CC++ & Algorithm

线性筛法求积性函数

极难3
语言版本：通用
概述：线性筛（欧拉筛）可以在线性时间内求出1到n所有数的某种积性函数值，是算法竞赛中批量计算积性函数的标准方法。

线性筛法求积性函数——如何快速算出1到100万的函数值？

想象一下，你有一堆数字，要计算每个数的“欧拉函数”（就是比它小且与它互质的数的个数）。比如 φ(12)=4（因为1,5,7,11）。如果直接每个数分解质因数，算100万个数会非常慢（可能电脑要跑几秒钟甚至更久）。而线性筛法就像一条流水线，只用扫一遍就能同时筛出所有质数，并顺便算出每个数的积性函数值，速度极快。

1. 为什么需要线性筛？——从生活例子理解“避免重复劳动”

假设你要给全班40个同学发作业本，每个同学都有一份名单，要算出自己的“好朋友个数”（类似欧拉函数）。最笨的办法是：每个人挨个检查全班同学是不是好朋友，那样每个同学要检查40次，总共1600次。但如果你先按“身高”排好队（相当于筛质数），然后让个子最小的同学先告诉别人“我是你最小的朋友”，这样每个人只用和自己前面的人比一次，总共40次就够了。

核心思想：

  • 每个合数（比如12=2×2×3）只会被它的最小质因子（2）筛掉一次。
  • 在筛的过程中，我们能知道这个合数的质因子信息，从而用积性函数的递推公式快速得到函数值。

2. 关键概念——你需要记住的“三件法宝”

（1）积性函数的“分拆”性质

如果函数 f(n) 满足：当 a 和 b 互质时，f(a×b)=f(a)×f(b)，那么就是积性函数。
比如欧拉函数 φ(n)、除数函数 d(n)（因子个数）、莫比乌斯函数 μ(n) 都是积性函数。

怎么用它？
对于任意数 n，把它写成最小质因子的幂次乘以剩下的互质部分：
n = p^e × m，其中 p 是 n 的最小质因子，且 p 不整除 m。
那么 f(n) = f(p^e) × f(m)。
所以只要知道 f(p^e)（质数幂上的值），就能算 f(n)。

（2）最小质因子的幂次——low 数组的作用

在线性筛中，我们需要记录每个数的最小质因子的幂次（即 p^e，其中 p 是最小质因子，e 是指数）。
用 low[n] 表示这个值。例如：

  • low[12]：12=2²×3，最小质因子2，指数2，所以 low[12]=2²=4。
  • low[18]=2¹×3²→最小质因子2，指数1，所以 low[18]=2。
  • 如果 n 本身就是质数，比如 n=7，则 low[7]=7。

为什么需要 low？
因为在递推时，我们经常遇到两种情况：

  • 如果当前数 i 和质数 p 互质（p不是i的因子），那么 i×p 的最小质因子就是 p，且 low[i×p]=p。
  • 如果 p 是 i 的因子，那么 i×p 的最小质因子还是 p，但指数增加了：low[i×p] = low[i] × p。

（3）两种递推情况——就像分水果

假设你已经知道了 f(i) 的值，现在要算 f(i×p)（p是质数）：

情况1：p 不是 i 的因子（即 i % p != 0）

  • 那么 p 和 i 互质，i×p 可以拆成 i 和 p 两部分，并且这两部分互质。
  • 直接相乘：f(i×p) = f(i) × f(p)。
  • 同时 low[i×p] = p（因为p是新的最小质因子）。

情况2：p 是 i 的因子（即 i % p == 0）

  • 此时 p 已经是 i 的最小质因子，所以 i×p 的最小质因子还是 p，但指数加1。
  • 不能直接相乘，因为 p 与 i 不互质。我们需要拆成最小质因子的幂次和互质部分：
    • 令 i = low[i] × rest，其中 low[i]=p^e，rest 与 p 互质。
    • 那么 i×p = (low[i]×p) × rest = low[i×p] × rest。
  • 根据积性：f(i×p) = f(low[i×p]) × f(rest)。
  • 而 f(rest) = f(i) / f(low[i])（因为积性，且 rest 与 low[i] 互质）。
  • 所以公式可以写成：f(i×p) = f(low[i]×p) × f(i) / f(low[i])。

简化技巧：
对于欧拉函数，有一个特别简单的递推公式：

  • 如果 p 是 i 的最小质因子，那么 φ(i×p) = p × φ(i)。
    这个公式的推导在后面的欧拉函数例子中会详细解释。

3. 以欧拉函数为例——最常用的积性函数

欧拉函数的性质（复习）

  • 若 p 是质数：φ(p) = p-1。
  • 若 n = p^e：φ(p^e) = p^e - p^{e-1} = p^{e-1} × (p-1)。
  • 积性：若 gcd(a,b)=1，则 φ(ab)=φ(a)φ(b)。

为什么当 p 是 i 的因子时 φ(i×p)=p×φ(i)？

设 i = p^e × m（m与p互质），则 i×p = p^{e+1} × m。
由于 φ 在质数幂上有公式：φ(p^{e+1}) = p^{e+1} - p^e = p × (p^e - p^{e-1}) = p × φ(p^e)。
而 φ 是积性函数，且 p^{e+1} 与 m 互质，所以：
φ(i×p) = φ(p^{e+1}) × φ(m) = [p × φ(p^e)] × φ(m) = p × [φ(p^e) × φ(m)] = p × φ(i)。
完美！

用生活例子理解这个递推

想象你有两种颜色的球：红色球（代表质数因子）和蓝色球（代表其他因子）。如果新加进来的红色球和已有的红色球颜色相同（即p是i的因子），那么红色球的数量会多一个，但它们的“贡献”是乘以红色球的个数（因为欧拉函数对相同质因子的处理是乘法）。如果新加的是新颜色（即p与i互质），那么两种颜色独立，直接相乘。

4. 新手容易犯的 5 个错误

  1. 忘记初始化 phi[1] 和 low[1]
    虽然1是特殊的，但很多积性函数定义 f(1)=1。如果不初始化，后面的计算会出错。

  2. 把 break 放错位置
    在情况2（i % p == 0）中，处理完 i×p 后必须立即 break，因为对于更大的质数，i×p 的最小质因子已经不是它们了（否则会重复计算）。如果忘记break，后面会错误地更新同一个合数多次。

  3. 混淆 low[x] == x 的条件
    “low[x] == x” 表示 x 是某个质数的幂（即 x = p^e）。但这个判断的前提是 i % p == 0 且我们已经更新了 low[x] = low[i] * p。只有当 i 本身也是质数幂时，low[x]==x 才成立。如果 i 有其他质因子，那么 low[x] 只是最小质因子的幂次，小于 x。

  4. 运算溢出
    比如 phi[i] * p 可能超过 int 范围，应该使用 long long。

  5. 不理解为什么情况2中直接 phi[x] = phi[i] * p 就是正确的
    很多新手会想：如果 i 不是质数幂（比如 i=12，含有因子2²×3），那么 i×p 的最小质因子幂次 low[12×2] = 2³ = 8，但 12×2=24 还有因子3，所以 low[24]=8，而 phi[24] 不能简单等于 phi[12] * 2。实际上，上一节我们已经证明了只有当 low[x]==x（即 x 是纯质数幂）时才能直接用 phi[i] * p。对于不是纯质数幂的情况，我们需要拆成 low 和剩余部分。代码已经通过 if (low[x]==x) 分支处理了，所以没问题。

5. 完整可运行的代码示例（C++ 和 Python）

C++ 实现（带详细注释）

#include <iostream>
#include <vector>
using namespace std;

const int MAXN = 1000000;
long long phi[MAXN + 1];      // 欧拉函数值
int low[MAXN + 1];            // 最小质因子的幂次，即 p^e
vector<int> primes;           // 质数列表
bool is_composite[MAXN + 1] = {false}; // 标记是否为合数

void linear_sieve_phi(int n) {
    phi[1] = 1;               // φ(1)=1
    low[1] = 1;               // low[1]=1
    for (int i = 2; i <= n; i++) {
        if (!is_composite[i]) {
            primes.push_back(i);
            // i 是质数
            phi[i] = i - 1;   // 质数的欧拉函数 = p-1
            low[i] = i;       // 最小质因子的幂次就是它本身 (p^1)
        }
        // 枚举所有小于等于i的最小质因子的质数
        for (size_t j = 0; j < primes.size() && i * primes[j] <= n; j++) {
            int p = primes[j];
            int x = i * p;          // 我们要计算的合数
            is_composite[x] = true; // 标记为合数

            if (i % p == 0) {
                // 情况2：p是i的最小质因子（因为p <= i的最小质因子）
                low[x] = low[i] * p;   // 更新最小质因子的幂次
                if (low[x] == x) {
                    // 如果low[x]==x，说明x是质数幂，比如 p^{e+1}
                    // 那么可以用递推公式：φ(p^{e+1}) = p * φ(p^e)
                    phi[x] = phi[i] * p;
                } else {
                    // 否则，x不是纯质数幂，需要拆成 low[x] 和 x/low[x] 两部分
                    // 因为 low[x] 与 x/low[x] 互质
                    phi[x] = phi[low[x]] * phi[x / low[x]];
                }
                break; // 注意：必须break，因为更大的质数p'会使得i*p'的最小质因子不是p'
            } else {
                // 情况1：p与i互质，p是x的最小质因子
                low[x] = p;
                phi[x] = phi[i] * phi[p]; // 积性直接相乘
            }
        }
    }
}

int main() {
    int N;
    cout << "请输入 N (建议不超过 1000000): ";
    cin >> N;
    if (N > MAXN) {
        cerr << "N 太大，请修改 MAXN" << endl;
        return 1;
    }
    linear_sieve_phi(N);
    cout << "前10个欧拉函数值:" << endl;
    for (int i = 1; i <= min(N, 10); i++) {
        cout << "φ(" << i << ") = " << phi[i] << endl;
    }
    return 0;
}

Python 实现（带详细注释）

MAXN = 1000000
phi = [0] * (MAXN + 1)          # 欧拉函数值列表
low = [0] * (MAXN + 1)          # 最小质因子幂次列表
primes = []                     # 质数列表
is_composite = [False] * (MAXN + 1)  # 合数标记

def linear_sieve_phi(n):
    phi[1] = 1
    low[1] = 1
    for i in range(2, n + 1):
        if not is_composite[i]:
            primes.append(i)
            phi[i] = i - 1
            low[i] = i           # 质数的最小质因子幂次就是本身
        # 枚举质数
        for p in primes:
            if i * p > n:
                break
            x = i * p
            is_composite[x] = True
            if i % p == 0:
                # 情况2：p是i的最小质因子
                low[x] = low[i] * p
                if low[x] == x:
                    # x是质数幂
                    phi[x] = phi[i] * p
                else:
                    # x不是纯质数幂，拆开计算
                    phi[x] = phi[low[x]] * phi[x // low[x]]
                break  # 重要：必须break
            else:
                # 情况1：p与i互质
                low[x] = p
                phi[x] = phi[i] * phi[p]

# 使用示例
N = int(input("请输入 N (建议不超过 1000000): "))
linear_sieve_phi(N)
print("前10个欧拉函数值:")
for i in range(1, min(N, 10) + 1):
    print(f"φ({i}) = {phi[i]}")

6. 推广到任意积性函数——只要知道“质数幂上的公式”就行

你可能会问：如果我要算的不是欧拉函数，而是其他积性函数，比如除数函数 d(n)（因子个数）或莫比乌斯函数 μ(n)，该怎么办？

核心思想不变：

  1. 对于质数幂 p^e，你需要知道 f(p^e) 的公式（或递推关系）。
  2. 在线性筛的两种情况下，依然可以用同样的分拆方法：
    • 情况1（互质）：f(i×p) = f(i) × f(p)
    • 情况2（p是i的因子）：f(i×p) = f(low[i]×p) × f(i) / f(low[i])
      或者更简单：如果 f 在质数幂上有递推公式，可以直接用 f(i×p) = f(i) * 某个倍数（比如欧拉函数就是乘以 p）。

例子1：除数函数 d(n)

  • d(p^e) = e+1（因为因子有 p^0, p^1, ..., p^e 共 e+1 个）。
  • 在筛法中，情况2（p是i的因子）时，我们需要知道 d(low[i]*p)，这可以通过记录指数来实现。
    具体做法：除了 low 数组，还可以维护一个 cnt 数组表示最小质因子的指数 e。
    那么 d(i×p) 的计算：如果 p 是 i 的因子，则 cnt[i×p] = cnt[i] + 1，然后 d(i×p) = d(i) / (cnt[i]+1) * (cnt[i×p]+1)。
    如果互质，则 d(i×p) = d(i) * d(p) = d(i) * 2。

例子2：莫比乌斯函数 μ(n)

  • μ(1) = 1
  • 对于质数 p：μ(p) = -1
  • 对于质数幂 p^e (e>=2)：μ(p^e) = 0
  • 积性性质：若 a 和 b 互质，则 μ(ab) = μ(a) μ(b)。

在线性筛中：

  • 情况1（互质）：μ(i×p) = μ(i) * (-1) = -μ(i)
  • 情况2（p是i的因子）：此时 i×p 含有至少 p^2，所以 μ(i×p) = 0（因为有平方因子）。

这样，用同样的框架，只需要修改 phi 部分为对应的函数名和规则即可。

7. 小结与练习

线性筛法求积性函数是数论编程的“神兵利器”，它能在 O(N) 的时间内算出所有 1~N 的函数值，是很多算法（比如莫比乌斯反演、杜教筛）的基础。理解它的关键是：

  • 每个合数只被最小质因子访问一次。
  • 利用 low 数组记录最小质因子幂次，实现递推。

练习题目：

  1. 用线性筛法求 1 到 N 的除数函数 d(n) （因子个数）的值。
  2. 用线性筛法求 1 到 N 的莫比乌斯函数 μ(n) 的值。

相关知识点指引：

  • 如果你对欧拉函数还不熟悉，可以回顾“欧拉函数及其性质”。
  • 学习完线性筛后，可以继续学习“莫比乌斯反演”，它大量使用线性筛预处理。
  • 如果想了解更高效的求积性函数方法，可以看看“杜教筛”。

现在，打开你的编程环境，亲手实现一下线性筛求欧拉函数，再试试求除数函数吧！

例题精讲

1单选题

线性筛法（欧拉筛）与普通埃氏筛的最大区别在于：

A线性筛法使用嵌套循环，埃氏筛使用单层循环
B线性筛法保证每个合数只被其最小质因子筛掉一次
C线性筛法只能用于质数判定，不能用于积性函数
D线性筛法的时间复杂度为O(n log log n)
2判断题

利用线性筛法可以批量计算任意积性函数f在1到n上的值，只需预先定义f(1)并给出f(p^k)与f(p^{k-1})的递推关系即可。

3填空题
下面是用线性筛法求欧拉函数φ(n)的代码片段（1≤n≤N），请补全空缺部分。
int phi[N], primes[N], cnt;
bool is_composite[N];
void linear_sieve(int n) {
    phi[1] = 1;
    for (int i = 2; i <= n; ++i) {
        if (!is_composite[i]) {
            primes[cnt++] = i;
            phi[i] = i - 1;
        }
        for (int j = 0; j < cnt && i * primes[j] <= n; ++j) {
            is_composite[i * primes[j]] = true;
            if (i % primes[j] == 0) {
                phi[i * primes[j]] = phi[i] * primes[j];
                ___;
            } else {
                phi[i * primes[j]] = phi[i] * (primes[j] - 1);
            }
        }
    }
}
4单选题

在用线性筛法求莫比乌斯函数μ(n)时，对于质数p，μ(p)的正确初始值是什么？

Aμ(p) = 1
Bμ(p) = -1
Cμ(p) = 0
Dμ(p) = p
5填空题
以下是用线性筛法求约数个数d(n)（n的正因子个数）的代码片段，已知d(1)=1，且若n = i * p（p为质数），记min_pow[i]表示i的最小质因子的幂次。请补全空缺处，使得当i % p == 0时能正确递推。
int d[N], min_pow[N];
void sieve(int n) {
    d[1] = 1;
    for (int i = 2; i <= n; ++i) {
        if (!is_composite[i]) {
            primes[cnt++] = i;
            d[i] = 2;
            min_pow[i] = 1;  // 质数的最小质因子幂次为1
        }
        for (int j = 0; j < cnt && i * primes[j] <= n; ++j) {
            int p = primes[j];
            is_composite[i * p] = true;
            if (i % p == 0) {
                min_pow[i * p] = min_pow[i] + 1;
                d[i * p] = d[i] / (min_pow[i] + 1) * (min_pow[i * p] + 1);
                ___;
            } else {
                min_pow[i * p] = 1;
                d[i * p] = d[i] * d[p];
            }
        }
    }
}