ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

从欧拉函数到线性筛:高效计算互质数的算法实践

从欧拉函数到线性筛:高效计算互质数的算法实践 1. 从一道经典算法题说起Relatives如果你刷过一些在线评测平台OJ的题目或者正在准备算法竞赛大概率见过这样一道题给定一个正整数n要求计算在1到n-1的范围内有多少个整数与n互质即最大公约数gcd(i, n) 1。这个问题的答案就是数论中大名鼎鼎的欧拉函数φ(n)的值。题目名字往往就叫 “Relatives” 或者 “Euler‘s Totient Function”。我第一次遇到它时觉得思路很直接遍历1到n-1逐个计算最大公约数统计互质的个数。写了个循环提交然后——毫不意外地收到了 “Time Limit Exceeded”超时的判决。当n的范围大到10^7甚至更高时这种O(n log n)级别的暴力算法在竞赛的严格时间限制下完全不可行。这道题就像一个分水岭它逼着你不能停留在暴力求解的舒适区必须去掌握更高效的工具素数筛法和欧拉函数的线性筛法。这不仅仅是解决一道题而是打开了一扇门门后是解决一大类数论问题的通用高效方案。今天我们就来彻底拆解 “Relatives” 这道题背后的核心如何将素数筛与欧拉函数的计算高效结合实现O(n)时间复杂度内求出1到n所有数的欧拉函数值。我会从最基础的原理讲起一步步推导到最终优化的代码实现并分享我在实现过程中踩过的坑和调试技巧。2. 欧拉函数定义、性质与暴力解法欧拉函数φ(n)对于正整数n表示小于等于n的正整数中与n互质的数的数目。例如φ(8) 4因为1, 3, 5, 7这四个数与8互质。2.1 核心性质与计算公式理解欧拉函数的几个关键性质是后续利用筛法高效计算的基础积性函数如果两个正整数a和b互质gcd(a, b) 1那么φ(a*b) φ(a) * φ(b)。这是筛法能够递推求解的核心。质数幂次对于一个质数p和正整数k有φ(p^k) p^k - p^(k-1) p^(k-1) * (p - 1)。直观理解在1到p^k中只有p的倍数共p^(k-1)个与p^k不互质。通用公式基于算术基本定理任何大于1的整数n都可以唯一分解为质因数的乘积n p1^k1 * p2^k2 * ... * pm^km。那么φ(n) n * (1 - 1/p1) * (1 - 1/p2) * ... * (1 - 1/pm)。这个公式直接给出了计算方法但需要质因数分解。2.2 暴力解法的局限性与复杂度分析最直观的解法就是根据定义或通用公式来实现。方法一遍历求gcdint phi_bruteforce(int n) { int count 0; for (int i 1; i n; i) { if (gcd(i, n) 1) count; } return count; }这种方法的时间复杂度是O(n log n)因为每次gcd计算是O(log n)。当n10^6时循环百万次在竞赛中已经非常吃力n10^7时必然超时。方法二利用通用公式单次计算优化int phi_formula(int n) { int result n; int temp n; for (int p 2; p * p temp; p) { if (temp % p 0) { while (temp % p 0) temp / p; result - result / p; // 等价于 result * (1 - 1/p) } } if (temp 1) { // 处理剩余的一个大于sqrt(n)的质因子 result - result / temp; } return result; }这种方法基于质因数分解时间复杂度为O(sqrt(n))。对于单次查询φ(n)这已经是很好的算法了。但是如果题目要求我们输出φ(1)到φ(n)的所有值这是很多题目的变体也是理解筛法的典型场景那么对每个数都执行一次O(sqrt(n))的分解总复杂度约为O(n sqrt(n))对于n10^6来说运算量达到10^9级别依然会超时。提示这里就引出了问题的关键。当需要批量、高效地计算大量欧拉函数值时我们必须寻找一种能够利用之前计算结果、避免重复分解质因数的算法。这就是筛法登场的时候。3. 素数筛法埃氏筛与线性筛的演进为了高效计算欧拉函数我们首先需要高效地获取素数信息。素数筛法是基础。3.1 埃拉托斯特尼筛法埃氏筛埃氏筛的思想非常简单从2开始将每个素数的所有倍数标记为合数。const int MAXN 10000000; bool is_prime[MAXN1]; vectorint primes; void eratosthenes_sieve(int n) { fill(is_prime, is_prime n 1, true); is_prime[0] is_prime[1] false; for (int i 2; i n; i) { if (is_prime[i]) { primes.push_back(i); // 从 i*i 开始标记因为 i*k (ki) 已经被更小的素数标记过了 if ((long long)i * i n) { for (int j i * i; j n; j i) { is_prime[j] false; } } } } }复杂度分析埃氏筛的时间复杂度为O(n log log n)空间复杂度O(n)。它已经足够快并且代码易于理解和记忆。对于n10^7它可以在毫秒级完成。埃氏筛的局限性每个合数会被它的所有质因子重复标记。例如30 2*3*5它会被2,3,5各标记一次。虽然不影响正确性但存在微小的效率损失。更重要的是在后续与欧拉函数结合进行线性筛时我们需要确保每个数只被其最小的质因子筛掉一次埃氏筛无法满足这个条件。3.2 欧拉筛线性筛线性筛的核心改进在于确保每个合数只被其最小的质因子筛除一次。这是实现O(n)复杂度的关键。const int MAXN 10000000; bool is_composite[MAXN1]; // 标记是否为合数 vectorint primes; void linear_sieve(int n) { fill(is_composite, is_composite n 1, false); for (int i 2; i n; i) { if (!is_composite[i]) { primes.push_back(i); } // 遍历当前已知的所有素数 for (int j 0; j primes.size(); j) { long long multiple (long long)i * primes[j]; if (multiple n) break; is_composite[multiple] true; // 关键步骤如果 primes[j] 是 i 的质因子则跳出循环 if (i % primes[j] 0) { break; } } } }原理解析外层循环i遍历所有数它既可能是素数也可能是合数。如果i是素数就加入素数表。内层循环用当前数i去乘上素数表里所有不大于i的最小质因子的素数。if (i % primes[j] 0) break;是灵魂语句。它保证了每个合数multiple i * primes[j]只会被其最小的质因子primes[j]筛掉。因为当primes[j]能整除i时i可以表示为primes[j] * k。那么对于下一个素数primes[j1]合数i * primes[j1] primes[j] * k * primes[j1]。这个数的最小质因子是primes[j]它本应该由(k * primes[j1])这个更大的数在将来乘上primes[j]时筛掉。如果现在不break就会用primes[j1]这个非最小质因子提前筛掉它导致重复。线性筛的复杂度严格是O(n)因为它每个合数只被标记一次。它为接下来线性求解欧拉函数提供了完美的框架。4. 线性筛法求欧拉函数原理与推导现在我们将线性筛和欧拉函数的计算融合在一起。目标是在筛出素数的同时利用欧拉函数的积性性质递推求出phi[1]到phi[n]的所有值。我们维护两个数组is_composite[MAXN]用于筛法phi[MAXN]存储欧拉函数值。显然phi[1] 1。我们需要考虑三种情况它们对应着线性筛中内层循环的break语句前后4.1 情况一i是素数如果一个数i是素数那么小于它的所有正整数1到i-1都与它互质。所以φ(i) i - 1。 在代码中当我们发现!is_composite[i]时执行phi[i] i - 1; primes.push_back(i);4.2 情况二i是合数且primes[j]不能整除i设multiple i * primes[j]且i % primes[j] ! 0。 这意味着primes[j]是multiple的一个新的、最小的质因子因为primes[j]比i的任何质因子都小不这里需要理解primes[j]是素数表中的素数按顺序遍历。由于i % primes[j] ! 0所以primes[j]不是i的因子。又因为primes[j]小于等于i的最小质因子不一定但关键点是primes[j]和i互质。 既然primes[j]与i互质根据欧拉函数的积性我们有φ(multiple) φ(i) * φ(primes[j]) φ(i) * (primes[j] - 1)。4.3 情况三i是合数且primes[j]能整除i设multiple i * primes[j]且i % primes[j] 0。 这意味着primes[j]是i的一个质因子。此时multiple与i的质因子集合完全相同只是primes[j]这个因子的指数增加了1。 设i primes[j]^k * m其中m不能被primes[j]整除。 那么multiple primes[j]^(k1) * m。 根据欧拉函数公式φ(i) φ(primes[j]^k) * φ(m) (primes[j]^k - primes[j]^(k-1)) * φ(m)φ(multiple) φ(primes[j]^(k1)) * φ(m) (primes[j]^(k1) - primes[j]^k) * φ(m)观察可得φ(multiple) φ(i) * primes[j]。 因为(primes[j]^(k1) - primes[j]^k) primes[j] * (primes[j]^k - primes[j]^(k-1))。更直观的理解当primes[j]是i的质因子时multiple的所有质因子和i完全一样。在计算φ(n)n * Π(1-1/p)时multiple比i多乘了一个primes[j]并且在连乘积部分Π(1-1/p)中因子(1-1/primes[j])是相同的。所以φ(multiple) φ(i) * primes[j]。这三种情况覆盖了线性筛过程中生成的所有合数multiple并且每个multiple只被其最小质因子primes[j]访问一次因此我们可以在O(n)时间内计算出所有phi值。5. 完整代码实现与逐行解析结合上述原理我们可以写出同时进行线性筛和欧拉函数计算的完整代码。这里以求解 “Relatives” 问题计算φ(n)以及其扩展问题输出前n个欧拉函数值为例。#include iostream #include vector #include cstring using namespace std; const int MAXN 10000000; // 根据题目要求调整 // 全局数组 bool is_composite[MAXN 1]; // 标记合数 int phi[MAXN 1]; // 存储欧拉函数值 vectorint primes; // 存储素数 // 线性筛法计算 1~n 的欧拉函数 void euler_sieve(int n) { // 初始化 fill(is_composite, is_composite n 1, false); phi[1] 1; // 定义 for (int i 2; i n; i) { // 情况一i是素数 if (!is_composite[i]) { primes.push_back(i); phi[i] i - 1; // 素数的欧拉函数值为 i-1 } // 用当前数 i 和素数表里的素数相乘筛去合数 for (int j 0; j primes.size(); j) { long long multiple (long long)i * primes[j]; if (multiple n) break; // 超过范围退出内层循环 is_composite[multiple] true; // 标记 multiple 为合数 // 关键判断primes[j] 是否是 i 的最小质因子 if (i % primes[j] 0) { // 情况三primes[j] 能整除 i phi[multiple] phi[i] * primes[j]; break; // 保证每个合数只被最小质因子筛一次 } else { // 情况二primes[j] 与 i 互质 phi[multiple] phi[i] * (primes[j] - 1); } } } } int main() { int n; // 假设题目要求计算 φ(n) // cin n; // euler_sieve(n); // 如果只求一个筛到n即可 // cout phi[n] endl; // 更常见的场景预处理 1~MAXN 的所有 phi 值 euler_sieve(MAXN); // 示例输出前 20 个欧拉函数值 for (int i 1; i 20; i) { cout phi( i ) phi[i] endl; } // 示例回答多次查询 int query; while (cin query query 0) { cout phi[query] endl; } return 0; }逐行解析与注意事项数组大小MAXN是预定义的最大范围数组大小要设为MAXN1因为我们的索引是从1到n。这是新手常犯的“差一错误”。数据类型multiple i * primes[j]这个计算可能溢出。i和primes[j]都是int但乘积可能超过int范围当n较大时。因此必须使用long long进行强制转换和比较。这是算法题中极其常见的坑点。初始化phi[1] 1是定义。记得在循环开始前设置好。内层循环条件multiple n时break这是必要的边界检查避免数组越界。标记合数is_composite[multiple] true;这行代码放在判断i % primes[j] 0之前还是之后答案是之前。因为无论属于情况二还是情况三multiple都是合数都需要被标记。先标记合数再根据互质关系计算phi[multiple]逻辑是清晰的。break的位置break语句只在i % primes[j] 0时执行。这保证了线性筛的正确性也是递推计算phi的分支点。6. 性能对比、应用场景与内存优化6.1 性能对比我们来对比一下三种方法在n10^7量级下的表现理论估算暴力遍历求gcdO(n log n) ≈ 10^7 * 24 ≈ 2.4e8次运算严重超时。单个数公式法求所有值O(n sqrt(n)) ≈ 10^7 * 3162 ≈ 3e10次运算不可能完成。线性筛法O(n) ≈ 10^7次运算在普通PC上可在0.1-0.3秒内完成完全满足竞赛要求。6.2 典型应用场景掌握了线性筛求欧拉函数你就能解决一大类问题原题“Relatives”直接输出phi[n]。区间互质个数统计有些题目会问[L, R]区间内与某个数M互质的数的个数。可以利用欧拉函数和容斥原理求解。欧拉函数前缀和求Σφ(i) (i1 to n)。预处理出所有phi[i]后前缀和可以在O(n)内完成之后每次查询O(1)。这在莫比乌斯反演等问题中非常常见。模运算相关欧拉定理a^φ(m) ≡ 1 (mod m)当gcd(a, m)1是RSA加密等算法的基础。快速计算φ(m)是关键步骤。6.3 内存优化技巧我们的代码使用了三个数组is_composite,phi,primes。对于n10^7bool is_composite[MAXN1]约 10 MB。int phi[MAXN1]约 40 MB。vectorint primes素数个数约为n / ln(n) ≈ 6e5存储约 2.4 MB。 总计约 52 MB这在大多数OJ的256MB内存限制下是可行的。如果内存非常紧张可以考虑以下优化省略is_composite数组我们可以在标记合数时利用phi数组的初始值来判断。例如初始化phi[i] i。在筛法过程中如果一个数phi[multiple]被计算过不再是初始值multiple那它就一定是合数。但这样会稍微增加判断逻辑的复杂度。使用bitset将is_composite改为bitsetMAXN1可以将标记数组的内存占用减少到约 1.25 MB。分块筛法当n极大如10^9时无法一次性分配O(n)内存。可以采用分段筛法每次只处理一段区间内存占用取决于区间大小但时间复杂度会略有增加。注意在竞赛中除非特别必要优先选择代码清晰、易于调试的实现。内存优化通常是在明确超限后再考虑。7. 调试技巧与常见错误排查实现线性筛求欧拉函数时以下几个错误点非常高发整数溢出前面提到的multiple i * primes[j]溢出。务必使用long long。// 错误 if (i * primes[j] n) break; // 正确 if ((long long)i * primes[j] n) break;数组越界确保循环边界是i n数组声明大小足够。在本地测试时可以使用valgrind或AddressSanitizer等工具检查内存错误。初始化遗漏忘记初始化phi[1] 1或is_composite数组导致结果错误。可以在函数开头显式地初始化所有相关数组。逻辑分支错误混淆情况二和情况三的计算公式。一个简单的记忆方法是如果primes[j]能整除i说明primes[j]已经在i的质因子集合里了那么φ(i*primes[j])就只是把i放大了primes[j]倍所以乘primes[j]。如果不能整除说明primes[j]是一个新的质因子引入了一个新的(1-1/p)项所以乘(primes[j]-1)。验证方法小数据暴力对拍写一个O(n sqrt(n))的单点求φ函数与筛法结果对比n从1到1000的所有值。这是最有效的调试手段。利用已知数列欧拉函数前几项为1, 1, 2, 2, 4, 2, 6, 4, 6, 4, 10, 4, 12, 6, 8, 8, 16, 6, 18, 8, ... 可以手动核对。性质检验对于质数p检查phi[p] p-1。对于n2检查phi[n]是否为偶数这是一个性质。计算Σφ(d) (d|n)应等于n。性能测试在本地用n1e7测试运行时间确保在合理范围内1秒左右。如果时间过长检查是否有不必要的递归、多层循环或低效操作。8. 举一反三线性筛法求解其他积性函数线性筛法的威力远不止于求欧拉函数。它实际上是一个框架适用于任何积性函数。只要你能定义好这个函数在素数p、素数幂p^k以及互质情况下的递推关系就能套用这个模板。以莫比乌斯函数μ(n)为例μ(1) 1若n有平方因子则μ(n) 0若n是k个不同质数的乘积则μ(n) (-1)^k在线性筛中当i为素数mu[i] -1。当i % primes[j] ! 0mu[multiple] -mu[i]因为增加了一个新的质因子。当i % primes[j] 0mu[multiple] 0因为出现了重复质因子即平方因子。以约数个数函数d(n)为例设n p1^k1 * p2^k2 * ...则d(n) (k11)*(k21)*...我们需要额外维护一个数组cnt[i]记录i的最小质因子的指数。 在线性筛中当i为素数d[i] 2,cnt[i] 1。当i % primes[j] ! 0multiple增加了一个新的质因子primes[j]指数为1。所以d[multiple] d[i] * 2,cnt[multiple] 1。当i % primes[j] 0multiple的最小质因子primes[j]的指数增加了1。所以cnt[multiple] cnt[i] 1并且d[multiple] d[i] / cnt[multiple] * (cnt[multiple] 1)。通过这两个例子你可以看到线性筛框架的强大与统一。掌握它你就掌握了一大批数论函数的高效批量计算方法。下次遇到需要预处理约数和、因子和等问题时不妨先想想能否用线性筛解决。
返回列表