从质因数分解到约数计算：数论核心算法与工程实践

1. 从一道经典面试题说起约数到底有多重要如果你刷过一些编程面试题或者参加过信息学竞赛大概率遇到过这样一类问题给定一个正整数 N求它的所有约数之和或者求 1 到 N 中所有数的约数个数之和。乍一看这似乎只是一个简单的数学计算暴力枚举不就行了但当你真正动手去写面对 N 可能高达 10^12 甚至更大的情况时你才会意识到问题的复杂性。暴力枚举的时间复杂度是 O(N)这在数据规模稍大时是完全不可接受的。这时数论中关于约数的系统性知识就从“课本上的公式”变成了解决实际问题的“救命稻草”。约数也叫因数是数论中最基础、最核心的概念之一。它描述了一个整数能被哪些整数整除。围绕约数有三个核心的计算问题约数个数、约数和以及最大公约数。这三个问题不仅是数学理论上的优美存在更是算法竞赛、密码学如 RSA 算法、数据压缩、乃至日常工程中优化资源分配如分页、任务调度的底层逻辑。很多人学习数论觉得抽象正是因为缺少了从“公式”到“代码”从“理论”到“应用”的桥梁。本文的目的就是拆掉这座桥的脚手架让你能一步步走过去真正掌握如何高效、优雅地处理约数问题。我会结合大量代码示例和实战场景让你不仅记住公式更理解其背后的“为什么”和“怎么做”。2. 基石算术基本定理与质因数分解在深入探讨约数的各种计算之前我们必须先夯实一个最基础、最重要的数论定理——算术基本定理。它是我们所有高效算法的理论源头。算术基本定理指出任何一个大于1的自然数 N都可以唯一地分解成有限个质数的乘积。这里的“唯一”是指如果不考虑质因数的排列顺序那么这种分解方式是唯一的。用公式表示就是N p1^a1 * p2^a2 * ... * pk^ak其中p1, p2, ..., pk是互不相同的质数a1, a2, ..., ak是正整数表示对应质因数的指数。为什么这个定理如此关键因为它将我们对一个庞大整数 N 的研究转化为了对若干个质数及其指数(pi, ai)的研究。质数是数的“原子”而算术基本定理告诉我们任何合数都是由这些“原子”以特定方式组合而成的。这极大地简化了问题。2.1 质因数分解的实战实现理论很美但如何用代码实现呢最常用的方法是试除法。核心思路从最小的质数2开始逐个尝试去除 N。如果能整除就记录这个质因数并用 N 不断除以它直到不能整除为止然后增加试除的质数。基础版本代码Pythondef prime_factors(n): factors [] i 2 # 处理因子2可以单独提出来加速因为2是唯一的偶质数 while i * i n: if n % i 0: cnt 0 while n % i 0: n // i cnt 1 factors.append((i, cnt)) # 记录质因数i和它的指数cnt i 1 if i 2 else 2 # 除了2之后只检查奇数这是一个小优化 # 循环结束后如果n还大于1那么它本身就是一个质数 if n 1: factors.append((n, 1)) return factors # 示例分解 360 print(prime_factors(360)) # 输出[(2, 3), (3, 2), (5, 1)]即 360 2^3 * 3^2 * 5^1代码解读与避坑点循环条件while i * i n这是关键优化。因为如果 N 有一个大于sqrt(N)的质因子那么它必然对应一个小于sqrt(N)的因子否则乘积就超过 N 了。所以当i大于sqrt(N)时N 要么是1要么是一个大于sqrt(N)的质数。这个优化将时间复杂度从 O(N) 降到了 O(√N)。单独处理2和奇数i 1 if i 2 else 2这行代码是一个有效的微优化。因为除了2以外所有质数都是奇数所以检查完2后我们可以只检查奇数减少一半的循环次数。最后的if n 1这个判断至关重要。当循环结束时如果n不等于1那么它一定是那个大于原数平方根的质因子。例如分解n13质数循环不会进入因为 2*2 13此时n仍为13需要单独加入结果。注意试除法对于单个大数的分解是有效的但对于极大整数如几百位则需要更复杂的算法如 Pollard Rho。在算法竞赛和大多数工程场景中试除法及其优化变种已经足够。有了质因数分解的结果[(p1, a1), (p2, a2), ...]我们就可以像搭积木一样推导出关于 N 的所有约数性质。3. 核心推导约数个数的公式与生成约数个数顾名思义就是整数 N 的正约数的总个数。我们记作d(N)或σ0(N)。3.1 公式推导为什么是乘积假设N p1^a1 * p2^a2 * ... * pk^ak。 N 的任何一个约数d必然是由这些质因数p1, p2, ..., pk组合而成的并且每个质因数pi的指数不能超过它在 N 中的指数ai但可以取 0即不包含该质因子。因此对于质因数p1在约数d中它的指数有(a1 1)种选择0, 1, 2, ..., a1。 同理对于p2有(a2 1)种选择。 ... 对于pk有(ak 1)种选择。由于每个质因数的选择是独立的根据乘法原理所有可能的约数总数就是各质因数指数加1的连乘积。约数个数公式d(N) (a1 1) * (a2 1) * ... * (ak 1)示例N 360 2^3 * 3^2 * 5^1d(360) (31) * (21) * (11) 4 * 3 * 2 24所以360 有 24 个正约数。3.2 实战如何枚举所有约数知道个数很重要但有时我们需要列出所有约数本身例如找一个数的所有真约数来判断它是否为完全数。暴力枚举从1到N显然太低效。我们可以利用质因数分解的结果通过递归或迭代生成所有组合。方法DFS深度优先搜索生成思路是遍历每个质因数决定其在当前约数中的指数从0到ai并累乘。def get_factors_from_pf(prime_factors): 根据质因数分解结果生成所有约数 factors [1] for p, a in prime_factors: current_factors [] # 对于当前质因数p考虑其指数从0到a for i in range(a 1): # 将现有所有约数乘以 p^i得到新的约数集合 for f in factors: current_factors.append(f * (p ** i)) factors current_factors # 更新约数列表 factors.sort() # 排序可选 return factors # 使用上一节的分解结果 pf [(2, 3), (3, 2), (5, 1)] all_factors get_factors_from_pf(pf) print(f“约数个数{len(all_factors)}”) # 输出24 print(f“所有约数{all_factors}”) # 输出[1, 2, 3, 4, 5, 6, 8, 9, 10, 12, 15, 18, 20, 24, 30, 36, 40, 45, 60, 72, 90, 120, 180, 360]复杂度分析生成的约数个数就是d(N)。这个算法的时间复杂度是 O(d(N))空间复杂度也是 O(d(N))。这比 O(N) 的暴力枚举要好得多因为d(N)通常远小于N。对于N360d(N)24而N360优势明显。实操心得在需要频繁获取一个数的约数时先进行质因数分解并缓存结果再动态生成约数列表是更高效的策略。特别是在一些数论动态规划问题中这个技巧非常有用。4. 进阶计算约数和的公式与应用约数和是指 N 的所有正约数相加之和。我们记作σ(N)。4.1 公式推导等比数列求和与乘法原理同样从质因数分解N p1^a1 * p2^a2 * ... * pk^ak出发。 考虑 N 的任意一个约数d它包含质因数p1的指数可以是 0, 1, ..., a1。那么在所有约数中p1这部分贡献的和是多少这相当于求一个等比数列的和p1^0 p1^1 p1^2 ... p1^a1。 根据等比数列求和公式这个和等于(p1^(a11) - 1) / (p1 - 1)。由于每个约数都是各质因数幂次的乘积并且不同质因数的贡献是独立的根据乘法原理所有约数的总和就等于每个质因数对应的等比数列和的乘积。约数和公式σ(N) [ (p1^(a11) - 1) / (p1 - 1) ] * [ (p2^(a21) - 1) / (p2 - 1) ] * ... * [ (pk^(ak1) - 1) / (pk - 1) ]示例计算σ(360)。360 2^3 * 3^2 * 5^1对于p12, a13:(2^(31)-1)/(2-1) (16-1)/1 15对于p23, a22:(3^(21)-1)/(3-1) (27-1)/2 13对于p35, a31:(5^(11)-1)/(5-1) (25-1)/4 6σ(360) 15 * 13 * 6 1170我们可以验证一下上一节列出的360的所有约数之和1234568910121518202430364045607290120180360 1170。完全正确。4.2 代码实现与模运算处理在编程中直接套用公式计算即可但需要注意两点大数幂运算当ai较大时计算pi^(ai1)可能溢出。需要使用快速幂算法并在必要时进行取模运算。除法与模逆元如果题目要求结果对一个质数MOD取模这在竞赛中很常见那么公式中的除法/(pi-1)不能直接进行。因为模运算下没有直接的除法。此时需要计算(pi-1)在模MOD下的乘法逆元将除法转化为乘法。代码示例含取模def fast_pow(base, exp, modNone): 快速幂算法支持取模 result 1 while exp 0: if exp 1: # 如果指数是奇数 result result * base if mod is not None: result % mod base base * base if mod is not None: base % mod exp 1 # 指数右移一位除以2 return result def mod_inv(x, mod): 求x在模mod下的乘法逆元mod是质数时可用费马小定理 return fast_pow(x, mod-2, mod) def sum_of_factors(n, modNone): 计算约数和可选取模 pf prime_factors(n) # 复用之前的质因数分解函数 ans 1 for p, a in pf: # 计算 (p^(a1) - 1) numerator fast_pow(p, a 1, mod) - 1 if mod is not None: numerator % mod # 计算 (p - 1) denominator p - 1 if mod is not None: # 模运算下需要乘以分母的逆元 inv_denom mod_inv(denominator, mod) term (numerator * inv_denom) % mod else: term numerator // denominator # 整数除法 ans ans * term if mod is not None: ans % mod return ans print(f“约数和 σ(360) {sum_of_factors(360)}”) # 输出1170重要提示上述mod_inv函数使用费马小定理仅当模数mod为质数时才成立。如果mod不是质数则需要使用扩展欧几里得算法求逆元。这是数论和组合数学中一个非常常见的坑点。5. 最大公约数算法、性质与实战最大公约数指两个或多个整数共有约数中最大的一个。记作gcd(a, b)。它是数论中应用最广泛的概念之一。5.1 欧几里得算法辗转相除法为什么它这么快这是计算最大公约数最著名、最高效的算法。其原理基于一个核心定理gcd(a, b) gcd(b, a % b)其中a % b表示 a 除以 b 的余数这个定理的直观理解假设d是a和b的公约数那么d一定能整除a和b。那么a可以表示为a k*b r其中r a % b。因为d能整除a和b所以它也一定能整除r a - k*b。因此d也是b和r的公约数。反之亦然。所以a, b的公约数集合与b, r的公约数集合完全相同自然它们的最大公约数也相同。递归实现def gcd_euclid(a, b): if b 0: return a return gcd_euclid(b, a % b)迭代实现更常用避免递归栈溢出def gcd(a, b): while b ! 0: a, b b, a % b return a复杂度分析欧几里得算法的时间复杂度是O(log(min(a, b)))。这是一个非常高效的算法即使对于非常大的整数如10^18级别也能瞬间得出结果。其复杂度与数字的位数成线性关系而不是与数字本身的大小成线性关系。5.2 扩展欧几里得算法解不定方程 ax by gcd(a, b)这是欧几里得算法的扩展它不仅能求出gcd(a, b)还能找到一对整数(x, y)使得a*x b*y gcd(a, b)成立。这个方程称为贝祖等式。算法原理递归理解 在辗转相除法的最后一步我们有b0,gcd(a,0)a。此时等式a*x 0*y a显然有一组解(x1, y0)。 然后我们“回溯”这个过程。假设我们已经知道了b和a%b对应的系数即b * x1 (a % b) * y1 gcd(b, a%b) gcd(a, b)我们知道a % b a - (a // b) * b代入上式b * x1 (a - (a//b)*b) * y1 gcd(a, b)整理得a * y1 b * (x1 - (a//b) * y1) gcd(a, b)所以对于(a, b)对应的系数(x, y)就是(y1, x1 - (a//b)*y1)。递归实现def ext_gcd(a, b): 返回 (gcd, x, y) 使得 a*x b*y gcd if b 0: return a, 1, 0 gcd_val, x1, y1 ext_gcd(b, a % b) x y1 y x1 - (a // b) * y1 return gcd_val, x, y # 示例解 35x 15y gcd(35,15) g, x, y ext_gcd(35, 15) print(f“gcd(35,15){g}, 系数x{x}, y{y}”) # 输出gcd(35,15)5, x1, y-2 # 验证35*1 15*(-2) 35 - 30 5正确。迭代实现更高效def ext_gcd_iter(a, b): x0, x1, y0, y1 1, 0, 0, 1 # 初始化系数矩阵 while b ! 0: q a // b a, b b, a % b # 同时更新系数 x0, x1 x1, x0 - q * x1 y0, y1 y1, y0 - q * y1 return a, x0, y0 # a 现在是 gcd扩展欧几里得算法的应用场景求解模线性方程求a*x ≡ c (mod m)的解。该方程有解的充要条件是gcd(a, m) | c。我们可以先解a*x m*y gcd(a, m)得到一组特解(x0, y0)然后进行变换得到原方程的解。求乘法逆元在模运算中若gcd(a, m) 1则a在模m下有逆元。这个逆元就是扩展欧几里得算法解a*x m*y 1得到的x模m后的正数。这比费马小定理更通用不要求m是质数。中国剩余定理的求解过程中也需要用到扩展欧几里得算法来合并同余方程组。5.3 多个数的最大公约数与最小公倍数多个数的最大公约数gcd(a, b, c) gcd(gcd(a, b), c)。可以依次计算。最小公倍数两个数a, b的最小公倍数lcm(a, b)与最大公约数有一个非常重要的关系lcm(a, b) a * b / gcd(a, b)这个公式非常高效因为我们有快速的gcd算法。 同样多个数的最小公倍数lcm(a, b, c) lcm(lcm(a, b), c)。代码示例def lcm(a, b): return a // gcd(a, b) * b # 先除后乘避免中间结果溢出 def gcd_list(nums): result nums[0] for num in nums[1:]: result gcd(result, num) return result def lcm_list(nums): result nums[0] for num in nums[1:]: result lcm(result, num) return result print(f“gcd(12, 18, 24) {gcd_list([12, 18, 24])}”) # 输出6 print(f“lcm(12, 18, 24) {lcm_list([12, 18, 24])}”) # 输出72避坑提醒计算lcm时务必使用a // gcd(a, b) * b的顺序而不是a * b // gcd(a, b)。虽然数学上等价但在编程中a * b可能会先溢出尤其是使用C/Java等语言的int类型时先进行除法可以减小中间值的大小。6. 综合实战从理论到解题的跨越掌握了基本公式和算法我们来看几个综合性的问题感受一下如何将这些知识串联起来解决实际问题。6.1 问题一求 1! 到 N! 每个数的约数个数这是一个经典的递推问题。直接对每个阶乘数进行质因数分解再计算d(N!)效率太低。我们需要利用阶乘的性质。核心观察N! 1 * 2 * 3 * ... * N。N!的质因数分解等于1, 2, ..., N每个数质因数分解的“总和”。更具体地说对于某个质数p它在N!中的指数a等于a floor(N/p) floor(N/p^2) floor(N/p^3) ...floor表示向下取整 这个公式的含义是1~N中有floor(N/p)个数至少包含一个因子p有floor(N/p^2)个数至少包含两个因子p这些数在第一次计数时已被算过一次所以这里再补一次... 以此类推。算法步骤用筛法如埃拉托斯特尼筛法预处理出所有不超过 N 的质数。对于每个质数p用上述公式计算它在N!中的指数a。根据约数个数公式d(N!) Π (a_i 1)。代码实现def count_divisors_of_factorial(n): # 1. 筛出所有质数 is_prime [True] * (n 1) primes [] for i in range(2, n 1): if is_prime[i]: primes.append(i) for j in range(i * i, n 1, i): is_prime[j] False # 2. 对每个质数计算在 n! 中的指数 result 1 MOD 10**9 7 # 通常结果很大需要取模 for p in primes: exp 0 power p while power n: exp n // power power * p result (result * (exp 1)) % MOD return result print(f“10! 的约数个数模 1e97: {count_divisors_of_factorial(10)}”)6.2 问题二求 1 到 N 中所有数的约数个数之和即求S(N) d(1) d(2) ... d(N)。N 可以很大比如 10^7。暴力法对每个i计算d(i)复杂度 O(N√N)不可接受。数论分块/贡献法转换视角。不考虑每个数有多少个约数而是考虑每个数k是多少个数的约数。 对于k在1~N中k的倍数有floor(N/k)个。也就是说k作为约数总共出现了floor(N/k)次。 因此S(N) Σ_{k1}^{N} floor(N/k)。直接计算这个求和式可以计算但复杂度仍是 O(N)。进一步优化观察floor(N/k)的值当k较大时这个值变化很慢。实际上floor(N/k)的取值只有大约2√N种。我们可以找到使floor(N/k)相同的k的区间[l, r]然后批量计算贡献。 对于给定的kfloor(N/k) v那么最大的k记为r满足floor(N/r) v的是r floor(N / v)。下一个区间从l r1开始。优化后算法def sum_of_divisor_counts(n): total 0 k 1 while k n: v n // k # 找到使得 n // r v 的最大 r r n // v # 区间 [k, r] 内每个数作为约数出现的次数都是 v # 这个区间对总和的贡献是区间长度 * v total (r - k 1) * v k r 1 # 跳到下一个区间起点 return total print(f“1到100的约数个数之和{sum_of_divisor_counts(100)}”) # 可以验证暴力计算 d(1)...d(100) 也是这个结果。这个算法的时间复杂度是O(√N)可以处理非常大的 N如 10^12。6.3 问题三判断两个数是否互质并求模逆元互质如果gcd(a, b) 1则a和b互质。模逆元如果gcd(a, m) 1则存在整数x使得a*x ≡ 1 (mod m)。x称为a模m的逆元记作a^{-1} (mod m)。我们可以用扩展欧几里得算法一次性解决这两个问题。def mod_inv_extgcd(a, m): 使用扩展欧几里得算法求逆元并判断是否互质 g, x, y ext_gcd_iter(a, m) if g ! 1: # gcd(a, m) ! 1说明a和m不互质逆元不存在 return None, False else: # 逆元存在将x调整到[0, m)范围内 inv x % m return inv, True a, m 17, 3120 inv, is_coprime mod_inv_extgcd(a, m) if is_coprime: print(f“{a} 和 {m} 互质{a} 模 {m} 的逆元是 {inv}”) print(f“验证{a} * {inv} % {m} {a * inv % m}”) else: print(f“{a} 和 {m} 不互质逆元不存在”)7. 性能优化与边界条件处理在实际编码和解题中除了掌握核心算法处理边界条件和进行优化同样重要。7.1 质因数分解的优化我们之前给出的试除法已经包含了i*i n和跳过偶数的优化。对于需要频繁分解很多数字的场景可以进一步优化预处理质数表先用线性筛法欧拉筛生成一个范围内的所有质数存放在列表中。分解时只用这些质数去试除避免了合数的无效判断。def linear_sieve(limit): is_prime [True] * (limit 1) primes [] for i in range(2, limit 1): if is_prime[i]: primes.append(i) for p in primes: if i * p limit: break is_prime[i * p] False if i % p 0: # 保证每个合数只被最小的质因子筛掉 break return primes # 假设我们只处理 n 10^6 的数 prime_list linear_sieve(10**6) def prime_factors_fast(n, primes): factors [] for p in primes: if p * p n: break if n % p 0: cnt 0 while n % p 0: n // p cnt 1 factors.append((p, cnt)) if n 1: factors.append((n, 1)) return factors当需要分解的数都在预处理范围内时这个方法非常快。记忆化缓存如果同一个数会被多次分解可以使用字典缓存结果。7.2 处理大数与溢出在计算约数和σ(N)时尤其是涉及取模运算要特别注意使用快速幂计算p^(a1)。模运算下的除法必须转换为乘以逆元。在计算lcm(a, b)时使用a // gcd(a, b) * b防止溢出。在 Python 中整数不会溢出但在 C/Java 中中间结果a * b可能超出long long范围需要特别注意或者使用__int128如果支持。7.3 特殊情况的考虑N 0 或 10 的约数定义在数学中不统一通常认为所有非零整数都是0的约数但这样约数个数无限在编程题中通常不会出现。1 的约数只有它自身质因数分解为空列表约数个数d(1)1约数和σ(1)1。在实现函数时要处理好边界。负整数通常数论讨论的是正整数。如果输入可能是负数可以先取绝对值处理。重复计算在循环中多次调用gcd,prime_factors等函数可能是性能瓶颈。考虑在外层预处理或缓存。8. 总结与个人心得走完这一趟从定义、公式推导到算法实现、实战应用的旅程你应该能感受到数论中的约数问题并不是一堆枯燥的公式而是一套强大的、有内在联系的工具箱。算术基本定理是这个工具箱的基石它将复杂的大数问题分解为质因数的组合问题。约数个数和约数和的公式是两块核心的模板让你能绕过暴力枚举直接通过质因数分解的结果进行高效计算。欧几里得算法及其扩展则是处理公约数、模运算和不定方程的瑞士军刀其简洁与高效令人赞叹。我个人在多年的算法竞赛和工程开发中处理约数相关问题的核心体会是“先分解后组合”。无论问题看起来多复杂只要它和整数的整除性有关第一步就应该想到质因数分解。分解之后问题的复杂度就从关于N降到了关于N的质因子个数和指数上这通常是一个巨大的优化。另一个深刻的教训是关于模运算中的逆元。早期我曾多次因为忘记处理模下的除法而得到错误结果。务必牢记(a / b) % MOD ≠ (a % MOD) / (b % MOD)。正确的做法是求b模MOD的逆元然后计算(a % MOD) * inv(b) % MOD。最后对于算法竞赛选手我强烈建议将以下内容作为模板熟练背诵并理解试除法质因数分解含优化。欧几里得算法求最大公约数。扩展欧几里得算法求逆元和解线性方程。快速幂算法。基于质因数分解计算约数个数、约数和的函数。数论分块求∑ floor(N/i)的代码。这些是解决大多数约数相关问题的“原子操作”。掌握了它们你就能像搭积木一样组合出解决更复杂问题的方案。数论的美在于其严谨的逻辑和广泛的应用希望这篇长文能帮你打开这扇门看到门后那片既深邃又实用的风景。下次当你再看到“求约数之和”这样的问题时希望你的第一反应不再是循环枚举而是会心一笑然后优雅地写下质因数分解的代码。