20% Project - Project Euler
Euler discovered the remarkable quadratic formula:
n^2 + n + 41
It turns out that the formula will produce 40 primes for the consecutive values n = 0 to 39. However, when n = 40, 40^2 + 40 + 41 = 40(40 + 1) + 41 is divisible by 41, and certainly when n = 41, 41^2 + 41 + 41 is clearly divisible by 41.
The incredible formula n^2 − 79n + 1601 was discovered, which produces 80 primes for the consecutive values n = 0 to 79. The product of the coefficients, −79 and 1,601, is −126,479.
Considering quadratics of the form:
n^2 + an + b, where |a| < 1,000 and |b| < 1,000 where |n| is the modulus/absolute value of n e.g. |11| = 11 and |−4| = 4
Find the product of the coefficients, a and b, for the quadratic expression that produces the maximum number of primes for consecutive values of n, starting with n = 0.
Initially I planned to tackle this problem from a strictly brute force perspective. I have found that working on a simple solution to a problem usually gives me insights to more elegant solutions.
Therefore I figured that the most straight forward and inelegant method was to produce a 1999 by 1999 array with individual elements of the array corresponding to the maximum number of consecutive primes generated by that combination of a and b.
| Thus the initial element of the array in position | 1,1 | would be the maximum number of consecutive primes generated by n2 + (-999)n + (-999), and the final position of | 1999,1999 | would be n2 + (999)n + (999). |
With this plan in mind I tried to think up a quick way to calculate the maximum number of consecutive primes for a given a and b.
I started by having Mathematica generate the first 101 values from the example given in the problem by using a=1 and b=41.
I then used PrimeQ to test whether each value was prime.
Then used Boole to convert these true false values to 0’s and 1’s.
And finally FirstPosition to find the first false in the list.
Note that the answer of 41 corresponds to the 39 primes listed in the example since this method doesn’t count the primes from 0-39 but from 1-40 and it found a non-Prime number in position 41.
So…
This method is horribly inelegant, however it does work so it’s time to see if I can use this method to generate the 1,999 by 1999 array I initially planned for.
I begin by trying my method on a smaller array and seeing if it would work at all.
Note that this entry should be much larger as it displays a 51×51 array of values.
To find the maximum element of this array I then use the Position and Max commands. Position gives the location in an array of a specified value, and Max gives the largest element of the array.
Since we are counting up from 0 we have now found that the Maximum number of consecutive primes is generated with a=1 and b=41 which matches the example.
One area of concern going forward is that it takes 0.655,204 secs for Mathematica to create the 51×51 array I used above. Since a 1,999×1,999 array is almost 1,600% larger than a 50×50 array, that means that generating a 1999×1999 array should take at least 17 minutes and 28.33 second. Most likely longer considering that the edges of the array will be dealing with a and b values greater than 900.
Still I have free time and I’m not paying for processing time so we might as well go forward and see how long it takes.
So 21 minutes and 14 seconds later we have our array and can try to find the Maximum number of consecutive primes.
And we have a max at |939,1971|, so 939-1000=-61 and 1971-1000=971 are are a and b values.
a*b=-61*971=-59231
And we thus have our answer for Problem 27.
Now we can focus on improving this process.
We could start by recognizing that b has to be a prime number. It was only after I had finished this problem that I realized that since n starts at 0 the first element will always be 0*0 + a*0 + b. Since there are 168 primes less than 1000 that means I only needed 168 rows instead of 1999.
A 1999×168} array is 8.4%} the size of a 1999×1999} array and only took Mathematica 1} minute and `56} seconds to calculate.
A larger time sink is that I am actually calculating out the first 100 n values for every combination of a and b tested. Thus even if the second value calculated is non-Prime I still waste time evaluating the next 98 values.
Using the following code I can just check if the n values are prime and stop when a non-Prime is found.
To better understand how much time this saves, I calculated a 50×50 array using the old method and the new method and had Mathematica generate an ArrayPlot to save space.
The second method took 6.1% of the time taken by the first method.
Using both methods together the entire process of generating the 1999×168 array takes 10 seconds.
I haven’t looked into it, but I suspect that I could similarly shave a couple of seconds off of the total running time by finding some trick for picking a values.
Conclusion
I started by creating a rather inelegant method for finding the number of consecutive primes generated on a 1999×1999 matrix. This method took 21 minutes and 14 seconds.
I next examined the problem more closely and found areas that could be refined so that I found that same solution in 10 seconds.
I’m confident that the majority of the improvement in performance came from optimizing the operation that found the number of consecutive primes.
Using the above code I found that 280,607 times the second number tested was non-Prime. Since there are only 336,000 entries in the array that means that 83.5% of the time while generating the original array I could have stopped after calculating the second value.
This was a longer write up than I originally expected, but it was a rather involved problem that benefited greatly from optimization.