Calculating nCr mod p efficiently when n is very large

Question

I need to calculate nCr mod p efficiently. Right now, I have written this piece of code, but it exceeds the time limit. Please suggest a more optimal solution.

For my case, p = 10^9 + 7 and 1 ≤ n ≤ 100000000

I have to also make sure that there is no overflow as nCr mod p is guaranteed to fit in 32 bit integer, however n! may exceed the limit.

def nCr(n,k):
    r = min(n-k,k)
    k = max(n-k,k)
    res = 1
    mod = 10**9 + 7

    for i in range(k+1,n+1):
        res = res * i
        if res > mod:
            res = res % mod

    res = res % mod
    for i in range(1,r+1):
        res = res/i
    return res

PS : Also I think my code may not be completely correct. However, it seems to work for small n correctly. If its wrong, please point it out !

Why are you worrying about overflow? Python's integer types don't have a fixed storage space; it will take as much storage space as it needs. — Ramchandra Apte, Nov 02 '13 at 06:33
@RamchandraApte: while that is true, it /is/ possible to cause an int overflow in python — inspectorG4dget, Nov 02 '13 at 06:34
@inspectorG4dget Well, the only way I know in Python 3 is to use certain C functions, such as those in the math module. Not sure about Python 2. — Ramchandra Apte, Nov 02 '13 at 06:36
Also, it exceeds the time limit. Can I optimize it further ? — OneMoreError, Nov 02 '13 at 06:36
@RamchandraApte: you’re right. `math.e(BIG_NUMBER)` will cause an overflow. Haven’t seen it happen “naturally”, though — inspectorG4dget, Nov 02 '13 at 06:43
@inspectorG4dget: this is a float ("double"). This overflows, but the integers don't. — Armin Rigo, Nov 02 '13 at 08:17

score 12 · Answer 1 · edited Jan 14 '15 at 16:09

From http://apps.topcoder.com/wiki/display/tc/SRM+467 :

long modPow(long a, long x, long p) {
    //calculates a^x mod p in logarithmic time.
    long res = 1;
    while(x > 0) {
        if( x % 2 != 0) {
            res = (res * a) % p;
        }
        a = (a * a) % p;
        x /= 2;
    }
    return res;
}

long modInverse(long a, long p) {
    //calculates the modular multiplicative of a mod m.
    //(assuming p is prime).
    return modPow(a, p-2, p);
}
long modBinomial(long n, long k, long p) {
// calculates C(n,k) mod p (assuming p is prime).

    long numerator = 1; // n * (n-1) * ... * (n-k+1)
    for (int i=0; i<k; i++) {
        numerator = (numerator * (n-i) ) % p;
    }

    long denominator = 1; // k!
    for (int i=1; i<=k; i++) {
        denominator = (denominator * i) % p;
    }

    // numerator / denominator mod p.
    return ( numerator* modInverse(denominator,p) ) % p;
}

Notice that we use modpow(a, p-2, p) to compute the mod inverse. This is in accordance to Fermat's Little Theorem which states that (a^(p-1) is congruent to 1 modulo p) where p is prime. It thus implies that (a^(p-2) is congruent to a^(-1) modulo p).

C++ to Python conversion should be easy :)

Also, note that modPow is already available in the form of `pow()` in Python. — Ramchandra Apte, Nov 02 '13 at 07:18
Helped me out a ton. Was hard to find this implementation elsewhere. — ryan1234, Sep 07 '14 at 04:10

score 2 · Answer 2 · answered Nov 02 '13 at 08:20

About the last question: I think that the mistake in your code is to compute the product, reduce it modulo k, and then divide the result by r!. This is not the same as dividing before reducing modulo k. For example, 3*4 / 2 (mod 10) != 3*4 (mod 10) / 2.

Calculating nCr mod p efficiently when n is very large

2 Answers2