算法流程 由于 � m 不是质数，我们考虑以下的这三个步骤。它们组成了扩展 Lucas 的思路框架。 我们将模数 � m 分解质因数： � = � 1 � 1 � 2 � 2 ⋯ � � � � m=p 1 a 1p 2 a 2⋯p n a n； 我们构造方程组： { � ≡ � 1 ( m o d � 1 � 1 ) � ≡ � 2 ( m o d � 2 � 2 ) ⋯ ⋯ � ≡ � � ( m o d � � � � ) ⎩ ⎨ ⎧x≡r 1(modp 1 a 1) x≡r 2(modp 2 a 2) ⋯⋯ x≡r n(modp n a n)注意到这些方程里的模数两两互质。我们使用 CRT 合并求解 � x ，即为我们想要的答案 ( � � ) m o d � ( r n)modm 。 接下来的问题就完全聚焦到了 � � r i的构造上来。我们考虑如下流程： 明确我们要构造 � � ≡ ( � � ) ( m o d � � � � ) r i≡( r n)(modp i a i) ； 转化问题。我们需要求解 ( � � ) m o d � � ( r n)modp a ，其中 � p 是质数。即我们要求解如下式子： � ! � ! ( � − � ) ! m o d � � n!(n−r)! n!modp a 注意到分母 � ! ( � − � ) ! r!(n−r)! 与 � � p a 不一定互质。那么就意味着我们没有办法求分母的逆元，进一步我们无法求解这个分式。考虑我们该如何改写它。 我们目标是使分母和 � � p a 互质。为此我们将式子改写为： � � ( � ! ) � ( � ! ) � ( ( � − � ) ! ) � m o d � � p v (r!) p((n−r)!) p(n!) pmodp a 其中 � = � � ( � ! ) − � � ( � ! ) − � � ( ( � − � ) ! ) v=v p(n!)−v p(r!)−v p((n−r)!) 。注意到若 � ≥ � v≥a ，则结果为 0 0 。 计算 � � ( � ! ) v p(n!) 。这可以利用 Legendre 公式求出： � � ( � ! ) = ∑ � = 1 ∞ ⌊ � � � ⌋ v p(n!)= i=1 ∑ ∞⌊ p i n⌋ 计算 ( � ! ) � (n!) p。利用前置知识中提到的递归公式： ( � ! ) � ≡ ( ± 1 ) ⌊ � / � � ⌋ ( ⌊ � / � ⌋ ! ) � ( ∏ 1 ≤ � ≤ ( � m o d � � ) 且 � ∤ � � ) ( m o d � � ) (n!) p≡(±1) ⌊n/p a ⌋ (⌊n/p⌋!) p1≤j≤(nmodp a )且j∤p ∏j(modp a ) 一些注意事项： i. ( ⌊ � / � ⌋ ! ) � (⌊n/p⌋!) p可以利用循环来展开递归，加快速度； ii. ∏ � ∏j 项，等价于在 ( � m o d � � ) ! (nmodp a )! 中去掉 � p 因子，故这里可以直接预打表，用 f[i] 来表示 � ! i! 中去掉 � p 因子的值。这个表只用打到 � < � � i
#include #define fastio ios::sync_with_stdio(0);cin.tie(0);cout.tie(0) #define int long long #define MIKU 0 using namespace std; int n, m, p; //求解 C(n,m) mod p。 //扩展欧几里得求逆元。 int exgcd(int a, int b, int &x, int &y) { if(b == 0) {x = 1; y = 0; return a;} int d = exgcd(b, a%b, y, x); y -= a / b * x; return d; } int inv(int a, int m) { int x, y; exgcd(a, m, x, y); return (x % m + m) % m; } //对每个质因数构造 r_i struct BPP { int p, a, pa; //p:模数 a:次数 pa:p^a。 vector f; //预处理 f[i] ，表示 i! 去掉 p 因子。 //构造函数，计算出 pa 和 f。 BPP(int prime, int power) : p(prime), a(power) { pa = 1; for(int i=1; i<=a; i++) pa *= p; f.resize(pa); f[0] = 1; for(int i=1; i1; n/=p) { if((n / pa) & neg) res = pa - res; //这一步完成符号修正，即公式中第一个因子。 res = res * f[n%pa] % pa; } return res; } //利用 vp(n!) 和 (n!)_p 算 C(n,r) mod p^a。 int C(int n, int r) { int v = vp(n) - vp(r) - vp(n-r); if(v >= a) return 0; int res = fac(n) * inv(fac(r) * fac(n-r) % pa, pa) % pa; for(; v; v--) res = res * p % pa; return res; } }; //扩展 Lucas ，负责完成质因数分解和 CRT 合并。 int exLucas(int n, int r, int m, int res = 0) { int tmp = m; for(int p=2; p*p<=tmp; p++) { //试除法分解质因数。 if(tmp % p) continue; int a = 0, pa = 1; for(; tmp%p==0; tmp/=p) a ++, pa *= p; BPP bpp = BPP(p, a); int ri = bpp.C(n, r); int mi = m / pa; res = (res + mi * ri % m * inv(mi, pa)) % m; // CRT 合并公式。 } if(tmp > 1) { BPP bpp = BPP(tmp, 1); int ri = bpp.C(n, r); int mi = m / tmp; res = (res + mi * ri % m * inv(mi, tmp)) % m; } return res; } signed main() { fastio; cin>>n>>m>>p; cout< � s m>n 是无解。 主要难点是不保证模数为质数。这就需要用到扩展 Lucas 定理。 下面贴出代码。相比于上面的模板代码，只改动了主函数。 Show Me the Code
---------------------
著作权归作者所有。
商业转载请联系作者获得授权，非商业转载请注明出处。
作者：Yzu_EtherealYz
源地址：https://www.cnblogs.com/Yzu-EtherealYz/p/20334937
来源：博客园cnblogs
© 版权声明：本文为博主原创文章，转载请附上博文链接！