【JN大周报】组合数的求法

说到组合数，哇塞，那可是。。。

想必组合数是什么都应该知道，下面给出组合数的几种求法。

---

首先，公式法：\(C_{n}^{m}=\frac{n!}{m!*(n-m)!}\),

考虑如果有模数，怎么办，直接求逆元，那就有点慢了

有一种线性求组合数逆元的算法：

设f(x)函数为逆元函数

\(f(n!)=f((n-1)!*n)\)

\(= f((n-1)!)*f(n)\)

则

\(f((n-1)!)=f(n!)*f(f(n))\)

因为\(f(f(n))=n\)

所以\(f(n)=f(n+1)*(n+1)\)

---

那么第二种，又可以理解成杨辉三角，也可以理解乘动态规划：

\(C_{n}^{m}=C_{n-1}^{m}+C_{n-1}^{m-1}\)

动态规划理解：对于最后一个元素是否取。

---

考虑如果n,m超级大，p不大，且是质数怎么做？

那么进入正题——卢卡斯定理。

简单来讲：

设 求\(C_{n}^{m}\),\(n=k*p+r\),\(m=q*p+s\)

那么\(C_{n}^{m}=C_{k}^{q}*C_{r}^{s}\)

证明：

（前置知识：二项式定理，超弱同余

首先先证明一个结论:\({(1+x)^p} \equiv {(1+x^p)} \pmod {p}\)

根据二项式定理：

\((1+x)^p=\sum_{i=0}^{p}C_{p}^{i}*x^i*1^{p-i}\)

那么很明显除了\(C_{p}^{p}\)和\(C_{p}^{0}\)外，其他都是可以被p整除的，那么%出来都是0

则原式化为：

\((1+x)^p \equiv C_{p}^{0}*x^0*i^p+C_{p}^{p}*x^p*1^0 \pmod {p}\)

​ $ \equiv 1+x^p \pmod {p}$

有了这个性质我们来考虑一个式子：

\((1+x)^n=((1+x)^{p})^k*(1+x)^r\)

​ \(=(1+x^p)^k*(1+x)^r\)

​ \(=(\sum_{i=0}^{k}C_{k}^{i}*x^{p*i})\)\(*\)$(\sum_{j=0}^{r}C_{j=0} *x^{j}) $

在展开左边：

\((1+x)^n=\sum_{i=0}^{n}C_{n}^{i}*x^{i}\)

考虑到我们要求的是\(C_{n}^{m }\)

带入i=m左边系数：\(C_{n}^{m}\)

要让右边x的次数为m

则要让\(p*i+j=m\)

很明显，\(i=q,j=s\)

那么在带入C

\(C_{n}^{m }\equiv C_{k}^{q}*C_{r}^{s} \pmod {p}\)

---

我同余等于傻傻分不清。大佬轻喷 柯w柯。

那么奉上卢卡斯代码：

#include<bits/stdc++.h>
#define int long long
typedef long long LL;
using namespace std;
int n,m,p,T;
const int MAXP=1e5;
int inv[MAXP+5],fac[MAXP+5];
inline int power(int a,int b,int p){
    int ret=1;
    for(;b;a=(a*a)%p,b>>=1)if(b&1) ret=(ret*a)%p;
    return ret;
}
inline int C(int n,int m,int p){
    if(m>n) return 0;
    if(m==n) return 1;
    return ((fac[n]*power(fac[m],p-2,p)%p)%p*power(fac[n-m],p-2,p)%p)%p;
}
inline int Lucas(int n,int m,int p){
    if(m==0) return 1;
    return Lucas(n/p,m/p,p)*C(n%p,m%p,p)%p;
}
signed main(){
    scanf("%lld",&T);
    while(T--){
        scanf("%lld%lld%lld",&n,&m,&p);
        fac[0]=fac[1]=1;
        for(int i=2;i<=p;++i) fac[i]=(fac[i-1]*i)%p;
        printf("%lld\n",Lucas(n+m,m,p));
    }
    return 0;
}