线性筛还能这么用？一个‘球盒问题’带你玩转因子个数统计与模数玄机

线性筛的魔法改造：用因子个数统计破解球盒难题

在算法竞赛中，我们常常会遇到一些看似是组合数学问题，实则暗藏数论玄机的题目。今天要探讨的这个"球盒问题"就是典型代表——将n个球放入n个盒子，要求每个盒子里的球与其编号的因子个数相同。初看像是排列组合，实则核心在于高效计算因子个数和巧妙利用模数特性。而解决这个问题的关键钥匙，竟是我们熟悉的线性筛法。

1. 问题本质与数学洞察

让我们先抛开代码，从数学角度理解这个问题。题目要求每个盒子的编号和放入球的编号具有相同数量的因子。这意味着我们需要：

1. 对所有1到n的数字计算其因子个数
2. 将相同因子个数的数字分为一组
3. 每组内部的排列方案数为该组大小的阶乘
4. 最终答案为各组阶乘的乘积

关键观察点：

• 因子个数相同的数字之间可以任意排列
• 当模数较小时（本题为500009），存在临界点使得某个因子个数的出现次数超过模数，此时整个乘积必然为0
• 通过打表发现，当n≥2,250,000时必然存在这样的因子个数

# 简单示例：n=3时的计算过程 因子个数统计： 1: 1 (因子数1) 2: 2 (因子数2) 3: 2 (因子数2) 分组情况： 因子数为1的数字：[1] 因子数为2的数字：[2,3] 方案计算： 1! * 2! = 1 * 2 = 2

2. 线性筛法的深度改造

标准的线性筛法用于高效找出素数，但我们可以改造它，使其在筛素数的同时计算每个数的因子个数。这需要深入理解线性筛的工作原理并进行巧妙扩展。

2.1 标准线性筛回顾

标准线性筛的核心思想是：

• 维护一个素数列表
• 对每个数i，用已知素数p去筛i*p
• 当i能被p整除时停止，确保每个合数只被最小素因子筛掉

// 标准线性筛伪代码 vector<bool> is_prime(n+1, true); vector<int> primes; for (int i = 2; i <= n; i++) { if (is_prime[i]) primes.push_back(i); for (int p : primes) { if (i*p > n) break; is_prime[i*p] = false; if (i % p == 0) break; } }

2.2 扩展为因子个数筛

要计算因子个数，我们需要记录：

1. 每个数的最小素因子及其幂次
2. 利用数论公式：若n = p₁^a₁ * p₂^a₂...pₖ^aₖ，则因子个数d(n) = (a₁+1)(a₂+1)...*(aₖ+1)

改造后的筛法分为三个阶段处理：

阶段一：处理p³ ≤ n的素数

• 对每个素数p，处理其所有幂次p^k
• 对每个p^k的倍数，乘以(k+1)的贡献

阶段二：处理p² ≤ n的剩余素数

• 类似阶段一，但只需处理k=1和k=2的情况

阶段三：处理剩余的大素数

• 这些素数在平方以上，每个只能贡献一个素因子，即乘以2

vector<int> ndivisors(MAXN + 1, 1); // 初始化为1，因为1没有素因子 // 阶段一：处理p³ ≤ MAXN的素数 for (int p = 2; p*p*p <= MAXN; ++p) { if (ndivisors[p] == 1) { // p是素数 int pk = p; for (int k = 1; ; ++k) { int mm = 0; for (int q = pk; q <= MAXN; q += pk) { ++mm; if (mm == p) mm = 0; // 避免重复计算 else ndivisors[q] *= (k + 1); } LL pk2 = LL(pk) * p; if (pk2 > MAXN) break; pk = int(pk2); } } }

3. 模数特性的巧妙利用

本题模数为500009，这个"小模数"带来了重要的优化机会。根据威尔逊定理的扩展思想，当某个因子个数的出现次数超过模数时，其阶乘必然包含模数的倍数，因此整个乘积模500009的结果为0。

关键优化步骤：

1. 预处理阶段设定上限MAXN=2,250,000
2. 当n ≥ MAXN时直接返回0
3. 否则使用预处理的因子个数统计结果

vector<unsigned> cnt(MAXN + 1, 0); vector<unsigned> res(MAXN + 1, 0); res[0] = 1; for (int n = 1; n <= MAXN; ++n) { int d = ndivisors[n]; ++cnt[d]; res[n] = res[n-1] * (cnt[d] % 1024u) % unsigned(MOD); if (cnt[d] >= 1024u) { unsigned a = res[n-1] * (cnt[d] >> 10) % unsigned(MOD); res[n] = (res[n] + (a << 10)) % unsigned(MOD); } if (res[n] == 0) break; // 提前终止 }

4. 性能优化与实现细节

为了处理T≤1e5的大规模查询，我们需要：

1. 预处理所有结果：预先计算1到MAXN的所有答案
2. 高效查询：O(1)时间响应每个查询
3. 内存优化：使用紧凑的数据结构存储中间结果

阶乘计算的优化技巧：

• 将计数器cnt[d]拆分为低10位和高位部分
• 分别计算贡献后再合并，避免大数运算
• 一旦结果变为0即可提前终止后续计算

注意：在实际编码中，位运算技巧可以显著提升模运算效率，但需要确保不会导致数值溢出。

5. 从数学到代码的完整实现

将上述洞察转化为完整解决方案，我们需要：

1. 初始化因子个数数组ndivisors
2. 三阶段筛法计算每个数的因子个数
3. 统计每个因子个数的出现次数cnt[d]
4. 预处理答案数组res[n]
5. 处理查询

完整解决方案框架：

#include <vector> using namespace std; const int MOD = 500009; const int MAXN = 2250000; vector<int> precompute() { vector<int> ndivisors(MAXN + 1, 1); // 三阶段筛法计算因子个数 // ... (省略具体实现) vector<unsigned> cnt(MAXN + 1, 0); vector<unsigned> res(MAXN + 1, 0); res[0] = 1; for (int n = 1; n <= MAXN; ++n) { int d = ndivisors[n]; ++cnt[d]; // ... (阶乘计算优化) if (res[n] == 0) break; } return res; } int main() { auto res = precompute(); // 处理查询 int T, n; scanf("%d", &T); while (T--) { scanf("%d", &n); printf("%d\n", n <= MAXN ? res[n] : 0); } return 0; }

6. 同类问题的扩展思考

这种基于线性筛的扩展技术可以应用于多种数论问题：

1. 因子和计算：类似因子个数，但改为计算(p^(a+1)-1)/(p-1)的乘积
2. 欧拉函数计算：记录每个数的最小素因子来递推计算
3. 莫比乌斯函数：在筛法过程中标记平方因子

性能对比表：

在实际比赛中，这种对经典算法的深度理解和灵活改造能力，往往是解决难题的关键。当遇到模数特殊的问题时，不妨多思考模数本身的性质可能带来的优化机会。