求C(n,m)%p(杨辉三角，逆元，Lucas定理)

注意： 本文主要探讨以下4种情况

1. 杨辉三角预处理组合数
2. m<p且p是素数(m≤1e5~1e6,n≤1e9,p≤1e9)：
3. 只需p是素数(m≤1e5~1e6,n≤1e9,p≤1e9)：
4. 只需p是素数(m≤n≤1e18,p≤1e5)：Lucas定理

①杨辉三角$O(n^2)$预处理组合数

例题：

代码：

#include <iostream>
#include <algorithm>
using namespace std;
const int N = 2010, mod = 1e9 + 7;
int c[N][N];
//循环i从0开始,j只循环一次:c[0][0] = 1,然后就进入i=1,不会RE
void init(){
    //从0开始方便
    for (int i = 0; i < N; i ++ )
        for (int j = 0; j <= i; j ++ )
            if (!j) c[i][j] = 1;//规定
            //组合数公式
            else c[i][j] = (c[i - 1][j] + c[i - 1][j - 1]) % mod;
}
int main(){
    int n;
    init();
    scanf("%d", &n);
    while (n -- ){
        int a, b;
        scanf("%d%d", &a, &b);
        printf("%d\n", c[a][b]);
    }
    return 0;
}

②m < p且p是素数(m≤1e5~1e6,n≤1e9,p≤1e9)

组合数求解公式可写为：
$$C_n^m=\frac{\left(n-m+1\right)\times \left(n-m+2\right)\times \dots \times \left(n-m+m\right)}{1\times 2\times \dots \times m}$$
经变形得到
$$\frac{\frac{\frac{\frac{\left(n-m+1\right)}{1}\times \left(n-m+2\right)}{2}\times \dots }{\dots }\times \left(n-m+m\right)}{m}$$

由此得出传统$O(m)$求法：

#define ll long long
ll C(ll n,ll m){
    ll ans=1;
    for(ll i=1;i<=m;i++){
        ans=ans*(n-m+i)/i;//注意先乘再除 
    }
}
ps1.易证每次除法都是整除
ps2.此方法(无取余)在n=62,m=31时开始溢出 

逆元定义：

假设 $a,b,m$ 是整数，$m>1$，且有 $a\ast b\equiv 1(modm)$ 成立，那么就说 $a,b$ 互为模 $m$ 的逆元，一般也记作 $a\equiv 1/b(modm)$ 或 $b\equiv 1/a(modm)$
通俗地说，如果两个整数的乘积模 $m$ 后等于 $1$ ，就称 $a,b$ 互为 $m$ 的逆元

逆元含义：

模$n$意义下，一个数$a$如果有逆元$x$，那么除以$a$相当于乘以$x$

逆元用处：

对乘法来说有$$（b\times a）\%m=（（b\%m）（a\%m）\%m$$
成立，但是对除法来说，$$（b/a）\%m=（（b\%m）/（a\%m））\%m$$
却不成立，$$（b/a）\%m=（（b\%m）/a）\%m$$也不成立，这时就需要逆元来计算$$（b/a）\%m$$
通过找到 $a$ 模 $m$ 的逆元，就有$$（b/a）\%m=（b\ast x）\%m$$成立 （证明过程略）

由于不能在除法时直接模上p，因此如果p是素数，可以使用扩展欧几里得算法或者费马小定理求出i模p的逆元，然后将除法取模转化为乘法取模来解决。需要注意的是，此时必须满足m<p，否则中间过程求逆元可能失效（即i是p倍数的情况）

显然这种做法的时间复杂度为$O(mlogm)$，其中$O(logm)$是计算逆元的复杂度。因此，这种做法能支持m≤1e5的情况（若时间方面较为宽松，则m≤1e6问题也不大），且对n和p的范围限制不大（例如n≤1e9,p≤1e9是可行的）,但是p必须是素数。

引入逆元后$O(mlogm)$的求法：

ps.m<p且p是素数(m≤1e5~1e6,n≤1e9,p≤1e9) 
C(int n,int m){
    int ans=1;
    for(int i=1;i<=m;i++){
        ans=1ll*ans*(n-m+i)%mod;//1ll*不写会爆掉 
        ans=1ll*ans*inverse(i)%mod;//求i模p的逆元
    }
    return ans;
}

补充1：费马小定理求逆元

费马小定理： 设$m$是素数，$a$是任意整数且$a\not\equiv 0\left(modm\right)$,则$a^{m-1}\equiv 1\left(modm\right)$

此时有$a\ast a^{m-2}\equiv 1\left(modm\right)$，因此$a$的逆元为$a^{m-2}(modm)$

求解逆元：

const int mod=1e9+7;
int qpow(int a,int k){
    int ans=1;
    while(k){
        if(k&1)
            ans=1ll*ans*a%mod;
        a=1ll*a*a%mod;
        k>>=1;
    }//k>>1：编译器检查不出来的错误
    return ans;
}
int inverse(int a,int mod){
    return qpow(a,mod-2);
}

整体代码：

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
#define f(i,a,b) for(int i=a;i<b;i++)
#define ff(i,a,b) for(int i=a;i<=b;i++)
const ll mod=1e9+7;
ll qpow(ll a,ll k){
    ll ans=1;
    while(k){
        if(k&1)
            ans=ans*a%mod;
        a=a*a%mod;
        k>>=1;
    }
    return ans;
}
ll inverse(ll a){
    return qpow(a,mod-2);
}
ll C(int n,int m){
    ll ans=1;
    for(int i=1;i<=m;i++){
        ans=ans*(n-m+i)%mod;//1LL*不写会爆掉
        ans=ans*inverse(i)%mod;//求i模p的逆元
    }
    return ans;
}
int main(){
    int n;
    cin>>n;
    while(n--){
        int a,b;
        cin>>a>>b;
        cout<<C(a,b)<<endl;
    }
    return 0;
}

补充2：讨论逆元失效条件

由定义知，求$a$模$m$的逆元，就是求解同余式$ax\equiv 1(modm)$，并且在实际使用中，一般把$x$的最小正整数解称为$a$模$m$的逆元，因此逆元经常被默认为$x$的最小正整数解。

显然，同余式$ax\equiv 1(modm)$是否有解取决于$1\%gcd(a,m)$是否为0，而这等价于$gcd(a,m)$是否为1：
①如果$gcd(a,m)\ne 1$，那么同余式$ax\equiv 1(modm)$无解，$a$不存在模$m$的逆元。
②如果$gcd(a,m)=1$，那么同余式$ax\equiv 1(modm)$在$(0,m)$上有唯一解，可以通过求解$ax+my=1$得到。

线性逆元($O\left(n\log n\right)\to O\left(n\right)$)

粗暴引入公式：
$$inv\left[x\right]\equiv \left(mod-\lfloor \frac{mod}{x}\rfloor \right)\ast inv\left[mod\%x\right]\%mod$$
（推导过程可看逆元详解）
这就意味着我们以一种DP的思路从小到大更新inv数组，可以以O(n)的复杂度更新出来

typedef long long ll;
inv[1]=1;
for(int i=2;i<=n;i++)
    inv[i]=((ll)(MOD-MOD/i)*inv[MOD%i])%MOD;

小技巧：

inv[2]=(mod+1)/2=mod-mod/2;

逆元预处理(组合数查询为$O\left(1\right)$级别)

typedef long long ll;
const int maxn = 2e6, mod = 911451407;
int re[maxn + 10], inv[maxn + 10], fac[maxn + 10];
inline void init(int n) {
    re[0] = inv[1] = fac[0] = re[1] = fac[1] = 1;
    for (int i = 2; i <= n; ++i) {
        fac[i] = (ll)fac[i - 1] * i % mod;
        inv[i] = (ll)(mod - mod / i) * inv[mod % i] % mod;//线性逆元 
        re[i] = (ll)re[i - 1] * inv[i] % mod;//前缀积    
    }
}
inline int C(int a,int b) {
    if (a < 0) return 0;
    return (ll)fac[a] * re[b] % mod * re[a - b] % mod;
}

③只需p是素数(m≤1e5~1e6,n≤1e9,p≤1e9)

事实上，由于$$C_n^m=\frac{\left(n-m+1\right)\times \left(n-m+2\right)\times \dots \times \left(n-m+m\right)}{1\times 2\times \dots \times m}$$
其中分子一定能够被分母整除，因此分子中含质因子p的个数必须不少于分母中含质因子p的个数。于是可以在for循环的过程中单独处理质因子p，即统计分子比分母中多含质因子p多少个，并在此过程中去除分子和分母中的质因子p，而对其余部分正常计算逆元；统计完毕后，如果分子中的质因子p的个数大于分母中的质因子p的个数，那么直接返回0；否则返回计算的结果。

int qpow(int a,int k){
    int ans=1;
    while(k){
        if(k&1)
            ans=1ll*ans*a%mod;
        a=1ll*a*a%mod;
        k>>=1;
    }//k>>1：编译器检查不出来的错误
    return ans;
}
int inverse(int a,int mod){
    return qpow(a,mod-2);
}


//m任意,且p是素数(m≤1e5~1e6,n≤1e9,p≤1e9)
int C(int n,int m,int p){
    int ans=1,numP=0;//ans存放结果，numP统计分子中的p比分母中的p多几个 
    for(int i=1;i<=m;i++){
        int temp=n-m+i;//分子
        while(temp %p==0){//去除分子中的所有p，同时累计numP 
            numP++;
            temp/=p;
        }
        ans=ans*temp%p;//乘以分子中除了p以外的部分
        temp=i;//分母
        while(temp%p==0){//去除分母中的所有p，同时减少numP 
            numP--;
            temp/=p;
        }
        ans=ans*inverse(temp,p)%p;//除以分母中除了p以外的部分 
    }
    if(numP>0)return 0;//分子中p的个数多于分母，直接返回0
    else return ans;//分子中p的个数等于分母，返回计算的结果
}

由于引入了numP的计算过程，这种做法的时间复杂度为$O(mlogn)$，但还是能支持m≤1e5的情况（同样如果设备符合要求，则m≤1e6问题也不大），且对n和p的范围限制不大（例如n≤1e9,p≤1e9是可行的）,但是p必须是素数。

④只需p是素数(m≤n≤1e18,p≤1e5)：Lucas定理

图片转载自：（Lucas定理）
代码：

int lucas(int n,int m){
    if(m==0)return 1;
    return C(n%p,m%p)*lucas(n/p,m/p)%p;
}

注意： 之后求C（…）用①的方法即可

例题1：莫的难题

莫的难题

简介：

埋和莫曾经是好朋友。埋是文科学霸，而莫却只是一个 OI 蒟蒻。一天，埋碰到一道难题跑来问莫。题目是这样的：有五个数字，分别是 5、2、1、3、9.莫可以取任意数字，每个数字可以取无限次。如：取两个 5，则组合为：55；取 2 与 1，则组合为：21。现在要问你所有组合中第 $C(n,m)\%1e9+7(n>=m)$ 个数有多大？

思路：

找规律发现思路为：将求得的组合数转化成五进制数，然后再把字符串转化回1，2，3，5，9和0，1，2，3，4对应的字符串。

解法一

//注意如果用%lld来读入会爆炸 
#include <bits/stdc++.h>
#define sc(x) scanf("%d", &(x))
using namespace std;
typedef long long ll;
const int mod=1e9+7;
int qkpow(int a, int b) {
    int ans = 1;
    while (b) {
        if (b & 1) ans = 1LL*ans * a % mod;
        a = 1LL*a * a % mod;
        b >>= 1;
    }
    return ans;
}
int getInv(int a) { return qkpow(a, mod - 2); }  //求一个数的逆元
int C(int n,int m){
    int ans=1;
    for(int i=1;i<=m;i++){
        ans=1LL*ans*(n-m+i)%mod;
        ans=1LL*ans*getInv(i)%mod;//求i模p的逆元
    }
    return ans;
}
char mp[5]={'1','2','3','5','9'};
int main() {
    int t;sc(t);
    while (t--) {
        int n, m; 
        sc(n), sc(m);
        int a = C(n, m);
        string ans="";
        while(a>0){
            --a;//精髓
            ans+=mp[a%5];
            a/=5;
        }
        reverse(ans.begin(),ans.end());
        cout<<ans<<endl;
    }
    return 0;
}

解法二(杨辉三角优化)

（参考题解：【题解】牛客IOI周赛17-普及组）

因为有多次询问，所以如果在数据极限的情况下可能还要用逆元预处理。这种方法太麻烦了，这里建议使用杨辉三角。首先，杨辉三角若查询$C(n,m)$，只需要$O(1)$的复杂度( 预处理复杂度为$O(n^2)$)。其次，杨辉三角递推式中仅含有加法，所以可以直接取模，省事很多（注意题目：$1<=m<=n<=100$）

#include<bits/stdc++.h>
using namespace std;
const int maxn=1e7+10;
const int mod = 1e9+7;
long long tgl[110][110], t, k;
int a[5] = {1, 2, 3, 5, 9};
int pt[maxn];
int main(){
    tgl[1][1] = 1;
    for (int i=2; i<=110; i++){
        for (int j=1; j<=i; j++) tgl[i][j] = (tgl[i-1][j] + tgl[i-1][j-1])%mod;
    }
    cin>>t;
    while (t --){
        int n, m, in=0; cin>>n>>m;
        k = tgl[n+1][m+1];
        while (k){
            k -= 1;
            pt[in++] = a[(k)%5];
            k /= 5;
        }
        for (int i=in-1; i>=0; i--) printf("%d", pt[i]);
        printf("\n");
    }
}

例题2：数列统计(逆元预处理)

数列统计

方法：

找规律发现原问题公式为$C_{L+x-2}^{x-1}$ (推导过程可看【题解】牛客IOI周赛17-普及组)

#include <cstdio>
#include <cstdlib>
#include <iostream>
using namespace std;
typedef long long ll;
const int maxn = 2e6, mod = 911451407;
int re[maxn + 10], inv[maxn + 10], fac[maxn + 10];
inline void init(int n) {
    re[0] = inv[1] = fac[0] = re[1] = fac[1] = 1;
    for (int i = 2; i <= n; ++i) {
        fac[i] = (ll)fac[i - 1] * i % mod;
        inv[i] = (ll)(mod - mod / i) * inv[mod % i] % mod;//线性逆元 
        re[i] = (ll)re[i - 1] * inv[i] % mod;//类比前缀和    
    }
}
inline int C(int a,int b) {
    if (a < 0) return 0;
    return (ll)fac[a] * re[b] % mod * re[a - b] % mod;
}
int main() {
    ios::sync_with_stdio(false), cin.tie(NULL), cout.tie(NULL);
    init(2e6);
    int T;
    cin >> T;
    int l, x;
    do {
        cin >> l >> x;
        --l;
        cout << C(l + x - 1, x - 1) << '\n';
    } while (--T);
    return 0;
}