Lucas定理

求C((m, n) % p) (p为素数)

表达式：C(n,m)%p=C(n/p,m/p)*C(n%p,m%p)%p

当p<1e6时

#include<algorithm>
#include<cstdio>
#include<cmath>
using namespace std;
typedef long long ll;
const int N = 100000;
ll f[N];
void init(int p) {
    f[0] = 1;
    for (int i=1; i<=p; ++i) f[i] = f[i-1] * i % p;
}

ll pow_mod(ll a, ll x, int p)   {
    ll ret = 1;
    while (x)   {
        if (x & 1)  ret = ret * a % p;
        a = a * a % p;
        x >>= 1;
    }
    return ret;
}

ll Lucas(ll n, ll k, int p) {       //C (n, k) % p
     ll ret = 1;
     while (n && k) {
        ll nn = n % p, kk = k % p;
        if (nn < kk) return 0;                   //inv (f[kk]) = f[kk] ^ (p - 2) % p
        ret = ret * f[nn] * pow_mod (f[kk] * f[nn-kk] % p, p - 2, p) % p;
        n /= p, k /= p;
     }
     return ret;
}

int main(void)  {
    int p = 10007, n , m;
    while(~scanf("%d %d", &n, &m)) {
        init (p);
        printf ("%lld\n", Lucas (n, m, p));
    }
    return 0;
}

扩展：当p较大（1e12左右）或者p不是素数时

#include<algorithm>
#include<cstdio>
#include<cmath>
using namespace std;
typedef long long ll;
ll POW(ll a, ll b, ll mod)
{
    ll ans=1;
    while(b)
    {
        if(b & 1) ans = ans * a % mod;
        a = a *a % mod;
        b >>= 1;
    }
    return ans;
}

ll POW(ll a, ll b)
{
    ll ans = 1;
    while(b)
    {
        if(b & 1) ans = ans * a;
        a = a * a;
        b >>= 1;
    }
    return ans;
}

ll exGcd(ll a, ll b, ll &x, ll &y)
{
    ll t, d;
    if(!b)
    {
        x=1;
        y=0;
        return a;
    }
    d = exGcd(b, a % b, x, y);
    t = x;
    x = y;
    y = t - a / b * y;
    return d;
}

bool modular(ll *a, ll *m, ll k)
{
    ll d, t, c, x, y, i;
    for(i = 2; i <= k; i++)
    {
        d = exGcd(m[1], m[i], x, y);
        c = a[i]-a[1];
        if(c % d) return false;
        t = m[i] / d;
        x = (c / d * x % t + t) % t;
        a [1] =m[1] * x + a[1];
        m [1] =m[1] * m[i] / d;
    }
    return true;
}

ll reverse(ll a, ll b)
{
    ll x, y;
    exGcd(a, b, x, y);
    return (x % b + b) % b;
}

ll C(ll n, ll m, ll mod)
{
    if(m > n) return 0;
    ll ans = 1, i, a, b;
    for(i = 1; i <= m; i++)
    {
        a = (n + 1 - i) % mod;
        b = reverse(i % mod, mod);
        ans = ans * a % mod * b % mod;
    }
    return ans;
}

ll C1(ll n, ll m, ll mod)
{
    if(m == 0) return 1;
    return C(n % mod, m % mod, mod) * C1(n / mod, m / mod, mod) % mod;
}

ll cal(ll n, ll p, ll t)
{
    if(!n) return 1;
    ll x = POW(p, t), i, y = n / x, temp = 1;
    for(i = 1; i <= x; i++) if(i % p) temp = temp * i%x;
    ll ans=POW(temp,y,x);
    for(i=y*x+1;i<=n;i++) if(i%p) ans=ans*i%x;
    return ans*cal(n/p,p,t)%x;
}

ll C2(ll n, ll m, ll p, ll t)
{
    ll x=POW(p, t);
    ll a, b, c, ap = 0, bp = 0, cp = 0, temp;
    for(temp = n; temp ;temp /= p) ap += temp / p;
    for(temp = m; temp; temp /= p) bp += temp / p;
    for(temp = n - m; temp; temp /= p) cp += temp / p;
    ap = ap - bp - cp;
    ll ans = POW(p, ap, x);
    a = cal(n, p, t);
    b = cal(m, p, t);
    c = cal(n - m, p, t);
    ans = ans * a % x * reverse(b, x) % x * reverse(c, x) % x;
    return ans;
}

//计算C(n,m)%mod
ll Lucas(ll n,ll m,ll mod)
{
    ll i, t, cnt = 0;
    ll A[205], M[205];
    for(i = 2; i * i <= mod; i++)
    if(mod % i == 0)
    {
        t = 0;
        while(mod % i == 0)
        {
            t++;
            mod /= i;
        }
        M[++cnt] = POW(i,t);
        if(t==1) A[cnt] = C1(n, m, i);
        else A[cnt] = C2(n, m, i, t);
    }
    if(mod > 1)
    {
        M[++cnt] = mod;
        A[cnt] = C1(n, m, mod);
    }
    modular(A, M, cnt);
    return A[1];
}
int main(void)  {
    ll mod = 10000006, n , m;
    while(~scanf("%lld %lld", &n, &m)) {

        printf ("%lld\n", Lucas (n, m, mod));
    }
    return 0;
}