组合数学(推导+模板)

组合数

$$C_n^m=C_{n-1}^m+C_{n-1}^{m-1}$$

$$C_n^m=C_n^{n-m}$$

杨辉三角递推式

O(n2)

	0	1	2	3	4	5	6	7	8	9	10
0	1
1	1	1
2	1	2	1
3	1	3	3	1
4	1	4	6	4	1
5	1	5	10	10	5	1
6	1	6	15	20	15	6	1
7	1	7	21	35	35	21	7	1
8	1	8	28	56	70	56	28	8	1
9	1	9	36	84	126	126	84	36	9	1
10	1	10	45	120	210	252	210	120	45	10	1

//范围小数大
int c[1001][1001];
void solve()
{
    for (int i = 0; i <= 10;i++){
        for (int j = 0; j <= i;j++){
            if(j==0){
                c[i][j] = 1;
            }
            else{
                c[i][j] = (c[i - 1][j - 1] + c[i - 1][j])%mod;
            }
        }
    }
}

快速幂&逆元

O(nlogp)

int fac[N], pr[N];//fac阶乘数组,pr逆元数组
int ksm(int d,int c,int mod){
    int ans = 1;
    while(c){
        if(c&1){
            ans = ans * d % mod;
        }
        c >>= 1;
        d = d * d % mod;
    }
    return ans;
}
void init(){
    // n<=1e7;
    fac[0] = pr[0] = 1;
    for (int i = 1; i <= ;i++){
        fac[i] = fac[i - 1] * i % mod;
        pr[i] = pr[i - 1] * ksm(i, mod - 2, mod)%mod;
    }
}
int c(int n,int k){
    return n < k ? 0 : ((fac[n] * pr[k]) % mod) * pr[n - k] % mod;
}

递推法逆元

O(n)
$$p\equiv 0(\mathrm{mod}p)$$

$$设p=k\ast i+r(1<r<i<p)$$

$$k\ast i+r\equiv 0(\mathrm{mod}p)$$

$$两边同乘i^{-1},r^{-1}得$$

$$k\ast r^{-1}+i^{-1}\equiv 0(\mathrm{mod}p)$$

$$i^{-1}\equiv -k\ast r^{-1}(\mathrm{mod}p)$$

$$i^{-1}\equiv p\ast r^{-1}-k\ast r^{-1}(\mathrm{mod}p)$$

$$i^{-1}\equiv (p-k)\ast r^{-1}(\mathrm{mod}p)$$

$$inv[i]=(p-p/i)\ast inv[p\%i]\%p$$

int fac[N],invfac[N],pr[N];
void init(){
    fac[1] = invfac[1] = invfac[0] = pr[1] = fac[0] = 1;
    //<=1e7
    for (int i = 2; i <= ;i++){
        fac[i] = fac[i - 1] * i%mod;
        pr[i] = pr[mod%i] * (mod - mod / i) % mod;
        invfac[i] = invfac[i - 1] * pr[i] % mod;
    }
}
int c(int n,int k){
    return n < k ? 0 : ((fac[n] * invfac[k]) % mod) * invfac[n - k] % mod;
}
/*
pr数组是一个用于计算模逆元的数组。模逆元是指在模运算下的乘法逆元，即对于给定的模数mod，对于任意非零整数a，存在一个整数b，使得a * b ≡ 1 (mod mod)。在这段代码中，pr数组用于存储每个数的模逆元。
具体来说，pr[i] 存储了 (mod - mod/i) * pr[mod%i] % mod 的结果。其中，pr[1] 被初始化为 1，而其余位置 i 的值通过递推计算得到。这样，通过查询 pr 数组可以获得任意数的模逆元。
*/

卢卡斯lucas定理

O(logpN⋅(p+log2p))

此方法用于非常大的 n 和 m ，但是 p 却比较小的计算。
$$C_n^m=C_{\frac{p}{n}}^{\frac{p}{m}}\times C_{n\%p}^{m\%p}(\mathrm{mod}p)$$

引理1：
$$C_p^x\equiv 0(\mathrm{mod}p)(0<x<p)$$

$$因C_p^x=\frac{P!}{x!(p-x)!}=\frac{p(p-1)!}{x(x-1)!(p-x)!}=\frac{p}{x}{C}_p^i{x}^i$$

$$故C_p^x\equiv p\times inv(x)C_{p-1}^{x-1}\equiv 0(\mathrm{mod}p)$$

引理2:

$$(1+x{)}^p\equiv 1+x^p(\mathrm{mod}p)$$

$$由二项式定理(1+x{)}^p\equiv \sum _{i=0}^p{C}_p^i{x}^i$$

$$由引理1可知，只剩i=0(C_p^0),p(C_p^p)两项，得证$$

推导证明
$$令n=ap+b,m=cp+d$$

$$(1+x{)}^n=C_n^0{x}^0+c_n^1{x}^1+...+C_n^k{x}^k+C_n^n{x}^n(\mathrm{mod}p)$$

$$(1+x{)}^n\equiv \sum _{i=0}^n{C}_n^i{x}^i(\mathrm{mod}p)$$

$$\equiv (1+x{)}^{ap+b}$$

$$\equiv ((1+x{)}^p{)}^a.(1+x{)}^b$$

$$\equiv \sum _{i=0}^a{C}_a^i{x}^{ip}.\sum _{j=0}^b{C}_b^j{x}^j(\mathrm{mod}p)$$

$$(1)式中x^m的系数是C_n^m$$

$$(2)式中x^m=x^{cp}.x^d的系数为C_a^c.C_b^d$$

$$故C_n^m\equiv C_a^c{C}_b^d(\mathrm{mod}p)$$

$$a=n/p,b=m/p,c=n\%p,d=m\%p$$

$$即C_n^m\equiv C_{n/p}^{m/p}{C}_{n\%p}^{m\%p}$$

int fac[N], pr[N];//fac阶乘数组,pr逆元数组
int ksm(int d,int c,int mod){
    int ans = 1;
    while(c){
        if(c&1){
            ans = ans * di % mod;
        }
        c >>= 1;
        d = d * d % mod;
    }
    return ans;
}
void init(){
    // n,m<=1e18;
    // mod<=1e6;
    fac[0] = pr[0] = 1;
    for (int i = 1; i <= mod;i++){
        fac[i] = fac[i - 1] * i % mod;
        pr[i] = pr[i - 1] * ksm(i, mod - 2, mod)%mod;
    }
}
int C(int n,int m){
    return fac[n] * pr[n - m] * pr[m]%mod;
}
int lucas(int n, int m)
{
    if (m == 0)
    {
        return 1;
    }
    return lucas(n / mod, m / mod) * C(n % mod, m % mod) % mod;
}

高精度组合数（不取模）

O(n)

C(5000,10000)是3009位

首先筛出5000以内的所有质数

再预 处 理 计 算 a 和 b 的 每 个 质 因 子 的 次 数

最后高精度计算结果
$$高精度计算C_n^m=P_1^{a_1}{P}_2^{a_2}...P_k^{a_k}$$

$$故a_1个数等于n!中的个数减去(n-m)!+m!中的个数\frac{p_1^x}{p_1^y{p}_1^z}$$

$$即a_1=x-y-z$$

$$a_1在x!中的个数su{m}_1=\frac{x}{p_1}+\frac{x}{p_1^2}+\frac{x}{p_1^3}+...+\frac{x}{p_1^k}$$

#include <bits/stdc++.h>
#define int long long
using namespace std;
const int N = 1e5 + 10;
int prime[N], cnt, vis[N];
int c[N];
void init()
{
    for (int i = 2; i < N; ++i)
    {
        if (!vis[i])
            prime[cnt++] = i;
        for (int j = 0; i * prime[j] < N; ++j)
        {
            vis[i * prime[j]] = 1;
            if (i % prime[j] == 0)
                break;
        }
    }
}
int get(int n, int p)
{ // n!中p的个数
    int tot = 0;
    while (n)
        tot += n / p, n /= p;
    return tot;
}
int getc(int n, int m, int p)
{ // c中p的个数
    return get(n, p) - get(n - m, p) - get(m, p);
}
void mul(int c[], int p, int &len)
{ // 高精度乘法
    int t = 0;
    for (int i = 0; i < len; i++)
    {
        t += c[i] * p;
        c[i] = t % 10;
        t /= 10;
    }
    while (t)
    {
        c[len++] = t % 10;
        t /= 10;
    }
}
signed main()
{
    ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
    init();
    int n, m;
    cin >> n >> m;
    int len = 1;
    c[0] = 1;
    for (int i = 0; i < cnt; ++i)
    {
        int p = prime[i];
        int s = getc(n, m, p);
        while (s--)
            mul(c, p, len);
    }
    for (int i = len - 1; i >= 0;i--)
        cout << c[i];
    cout << '\n';
    return 0;
}