Digit sum of central binomial coefficients The task is simply to see how much faster you can calculate n choose n/2 (for even n) than the builtin function in python.  Of course for large n this is a rather large number so rather than output the whole number you should output the sum of the digits. For example, for n = 100000, the answer is 135702. For n=1000000 it is 1354815.

Here is the python code:

from scipy.misc import comb
def sum_digits(n):
   r = 0
   while n:
       r, n = r + n % 10, n / 10
   return r
sum_digits(comb(n,n/2,exact=True))


Your score is (highest n on your machine using your code)/(highest n on your machine using my code).  Your code must terminate in 60 seconds or less.

Your program must give the correct output for all even n: 2 &lt;= n &lt;= (your highest n)

You can't use any builtin code or libraries which calculate binomial coefficients or values which can quickly be transformed into binomial coefficients.

You can use any language of your choice.



Leading answer   The current leading answer with an amazing 680.09 is by justhalf.
C++ (GMP) - (287,000,000 / 422,000) = 680.09

Shamelessly combine Kummer's Theorem by xnor and GMP by qwr.
Still not even close to the Go solution, not sure why.

Edit: Thanks Keith Randall for the reminder that multiplication is faster if the number is similar in size. I implemented multi-level multiplication, similar to memory coalescing concept on memory management. And the result is impressive. What used to take 51s, now takes only 0.5s (i.e., 100-fold improvement!!)


OLD CODE (n=14,000,000)
Done sieving in 0.343s
Done calculating binom in 51.929s
Done summing in 0.901s
14000000: 18954729

real    0m53.194s
user    0m53.116s
sys 0m0.060s

NEW CODE (n=14,000,000)
Done sieving in 0.343s
Done calculating binom in 0.552s
Done summing in 0.902s
14000000: 18954729

real    0m1.804s
user    0m1.776s
sys 0m0.023s


The run for n=287,000,000


Done sieving in 4.211s
Done calculating binom in 17.934s
Done summing in 37.677s
287000000: 388788354

real    0m59.928s
user    0m58.759s
sys 0m1.116s


The code. Compile with -lgmp -lgmpxx -O3

#include &lt;gmpxx.h&gt;
#include &lt;iostream&gt;
#include &lt;time.h&gt;
#include &lt;cstdio&gt;

const int MAX=287000000;
const int PRIME_COUNT=15700000;

int primes[PRIME_COUNT], factors[PRIME_COUNT], count;
bool sieve[MAX];
int max_idx=0;

void run_sieve(){
    sieve[2] = true;
    primes[0] = 2;
    count = 1;
    for(int i=3; i&lt;MAX; i+=2){
        sieve[i] = true;
    }
    for(int i=3; i&lt;17000; i+=2){
        if(!sieve[i]) continue;
        for(int j = i*i; j&lt;MAX; j+=i){
            sieve[j] = false;
        }
    }
    for(int i=3; i&lt;MAX; i+=2){
        if(sieve[i]) primes[count++] = i;
    }
}

mpz_class sum_digits(mpz_class n){
    clock_t t = clock();
    char* str = mpz_get_str(NULL, 10, n.get_mpz_t());
    int result = 0;
    for(int i=0;str[i]&gt;0;i++){
        result+=str[i]-48;
    }
    printf("Done summing in %.3fs\n", ((float)(clock()-t))/CLOCKS_PER_SEC);
    return result;
}

mpz_class nc2_fast(const mpz_class &amp;x){
    clock_t t = clock();
    int prime;
    const unsigned int n = mpz_get_ui(x.get_mpz_t());
    const unsigned int n2 = n/2;
    unsigned int m;
    unsigned int digit;
    unsigned int carry=0;
    unsigned int carries=0;
    mpz_class result = 1;
    mpz_class prime_prods = 1;
    mpz_class tmp;
    mpz_class tmp_prods[32], tmp_prime_prods[32];
    for(int i=0; i&lt;32; i++){
        tmp_prods[i] = (mpz_class)NULL;
        tmp_prime_prods[i] = (mpz_class)NULL;
    }
    for(int i=0; i&lt; count; i++){
        prime = primes[i];
        carry=0;
        carries=0;
        if(prime &gt; n) break;
        if(prime &gt; n2){
            tmp = prime;
            for(int j=0; j&lt;32; j++){
                if(tmp_prime_prods[j] == NULL){
                    tmp_prime_prods[j] = tmp;
                    break;
                } else {
                    mpz_mul(tmp.get_mpz_t(), tmp.get_mpz_t(), tmp_prime_prods[j].get_mpz_t());
                    tmp_prime_prods[j] = (mpz_class)NULL;
                }
            }
            continue;
        }
        m=n2;
        while(m&gt;0){
            digit = m%prime;
            carry = (2*digit + carry &gt;= prime) ? 1 : 0;
            carries += carry;
            m/=prime;
        }
        if(carries&gt;0){
            tmp = 0;
            mpz_ui_pow_ui(tmp.get_mpz_t(), prime, carries);
            for(int j=0; j&lt;32; j++){
                if(tmp_prods[j] == NULL){
                    tmp_prods[j] = tmp;
                    break;
                } else {
                    mpz_mul(tmp.get_mpz_t(), tmp.get_mpz_t(), tmp_prods[j].get_mpz_t());
                    tmp_prods[j] = (mpz_class)NULL;
                }
            }
        }
    }
    result = 1;
    prime_prods = 1;
    for(int j=0; j&lt;32; j++){
        if(tmp_prods[j] != NULL){
            mpz_mul(result.get_mpz_t(), result.get_mpz_t(), tmp_prods[j].get_mpz_t());
        }
        if(tmp_prime_prods[j] != NULL){
            mpz_mul(prime_prods.get_mpz_t(), prime_prods.get_mpz_t(), tmp_prime_prods[j].get_mpz_t());
        }
    }
    mpz_mul(result.get_mpz_t(), result.get_mpz_t(), prime_prods.get_mpz_t());
    printf("Done calculating binom in %.3fs\n", ((float)(clock()-t))/CLOCKS_PER_SEC);
    return result;
}

int main(int argc, char* argv[]){
    const mpz_class n = atoi(argv[1]);
    clock_t t = clock();
    run_sieve();
    printf("Done sieving in %.3fs\n", ((float)(clock()-t))/CLOCKS_PER_SEC);
    std::cout &lt;&lt; n &lt;&lt; ": " &lt;&lt; sum_digits(nc2_fast(n)) &lt;&lt; std::endl;
    return 0;
}