专题·Lucas定理【including Lucas定理，扩展Lucas

初见安~这里是数论专题（6）【详见数论专栏

本篇有前置知识点需要掌握，建议先了解下：费马小定理，中国剩余定理，乘法逆元

一、Lucas定理

Lucas定理用于求解 $C_{n}^{m}\textrm{} modp$ 的组合数取模的问题。其中p为质数。

组合数取模我们可以用分解质因数再分别相乘取模相加，但有时也很容易超时。

Lucas定理中，将n和m转化为p进制，即有：

$n = n_k * p^k + n_{k - 1} * p^{k -1} +...+n_1 * p + n_0$

$m = m_k * p^k + m_{k - 1} * p^{k -1} +...+m_1 * p + m_0$

则有：

$C_{n}^{m} \equiv C_{n_k}^{m_k} * ...*C_{n_0}^{m_0} \equiv \prod_{i = 0}^{k}C_{n_i}^{m_i}(modp)$

这里问题就来了——组合数取模，最后还是转化成了一堆组合数的累乘，有什么有吗？？？

其实就相当于把一个很大的乘法分解了，使其不容易越界，不容易超时。

而我们把C（n， m）化小了过后，可以选择直接暴力用组合公式： $C_{n}^{m} = \frac{n!}{m!(n - m)!}$

当然分子和分母是有可以化简的部分的，可以化简为 $\frac{(n - m + 1) * (n - m) *...*n}{m!}$

其中要除以m！，我们可以直接乘m的逆元 $m^{-1}$ 。而根据费马小定理（传送门同逆元），p为质数，我们可以得出：设x = m！，则有 $x^{p - 1} \equiv1(modp)$ 。而根据逆元的性质，我们设x的逆元为a，必有 $a * x \equiv 1(modp)$ 。两式联立便可以得到：

$a = x^{p - 2}$ 。所以我们用快速幂就可以求得x的逆元了。

核心代码实现如下：

typedef long long ll;//求其逆元
ll mod_power(ll a, ll b, ll p)
{ll ans = 1;while(b){if(b & 1) ans = ans * a % p;a = a * a % p;b >>= 1;}return ans;
}//求组合数。因为分解得较小了，所以可以用暴力
ll C(ll n, ll m, ll p)
{if(m > n) return 0;ll c1 = 1, c2 = 1;for(int i = n - m + 1; i <= n; i++)//由组合数公式化简得c1 = c1 * i % p;for(int i = 2; i <= m; i++)c2 = c2 * i % p;return c1 * mod_power(c2, p - 2, p) % p;//乘其逆元，是一种优化
}//Lucas定理
ll Lucas(ll n, ll m, ll p)
{if(!m) return 1;return C(n % p, m % p, p) * Lucas(n / p, m / p, p) % p;//累乘，递归实现
}

二、扩展Lucas定理（exLucas）

各种的所谓扩展，都是在原有的基础上去掉一些限制条件。比如——扩展Lucas定理，就是求解 $C_{n}^{m} modp$ 的问题，且不保证p为质数。

很明显，如果p过大，会很容易导致我们求组合数的时候越界，而且不是质数，不能单纯地用Lucas定理。所以——既然不是质数，变成质数就行了嘛。怎么变？分解质因数。

（1）我们将p质因数分解为t个质数的k次方的乘积： $p = p_1^{k_1} * p_2^{k_2} *...*p_t^{k_t}$ 。

然后运用到一个显而易见的结论—— $x \equiv x mod p(mod p)$ 来建立同余方程：

$\left\{\begin{matrix} C_{n}^{m} \equiv C_{n}^{m} mod p_1^{k_1}(modp_1^{k_1})& & & & & & & & & \\ C_{n}^{m} \equiv C_{n}^{m} mod p_2^{k_2}(modp_2^{k_2})& & & & & & & & & \\ ......& & & & & & & & & \\ C_{n}^{m} \equiv C_{n}^{m} mod p_t^{k_t}(modp_t^{k_t})& & & & & & & & & \end{matrix}\right.$

又因为满足mod的数为分解后的质因数的幂，所以互质，可以用中国剩余定理求解。看起来你还不是要求 $C_{n}^{m}$ ，而且求很多次。事实上，因为mod的数已经分解了，所以已经是比暴力求阶乘要优很多了。所以这就是我们的第一步。

（2）接下来我们要求 $C_{n}^{m} mod p_i^{k_i}$ 了。当然，入口可以直接设立在分解质因数的时候。参考Lucas定理的思路，我们这里用排列组合公式来化简。同样的要运用到逆元。但是——这里要mod的数并不是质数，如果我们为了求逆元而用费马小定理去套用欧拉函数，那是一定复杂化了的。不如直接点，就像把p变成质数一样——我们也把n!、m！和（n - m）！变成和pi的ki次方互质的数，也就是除以其最大包含有pi^ki的约数。（比如100和5不互质，但100/（5^2）= 4就可以了）就可以转化为：（这里的k和第一步中质因数分解的k没有任何关系）

$\frac{n!}{m!(n - m)!} = \frac{\frac{n!}{p_i^{k_1}}}{\frac{m! * (n - m)!}{p_i^{k_2} * p_i^{k_3}}} * p^{k_1 - k_2 - k_3}$ 。这样一来我们就可以用到逆元了——求 $\frac{m!}{p^{k_2}}$ 和 $\frac{(n - m)!}{p^{k_3}}$ 的逆元。当然，因为会用到多次，如果调用快速幂就不方便了，所以我们选用原始一点的扩欧定理，将其转化为线性同余方程来求解。【可能用快速幂也行的通】

（3）在第二步的基础上，我们又将问题转化为了——如何求解（这里就简写p和k了，k是（1）中的意义，a是（2）中的k，是那么个意思就行了） $\frac{n!}{p^a} mod p^k$ 。很明显，p^a我们可以算，p^k我们可以传过来，所以关键就在于n！。暴力？绝对不行。不如找找有没有能简化的步骤。

举个例子来看看吧——19！ % 3^2。

其中， $19! = 1 * 2 * 3 * 4 * 5 * .. *19%u3002$ 。

不难发现，可以把3的倍数提出来打包一下变成： $3 * (1 * 2 * 3 * 4 * 5 * 6) * 1 * 2 * 4 * ... *17 * 19$

也就是 $3 * 6! * 1 * 2 * 4 * ... * 17 * 19$

至于后面那一串，我们展开来看（建议手算），可以发现：每3^2长度mod 3^2的值都是一样的。也就是说：

$(1 * 2 * 4 * 5 * 7 * 8) \equiv(10 * 11 * 13 * 14 * 16 * 17) \equiv 8(mod9)$

所以至此，我们可以把求阶乘的这一步分解为三个部分——p幂次的阶乘*p，n / （p^k）个同余的式子和最后剩下的n % （p^k）的部分。其中求p幂次的阶乘，我们亦可以递归调用求阶乘的fac函数。所以到这里——这个问题就解决啦！！！其中有不少为了优化而略显复杂的操作，但逻辑上还是十分清晰的。

下面上代码——【过程太长，就上完整代码了】

#include<bits/stdc++.h>
using namespace std;
long long c[1000005], a[1000005];
long long power(long long a, long long b, long long p)//快速幂
{long long ans = 1;while(b){if(b & 1) ans = ans * a % p;a = a * a % p;b >>= 1;}return ans;
}long long fac(long long n, long long p, long long pk)//求阶乘
{if(!n) return 1;long long ans = 1;for(int i = 1; i < pk; i++)if(i % p) ans = ans * i % pk; //同余部分 ans =  power(ans, n / pk, pk);for(int i = 1; i <= n % pk; i++)//剩余无法凑同余的部分if(i % p) ans = ans * i % pk; return ans * fac(n / p, p, pk) % pk;
}long long exgcd(long long a, long long b, long long &x, long long &y)
{if(!b){x = 1, y = 0;return a;}long long xx, yy, g = exgcd(b, a % b, xx, yy);x = yy;y = xx - a / b * yy;return g;
}long long inv(long long a, long long p)//求逆元
{long long y, x;exgcd(a, p, x, y);//扩欧定理求解return (x % p + p) % p;
} long long C(long long n, long long m, long long p, long long pk)//求组合数
{if(m > n) return 0;long long f1 = fac(n, p, pk), f2 = fac(m, p, pk), f3 = fac(n - m, p, pk), cnt = 0;//因为要求逆元，所以三个的阶乘都要老老实实算出来for(long long i = n; i; i /= p)cnt += i / p;       for(long long i = m; i; i /= p)cnt -= i / p;for(long long i = n - m; i; i /= p)cnt -= i / p;//这里不能用费马，p不为素数。//cnt计算的是step2中的(k1 - k2 - k3)return f1 * inv(f2, pk) % pk * inv(f3, pk) % pk * power(p, cnt, pk) % pk;
}long long CRT(long long cnt)//中国剩余定理
{long long M = 1, ans = 0;for(int i = 1; i <= cnt; i++)M *= c[i];//p的值变了，所以要重新计算for(int i = 1; i <= cnt; i++){ans = (ans + a[i] * (M / c[i]) % M * inv(M / c[i], c[i]) % M) % M;} return ans;
}long long exlucas(long long n, long long m, long long p)
{long long temp, cnt= 0;for(int i = 2; p > 1 && i <= p / i; i++)//分解质因数{long long tmp = 1; while (p % i == 0)p /= i, tmp *= i;if (tmp > 1)//有这个p a[++cnt] = C(n, m, i, tmp), c[cnt] = tmp;//c里存pk }if(p > 1) c[++cnt] = p, a[cnt] = C(n, m, p, p);return CRT(cnt);
}int main()
{long long m, n, p;scanf("%lld%lld%lld", &n, &m, &p);printf("%lld", exlucas(n, m, p));
}

看起来代码真的很长。其实只要理解了推导过程就还是很好记的~

以上就是关于Lucas定理的内容啦~

迎评：）
——End——