LightOJ 1067 Combinations

Posted 2021-03-10 chinesepikaync

求 (C_n^k\%p)，(C) 为组合数，(p=1000003)，是个质数

---

前置芝士：

费马小定理

(a^{p-1} equiv a mod p) , (p) 为质数（费马小定理）

二项式定理：
[ (a+b)^n=sum_{i=0}^n{C_n^i}*a^i*b^{n-i} ]
当然，因为组合数是对称的，所以也可以写成
[ (a+b)^n=sum_{i=0}^n{C_n^i}*a^{n-i}*b^i ]
柿子1

((a+b)^p equiv a^p+b^p mod p)，(p) 为质数

证明：

由把 (a+b) 看成一个整体，由费马小定理得 ((a+b)^p equiv (a+b) mod p)

由费马小定理得 (a equiv a^p mod p) , (b equiv b^p mod p)

所以 ((a+b)^p equiv (a+b) equiv (a^p+b^p) mod p)

柿子2

((1+x)^p equiv (1+x^p) mod p)

直接代入柿子1即可

卢卡斯定理
[ C_n^m\%p=C_{lfloor frac np floor}^{lfloor frac mp floor}*C_{lfloor n\%p floor}^{lfloor m\%p floor}\%p ]
证明：

设 (b=n\%p)，

有
[ (1+x)^n=(1+x)^{lfloor frac np floor*p}*(1+x)^b ]
由柿子2得
[ equiv (1+x^p)^{lfloor frac np floor}*(1+x)^b mod p ]
由二项式定理得
[ equiv (sum_{i=0}^{{lfloor frac np floor}}{C_{lfloor frac np floor}^i}*x^{pi})*(sum_{j=0}^bC_b^j*x^j) mod p ]
我们把这个柿子叫做 (tanao) 柿 ，把它暴力展开一下
[ tanao=(C_{lfloor frac np floor}^0*x^{0p}+C_{lfloor frac np floor}^1*x^{1p}+...+C_{lfloor frac np floor}^{lfloor frac np floor}*x^{{lfloor frac np floor}*p})*(C_b^0*x^0+C_b^1*x^1+...+C_b^b*x^b) ]
我们把 (tanao) 柿 里面所有的 (x^m) 项都拿出来加在一起记为 (Q) ，那么这个 (Q=A*x^m)

其中，(A=sum C_{lfloor frac np floor}^i*C_b^j) , ((pi+j=m))，因为 (b=n\%p,jle b) 所以 (j<p)

我们再把原来那个 ((1+x)^n) 用二项式定理展开
[ (1+x)^n=sum_{i=0}^n{C_n^i}*x^i ]
我们发现 (x^i) 项对应的系数恰好是 (C_n^i)

也就是说，我们上面 (x^m) 项对应的系数恰好是 (C_n^m)

即，(C_n^m equiv A equiv sum C_{lfloor frac np floor}^i*C_b^j mod p) , ((pi+j=m,j<p))

所以我们现在有，(pi+j=m,j<p,b=n\%p)

可以得出
[ left{ egin{aligned} i=lfloor frac mp floor,j=m\%p, =n\%p, end{aligned} ight. ]
然后你定睛一看，(i,j,b) 都是唯一的。。。。

那就把可以把 (Sigma) 扔掉了，然后回代
[ C_n^m equiv C_{lfloor frac np floor}^{lfloor frac mp floor}*C_{lfloor n\%p floor}^{lfloor m\%p floor} mod p ]
即
[ C_n^m \%p = C_{lfloor frac np floor}^{lfloor frac mp floor}*C_{lfloor n\%p floor}^{lfloor m\%p floor} \%p ]
(wonderful!)

好，我们现在可以开始做题了

先把 (C_n^m) 里的 (n,m) 用卢卡斯递归降到 (p) 以下，然后再计算，递归回来

(n,m) 小于 (p) 的时候就好做了，先预处理出 (p) 以内的阶乘 (f[i]=1*2*...*n\%p) (O(n)) 预处理出来，然后用阶乘公式 (C_n^m=frac {n!}{m!(n-m)!}) 计算，注意这里是膜 (p) 意义下，所以算除法要用乘法逆元，直接费马小定理即可

时间复杂度 (O()很快())

// This code Write By chtholly_micromaker(MicroMaker)
#include <cstdio>
#include <cctype>
#define reg register
#define int long long
using namespace std;
const int p=1e6+3;
template <class t> inline void rd(t &s)
{
    s=0;
    reg char c=getchar();
    while(!isdigit(c))
        c=getchar();
    while(isdigit(c))
        s=(s<<3)+(s<<1)+(c^48),c=getchar();
    return;
}
int f[p+100];
inline void Init()
{
    f[0]=1;
    for(int i=1;i<=p;++i)
        f[i]=f[i-1]*i%p;
    return;
}
inline int fastpow(int a,int b)
{
    reg int res=1;
    a%=p;
    for(;b;b>>=1,a=a*a%p)
        if(b&1)
            res=res*a%p;
    return res;
}
inline int inv(int x)
{
    return fastpow(x,p-2);
}
inline int C(int n,int m)
{
    if(m>n)
        return 0;
    return f[n]*(fastpow(f[m]*f[n-m]%p,p-2))%p;
}
inline int Lucas(int n,int m)
{
    if(!m)
        return 1;
    return C(n%p,m%p)*Lucas(n/p,m/p)%p;
}
inline void work()
{
    int n,m;
    rd(n);rd(m);
    printf("%lld
",Lucas(n,m));
}
signed main(void)
{
    Init();
    int t;rd(t);
    for(int i=1;i<=t;++i)
        printf("Case %lld: ",i),work();
    return 0;
}