While browsing the Project Euler problems, I found problem 97, which is an incredibly straightforward problem that is buried within a sea of much more challenging problems. It’s called “Large non-Mersenne Prime” and it states the following.

The first known prime found to exceed one million digits was discovered in 1999, and is a Mersenne prime of the form $2^{6972593}āˆ’1$; it contains exactly 2,098,960 digits. Subsequently other Mersenne primes, of the form $2^pāˆ’1$, have been found which contain more digits.

However, in 2004 there was found a massive non-Mersenne prime which contains 2,357,207 digits: $28433\times 2^{7830457}+1$.

Find the last ten digits of this prime number.

Solution in Python

Since we’re only looking for the final ten digits of the number, we can form a solution modulo a sufficiently large number, effectively forgetting about any larger power digits. We’ll do all of the powers of two first, then multiply by the constant, then add one. We’ll take the first ten digits of the result and that should be the answer we’re looking for. I’d say that this one is pretty quick and painless.

#!/usr/bin/env python
 
import time
 
start = time.time()
 
n = 2
for i in range(7830456):
    n = (2 * n) % 10000000000
 
n *= 28433
n += 1
 
n = n % 10000000000
 
elapsed = time.time() - start
 
print "Result %s found in %s seconds" % (n, elapsed)

Then, executing, we have:

$ python problem-97.py 
Result 8739992577 found in 1.08037400246 seconds

Solution in C

This is a situation where the C solution is nearly as fast to write as the Python solution, and the code really follows the same structure.

#include <stdio.h>
#include <stdlib.h>
 
int main(int argc, char **argv) {
    unsigned long       i;
    unsigned long long  n = 2;
 
    for (i=0;i<7830456;i++) {
        n = (n << 1) % 10000000000;
    }
 
    n *= 28433;
    n += 1;
 
    n = n % 10000000000;
 
    printf("Result %lld\n", n);
 
    return 0;
}

When we execute, the result is returned very quickly. (I’ve optimized the machine code using the -O3 flag as you can see below.)

$ gcc -O3 problem-97.c -o problem-97
$ time ./problem-97
Result 8739992577
 
real	0m0.030s
user	0m0.028s
sys	0m0.004s

That feels pretty fast, but I think a bottleneck is popping up when we bitshift (multiply by 2) AND reduce the result modolo 10000000000 at every single iteration of the for loop. The bitshifting is fast, but the modular arithmetic isn’t so much. Let’s rewrite the code so that we’re bitshifting at every iteration, but only using modular arithmetic every several iterations.

#include <stdio.h>
#include <stdlib.h>
 
 
int main(int argc, char **argv) {
    unsigned long       i,j;
    unsigned long long  n = 2;
 
    j = 0;
    for (i=0;i<7830456;i++) {
        //n = (n << 1) % 10000000000;
        n <<= 1;
        j++;
        if (j%10==0) {
            j = 0;
            n = n % 10000000000;
        }
    }
 
    n *= 28433;
    n += 1;
 
    n = n % 10000000000;
 
    printf("Result %lld\n", n);
 
    return 0;
}

That runs in roughly 1/3 the time of the original C code.

$ gcc -O3 problem-97.c -o problem-97
$ time ./problem-97
Result 8739992577
 
real	0m0.009s
user	0m0.008s
sys	0m0.000s