组合数计算（Lucas 逆元 中国剩余定理）

Problem 1. treasure
Input file: treasure.in
Output file: treasure.out
Time limit: 1 second
Memory limit: 256 MB
最近Mr. H 得到了一个宝箱, 但宝箱被上了锁, 需要解决一个问题才能将其打开, 问题是这样的:
你有n 个不同的苹果, 你想从里面选出m 个来, 问方案数, 结果可能很大, 输出模M 后的结果输出, 其中
M 是k 个不同的素数p1; p2;…; pk 的乘积.
Mr. H 一眼就解决了这个问题, 但他知道你们最近才学了这方面的知识, 于是就把这个问题交给了你们.
Input
第1 行，3 个整数n m k, 意义如上。
第2 行k 个不同的素数:p1 p2 … pk.
Output
输出一个整数表示答案.
Sample
treasure.in
9 5 2
3 5
treasure.out
6
Note
. 对于10% 的数据，1 <= n <= 10^3,k = 1;
. 对于30% 的数据，1 <= n <= 10^5,k = 1,n < p1;
. 对于60% 的数据，1 <= n <= 10^18,k = 1;
. 对于100% 的数据，1 <= n <= 10^18, 1 <= k <= 10
. 对于所有数据，0 <= m <= n <= 10^18, 2 <= pi <= 10^5 且pi 互不相同.

思路：
求较大的组合数取模，需要逆元
再大一点，就要用Lucas了
这道题的mod不是质数，所以不能直接Lucas
所以对于每个质因数算出答案，得到一个形如一个ai=xi（mod mi）的式子
最后用中国剩余定理把所有式子合并起来就行了

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <iostream>
#define LL long long
using namespace std;

int k;
LL n, m, mud[100010];

LL mpow(LL a, LL b, LL mod) {
    LL ans;
    for( ans = 1; b; a=a*a%mod, b>>=1 ) 
        if(b & 1) ans = ans * a % mod;
    return ans;
}

LL Lucas(LL n, LL m, LL mod) {
    if(m > n) return 0;
    if(n > mod) return Lucas(n%mod, m%mod, mod) * Lucas(n/mod, m/mod, mod) % mod;
    else {
        LL ans = mud[n] * mpow(mud[m],mod-2,mod) % mod * mpow(mud[n-m],mod-2,mod);
        return ans % mod;
    }
}

void exgcd(LL a, LL b, LL &d, LL &x, LL &y) {
    if(b == 0) d = a, x = 1, y = 0;
    else {
        exgcd(b, a%b, d, y, x);
        y -= a / b * x;
    }
}

//ans = a0+t0*m0 = a1+t1*m1
//a0+t0*m0 = a1+t1*m1
//t0*m0 - t1*m1 = a1 - a0
//形同ax + by = c
//a = m0; b = m1; c = a1-a0; x = t0;
//解出ax + by = d(gcd(a,b))的x,y 
//ax + by = c中x=x*c/d; mod = b/d 
//ans = a0+x*m0(mod (m0+m1)/d)
//ans = (a0+x*m0) + t0*(m0+m1)/d
//ans = a0+t0*m0

void crt(LL &a0, LL &m0, LL a1, LL m1) {
    LL x, y, d;
    exgcd(m0, m1, d, x, y);
    LL p = m1 / d;
    x = ((a1 - a0) / d * x % p + p) % p;
    a0 += m0 * x;
    m0 = m0 * m1 / d;
}

int main() {
    freopen("treasure.in", "r", stdin);
    freopen("treasure.out", "w", stdout);
    scanf("%I64d%I64d%d", &n, &m, &k);
    LL a0 = 0, m0 = 1, a1, m1;
    for( int i = 1; i <= k; i++ ) {//把每个质因数拆开考虑，质因数mi可以得到一个ai=xi（mod mi） 
        scanf("%I64d", &m1);
        mud[0] = 1;
        for( int i = 1; i <= m1; i++ ) mud[i] = mud[i-1] * i % m1;//
        a1 = Lucas(n, m, m1);
        crt(a0, m0, a1, m1);//中国剩余定理把两个形同ai=xi（mod mi）的式子合并起来 
    }
    printf("%I64d\n", a0);
    return 0;
}