本文只解决一个问题：求$C_a^b\mathrm{\%}p$.

递归公式（杨辉三角)

原理

众所周知， $C_a^b=C_{a-1}^b+C_{a-1}^{b-1}$.

可以从组合数的现实意义上证明其正确性：

$C_a^b$表示从a个不同物体中选出b个的所有方案数。

设a个物体中有一个物体x，则$C_a^b$可以分为 包含x的方案 和不包含x的方案。

• 包含x， 则还需要从剩余的a-1个物体中再选b-1个， 即为$C_{a-1}^{b-1}$
• 不包含x， 则需要从另外的a-1个物体中选b个， 即为$C_{a-1}^b$

因此， $C_a^b=C_{a-1}^b+C_{a-1}^{b-1}$.

使用此公式可以以$O(ab)$的时间复杂度预处理出$C_0^0\mathrm{到}{C}_a^b$范围中所有的组合数。

此方法适用于a*b在$1e8$范围内且需要频繁查询的问题。

代码实现

c++

for (int i = 0; i <= n; i ++) {
 C[i][0] = 1;
 for (int j = 0; j <= i; j ++) {
  C[i][j] = (C[i - 1][j] + C[i - 1][j - 1]) % p;
 }
}

阶乘公式

原理

众所周知， $C_a^b=\frac{a!}{b!(a-b)!}$.

可以从排列数的角度证明其正确性：

$A_a^b$表示有先后顺序的从a个不同物体中选出b个的所有方案数，$A_a^b=\frac{a!}{(a-b)!}$。

若只选取固定的b个物品， 则选取顺序显然有$b!$种。

即$A_a^b=b!C_a^b$, $C_a^b=\frac{a!}{b!(a-b)!}$

通常的， 可以预处理 1~a的阶乘以及阶乘的逆元， 如此一来使用$O(n)$的时间复杂度预处理， 即可$O(1)$的完成每次查询。

但有一个值得注意的问题是， 既然需要计算逆元， 就得保证各阶乘的逆元存在。当p大于a时， 1~a的阶乘必然是存在的， 因为x与p互质、y与p互质， 则x*y同样也与p互质。

所以此方法适用于a小于$1e8$且a小于p的问题。

代码实现

fact[0] = inv[0] = 1;
for (int i = 1; i <= n; i ++) {
 fact[i] = (fact[i - 1] * i) % p;
 inv[i] = (inv[i-1] * modPow(i, p - 2, p)) % p;
}

Lucas定理

原理

Lucas定理:

对于素数p， 有

$$C_a^b\equiv C_{a\mathrm{\%}p}^{b\mathrm{\%}p}\times C_{\lfloor \frac{a}{p}\rfloor }^{\lfloor \frac{b}{p}\rfloor }(modp)$$

Lucas定理常见的应用场景是: a非常大， 而p较小。

在实现时， 因为p比较小， 我们一般直接计算$C_{a\mathrm{\%}p}^{b\mathrm{\%}p}$. 而$C_{\lfloor \frac{a}{p}\rfloor }^{\lfloor \frac{b}{p}\rfloor }$则递归的使用Lucas定理计算。

代码实现

ll lucas(ll a, ll b, ll p) {
    return b == 0 ? 1 % p : C(a % p, b % p, p) * lucas(a / p, b / p, p) % p;
}

模板题

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

输入格式

第一行包含整数 $n$。

接下来 $n$行，每行包含一组 $a,b,p$。

输出格式

共 $n$行，每行输出一个询问的解。

数据范围

$1\le n\le 20$$1\le b\le a\le {10}^{18}$, $1\le p\le {10}^5$,

输入样例

3
5 3 7
3 1 5
6 4 13

输出样例

3
3
2

参考代码

# include <bits/stdc++.h>
using namespace std;
using ll = long long;
const int N = 1e5;
ll fact[N], inv[N];

ll modPow(ll a, ll n, ll p) {
    ll res = 1 % p;
    while (n) {
        if (n & 1) res = res * a % p;
        a = a * a % p, n >>= 1;
    }
    return res;
}

ll C(ll a, ll b, ll p) {
    return fact[a] * inv[b] % p * inv[a - b] % p;
}

ll lucas(ll a, ll b, ll p) {
    return b == 0 ? 1 % p : C(a % p, b % p, p) * lucas(a / p, b / p, p) % p;
}

int main () {
    int n; cin >> n;
    while (n -- ) {
        ll a, b, p; cin >> a >> b >> p;
        fact[0] = inv[0] = 1 % p;
        for (ll i = 1; i < p; i ++) {
            fact[i] = fact[i - 1] * i % p;
            inv[i] = inv[i - 1] * modPow(i, p - 2, p) % p;
        }
        cout << lucas(a, b, p) << endl;
    }
    return 0;
}

时间复杂度

$lucas$函数的执行次数显然是${\log }_{p}a$左右次， 而每次执行的时间复杂度取决于$C$函数。

• 如果像如上代码一样， 预处理p内的阶乘fact和逆元inv， 则C函数的时间复杂度为$O(1)$, lucas算法的时间复杂度为$O(p+{\log }_{p}a)$
• 如果不予记录， 而是每次都在C函数内递推一遍， 则C函数的时间复杂度为$O(p)$, lucas算法的时间复杂度为$O(p{\log }_{p}a)$.

参考

二项式系数 - 维基百科，自由的百科全书 (wikipedia.org)

排列组合 - OI Wiki (oi-wiki.org)

卢卡斯定理 - OI Wiki (oi-wiki.org)

卢卡斯定理 - 维基百科，自由的百科全书 (wikipedia.org)

算法学习笔记(25): 卢卡斯定理 - 知乎 (zhihu.com)