简单的数论函数及其求和方法

目录

• 数论函数
  • 记号与约定
  • 线性筛
  • 两例变量分离
  • 关于 $\mu^2(i)$
• 狄利克雷卷积
  • 狄利克雷卷积
  • 基本性质
  • 几个简单的狄利克雷卷积
  • 狄利克雷卷积的求法 (特别鸣谢: pp_orange)
    • 当 $f,g$ 是普通函数时
    • 当 $f$ 为积性函数时
    • 当 $f,g$ 为积性函数时
• 杜教筛
  • 前置：$\sum \limits_{i = 1}^{n} i^k \ \ (k\leq 4)$ 公式
  • 流程
  • 两例简单的杜教筛
  • 点乘
  • 简单运用
    • 例 1
    • 例 2 (Luogu P3768 简单的数学题)
    • 例 3
  • [NOI2016] 循环之美
  • BZOJ3512: DZY Loves Math IV
• Min_25 筛
  • 综述
  • 流程
    • 第一步
    • 复杂度分析 · 一
    • 第二步
      • 综述
      • 方法一
      • 代码 (For Luogu P4213)
      • 方法二
  • 代码 (For LibreOJ 6784)
  • 小试牛刀
    • 例 1

数论函数#

记号与约定#

称函数 $f$ 为数论函数，当且仅当其定义域为 $\mathbb Z$，陪域为 $\mathbb C$。

称数论函数 $f$ 为积性函数当且仅当对于任意 $a,b \in \mathbb N_+,\gcd(a,b) = 1$，有 $f(ab) = f(a)f(b)$。

称数论函数 $f$ 为完全积性函数当且仅当对于任意 $a,b \in \mathbb N_+$，有 $f(ab) = f(a)f(b)$。

按照如下记号记录部分数论函数：

• $I(n) = 1$
• $\epsilon(n) = [n = 1]$
• $id(n) = n$
• $id_k(n) = n^k$
• $d(n) = \sum \limits_{d \mid n} 1$
• $\sigma(n) = \sum \limits_{d \mid n} d$
• $\mu,\varphi$ 定义照常，不再赘述。

其中 $I,\epsilon,id,id_k$ 完全积性。

另外，文中部分 $a/b,\frac{a}{b}$ 型分数根据上下文可能有向下取整（$\left \lfloor \frac ab \right \rfloor$）的含义，请注意甄别。

线性筛#

对于任意积性函数，我们可以线性筛求出。

我们只需要关心积性函数在素数幂 $p^k$ 处的取值。为了检验一个数 $i$ 是否是 $p^k$ 型的，我们可以维护 $low(i)$ 表示 $i$ 中最小质因子的部分。例如 $low(6) = 2,low(49) = 49,low(12)=4$。当且仅当 $i = low(i)$ 时 $i$ 是 $p^k$ 型的。

具体的 $k$ 和 $p$ 也是不难维护的。

下面的模板来自 Alex-Wei。

for(int i = 2; i < N; i++) {
  if(!vis[i]) pr[++pcnt] = i, f[i] = ..., low[i] = i; // 单独算 f(p)
  for(int j = 1; j <= pcnt && i * pr[j] < N; j++) {
    vis[i * pr[j]] = 1;
    if(i % pr[j] == 0) { // i 与 p 不互质
      low[i * pr[j]] = low[i] * pr[j];
      if(i == low[i]) f[i * pr[j]] = ...; // i = p ^ k，单独算 f(p ^ {k + 1})
      else f[i * pr[j]] = f[i / low[i]] * f[low[i * pr[j]]];
      break;
    }
   	low[i * pr[j]] = pr[j];
    f[i * pr[j]] = f[i] * f[pr[j]]; // i 与 p 互质，f(ip) = f(i)f(p) 
  }
}

两例变量分离#

• $d(ij)$

  $$d(ij) = \sum \limits_{x \mid i} \sum \limits_{y \mid j} [\gcd(x,y) = 1]$$

  构造双射容易证明。

• $\varphi(ij)$

  $$\varphi(ij) =\frac{\varphi(i)\varphi(j)\gcd(i,j)}{\varphi(\gcd(i,j))}$$

  证明：

  $$\varphi(i)\varphi(j)\gcd(i,j) = \gcd(i,j) \left(i \prod \limits_{p \mid i,p \in \mathbb P} \frac{p-1}{p}\right) \left(j \prod \limits_{q \mid j,q \in \mathbb P} \frac{q-1}{q}\right) \\ = ij\gcd(i,j) \left(\prod \limits_{p \mid ij,p \in \mathbb P} \frac{p-1}{p}\right) \left(\prod \limits_{q \mid \gcd(i,j),q \in \mathbb P} \frac{q-1}{q}\right) \\ = \varphi(ij)\varphi(\gcd(i,j))$$

关于 $\mu^2(i)$#

有性质：

$$\mu^2(i) = \sum \limits_{d^2 \mid i} \mu(d)$$

当 $\mu(i) \neq 0$ 时，必然不存在平方因子，因此恰有 $d = 1$ 造成 $1$ 的贡献。

当 $\mu(i) = 0$ 时，设 $p_1,p_2,\cdots,p_k$ 这些质因子在 $i$ 中的次数至少为 $2$。$d$ 要在右侧造成贡献，只能是这些质因子构成的集合的某个子集的乘积，根据二项式定理，贡献为 $0$。

狄利克雷卷积#

狄利克雷卷积#

狄利克雷卷积：是一个算子。对于数论函数 $f,g$，我们有

$$(f*g)(n) = \sum \limits_{d \mid n}f(n)g(\frac nd)$$

狄利克雷逆：若 $f* g = \epsilon$ 则称 $f,g$ 互为狄利克雷逆。

对于 $f(1) \neq 0$ 的 $f$，$f$ 的狄利克雷逆一定存在。显然这条件已经必要。如下推导得到 $g$：

$$\epsilon(n) = \sum \limits_{d | n} f(d) g(\frac n d) \\ f(1)g(n) = \epsilon(n) - \sum \limits_{d \mid n,d > 1}f(d)g(\frac{n}{d}) \\ g(n) = \frac{1}{f(1)}\left(\epsilon(n) - \sum \limits_{d \mid n,d > 1}f(d)g(\frac{n}{d})\right)$$

基本性质#

• 交换律

  $$f * g= g* f$$

• 结合律

  $$(f * g) * h = f* (g* h)$$

• 分配律

  $$(f\pm g) * h = f * h \pm g * h$$

• 消去律

  $$f * h = g * h \Leftrightarrow f = g \ \ (h(1) \neq 0)$$

• 两个积性函数的狄利克雷卷积是积性函数

• 积性函数的狄利克雷逆是积性函数

几个简单的狄利克雷卷积#

• $\mu * I = \epsilon$

  这也是为什么我们可以把 $[\gcd(x,y) = 1]$ 化开变成 $\sum \limits_{d \mid \gcd(x,y)} \mu(d)$。

• $\varphi * I = \text{id}$

  这也是为什么我们可以把 $\gcd(x,y)$ 化开变成 $\sum \limits_{d \mid \gcd(x,y)} \varphi(d)$。

• $I * I = d$

  根据 $d$ 的定义。

• $\text{id} * I = \sigma$

  根据 $\sigma$ 的定义。

狄利克雷卷积的求法 (特别鸣谢: pp_orange)#

$f,g$ 已知，每一项可以 $O(1)$ 求出。设我们要计算 $h = f*g$ 的前 $n$ 项。

当 $f,g$ 是普通函数时#

采用定义式 $O(n \log n)$ 计算。

Code：[V 我 50 我帮你写]。

当 $f$ 为积性函数时#

考察 $f = I$ 这一特例。此时等价于狄利克雷前缀和，复杂度为 $O(n \log \log n)$。

当 $f$ 为任意积性函数时，我们可以执行类似的操作。考虑枚举质数 $p$。初始我们的卷积中只有 $h'(n) = f(1)g(n)$ 这一项。每枚举一个质数 $p$，我们就枚举不含因子 $p$ 的正整数 $n$，把 $f(p)h(n),f(p^2)h(n),f(p^3)h(n),\cdots\ (p \nmid n)$ 这样的项贡献到对应的新的 $h'(np),h'(np^2),h'(np^3),\cdots$ 上。这样我们就相当于在 $p$ 这个维度做了一遍前缀和，现在 $h'(n)$ 里面包含了所有 $f(x)g(n/x)$ 的项，其中 $x$ 整除 $n$，且 $x$ 的质因子集合中只有我们枚举过的质数。

因为高维前缀和是对的，那么这个操作也应当是对的。此时复杂度为 $\sum \limits_{1 \leq p^k \leq n,k \ge 1,p \in \mathbb P} \left \lfloor \frac{n}{p^k} \right \rfloor = O(n \log \log n)$。这是因为对于每个 $p$，后面的部分都可以看作余项（等比数列求和不超过常数倍的第一项），因此其复杂度和狄利克雷前缀和 / 埃氏筛法相同。

Code：[V 我 50 我帮你写]。

当 $f,g$ 为积性函数时#

时间复杂度为 $O(n)$。

我们只需要求出 $h$ 在素数幂处的取值。对于每个素数直接暴力，小于等于 $\sqrt n$ 的素数复杂度可以看成 $O\left(\dfrac {\sqrt n}{\log \sqrt n} \times \log^2 n\right) = O(\sqrt n \log n)$，而对于大于 $\sqrt n$ 的素数只会有一项。因此，计算素数幂处取值的复杂度可以看作低阶项，不影响用 $O()$ 计算的复杂度。

Code：[V 我 50 我帮你写]。

杜教筛#

前置：$\sum \limits_{i = 1}^{n} i^k \ \ (k\leq 4)$ 公式#

• $k = 1$：$n(n+1)/2$。
• $k = 2$：$n(n+1)(2n+1)/6$。
• $k = 3$：$n^2(n+1)^2/4$。
• $k = 4$：$n(n+1)(2n+1)(3n^2+3n-1)/30$。已经很少用到了。
• 对于 $k > 4$，使用连续点值拉格朗日插值可以 $O(k)$ 计算一次。

流程#

考虑求某数论函数 $f$ 的前缀和 $S(n) = \sum \limits_{i = 1}^n f(i)$。杜教筛的想法是通过一个比较好的数论函数 $g$ 和比较好的 $f*g$ 把 $f$ 的前缀和化简。这和 $f,g,f*g$ 是否是积性函数无关。

我们知道

$$\sum \limits_{i = 1}^{n} (f*g)(i) = \sum \limits_{i = 1}^{n} g(i) S(\left \lfloor \frac ni \right \rfloor) \\ g(1)S(n ) = \sum \limits_{i = 1}^{n} (f*g)(i) - \sum \limits_{i = 2}^{n} g(i) S(\left \lfloor \frac ni \right \rfloor)$$

如果 $f*g$ 和 $g$ 的前缀和比较简单，我们就可以快速求出 $f$ 的前缀和了。

考虑如下伪代码：

inline int SumF(int n){ // 注意：这里计算的是 g(1)S(n). 但大部分时候 g 是积性函数或其点乘，因此 g(1)=1.
	int ans=SumFG(n);
	for(int l=2,r;l<=n;l=r+1){
		r=n/(n/l);
		ans-=(SumG(r)-SumG(L-1))*SumF(n/l);
	}
	return ans;
}

当使用记忆化的时候，根据整除分块经典结论，我们只会用到不超过 $2 \sqrt n$ 个值。

那么 $T(n) = O(\sqrt n) + \sum \limits_{i = 1}^{\sqrt n}\left(O(\sqrt i) +O\left(\sqrt {\dfrac{n}{i}}\right)\right)$。求和改积分容易证得 $T(n) = O(n^{3/4})$。

否则时间复杂度不详（OI-Wiki 对不加以记忆化的时间复杂度证明作出了证伪：不能直接将第二层之后的 $T\left(\sqrt{\dfrac{n}{ij}}\right)$ 看作余项）。

如果我们预处理部分 $S(i) (i = 1,2,\cdots,m) \ \ (m \geq \sqrt n)$，设预处理复杂度为 $T_0(m)$，则总复杂度为

$$T(n) = T_0(m) + \sum \limits_{k \in \{\left\lfloor n/x\right\rfloor\},k > m} T(k) \\ T_0(m) + \sum \limits_{k = 1}^{n/m} O(\sqrt{\frac nk}) \\ = O\left(T_0(m) +\frac{n}{\sqrt m}\right)$$

那么一般 $T_0(m) = O(m)$，此时 $m = \Theta(n^{2/3})$ 最优。

详见 OI-Wiki。

两例简单的杜教筛#

• $f = \mu$

  根据 $\mu * I = \epsilon$，取 $f = \mu,g = I, f*g = \epsilon$。

• $f = \varphi$

  根据 $\varphi * I = id$，取 $f = \varphi,g = I,f * g = id$。

点乘#

定义 $f,g$ 点乘的结果 $f\cdot g$ 是一个函数，满足 $(f\cdot g)(n) = f(n)g(n)$。

对于数论函数 $A,B,C$，当 $C$ 是完全积性函数的时候，我们有

$$(A\cdot C)*(B \cdot C) = (A*B) \cdot C$$

那么当我们要求前缀和的函数 $f$ 为两个简单函数点乘的时候，我们可以运用这个等式，从而寻求形式较好的 $g$ 和 $f*g$，然后杜教筛解决。

下面是两例比较基本的点乘处理方法。

• $f = \mu \cdot id_k$：取 $g = id_k$。

  $$h = (\mu \cdot id_k )* id_k = (\mu \cdot id_k) * (I * id_k) = (\mu * I) * id_k = \epsilon$$

• $f = \varphi \cdot id_k$：取 $g = id_k$。

  $$h = (\varphi \cdot id_k )* id_k = (\varphi \cdot id_k) * (I * id_k) = (\varphi * I) * id_k = id_{k+1}$$

下面是抄的一例点乘：

• $f = \mu^2 * (\mu \cdot id)$：取 $g = id$。
  $$h = \mu^2 * (\mu \cdot id) * id = \mu^2 * ((\mu * I) \cdot id) = \mu ^2 * \epsilon = \mu ^2$$
  而 $\mu^2$ 的前缀和可以通过公式 $\mu^2(i) = \sum \limits_{d^2 \mid i} \mu(d)$ 计算：简单地交换求和号，前缀和变为 $\sum \limits_{d = 1}^{\sqrt n} \mu(d) \left \lfloor \frac{n}{d^2} \right \rfloor$。直接计算即可做到 $O(\sqrt n)$，求出 $\mu(i)$ 的前缀和可以进一步做到 $O(n^{1/3})$。但是在低于 $O(\sqrt n)$ 时间复杂度内求出 $\mu(i)$ 的前缀和需要再写筛法，其实就不如直接 $O(\sqrt n)$ 求了。

另一例略显奔放的点乘（出自某场联考的最后一部分）：

• $f = id_k \cdot ((\mu \cdot id_k) * I)$：取 $g = id_{2k}$。
  $$h = (id_k\cdot ((\mu \cdot id_k) * I) ) * (id_k \cdot id_k)= ({\color{red}(\mu \cdot id_k) * id_k}* I) \cdot id_k \\ = (\epsilon * I) \cdot id_k \\ = id_k$$
  红色部分是第一例点乘中出现的结果，它等于 $\epsilon$。

小结：欲求 $f = A\cdot C$ 形的前缀和，我们希望 $g$ 形如 $B \cdot C$，且 $g$ 的前缀和易求。然后，运用上面的等式将 $f * g$ 化作 $(A*B) \cdot C$。适当地选取 $B$ 可以将 $A*B$ 化作简单的函数，从而 $f*g$ 的前缀和易求。

简单运用#

例 1#

求

$$\sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^{m} \varphi(\gcd(i,j))$$

对某个大质数取模

$1 \leq n \leq m \leq 10^9$

经过平凡的推导可以得出原式等于

$$\sum \limits_{T = 1}^{n} \left \lfloor \frac nT \right \rfloor \left \lfloor \frac mT \right \rfloor \sum \limits_{i \mid T} \mu(i) \varphi\left(\frac Ti\right)$$

想要整除分块，则需要求出 $\mu * \varphi$ 的前缀和。

令 $f = \mu * \varphi,g = I,f*g = \varphi$。利用杜教筛求出 $\varphi$ 和 $\mu * \varphi$ 的前缀和即可。注意到复杂度仍然是 $O(n^{2/3})$ 的，因为我们可以线性筛预处理 $f$ 的前 $O(n^{2/3})$ 项。

例 2 (Luogu P3768 简单的数学题)#

求

$$\sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^{n} ij\gcd(i,j)$$

对某个大质数取模

$1 \leq n \leq 10^9$

经过平凡的推导可以得出原式等于

$$\sum \limits_{d = 1}^{n} \varphi(d) d^2 S\left(\left \lfloor \frac ni\right \rfloor\right)^2, \text{where } S(x) = \frac{x(x+1)}{2}$$

我们知道了怎么处理点乘，因此杜教筛也是容易的：$f = \varphi \cdot id_2,g = id_2,(f*g) = id_3$。

例 3#

求

$$\sum \limits_{i = 1}^{n} \sum \limits_{j =1}^{n} \operatorname{lcm}(i,j)^k$$

$1 \leq n \leq 10^9,1 \leq k \leq 1000$

这碟醋终于出现了！

平凡推导可得原式等于

$$\sum \limits_{T=1}^n S_k\left(\dfrac nT\right)^2 T^k \sum \limits_{t \mid T} \mu(t) \cdot t^k,\\ \text{where }S_k(n) = \sum \limits_{i = 1}^{n} i^k$$

按照上面点乘处的推导，令 $f(n) = n^k \sum \limits_{d \mid n} \mu(d) \cdot d^k$，请出 $g(n) = n^{2k}$，$f * g = id_k$。

一个粗劣的复杂度上界是 $O((\sqrt n + n^{2/3})k)$，但事实上整个杜教筛过程中也只会用到 $O(\sqrt n)$ 个 $S_k$ 的值，对它们记忆化可以分析出更好的复杂度 $O(\sqrt n k+n^{2/3})$。

[NOI2016] 循环之美#

简要题面：

$$\sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^m [\gcd(i,j) = 1][\gcd(j,k) = 1]$$

$1 \leq n,m \leq 10^9,1 \leq k \leq 10^3$

这题处理技巧比较奇怪，我们不太能直接把两个艾弗森括号都拆成 $\mu$，而是先保留下来 $[\gcd(j,k) = 1]$ 的那部分。在这里，和一个给定常数互质这个条件是比较好的。

$$\sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^m [\gcd(i,j) = 1][\gcd(j,k) = 1] \\ = \sum \limits_{j = 1}^m [\gcd(j,k) = 1] \sum \limits_{i = 1}^{n} [\gcd(i,j) = 1] \\ = \sum \limits_{j = 1}^m [\gcd(j,k) = 1] \sum \limits_{i = 1}^{n} \sum \limits_{d \mid \gcd(i,j)} \mu(d) \\\\ = \sum \limits_{i = 1} ^{\min(n,m)} \mu(d) [\gcd(d,k) = 1] \left \lfloor \frac nd \right \rfloor \sum \limits_{j = 1}^{m/d} [\gcd(j,k) = 1]$$

然后要求两个东西：

$$f(n) = \sum \limits_{i = 1}^{n} [\gcd(i,k) = 1] \\ g(n,k) = \sum \limits_{i = 1}^{n} [\gcd(i,k) = 1] \mu(i)$$

对于 $f(n)$，根据 $\gcd(i,k) = \gcd(i+k,k)$，艾弗森括号的取值是 $k$ - 循环的，我们只需要预处理 $f(1) \sim f(k)$ 即可。

对于 $g(n,k)$ 则比较需要有想象力。

$$g(n,k) = \sum \limits_{i = 1}^{n} \mu(i) \sum \limits_{d \mid \gcd(i,k)} \mu(d) \\ = \sum \limits_{1 \leq d \leq n, d \mid k} \mu(d) \sum \limits_{i=1}^{n/d} \mu(id) \\ = \sum \limits_{1 \leq d \leq n, d \mid k} \mu^2(d) \sum \limits_{i = 1}^{n/d} [\gcd(i,d) = 1]\mu(i) \\ = \sum \limits_{1 \leq d \leq n,d \mid k} \mu^2(d) g(n/d,d)$$

递归计算，直到递归到 $d = 1$ 时用杜教筛求 $\mu$ 前缀和。

注意到最终要算 $O(\sqrt n)$ 个 $g(n,...)$，$...$ 这一维是 $k$ 的因子，因此总复杂度为 $O(d(k)\sqrt n+ n^{2/3} + k)$。

BZOJ3512: DZY Loves Math IV#

求

$$\left(\sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^{m} \varphi(ij) \right)\bmod 998244353$$

$1 \leq n \leq 10^5,1 \leq m \leq 10^9,n \leq m$

根据《毒瘤之神的考验》一题快进到

$$\sum \limits_{T =1}^{n} \left(\sum \limits_{i = 1}^{n/T} \varphi(iT) \sum \limits_{j = 1}^{m/T} \varphi(jT) \right)\sum \limits_{d \mid T} \frac{d\mu\left(\frac Td\right)}{\varphi(d)}$$

记 $S(a,b) = \sum \limits_{i = 1}^{a} \varphi(ib),f(n) = \sum \limits_{d \mid n} \frac{d\mu\left(\frac Td\right)}{\varphi(d)}$，则上式可写作

$$\sum \limits_{i = 1}^{n} f(i)S(n/i,i)S(m/i,i)$$

在《毒瘤之神的考验》一题中，要用到的全体 $S(x,y)$ 可 $O(n \log n)$ 预处理，从而最终复杂度也是 $O(n \log n)$ 的。

本题需要神秘的观察。我们考虑类似扰动法的策略，

$$\text{let } A(n,m) = \sum \limits_{i = 1}^{n} \sum \limits_{j = 1}^{m} \varphi(ij) = \sum \limits_{i = 1}^{n} f(i)S(n/i,i)S(m/i,i), \text{then} \\ S(n,k) = \sum \limits_{i = 1}^{n} \varphi(ik) \\ = A(n,k) - A(n,k-1) \\ = \sum \limits_{i = 1}^{n} f(i)S(n/i,i)\left(S(k/i,i) - S((k-1)/i,i)\right) \\ = \varphi(k) \sum \limits_{i \mid k} f(i) S(n/i,i)$$

注意到 $k$ 这一维是小于等于给定的 $n$ 的（和此处的变量 $n$ 作区分），于是对于 $S(m/i,i)$ 这样的状态来说，计算它们的复杂度不太大。我们可以采用记忆化 + 预处理小的部分来减小递归的复杂度。

Min_25 筛#

综述#

考虑求一类积性函数的前缀和，其中 $f$ 是积性函数，满足：

• 对于任意素数 $p$，$f(p)$ 处取值是关于 $p$ 的低次多项式。
• $f(p^c)$ 处取值容易求得。

则 Min_25 筛可以在 $O(\frac{n^{3/4}}{\log n})$ 或 $O(n^{1 - \epsilon})$ 的复杂度内求出 $f$ 的前缀和。

流程#

Min_25 筛的求解分为两步：

1. 对于 $W(n) = \left\{\left\lfloor \dfrac{n}{x}\right\rfloor \mid 1 \leq x \leq n\right\}$，考虑其中每一个元素 $m$，求出 $\sum \limits_{1 \leq i \leq m,i \in \mathbb P} f(i)$，即前缀素数点处取值之和；
2. 根据第一步求得的 $O(\sqrt n)$ 个答案进而求解出整个前缀和。

第一步#

因为 $f(p)$ 的取值是关于 $p$ 的低次多项式，因此，在这一步中，我们只需要考虑 $f(p) = p^k$ 的情形如何求解（对于项数更多的情形，分别求解后求和即可）。

我们通常找一个 $f'$ 来替代 $f$。这样的 $f'$ 具有完全积性，且在素数 $p$ 处和 $f$ 相同。这里考虑函数 $f'(x) = x^k$。

Min_25 筛的第一步，我们要对于每个 $m$，求出数组 $g$。其中 $g$ 的定义如下：

$$g(m,j) = \sum \limits_{1 \leq i \leq m,(i \in \mathbb P \operatorname {or} L(i) >P_j)} f'(i)$$

其中 $L(i)$ 表示 $i$ 的最小质因子，$P_j$ 表示第 $j$ 个质数，$P_0 = 0$。上式即前缀素数点处或最小质因子大于 $P_j$ 处 $f'$ 取值之和（注意是 $f'$ 而非 $f$。我们用 $f'$ 暂时替代了 $f$）。

对于 $1$，其不存在最小质因子，故 $[L(1) > P_j]$ 总是 $0$。此处可改作 $2 \leq i \leq m$，但为了推导方便保留可能的 $i = 1$ 取值，虽然 $f'(1)$ 永远不能被统计进贡献。

设 $q(m)$ 表示最大的 $j$ 使得 $P_j^2 \leq m$，那么 $g(m,q(m))$ 即为 $\sum \limits_{1 \leq i \leq m,i \in \mathbb P} f(i)$。

考虑递推求出 $g(m,j)$。初始有 $g(m,0) = \sum \limits_{i = 2}^{m} i^k$。对于 $j \geq 1$，我们会在 $g(m,j-1)$ 的基础上删除若干个点。质数是不会被删除的，因此我们会恰好删除最小质因子是 $P_j$ 的非质数点。

若 $P_j^2 > m$，$g(m,j)$ 不会有变化。

若 $P_j^2 \leq m$，形如 $P_jx(P_jx \leq m)$，且 $P_jx$ 的最小质因子恰好是 $P_j$ 的点会造成贡献。这可以表示为下面的式子：

$$f'(P_j)\left( g\left(\left \lfloor \dfrac{m}{P_j} \right \rfloor,j-1\right)-g(P_{j-1},j-1)\right)$$

来解释一下这个式子：首先，我们可以根据完全积性把 $f'(P_jx)$ 拆开。现在考虑 $x$ 的取值，很显然是来源于 $g\left(\left \lfloor \dfrac{m}{P_j} \right \rfloor,j-1\right)$。但是既小于等于 $\left \lfloor \dfrac{m}{P_j} \right \rfloor$ 又小于等于 $P_{j-1}$ 的质数是不能作为 $x$ 的，因为我们有 $P_{j-1} < P_j \leq \sqrt m$，所以 $\left \lfloor \dfrac{m}{P_j} \right \rfloor \geq P_{j-1}$，只需要考虑 $g(P_{j-1},j-1)$ 就可以了。

因此，我们得到了下式：

$$g(m,j) = \begin{cases} g(m,j-1) & P_j^2 > m \\ g(m,j-1)-f'(P_j) \times \left( g\left(\left \lfloor \dfrac{m}{P_j} \right \rfloor,j-1\right)-g(P_{j-1},j-1)\right) & P_j^2 \leq m\end{cases}$$

复杂度分析 · 一#

现在考虑复杂度。根据经典结论 $\left \lfloor \dfrac{n}{ab} \right \rfloor=\left \lfloor \dfrac{\frac{n}{a}}{b} \right \rfloor$，我们对于 $m \in W(n)$ 有 $W(m) \subseteq W(n)$。同时可以证明 $\forall 1 \leq i \leq \sqrt n$，$i \in W(n)$，因此对于任意 $m \in W(m)$ 我们求 $g(m,j)$ 的时候，只会遍历 $g$ 数组中第一维在 $W(n)$ 中的元素。同时，求出单个状态的复杂度是 $O(1)$ 的，因此状态数即为复杂度。

根据范围素数近似估计，复杂度可写作

$$T(n) = \sum \limits_{i \leq \sqrt n} O(\pi(\sqrt i)) + O\left(\pi\left (\sqrt \frac ni\right) \right) \\ = \sum \limits_{i \leq \sqrt n} O\left(\pi\left (\sqrt \frac ni\right) \right) \\ = \sum \limits_{i \leq \sqrt n} O\left( \dfrac{\sqrt \frac ni }{\ln \sqrt {\frac ni}}\right) \\ $$

注意到 $i \leq \sqrt n$，因此 $\ln \sqrt \frac ni$ 只和 $\ln n$ 差常数倍。直接将这部分提出来，上面的部分积分，得到

$$T(n ) =O\left(\frac{n^{3/4}}{\log n}\right)$$

第二步#

综述

现在我们已经知道所有 $\sum \limits_{1 \leq i \leq m,i \in \mathbb P} f(i),\ \text{where}\ m \in W(n)$ 了。将其简记为 $g(m)$。考虑使用这些 $g(m)$ 求出原来的解。

设 $S(m,j)$ 表示 $\sum \limits_{1 \leq i \leq m}[L(i) > P_j] f(i)$，即最小质因子大于 $P_j$ 的 $f(i)$ 之和（我们现在抛弃 $f'$ 了，同时也不再限定 $i$ 是质数；只不过 $1$ 仍然被抛弃在外）。这样再调用 $S(n,0)$ 就得到答案了。

方法一

这种方法复杂度为 $O(n^{1 - \epsilon})$，在 $n \leq 10^{13}$ 的时候也可以看作 $O\left(\frac{n^{3/4}}{\log n}\right)$，够用。实际测试略快于第二种。这里着重介绍。

仍然考虑如何求解。分 $i$ 是质数或合数两种情况讨论：

• $i$ 是质数

  根据 $g$ 的定义，这部分的贡献是 $g(n) - g(P_j)$。

• $i$ 是合数

  枚举最小质因子 $P_k (k > j)$，枚举指数 $e$，利用积性函数的性质进行变量分离，得到

  $$\sum \limits_{k > j} \sum \limits_{e \geq 1,P_{k}^e \leq m} f(P_{k}^e)\left(S\left(\left\lfloor\frac{m}{P_k^e}\right\rfloor,k\right)+ [e>1]\right)$$

  这里 $[e>1]$ 有些耐人寻味。考虑前面的 $S\left(\left\lfloor\frac{m}{P_k^e}\right\rfloor,k\right)$ 没有考虑 $P_{k}^e$ 自己的贡献。如果 $e = 1$，那么 $P_k^e$ 就是 $P_k$，应该归到 $i$ 是质数的分类里，没有贡献。否则，$P_{k}^e$ 会造成 $f(P_{k}^e)$ 的贡献。

则 $S(n,0) + 1$ 就是答案（我们之前没有考虑 $f(1)$ 的贡献）。

代码 (For Luogu P4213)

略好于之前的版本。

# include <bits/stdc++.h>

const int N=1000010,INF=0x3f3f3f3f,mod=1e9+7,inv6=166666668,inv2=500000004;

typedef long long ll;
int pr[N],pc;
bool vis[N];

int g1[N],g2[N];

ll w[N],idx1[N],idx2[N];

int wtot;

ll MAXN;
ll n;

inline int read(void){
	int res,f=1;
	char c;
	while((c=getchar())<'0'||c>'9')
		if(c=='-')f=-1;
	res=c-48;
	while((c=getchar())>='0'&&c<='9')
		res=res*10+c-48;
	return res*f;
}
inline void init(void){
	for(int i=2;i<=MAXN;++i){
		if(!vis[i]) pr[++pc]=i;
		for(int j=1;i*pr[j]<=MAXN&&j<=pc;++j){
			vis[i*pr[j]]=true;
			if(i%pr[j]==0) break;
		}
	}
	return;
}
inline int sum1(ll x){
	x%=mod;
	return x*(x+1)%mod*inv2%mod;
}
inline int sum2(ll x){
	x%=mod;
	return x*(x+1)%mod*(2*x+1)%mod*inv6%mod;
}
inline int idx(ll x){
	return (x<=MAXN)?idx1[x]:idx2[n/x];
}

inline int adc(int a,int b){
	return (a+b>=mod)?(a+b-mod):(a+b); 
}
inline int dec(int a,int b){
	return (a<b)?(a-b+mod):(a-b);
}
inline void add(int &a,int b){
	a=adc(a,b);
}
inline void del(int &a,int b){
	a=dec(a,b);
}
inline int mul(int a,int b){
	return 1ll*a*b%mod;
}

ll S(ll x,int i){
	if(pr[i]>=x) return 0;
	
	int xid=idx(x),pid=idx(pr[i]);
	ll res=dec(dec(g2[xid],g1[xid]),dec(g2[pid],g1[pid]));
	
	for(int j=i+1;j<=pc&&1ll*pr[j]*pr[j]<=x;++j){
		for(ll e=1,p=pr[j],_p;p<=x;++e,p*=pr[j]){
			_p=p%mod;
			res=adc(res,mul(mul(dec(_p,1),_p),S(x/p,j)+(e>1)));
		}
	}
	return res;
}
int main(void){
	scanf("%lld",&n);
	MAXN=1ll*sqrt(n);
	init();
	for(ll l=1,r;l<=n;l=r+1){
		ll val=n/l;
		r=n/val,w[++wtot]=val;
		g1[wtot]=dec(sum1(val),1),g2[wtot]=dec(sum2(val),1);
		if(val<=MAXN) idx1[val]=wtot;
		else idx2[n/val]=wtot;
	}
	for(int i=1;i<=pc;++i){
		for(int j=1;1ll*pr[i]*pr[i]<=w[j]&&j<=wtot;++j){
			int mid=idx(w[j]/pr[i]),pid=idx(pr[i-1]);
			del(g1[j],mul(pr[i],dec(g1[mid],g1[pid])));
			del(g2[j],mul(mul(pr[i],pr[i]),dec(g2[mid],g2[pid])));
		}
	}
	printf("%d",adc(S(n,0),1));
	return 0;
}

方法二

再度审视 $S$ 的定义：$S(m,j)$ 表示 $\sum \limits_{1 \leq i \leq m}[L(i) > P_j] f(i)$，即最小质因子大于 $P_j$ 的 $f(i)$ 之和。

我们直接套用第一步中的递推方式。

$$S(m,j) = S(m,j+1) + \\ [P_{j+1} \leq \sqrt m]\sum \limits_{P_{j+1}^{e} \leq m} f(P_{j+1}^e)\left( S\left( \left \lfloor\frac{m}{P_{j+1}^e} \right \rfloor,j+1 \right) - g\left (\min\left(\left \lfloor\frac{m}{P_{j+1}^e} \right \rfloor,P_{j+1}\right)\right) + [e>1] \right)$$

这里的取 $\min$ 是不好的。我们注意到当 $P_{j+1}^{e+1} > m$ 的时候，$S() - g()$ 没有贡献，因为这意味着 $P_{j+1} > \dfrac{m}{P_{j+1}^e}$。于是可以改为

$$S(m,j) = S(m,j+1) + \\ [P_{j+1} \leq \sqrt m]\sum \limits_{P_{j+1}^{e+1} \leq m} f(P_{j+1}^e)\left( S\left( \left \lfloor\frac{m}{P_{j+1}^e} \right \rfloor,j+1 \right) - g\left (P_{j+1}\right) \right) + f(P_{j+1}^{e+1})$$

使用和第一步相同的方法进行递推即可。

至此可以使用与第一步类似的方法分析出 $O\left(\frac{n^{3/4}}{\log n}\right)$。

方法二的优势是，我们可以求出每个 $m \in W(n)$ 的前缀 $m$ 的答案。这无疑在整除分块需要计算 $O(\sqrt n)$ 个前缀的函数值之和的时候给我们提供了便利。

代码 (For LibreOJ 6784)#

# include <bits/stdc++.h>

const int N=2000010,INF=0x3f3f3f3f,mod=1e9+7,inv6=166666668,inv2=500000004;

typedef long long ll;
int pr[N],pc;
bool vis[N];

int g1[N],g0[N];

ll w[N],idx1[N],idx2[N];

int s[N];

int wtot;

ll MAXN;
ll n;

inline int read(void){
	int res,f=1;
	char c;
	while((c=getchar())<'0'||c>'9')
		if(c=='-')f=-1;
	res=c-48;
	while((c=getchar())>='0'&&c<='9')
		res=res*10+c-48;
	return res*f;
}
inline void init(void){
	for(int i=2;i<=MAXN;++i){
		if(!vis[i]) pr[++pc]=i;
		for(int j=1;i*pr[j]<=MAXN&&j<=pc;++j){
			vis[i*pr[j]]=true;
			if(i%pr[j]==0) break;
		}
	}
	return;
}
inline int sum1(ll x){
	x%=mod;
	return x*(x+1)%mod*inv2%mod;
}
inline int sum0(ll x){
	return x%mod;
}
inline int idx(ll x){
	return (x<=MAXN)?idx1[x]:idx2[n/x];
}

inline int adc(int a,int b){
	return (a+b>=mod)?(a+b-mod):(a+b); 
}
inline int dec(int a,int b){
	return (a<b)?(a-b+mod):(a-b);
}
inline void add(int &a,int b){
	a=adc(a,b);
}
inline void del(int &a,int b){
	a=dec(a,b);
}
inline int mul(int a,int b){
	return 1ll*a*b%mod;
}
inline int f(int x,int k){
	return x^k;
}

std::vector <int> vec;

int main(void){
	scanf("%lld",&n);
	MAXN=1ll*sqrt(n);
	init();
	for(ll l=1,r;l<=n;l=r+1){
		ll val=n/l;
		r=n/val,w[++wtot]=val;
		g1[wtot]=dec(sum1(val),1),g0[wtot]=dec(sum0(val),1);
		if(val<=MAXN) idx1[val]=wtot;
		else idx2[n/val]=wtot;
	}
	for(int i=1;i<=pc;++i){
		for(int j=1;1ll*pr[i]*pr[i]<=w[j]&&j<=wtot;++j){
			int mid=idx(w[j]/pr[i]),pid=idx(pr[i-1]);
			del(g1[j],mul(pr[i],dec(g1[mid],g1[pid])));
			del(g0[j],dec(g0[mid],g0[pid]));
		}
	}
	
	for(int i=1;i<=wtot;++i) s[i]=g1[i]=adc(dec(g1[i],g0[i]),(w[i]>=2)*2);
	for(int i=pc;i;--i){
		for(int j=1;1ll*pr[i]*pr[i]<=w[j]&&j<=wtot;++j){
			for(ll e=1,p=pr[i];p<=w[j]/pr[i];++e,p*=pr[i]){
				add(s[j],mul(f(pr[i],e),dec(s[idx(w[j]/p)],g1[idx(pr[i])])));
				add(s[j],f(pr[i],e+1));
			}
		}
	}
	vec.resize(wtot);
	for(int i=0;i<wtot;++i) vec[i]=adc(s[i+1],1);
	std::sort(vec.begin(),vec.end());
	vec.erase(std::unique(vec.begin(),vec.end()),vec.end());
	int ans=0;
	for(auto v:vec) ans^=v;
	printf("%d",ans);
	
	return 0;
}

小试牛刀#

例 1#

(杜教筛 例 3)

求

$$\sum \limits_{i = 1}^{n} \sum \limits_{j =1}^{n} \operatorname{lcm}(i,j)^k$$

$1 \leq n \leq 10^9,1 \leq k \leq 1000$

$$\sum \limits_{T=1}^n S_k\left(\dfrac nT\right)^2 T^k \sum \limits_{t \mid T} \mu(t) \cdot t^k$$

快进到这里。

那么我们只要瞪眼法看出后面是积性函数，直接上 Min_25 筛就完事了，不需要观察到可以卷 $id_{2k}$。

Typora 崩掉了，不写了。