Lucas 定理

组合意义天地灭。

Lucas 定理

问题 $1$：给定 $n,m\in \mathbb{N}$ 与 $p\in \mathbb{P}$，其中 $n$ 与 $m$ 相当大，而 $p$ 则相对较小，要求计算 $(\genfrac{}{}{0ex}{}{n}{m})modp$ 的值。

一般的预处理逆元以及递推的方法在 $n,m$ 充分大时均会失效，我们需要新的工具来解决这类问题。

对于 $n,m\in \mathbb{N}$ 与 $p\in \mathbb{P}$，有 $(\genfrac{}{}{0ex}{}{n}{m})\equiv \prod (\genfrac{}{}{0ex}{}{n_i}{m_i})(\mathrm{mod}p)$，其中 $n=\sum _{i\in \mathbb{N}}{n}_i{p}^i$，$m=\sum _{i\in \mathbb{N}}{m}_i{p}^i$，也即 $n$ 与 $m$ 的 $p$ 进制表示。

上边的这一坨就是 Lucas 定理了。不过这貌似跟我们平常使用的递归形式不太一样？

对于 $n,m\in \mathbb{N}$ 与 $p\in \mathbb{P}$，有 $(\genfrac{}{}{0ex}{}{n}{m})\equiv (\genfrac{}{}{0ex}{}{nmodp}{mmodp})(\genfrac{}{}{0ex}{}{\lfloor \frac{n}{p}\rfloor }{\lfloor \frac{m}{p}\rfloor })(\mathrm{mod}p)$。

这是我们平常见到的形式。不过，稍作观察就可以发现，这两种形式是等价的，在此不做证明。

下文中，如无特殊说明，均有 $n,m\in \mathbb{N}$ 且 $p\in \mathbb{P}$。

下面我们尝试证明 Lucas 定理。

引理 $1$：对于任意 $x\in \mathbb{N}$ 且 $x<p$，$(\genfrac{}{}{0ex}{}{p}{x})\equiv 0(\mathrm{mod}p)$。

这是好证明的，我们直接从定义下手即可。由于 $(\genfrac{}{}{0ex}{}{p}{x})=\frac{p!}{x!(p-x)!}$，又 $x<p$，因此 $x!$ 中一定没有质因子 $p$，也即 $x!$ 在 $modp$ 意义下存在逆元 $(x!{)}^{-1}$。而 $(p-x)!$ 同理，因此有：

$$(\genfrac{}{}{0ex}{}{p}{x})\equiv p!(x!{)}^{-1}((p-x)!{)}^{-1}(\mathrm{mod}p)$$

右边显然为 $p$ 的倍数。$◻$

引理 $2$：对于任意 $x\in \mathbb{N}$，有 $(1+x{)}^p\equiv 1+x^p(\mathrm{mod}p)$。

通过二项式定理，我们得到 $(1+x{)}^p=\sum _{0\le i\le p}(\genfrac{}{}{0ex}{}{p}{i})x^i$，而根据引理 $1$，我们立即得到上式中除去首项 $(\genfrac{}{}{0ex}{}{p}{0})x^0=1$ 与末项 $(\genfrac{}{}{0ex}{}{p}{p})x^p=x^p$ 外，其余所有项在 $modp$ 意义下均为 $0$，因此引理得证。$◻$

现在来做一些推导：

$$\begin{array}{rl}(1+x{)}^m& =(1+x{)}^{\lfloor \frac{m}{p}\rfloor p+mmodp}\\ & =((1+x{)}^p{)}^{\lfloor \frac{m}{p}\rfloor }(1+x{)}^{mmodp}\\ & \equiv (1+x^p{)}^{\lfloor \frac{m}{p}\rfloor }(1+x{)}^{mmodp}\\ & =\sum _{0\le i\le \lfloor \frac{m}{p}\rfloor }(\genfrac{}{}{0ex}{}{\lfloor \frac{m}{p}\rfloor }{i})x^{ip}\sum _{0\le i\le mmodp}(\genfrac{}{}{0ex}{}{mmodp}{i})x^i\\ & =\sum _{0\le i\le \lfloor \frac{m}{p}\rfloor }\sum _{0\le j\le mmodp}(\genfrac{}{}{0ex}{}{\lfloor \frac{m}{p}\rfloor }{i})(\genfrac{}{}{0ex}{}{mmodp}{j})x^{ip+j}(\mathrm{mod}p)\end{array}$$

第 $3$ 步的转化使用了引理 $2$。容易看出 $ip+j$ 取遍了所有 $[0,m]$ 的整数，于是考虑枚举 $k=ip+j$，则有：

$$(1+x{)}^m\equiv \sum _{0\le k\le m}(\genfrac{}{}{0ex}{}{\lfloor \frac{m}{p}\rfloor }{\lfloor \frac{k}{p}\rfloor })(\genfrac{}{}{0ex}{}{mmodp}{kmodp})x^k(\mathrm{mod}p)$$

又由二项式定理，我们立即得到：

$$\sum _{0\le k\le m}(\genfrac{}{}{0ex}{}{m}{k})x^k\equiv \sum _{0\le k\le m}(\genfrac{}{}{0ex}{}{\lfloor \frac{m}{p}\rfloor }{\lfloor \frac{k}{p}\rfloor })(\genfrac{}{}{0ex}{}{mmodp}{kmodp})x^k(\mathrm{mod}p)$$

那么 Lucas 定理是显而易见的。$◻$

实际使用中，由于给出的 $p$ 一般不会超过 ${10}^6$，我们通常考虑直接计算 $(\genfrac{}{}{0ex}{}{nmodp}{mmodp})$ 的值，而 $(\genfrac{}{}{0ex}{}{\lfloor \frac{n}{p}\rfloor }{\lfloor \frac{m}{p}\rfloor })$ 的值则递归的使用 Lucas 定理进行计算，其边界条件为 $m=0$ 时返回 $1$。如果预处理所有 $i\in [0,p)$ 的阶乘与其逆元的话，我们便可以 $\mathcal{O}(1)$ 计算前边一部分的值，那么容易得到一个 $\mathcal{O}(p\log p)+\mathcal{O}({\log }_{p}n)$ 的算法。

扩展 Lucas 定理

这玩意虽然叫定理，但事实上只是个算法。

发现 Lucas 定理只能解决 $p\in \mathbb{P}$ 的问题，而某些时候题目给定的模数不一定是质数，这个时候我们该怎么办呢？

问题 $2$：给定 $n,m,p\in \mathbb{N}$，其中 $n$ 与 $m$ 相当大，要求计算 $(\genfrac{}{}{0ex}{}{n}{m})modp$ 的值。

由于未保证 $p$ 为质数，普通的 Lucas 定理对此也无能为力。我们考虑将 $p$ 质因数分解，对每一个因子求解后得到一组同余方程，再用 CRT 来合并答案。

对于 $p\in \mathbb{P}$，有：

$$\begin{array}{rl}(\genfrac{}{}{0ex}{}{n}{m})mod{p}^k& =\frac{n!}{m!(n-m)!}mod{p}^k\\ & =\frac{\frac{n!}{p^{g(n)}}}{\frac{m!}{p^{g(m)}}\frac{(n-m)!}{p^{g(n-m)}}}{p}^{g(n)-g(m)-g(n-m)}mod{p}^k\end{array}$$

其中 $g(x)$ 为 $x!$ 中质因子 $p$ 的次数。由于 $\frac{x!}{p^{g(x)}}\perp p^k$，那么 $\frac{x!}{p^{g(x)}}$ 存在 $mod{p}^k$ 的逆元，于是我们考虑求出这个式子的值：

$$\frac{n!}{p^{g(n)}}mod{p}^k$$

考虑记上边的式子为 $f(n)$，然后再来做一些推导，我们考虑将 $n!$ 中所有的 $p$ 的倍数提出来，容易看出一共有 $\lfloor \frac{n}{m}\rfloor$ 个这样的数：

$$\begin{array}{rl}f(n)& =\frac{n!}{p^{g(n)}}\\ & =\frac{(p^{\lfloor \frac{n}{p}\rfloor }\prod _{1\le i\le \lfloor \frac{n}{p}\rfloor }i)(\prod _{1\le i\le n,i\not\equiv 0(\mathrm{mod}p)}i)}{p^{g(n)}}\\ & =\frac{(p^{\lfloor \frac{n}{p}\rfloor }(\lfloor \frac{n}{p}\rfloor )!)(\prod _{1\le i\le n,i\not\equiv 0(\mathrm{mod}p)}i)}{p^{g(n)}}\end{array}$$

由于 $p$ 的倍数在 $1\sim n$ 中分布均匀，并且 $1\sim p^k$ 在 $mod{p}^k$ 下构成一个完全剩余系，那么后面的积式是有循环节的，那么有：

$$\begin{array}{rl}f(n)& =\frac{(p^{\lfloor \frac{n}{p}\rfloor }(\lfloor \frac{n}{p}\rfloor )!)(\prod _{1\le i\le n,i\not\equiv 0(\mathrm{mod}p)}i)}{p^{g(n)}}\\ & \equiv \frac{(p^{\lfloor \frac{n}{p}\rfloor }(\lfloor \frac{n}{p}\rfloor )!)(\prod _{1\le i<p^k,i\not\equiv 0(\mathrm{mod}p)}i{)}^{\lfloor \frac{n}{p^k}\rfloor }(\prod _{\lfloor \frac{n}{p^k}\rfloor p^k\le i\le n,i\not\equiv 0(\mathrm{mod}p)}i)}{p^{g(n)}}(\mathrm{mod}{p}^k)\end{array}$$

这步变换不难理解。$\prod _{1\le i\le n,i\not\equiv 0(\mathrm{mod}p)}i$ 根据 $mod{p}^k$ 的值可以分为 $\lfloor \frac{n}{p^k}\rfloor$ 块，其中前边的块在没有去掉 $p$ 的倍数前均构成 $p^k$ 的一个完系，可以用快速幂求得，而最后一块是散块直接暴力累加就行。

下一步的化简就需要用到 $g$ 的性质了。容易注意到 $g$ 的一种递归式：

$$g(n)=\lfloor \frac{n}{p}\rfloor +g(\lfloor \frac{n}{p}\rfloor )$$

这非常显而易见，不予证明。那么我们可以将式子作进一步化简：

$$\begin{array}{rl}f(n)& \equiv \frac{(p^{\lfloor \frac{n}{p}\rfloor }(\lfloor \frac{n}{p}\rfloor )!)(\prod _{1\le i<p^k,i\not\equiv 0(\mathrm{mod}p)}i{)}^{\lfloor \frac{n}{p^k}\rfloor }(\prod _{\lfloor \frac{n}{p^k}\rfloor p^k\le i\le n,i\not\equiv 0(\mathrm{mod}p)}i)}{p^{g(n)}}\\ & =\frac{p^{\lfloor \frac{n}{p}\rfloor }(\lfloor \frac{n}{p}\rfloor )!}{p^{g(n)}}(\prod _{1\le i<p^k,i\not\equiv 0(\mathrm{mod}p)}i{)}^{\lfloor \frac{n}{p^k}\rfloor }(\prod _{\lfloor \frac{n}{p^k}\rfloor p^k\le i\le n,i\not\equiv 0(\mathrm{mod}p)}i)\\ & =\frac{p^{\lfloor \frac{n}{p}\rfloor }{p}^{g(\lfloor \frac{n}{p}\rfloor )}f(\lfloor \frac{n}{p}\rfloor )}{p^{g(n)}}(\prod _{1\le i<p^k,i\not\equiv 0(\mathrm{mod}p)}i{)}^{\lfloor \frac{n}{p^k}\rfloor }(\prod _{\lfloor \frac{n}{p^k}\rfloor p^k\le i\le n,i\not\equiv 0(\mathrm{mod}p)}i)\\ & =f(\lfloor \frac{n}{p}\rfloor )(\prod _{1\le i<p^k,i\not\equiv 0(\mathrm{mod}p)}i{)}^{\lfloor \frac{n}{p^k}\rfloor }(\prod _{\lfloor \frac{n}{p^k}\rfloor p^k\le i\le n,i\not\equiv 0(\mathrm{mod}p)}i)(\mathrm{mod}{p}^k)\end{array}$$

于是就做完了。最后再用 CRT 合并一下就可以求得答案。

通过分析可以得到，时间复杂度为 $\mathcal{O}(p^k{\log }_{p}n)$。

constexpr int N = 1050;
i64 n, m, p;
i64 a[N], b[N];
i64 fastPow(i64 base, i64 exp, i64 mod) {
  i64 res = 1;
  for (; exp; exp >>= 1) {
    if (exp & 1) (res *= base) %= mod;
    (base *= base) %= mod;
  }
  return res;
}
i64 exgcd(i64 a, i64 b, i64 &x, i64 &y) {
  if (!b) return x = 1, y = 0, a;
  i64 d = exgcd(b, a % b, y, x);
  y -= a / b * x;
  return d;
}
i64 inv(i64 num, i64 mod) {
  i64 x, y;
  exgcd(num, mod, x, y);
  return (x + mod) % mod;
}
i64 f(i64 n, i64 p, i64 pk) {
  if (!n) return 1;
  i64 cir = 1, rst = 1;
  for (i64 i = 1; i <= pk; i++) {
    if (i % p) (cir *= i % pk) %= pk;
  }
  cir = fastPow(cir, n / pk, pk);
  for (i64 i = pk * (n / pk); i <= n; i++) {
    if (i % p) (rst *= i % pk) %= pk;
  }
  return f(n / p, p, pk) * cir % pk * rst % pk;
}
i64 g(i64 n, i64 p) {
  i64 res = 0, cpy = p;
  for (; p <= n; p *= cpy) res += n / p;
  return res;
}
i64 C(i64 n, i64 m, i64 p, i64 pk) {
  i64 a = f(n, p, pk), b = inv(f(m, p, pk), pk), c = inv(f(n - m, p, pk), pk);
  return a * b % pk * c % pk * 
         fastPow(p, g(n, p) - g(m, p) - g(n - m, p), pk) % pk;
}
i64 exLucas(i64 n, i64 m, i64 mod) {
  i64 cpy = mod, res = 0;
  int cnt = 0;
  for (i64 i = 2; i * i <= cpy; i++) {
    if (cpy % i == 0) {
      i64 pk = 1;
      while (cpy % i == 0) cpy /= i, pk *= i;
      a[++cnt] = pk, b[cnt] = C(n, m, i, pk);
    }
  }
  if (cpy != 1) a[++cnt] = cpy, b[cnt] = C(n, m, cpy, cpy);
  for (int i = 1; i <= cnt; i++) {
    i64 k = mod / a[i], t = inv(k, a[i]);
    (res += b[i] * k % mod * t % mod) %= mod;
  }
  return res;
}