0

I'm calculating nCr % MOD for large values of n

I'm using the relation (n+1)Cr=(nCr)*(n+1)/(n+1-r)

I have to iterate over a loop for different values of n keeping r constant.

llu fact=1;
/*The loop begins from i>=M+1 */
fact=(fact*(i-1)*modInverse(i-M,MOD))%MOD;  // Single statement executed during each iteration of loop


  Here I'm calculating (i-1)C(M-1)
  Here M and MOD are constant values 
  MOD=1000000009 and llu refers to unsigned long long

What I'm doing exactly is

(n+1)Cr % MOD = [ (nCr) * (n+1) * modInverse(n+1-r) ] % MOD

Here modInverse calulates modular multiplicative inverse as according to following definition:

llu modPow(llu a, llu x, llu p)
{
    //calculates a^x mod p
    llu res = 1;
    while(x > 0)
    {
        if( x % 2 != 0)
        {
            res = (res * a) % p;
        }
        a = (a * a) % p;
        x /= 2;
    }
    return (res%MOD);
}
llu modInverse(llu a, llu p)
{
    //calculates the modular multiplicative of a mod m assuming p is prime
    return modPow(a, p-2, p);
}

Now the problem is that I'm not getting correct values of nCr for large values of n (order 10^6). Is my approach of

(n+1)Cr % MOD = [ (nCr) * (n+1) * modInverse(n+1-r) ] % MOD

conceptually wrong?

Keith Thompson
  • 254,901
  • 44
  • 429
  • 631
CPPCoder
  • 155
  • 1
  • 1
  • 10
  • I'm working on the same idea as (a*b)%c= ((a%c)*(b%c))%c – CPPCoder Mar 13 '14 at 20:42
  • Have you tried adding the suffix ULL or LLU when you divide by 2? x /= 2LLU; Usually while coding and not specifying the compiler which type the literals are brings out this kind of truncation problems. – Salvador Medina Jan 20 '15 at 01:21

1 Answers1

0

Yes, the formula at the bottom is mathematically sound.

However, by doing 2 multiplications before taking the modulus, this increases the chance of an overflow problem.

For example, MOD is O(10^10) so modInverse is also O(10^10). If n is O(10^6), then the product is O(10^26) which is O(2^86) and would cause an overflow for a uint64, and giving you the wrong answer. Consider instead:

(n+1)Cr % MOD = [ ([(nCr) * (n+1)] % MOD) * modInverse(n+1-r) ] % MOD

You may also be interested in my answer to this suspiciously similar question.

Community
  • 1
  • 1
Jim Wood
  • 891
  • 8
  • 20