组合数的计算方法

#计算组合数的几种方法

前置知识 费马小定理(求乘法逆元) $a^{p-1}\equiv 1 (mod p)$ 若$x$为逆元 $\frac{a}{b} \equiv a \times x (mod p)$

可以得到$\Longrightarrow x = b^{p-2}$

第一种 暴力

根据定义

$C_a^b = \frac{a\times (a-1) \times (a-2) \times ... \times (a-b+1)}{1 \times 2 \times 3 \times ... \times b}$

在 $a < 60$ 的时候可以直接用 $long long$ 暴力循环求。

long long C(int a, int b)
{
    long long res = 1;
    for (int i=1,j=a;i<=b;i++,j--)
    {
        res = res * j / i;
    }
    return res;
}

第二种 递推公式

给定 $n$ 组询问，每组询问给定两个整数 $a$，$b$，请你输出 $C_a^b\mod(10^9+7)$ 的值。

数据范围 $1≤n≤10000$ $1≤b≤a≤2000$

不难发现，$a$ 和 $b$ 的范围不大，所有的组合有 $a \times b$ 共$4\times10^6$ 种，将所有情况预处理出来，对每次查询$O(1)$，可以在规定时间内通过。 利用递推公式：

 
for (int i=0;i<N;i++)
{
    for (int j=0;j<=i;j++)
    {
        if (j == 0) c[i][j] = 1;
        else c[i][j] = c[i-1][j] + c[i-1][j-1];
    }
}

第三种 快速幂

给定 $n$ 组询问，每组询问给定两个整数 $a$，$b$，请你输出 $C_a^b\mod(10^9+7)$ 的值。

数据范围 $1≤n≤10000$ $1≤b≤a≤10^5$

相比于第一种，第二种的 $a$ 与 $b$ 扩大了范围。这时可以用组合数的定义式直接求，但要用到快速幂算法求逆元。 $C_{a}^{b} \mod p= \frac{a!}{b!\times(a-b)!} \mod p$ $C_{a}^{b} \mod p= a! \times (b!)^{-1} \times ((a-b)!)^{-1} \mod p$

在计算阶乘的时候，可以进行预处理，阶乘的逆元同理。

 
LL quickPower(LL a,LL b,LL m)
{
    LL res = 1;
    while (b)
    {
        if (b & 1) res *= a;
        res %= m;
        b >>= 1;
        a *= a;
        a %= m;
    }
    return res;
}
 
fact[0] = infact[0] = 1;
for (int i=1;i<N;i++)
{
    fact[i] = fact[i-1] * i % mod;
    infact[i] = infact[i-1] * quickPower(i, mod-2, mod) % mod;
}

逆元阶乘可以直接相乘，是不是说 $a$ 的逆元乘上 $b$ 的逆元就是 $a\times b$ 的逆元？

是的，如果 $x$ 是 $a$ 的逆元，$y$ 是 $b$ 的逆元，那么 $ax≡1$,$by≡1$所以 $abxy≡1$，所以 $xy$ 是 $ab$ 的逆元。

第四种 $Lucas$ 定理

$Lucas定理$

给定 $n$ 组询问，每组询问给定三个整数 $a,b,p$，其中 $p$ 是质数，请你输出 $C_a^b\mod p$ 的值。

数据范围 $1≤n≤20$ $1≤b≤a≤10^{18}$ $1≤p≤10^5$

在这种情况，测试数据组数比较小，但是 $a$ 和 $b$ 的范围很大，这时要用到 $Lucas$ 定理。

LL qmi(LL a, LL k) //快速幂
{
    LL res = 1;
    while (k)
    {
        if (k & 1) res *= a;
        res %= p;
        k >>= 1;
        a *= a;
        a %= p;
    }
    return res;
}
 
LL C(LL a,LL b) //定义求组合数
{
    LL res = 1;
    for (LL i=1, j=a;i<=b;i++,j--)
    {
        res = res * j % p;
        res = res * qmi(i, p-2) % p;
    }
    return res;
}
 
LL lucas(LL a, LL b) //Lucas定理
{
    if (a < p) return C(a,b);
    return C(a % p, b % p) * lucas(a / p, b / p) % p;
}