盧卡斯定理

Lucas 定理

引入

Lucas 定理用於求解大組合數取模的問題，其中模數必須為素數。正常的組合數運算可以通過遞推公式求解（詳見排列組合），但當問題規模很大，而模數是一個不大的質數的時候，就不能簡單地通過遞推求解來得到答案，需要用到 Lucas 定理。

定義

Lucas 定理內容如下：對於質數 \(p\)，有

\[ \binom{n}{m}\bmod p = \binom{\left\lfloor n/p \right\rfloor}{\left\lfloor m/p\right\rfloor}\cdot\binom{n\bmod p}{m\bmod p}\bmod p \]

觀察上述表達式，可知 \(n\bmod p\) 和 \(m\bmod p\) 一定是小於 \(p\) 的數，可以直接求解，\(\displaystyle\binom{\left\lfloor n/p \right\rfloor}{\left\lfloor m/p\right\rfloor}\) 可以繼續用 Lucas 定理求解。這也就要求 \(p\) 的範圍不能夠太大，一般在 \(10^5\) 左右。邊界條件：當 \(m=0\) 的時候，返回 \(1\)。

時間複雜度為 \(O(f(p) + g(n)\log n)\)，其中 \(f(n)\) 為預處理組合數的複雜度，\(g(n)\) 為單次求組合數的複雜度。

實現

C++Python

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

def Lucas(n, m, p):
  if m == 0:
      return 1
  return (C(n % p, m % p, p) * Lucas(n // p, m // p, p)) % p

證明

考慮 \(\displaystyle\binom{p}{n} \bmod p\) 的取值，注意到 \(\displaystyle\binom{p}{n} = \frac{p!}{n!(p-n)!}\)，分子的質因子分解中 \(p\) 的次數恰好為 \(1\)，因此只有當 \(n = 0\) 或 \(n = p\) 的時候 \(n!(p-n)!\) 的質因子分解中含有 \(p\)，因此 \(\displaystyle\binom{p}{n} \bmod p = [n = 0 \vee n = p]\)。進而我們可以得出

\[ \begin{align} (a+b)^p &= \sum_{n=0}^p \binom pn a^n b^{p-n}\\ &\equiv \sum_{n=0}^p [n=0\vee n=p] a^n b^{p-n}\\ &\equiv a^p + b^p \pmod p \end{align} \]

注意過程中沒有用到費馬小定理，因此這一推導不僅適用於整數，亦適用於多項式。因此我們可以考慮二項式 \(f^p(x)=(ax^n + bx^m)^p \bmod p\) 的結果

\[ \begin{align} (ax^n + bx^m)^p &\equiv a^p x^{pn} + b^p x^{pm} \\ &\equiv ax^{pn} + bx^{pm}\\ &\equiv f(x^p) \end{align} \]

考慮二項式 \((1+x)^n \bmod p\)，那麼 \(\displaystyle\binom n m\) 就是求其在 \(x^m\) 次項的取值。使用上述引理，我們可以得到

\[ \begin{align} (1+x)^n &\equiv (1+x)^{p\lfloor n/p \rfloor} (1+x)^{n\bmod p}\\ &\equiv (1+x^p)^{\lfloor n/p \rfloor} (1+x)^{n\bmod p} \end{align} \]

注意前者只有在 \(p\) 的倍數位置才有取值，而後者最高次項為 \(n\bmod p \le p-1\)，因此這兩部分的卷積在任何一個位置只有最多一種方式貢獻取值，即在前者部分取 \(p\) 的倍數次項，後者部分取剩餘項，即 \(\displaystyle\binom{n}{m}\bmod p = \binom{\left\lfloor n/p \right\rfloor}{\left\lfloor m/p\right\rfloor}\cdot\binom{n\bmod p}{m\bmod p}\bmod p\)。

exLucas 定理

Lucas 定理中對於模數 \(p\) 要求必須為素數，那麼對於 \(p\) 不是素數的情況，就需要用到 exLucas 定理。

過程

第一部分：中國剩餘定理

要求計算二項式係數 \(\binom{n}{m}\bmod M\)，其中 \(M\) 可能為合數。

考慮利用中國剩餘定理合併答案，這種情況下我們只需求出 \(\binom{n}{m}\bmod p^\alpha\) 的值即可（其中 \(p\) 為素數且 \(\alpha\) 為正整數）。

根據 唯一分解定理，將 \(M\) 質因數分解：

\[ M={p_1}^{\alpha_1}\cdot{p_2}^{\alpha_2}\cdots{p_r}^{\alpha_r}=\prod_{i=1}^{r}{p_i}^{\alpha_i} \]

對於任意 \(i,j\)，有 \({p_i}^{\alpha_i}\) 與 \({p_j}^{\alpha_j}\) 互質，所以可以構造如下 \(r\) 個同餘方程：

\[ \left\{ \begin{aligned} a_1\equiv \displaystyle\binom{n}{m}&\pmod {{p_1}^{\alpha_1}}\\ a_2\equiv \displaystyle\binom{n}{m}&\pmod {{p_2}^{\alpha_2}}\\ &\cdots\\ a_r\equiv \displaystyle\binom{n}{m}&\pmod {{p_r}^{\alpha_r}}\\ \end{aligned} \right. \]

我們發現，在求出 \(a_i\) 後，就可以用中國剩餘定理求解出 \(\displaystyle\binom{n}{m}\)。

第二部分：移除分子分母中的素數

根據同餘的定義，\(\displaystyle a_i=\binom{n}{m}\bmod {p_i}^{\alpha_i}\)，問題轉化成，求 \(\displaystyle \binom{n}{m} \bmod p^\alpha\)（\(p\) 為質數）的值。

根據組合數定義 \(\displaystyle \binom{n}{m} = \frac{n!}{m! (n-m)!}\)，\(\displaystyle \binom{n}{m} \bmod p^\alpha = \frac{n!}{m! (n-m)!} \bmod p^\alpha\)。

由於式子是在模 \(p^\alpha\) 意義下，所以分母要算乘法逆元。

同餘方程 \(ax \equiv 1 \pmod p\)（即乘法逆元）有解的充要條件為 \(\gcd(a,p)=1\)（裴蜀定理），

然而 無法保證有解，發現無法直接求 \(\operatorname{inv}_{m!}\) 和 \(\operatorname{inv}_{(n-m)!}\)，

所以將原式轉化為：

\[ \frac{\frac{n!}{p^x}}{\frac{m!}{p^y}\frac{(n-m)!}{p^z}}p^{x-y-z} \bmod p^\alpha \]

\(x\) 表示 \(n!\) 中包含多少個 \(p\) 因子，\(y, z\) 同理。

第三部分：Wilson 定理的推論

問題轉化成，求形如：

\[ \frac{n!}{q^x}\bmod q^k \]

的值。這時可以利用 Wilson 定理的推論。如果難以理解，可以看看下面的解釋。

解釋

一個示例：22! mod 9

先考慮 \(n! \bmod q^k\)，

比如 \(n=22, q=3, k=2\) 時：

\(22!=1\times 2\times 3\times 4\times 5\times 6\times 7\times 8\times 9\times 10\times 11\times 12\)

\(\times 13\times 14\times 15\times 16\times 17\times 18\times 19\times20\times21\times22\)

將其中所有 \(q\) 的倍數提取，得到：

\(22!=3^7 \times (1\times 2\times 3\times 4\times 5\times 6\times 7) \times (1\times 2\times 4\times 5\times 7\times 8\times 10 \times 11\times 13\times 14\times 16\times 17\times 19 \times 20 \times 22 )\)

可以看到，式子分為三個整式的乘積：

是 \(3\) 的冪，次數是 \(\lfloor\frac{n}{q}\rfloor\)；
是 \(7!\)，即 \(\lfloor\frac{n}{q}\rfloor!\)，由於階乘中仍然可能有 \(q\) 的倍數，考慮遞歸求解；
是 \(n!\) 中與 \(q\) 互質的部分的乘積，具有如下性質：
\(1\times 2\times 4\times 5\times 7\times 8\equiv10 \times 11\times 13\times 14\times 16\times 17 \pmod{3^2}\)，
即：\(\displaystyle \prod_{i,(i,q)=1}^{q^k}i\equiv\prod_{i,(i,q)=1}^{q^k}(i+tq^k) \pmod{q^k}\)（\(t\) 是任意正整數）。
\(\displaystyle \prod_{i,(i,q)=1}^{q^k}i\) 一共循環了 \(\displaystyle \lfloor\frac{n}{q^k}\rfloor\) 次，暴力求出 \(\displaystyle \prod_{i,(i,q)=1}^{q^k}i\)，然後用快速冪求 \(\displaystyle \lfloor\frac{n}{q^k}\rfloor\) 次冪。
最後要乘上 \(\displaystyle \prod_{i,(i,q)=1}^{n \bmod q^k}i\)，即 \(19\times 20\times 22\)，顯然長度小於 \(q^k\)，暴力乘上去。

上述三部分乘積為 \(n!\)。最終要求的是 \(\frac{n!}{q^x}\bmod{q^k}\)。

所以有：

\[ n! = q^{\left\lfloor\frac{n}{q}\right\rfloor} \cdot \left(\left\lfloor\frac{n}{q}\right\rfloor\right)! \cdot {\left(\prod_{i,(i,q)=1}^{q^k}i\right)}^{\left\lfloor\frac{n}{q^k}\right\rfloor} \cdot \left(\prod_{i,(i,q)=1}^{n\bmod q^k}i\right) \]

於是：

\[ \frac{n!}{q^{\left\lfloor\frac{n}{q}\right\rfloor}} = \left(\left\lfloor\frac{n}{q}\right\rfloor\right)! \cdot {\left(\prod_{i,(i,q)=1}^{q^k}i\right)}^{\left\lfloor\frac{n}{q^k}\right\rfloor} \cdot \left(\prod_{i,(i,q)=1}^{n\bmod q^k}i\right) \]

\(\displaystyle \left(\left\lfloor\frac{n}{q}\right\rfloor\right)!\) 同樣是一個數的階乘，所以也可以分為上述三個部分，於是可以遞歸求解。

等式的右邊兩項不含素數 q。事實上，如果直接把 n 的階乘中所有 q 的冪都拿出來，等式右邊的階乘也照做，這個等式可以直接寫成：

\[ \frac{n!}{q^{x}} = \frac{\left(\left\lfloor\frac{n}{q}\right\rfloor\right)!}{q^{x'}} \cdot {\left(\prod_{i,(i,q)=1}^{q^k}i\right)}^{\left\lfloor\frac{n}{q^k}\right\rfloor} \cdot \left(\prod_{i,(i,q)=1}^{n\bmod q^k}i\right) \]

式中的 \(x\) 和 \(x'\) 都表示把分子中所有的素數 \(q\) 都拿出來。改寫成這樣，每一項就完全不含 \(q\) 了。

遞歸的結果，三個部分中，左邊部分隨着遞歸結束而自然消失，中間部分可以利用 Wilson 定理的推論 0，右邊部分就是推論 2 中的 \(\prod_{j\geq 0}(N_j!)_p\)。

下面這種寫法，擁有單次詢問 \(O(p\log p)\) 的時間複雜度。其中 int inverse(int x) 函數返回 \(x\) 在模 \(p\) 意義下的逆元。

實現

LL calc(LL n, LL x, LL P) {
  if (!n) return 1;
  LL s = 1;
  for (LL i = 1; i <= P; i++)
    if (i % x) s = s * i % P;
  s = Pow(s, n / P, P);
  for (LL i = n / P * P + 1; i <= n; i++)
    if (i % x) s = i % P * s % P;
  return s * calc(n / x, x, P) % P;
}

LL multilucas(LL m, LL n, LL x, LL P) {
  int cnt = 0;
  for (LL i = m; i; i /= x) cnt += i / x;
  for (LL i = n; i; i /= x) cnt -= i / x;
  for (LL i = m - n; i; i /= x) cnt -= i / x;
  return Pow(x, cnt, P) % P * calc(m, x, P) % P * inverse(calc(n, x, P), P) %
         P * inverse(calc(m - n, x, P), P) % P;
}

LL exlucas(LL m, LL n, LL P) {
  int cnt = 0;
  LL p[20], a[20];
  for (LL i = 2; i * i <= P; i++) {
    if (P % i == 0) {
      p[++cnt] = 1;
      while (P % i == 0) p[cnt] = p[cnt] * i, P /= i;
      a[cnt] = multilucas(m, n, i, p[cnt]);
    }
  }
  if (P > 1) p[++cnt] = P, a[cnt] = multilucas(m, n, P, P);
  return CRT(cnt, a, p);
}

若不考慮 excrt 的複雜度，通過預處理 \(\frac{n!}{n以內的p的所有倍數的乘積}\bmod{p}\)，可以使時間複雜度優化至單次 \(O(p + \log p)\)。而如果 p 是固定的，我們在一開始就可以對 p 進行分解，並進行預處理，可以達到總複雜度 \(O(p + T\log p)\)。

習題

本页面最近更新：，更新历史
发现错误？想一起完善？在 GitHub 上编辑此页！
本页面贡献者：OI-wiki
本页面的全部内容在 CC BY-SA 4.0 和 SATA 协议之条款下提供，附加条款亦可能应用