How to find Choose(m, n) % P for large inputs?

Question

P = (10^9 + 7)
Choose(m, n) = m! / (n! * (m - n)!)

I want to calculate the value of Choose(m, n) mod P for large m and n. How can I do that in C++ ?

Answer 1

This is what I use as it has a fairly good range without too much intermediate overflow. However C(n,k) gets big fast, it is O(n^k) after all.

size_t N_choose_K(size_t n, size_t k)
{
    size_t numer = 1;
    size_t denom = 1;
    if (k > n - k) {
        k = n - k;
    }
    while (k > 0) {
        numer *= n;
        denom *= k;
        --n; --k;
    }
    return numer / denom;
}

EDIT: Assumes you need integral results. You can move to floating point results and get more range if you need it and can stand to lose the accuracy.

Answer 2

You can use the fact that multiplication is closed under Z_p meaning that: a b mod p = (a mod p) (b mod p) mod p. Using this theorem, one can calculate a^b mod p effectively (for instance, this Python implementation.

Next we can make use of Euler's theorem saying: a^-1 mod p=a^(p-2) mod p.

Now that we know these facts, we can come up with an effective solution: first we multiply over all elements in the numerator, this is thus a range from k+1 (inclusive) to n, and since this is a multiplication, we can always perform a modulo:

long long numerator(int n, int k, int p) {
    long long l = 1;
    for(int j = k+1; j <= n; j++) {
        l = (l*j)%p;
    }
    return l;
}

Now we still need to divide it by (n-k)!. We can do this by first calculating (n-k)! mod p like we already did in the previous code fragment:

long long denominator(int n, int k, int p) {
    long l = 1;
    for(int j = 2; j <= n-k; j++) {
        l = (l*j)%p;
    }
    return l;
}

Now in order to divide it, we can use Euler's theorem on the result of denominator. Therefore we first implement the pow function with modulo:

long long pow(long long a, int k, int p) {
    if(k == 0) {
        return 1;
    }
    long long r = pow((a*a)%p,k>>0x01,p);
    if((k&0x01) == 0x01) {//odd number
        r = (r*a)%p;
    }
    return r;
}

Now we can merge these together like:

long long N_choose_K(int n, int k, int p) {
    long long num = numerator(n,k,p);
    long long den = denominator(n,k,p);
    return (num*pow(den,p-2,p))%p;
}

So what you basically do is you determine the numerator num in Z_p, the value of the denominator den in Z_p, and then you use Euler's theorem to find the inverse of the denominator in Z_p, such that you can multiply and perform the last modulo operation. Then you can return it.