1. 项目概述当组合数遇上非质数模数在算法竞赛和数论研究中计算组合数 C(n, m) 对一个大整数 P 取模的结果是一个经典且高频的需求。当模数 P 是一个质数时我们有成熟的卢卡斯定理Lucas Theorem可以高效解决其核心思想是将 n 和 m 按 P 进制分解递归求解。然而现实世界尤其是题目里往往没那么“理想”——模数 P 经常不是一个质数它可能是一个合数甚至是质数幂的乘积。这时标准的卢卡斯定理就失效了。这正是“扩展卢卡斯定理”exLucas要解决的痛点它允许我们在模数为任意正整数不一定是质数的情况下计算组合数取模。这个“模板”项目的核心就是理解和实现 exLucas 算法。它不是一个单一的公式而是一个精巧的算法框架巧妙地将中国剩余定理、质因数分解、阶乘模质数幂的计算等多个数论工具串联起来。对于需要处理模数非质数的组合数问题比如某些计数类动态规划、多项式系数计算或者直接就是一道要求计算 C(n, m) mod P 的题目掌握 exLucas 是必不可少的技能。2. 算法核心思路拆解化整为零分而治之exLucas 算法的核心思想是“分解与合并”。面对一个非质数模数 P我们无法直接在其上定义模逆元因为不是所有数都有逆元所以需要转换思路。2.1 第一步模数质因数分解假设我们要计算 C(n, m) mod P其中 P 是任意正整数。 首先对模数 P 进行质因数分解P p1^k1 * p2^k2 * ... * pt^kt其中 pi 是互不相同的质数ki 是正整数。算法的目标转化为分别求出C(n, m) mod pi^ki对于每一个 i 的结果然后再利用中国剩余定理将这些结果合并得到最终C(n, m) mod P。为什么这么做因为模一个质数幂pi^ki的环境虽然比模质数复杂但比模任意合数要规则得多。我们可以在每个pi^ki的“局部”环境下解决问题。2.2 第二步解决子问题 C(n, m) mod p^k这是整个算法最难也是最核心的部分。对于某个质数 p 和指数 k我们需要计算C(n, m) mod p^k。组合数公式为C(n, m) n! / (m! * (n-m)!)。 在模p^k下分母m!和(n-m)!可能含有因子 p因此它们可能没有模p^k下的逆元因为与模数不互素。我们不能直接计算除法。exLucas 的巧妙之处在于对阶乘进行“净化”处理将阶乘 n! 写成两部分乘积的形式n! p^a * b。 其中a是 n! 中质因子 p 的个数b是 n! 中剔除了所有因子 p 之后剩下的部分并且b与 p 互质。例如对于 p2, k3计算 10! mod 8。 10! 12345678910。 其中因子 2 出现在 2, 42^2, 623, 82^3, 1025。 总共有 12131 8 个因子 2所以 a8。 剔除所有因子 2 后剩下的数相乘1131537195。计算这个乘积 mod 8得到 b。这样组合数就可以表示为C(n, m) [n! / (m! * (n-m)!)] [p^{a_n} * b_n] / [p^{a_m} * b_m * p^{a_{n-m}} * b_{n-m}] p^{a_n - a_m - a_{n-m}} * (b_n * inv(b_m) * inv(b_{n-m}))其中inv(x)表示 x 在模p^k下的逆元。由于 b_n, b_m, b_{n-m} 都与 p 互质所以它们在模p^k下存在逆元可以安全计算。因此问题进一步转化为两个子问题计算指数e a_n - a_m - a_{n-m}。计算“净化”后的部分b_n * inv(b_m) * inv(b_{n-m}) mod p^k。如果e k那么p^e mod p^k 0整个组合数模p^k就是 0。 否则最终结果就是(p^e * (b部分的结果)) mod p^k。2.3 第三步递归计算净化阶乘 b如何高效计算F(n) n! 中剔除所有因子 p 后的乘积 mod p^k即上面的 b 这里采用递归思想。观察 n!我们可以把 1 到 n 的数按是否包含因子 p 来分组不包含因子 p 的数1, 2, ..., p-1, p1, ..., 2p-1, ...包含至少一个因子 p 的数p, 2p, 3p, ..., floor(n/p) * p。对于包含因子 p 的数我们可以提取出一个公因子 p剩下的部分就变成了1, 2, 3, ..., floor(n/p)。而这正好是F(floor(n/p))要处理的问题因此递归公式如下F(n) [乘积_{i1, i%p!0}^{p^k} i]^{floor(n / p^k)} * [乘积_{i1}^{n mod p^k} i (if i%p!0)] * F(floor(n/p)) mod p^k公式解读乘积_{i1, i%p!0}^{p^k} i这是一个周期乘积。在模p^k的意义下从 1 到p^k这p^k个数中所有不与 p 互质的数即 p 的倍数我们暂时不考虑它们被归入递归部分。剩下的p^{k-1}*(p-1)个数构成一个“周期块”它们的乘积记为pre。因为模运算的周期性n!中完整的周期块有floor(n / p^k)个所以贡献是pre^{floor(n / p^k)}。乘积_{i1}^{n mod p^k} i (if i%p!0)这是最后一个不完整的周期块中不与 p 互质的数的乘积。F(floor(n/p))这是递归部分处理所有被提取出来的因子 p 之后剩下的那个整数序列的净化阶乘。同时计算指数a即 n! 中因子 p 的个数有一个著名的公式勒让德定理a floor(n/p) floor(n/p^2) floor(n/p^3) ...这个计算可以在递归计算F(n)的过程中顺便完成。2.4 第四步中国剩余定理合并结果经过第二步我们得到了 t 组同余方程x ≡ C(n, m) (mod p1^k1)x ≡ C(n, m) (mod p2^k2)...x ≡ C(n, m) (mod pt^kt)由于p1^k1, p2^k2, ..., pt^kt两两互质因为 pi 是不同质数满足中国剩余定理的应用条件。我们可以使用 CRT 求出唯一的x mod P这个 x 就是C(n, m) mod P。CRT 的合并公式为 设M PMi M / pi^kiti是Mi在模pi^ki下的逆元因为 Mi 与 pi^ki 互质逆元存在。 则解为x ≡ Σ(ai * Mi * ti) (mod M)其中ai C(n, m) mod pi^ki。3. 算法实现细节与关键步骤理解了思路我们来看具体的实现。实现 exLucas 主要需要三个函数快速幂、扩展欧几里得求逆元、以及核心的calc函数。3.1 辅助函数快速幂与扩展欧几里得这些是基础数论工具。// 快速幂计算 base^exp % mod long long qpow(long long base, long long exp, long long mod) { long long res 1; while (exp) { if (exp 1) res res * base % mod; base base * base % mod; exp 1; } return res; } // 扩展欧几里得求解 ax by gcd(a, b)返回 gcd并通过引用返回 x, y long long exgcd(long long a, long long b, long long x, long long y) { if (b 0) { x 1; y 0; return a; } long long d exgcd(b, a % b, y, x); y - a / b * x; return d; } // 求 a 在模 mod 下的逆元要求 gcd(a, mod) 1 long long inv(long long a, long long mod) { long long x, y; exgcd(a, mod, x, y); return (x % mod mod) % mod; // 保证返回正数 }3.2 核心函数计算 F(n) 和指数 a这个函数对应思路中的第二步计算对于特定的(p, pk)pk p^kC(n, m) mod pk的值。 我们实现一个函数lucas_pk(long long n, long long m, long long p, long long pk)。// 计算净化阶乘 F(n) mod pk同时返回指数 a通过引用 long long factorial_p(long long n, long long p, long long pk, long long a) { if (n 0) { a 0; return 1; } // 计算一个完整周期块的乘积 pre long long res 1; for (long long i 1; i pk; i) { if (i % p) { // 只取与 p 互质的数 res res * i % pk; } } res qpow(res, n / pk, pk); // 完整周期块的贡献 // 计算最后一个不完整周期块 for (long long i 1; i n % pk; i) { if (i % p) { res res * i % pk; } } // 递归计算 F(floor(n/p))并累加指数 a long long next_res factorial_p(n / p, p, pk, a); res res * next_res % pk; a n / p; // 勒让德公式累加因子 p 的个数 return res; } // 计算 C(n, m) mod pk long long C_pk(long long n, long long m, long long p, long long pk) { if (m n) return 0; long long a_n 0, a_m 0, a_nm 0; // 计算三个净化阶乘及对应的指数 long long f_n factorial_p(n, p, pk, a_n); long long f_m factorial_p(m, p, pk, a_m); long long f_nm factorial_p(n - m, p, pk, a_nm); long long e a_n - a_m - a_nm; if (e 0) { // 计算净化部分f_n * inv(f_m) * inv(f_nm) mod pk long long res f_n * inv(f_m, pk) % pk * inv(f_nm, pk) % pk; // 乘上 p^e res res * qpow(p, e, pk) % pk; return res; } else { // 如果指数为负说明分母中 p 的因子比分子多整个数不是整数但在模 pk 下我们按 0 处理 // 实际上在组合数定义中 m n且为整数e 不可能为负。这里为了安全可以返回 0 或报错。 return 0; } }3.3 主函数exLucas 与 CRT 合并现在实现最终的exLucas函数它负责质因数分解 P对每个质因子幂调用C_pk最后用 CRT 合并。// 扩展卢卡斯定理主函数计算 C(n, m) mod P long long exLucas(long long n, long long m, long long P) { if (m n) return 0; long long tmp P; vectorpairlong long, long long factors; // 存储 (质数 p, 幂次 pk) // 质因数分解 P for (long long i 2; i * i tmp; i) { if (tmp % i 0) { long long pk 1; while (tmp % i 0) { tmp / i; pk * i; } factors.emplace_back(i, pk); } } if (tmp 1) { factors.emplace_back(tmp, tmp); } // 分别计算 C(n, m) mod pk vectorlong long a(factors.size()), mod(factors.size()); for (size_t i 0; i factors.size(); i) { long long p factors[i].first; long long pk factors[i].second; a[i] C_pk(n, m, p, pk); mod[i] pk; } // 中国剩余定理合并 long long res 0, M P; for (size_t i 0; i factors.size(); i) { long long Mi M / mod[i]; long long ti inv(Mi, mod[i]); // 求 Mi 在模 mod[i] 下的逆元 res (res a[i] * Mi % M * ti % M) % M; } return res; }4. 边界处理、优化与常见问题4.1 边界情况与细节处理m n或m 0组合数定义为 0应在函数入口处判断并返回 0。P 1任何数模 1 都为 0可以直接返回 0。p^k可能很大在计算周期块乘积pre时pk可能达到1e9甚至更大虽然通常题目会控制循环for (i1; ipk; i)会超时。这是实现中的一个关键优化点。优化方法我们不需要真的计算1..pk的完整乘积。注意到我们只关心i % p ! 0的i并且是在模pk下计算。我们可以预处理一个数组pre[pk]但pk太大时不行。实际上由于我们最终要计算pre^{floor(n/pk)}而floor(n/pk)通常很小因为 n 有限我们可以不显式计算出整个pre而是在递归时直接计算(n % pk)!的部分并通过递归关系处理完整周期。上面给出的factorial_p实现是一种简化描述在pk较大时计算res qpow(res, n / pk, pk)中的res即pre会成为瓶颈。更高效的实现是避免直接计算pre而是通过递归将n缩小。更优的factorial_p实现思路long long factorial_p(long long n, long long p, long long pk) { if (n 0) return 1; // 递归计算 F(n) F(n/p) * (周期块乘积)^{n/pk} * (剩余块乘积) long long res factorial_p(n / p, p, pk); // 递归处理提取p后的部分 // 计算周期块和剩余块 long long period_prod 1; // 这里需要高效计算 1..pk 中非p倍数的乘积 mod pk long long rem_prod 1; // 计算 1..(n%pk) 中非p倍数的乘积 mod pk // ... 计算 period_prod 和 rem_prod ... res res * qpow(period_prod, n / pk, pk) % pk; res res * rem_prod % pk; return res; }计算period_prod和rem_prod如果pk不大比如 1e6可以预处理。如果pk很大可能需要更数学的方法或者题目通常会保证pk即 p^k不会太大使得O(pk)的循环可接受。这是算法的时间复杂度瓶颈之一通常认为p^k在1e6量级内是可接受的。求逆元inv(f_m, pk)务必确保f_m与pk互质即不含因子 p这是我们设计算法时保证的。使用扩展欧几里得求逆元是安全的。4.2 时间复杂度分析假设模数 P 分解为Π pi^ki。质因数分解 PO(√P)但 P 通常不大如 1e9 以内。对于每个质因子幂(p, pk)factorial_p递归深度为O(log_p n)。每次递归需要计算周期块乘积若pk不大可O(pk)预处理若直接循环则每次O(pk)和剩余块乘积O(n%pk)。粗略估计处理一个(p, pk)的时间复杂度约为O(pk log_p n)。CRT 合并O(t)t 是质因子个数很小。因此总时间复杂度主要取决于最大的pk值。如果最大的pk在百万级别算法可以在合理时间内运行。4.3 常见问题与调试技巧结果错误为 0检查是否m n导致直接返回 0。检查C_pk函数中当指数e k时返回p^e * ... mod pk如果e kp^e mod pk确实为 0这是正确的。但需确认e计算是否正确factorial_p中的a累加是否正确。检查模数pk是否在运算过程中因为溢出变成了负数或奇怪的值。确保使用long long并在乘法后及时取模。运行超时最大的可能是pk太大导致计算周期块乘积pre的循环O(pk)太慢。确认题目约束pk是否真的可能很大如 1e7。如果很大需要实现上文提到的优化版factorial_p避免直接循环计算pre。递归深度过大log_p n通常很小比如 n1e18, p2深度也就 60左右不是问题。逆元计算失败在C_pk中计算inv(f_m, pk)前理论上f_m与pk互质。如果求逆失败扩展欧几里得返回的 gcd 不为1说明factorial_p函数有 bug没有正确剔除因子 p。CRT 合并结果错误检查inv(Mi, mod[i])是否计算正确即Mi与mod[i]即pk是否互质。由于Mi P / pk而pk是p^kMi包含其他质因子肯定与pk互质所以逆元存在。检查合并公式res (res a * Mi % M * ti % M) % M;中的取模是否正确确保每次加法乘法后都取模M防止溢出。一个实用的调试方法用小的、手算可验证的样例进行测试。 例如计算C(5, 2) mod 6。分解6 2 * 3。计算C(5,2)10。10 mod 2 0。10 mod 3 1。解同余方程组x ≡ 0 (mod 2), x ≡ 1 (mod 3)。解得x ≡ 4 (mod 6)。所以C(5,2) mod 6 4。 用你的 exLucas 程序计算看结果是否为 4。5. 实战应用与扩展思考exLucas 算法虽然原理和实现略显复杂但它解决了模数非质数时组合数计算的根本问题。在算法竞赛中它通常以“模板题”的形式出现要求你实现它来计算C(n, m) mod P。理解其每一步的数学原理对于应对可能的变化至关重要。可能的变化包括模数 P 很大但质因子幂p^k较小这是 exLucas 发挥作用的典型场景。算法效率取决于最大的pk。需要计算多次组合数可以对每个质因子幂pk预处理周期块乘积pre这样在多次调用C_pk时可以避免重复计算pre提升效率。与其它数论定理结合有时题目可能要求计算C(n, m) mod P但 P 不是固定的或者需要处理更复杂的求和式。exLucas 可以作为其中一个模块。个人实现心得理解优于记忆尝试自己推导一遍factorial_p的递归公式比死记硬背代码更有用。理解了“提取因子p”和“周期块”的概念就能应对各种变体。重视边界仔细处理n0, m0, p^k1等情况并确保在C_pk中正确处理指数e为 0 或大于等于k的情况。测试驱动实现后务必用多个小数据包括边界数据和暴力计算对于小的 n, m, P进行对比测试确保正确性。再找一些标准题目如洛谷 P4720的样例进行测试。复杂度心里有数明确算法的瓶颈在于pk的大小。如果题目中pk可能很大比如由大质数组成需要意识到可能存在的性能问题或者考虑题目是否保证了pk不会太大。最后exLucas 是数论工具链中重要的一环。它将组合数计算、阶乘模运算、中国剩余定理、质因数分解等知识串联起来。掌握它不仅能解决一类特定问题更能加深你对模运算、同余方程和递归思想的理解。在遇到模数非质数的组合数问题时你便可以有条不紊地将其分解为质因子幂上的子问题然后逐个击破最后合并得到答案。这种“化整为零”的思路在算法设计中也是一种强大的策略。