预备
0.1 渐进符号
其实不少高等数学 / 数学分析教材在讲解无穷小的比较时已经相当严谨地介绍过大 O、小 O 记号，然而各种历史习惯记法的符号滥用（abuse of notation）[1] 直到现在都让笔者头疼. These notations seem to be innocent, but can be catastrophic without careful manipulation. For example,
n=O(n2)∧n2=O(n2)⟹n=n2
Knuth 在《具体数学》里举出的例子[2]. “=” 隐含的对称性使其在 g(x)=O(f(x)) 中格格不入. 事实上，将 O(f(x)) 看作“阶不高于 f(x) 的所有函数的集合”是比“某个阶不高于 f(x) 的函数”更严谨的理解. 因此，本文将使用 f(x)∈O(g(x)) （有时也记为 O(f(x))⊂O(g(x))）的集合论符号代替传统的 f(x)=O(g(x)) 记法.
n2sinn∈O(n2)⟹∑i=1ni2sini∈∑i=1nO(i2)⊂O(∑i=1ni2)⊂O(n3) 或更一般的， g(x)∈O(f(x))⟹∑P(n,i)g(i)∈∑P(n,i)O(f(i))⊂O(∑P(n,i)f(i))
没看出有啥问题，对吧？笔者在写作此文时犯了同样的错误. 请注意，大 O 记号的作用对象是函数，f(i) 是什么？它只是个函数值，是确定的数——这是因为 i 也是求和枚举中确定的数，而不是 n 这种真正代表变元的记号. 所以 O(f(i)) 是什么？它什么也不是.
这种错误的出现是在所难免的，我们太习惯用 x、x3+5x2+x 这种变元都不明确的记号来表示函数了[1] . 写成 f(x) 也不严谨，因为只有 f 才应代表函数本身，f(x) 只能是函数值. 这样我们就可以放心地写下 O(f)，不用担心把变元与确定值弄混了.
然而大家还是喜欢写 O(n2) 和 O(en2)，而不是奇怪的 O(id2) 和 O(exp∘id2). 所以，我们大概只能沿用这种不太严谨的记号，并时刻提醒自己加倍小心了. （形如 x↦ex2 的 λ 风格“匿名函数”记号可能更好？）
但上述命题从结论上是正确的. 正确的推导过程应为 ∑P(n,i)g(i)≤∑P(n,i)Cf(i)≤C∑P(n,i)f(i)∈O(∑P(n,i)f(i))
第一步是直接由大 O 记号的定义得到的结果.
Wikipedia[3] 中有一张详尽的表格介绍了各种渐进符号的定义，OI Wiki[4] 上也有极好的讲解，尚不熟练的读者可以参考. 有兴趣仔细研究的读者可以参考《具体数学》第九章[2] 、Wikipedia 及其 reference（个人推荐 Knuth 关于 O、Ω、Θ 的短文[5] ）. 本文除用 “∈” 和“⊂”替代 “=” 外，完全使用 Knuth 提议的记号体系.
0.2 调和数 H(n) / 调和级数
调和级数的部分和 H(n) 定义为 H(n)=∑i=1n1i 通过一些与 e 有关的数列放缩可以证明 limn→∞(H(n)−logn)=c，其中 c≈0.577 是 Euler 常数. 因此 nH(n)∼nlogn∈Θ(logn).
0.3 自然数等幂和 Pp(n) / p - 级数
p - 级数可视为调和级数的推广. 其部分和定义为 Pp(n)=∑i=1ni−p
p - 级数具有如下性质：
当 p>1 时，p - 级数收敛；
当 p=1 时，p - 级数是调和级数；
当 −∞<p<1 时，我们指出 Pp(n)∼11−pn1−p∈Θ(n1−p)
−∞<p<1 时 p - 级数的渐进估计可以从连续幂函数积分的角度理解. 证明这渐进性，离散情况下，可对 np 差分后前缀和 + 二项式定理得到高次项系数，或可用离散微积分理论得到精确表示（参见《具体数学》[6] ）；连续情况下，Lagrange 中值定理应为较简单的估计方法. 这里从略. 总之，我们得到： Pp(n)∈{Θ(n1−p)p<1Θ(nlogn)p=1Θ(1)p>1
1 约数函数 σz(n)
约数函数（Divisor Function，也可称为除数函数、因数函数）是与 n 的因子有关的一类函数，定义如下：
Definition 1 (约数函数) σz(n)=∑d∣ndz
当 z=0 时，σ0(n) 被称为约数个数函数（number-of-divisors function），常被记为 d(n) 或 τ(n). 当 z=1 时，σ1(n) 被称为约数和函数（sum-of-divisors function），常直接记为 σ(n).
Example 1 估计 σ0(n) 的渐进上界.
也就是估计 n 的因子的数量. 一个广为人知的上界是 2n，因为 n 的所有小于 n 的因子 d 均与另一因子 nd 一一对应.
事实上进一步可以证明 σ0(n)∈o(nϵ)∀ϵ>0[7] ，虽然这在 OI 中并不实用.
Example 2 估计 σ0^(n)=∑i=1nσ0(i) 的渐进上界.
即估计 1 到 n 中所有数因子个数的和. 这是一个形式上鲜为人知但其应用广为人知的例子. 变换求和顺序，容易得到
σ0^(n)=∑i=1nσ0(i)=∑i=1n∑d∣i1=∑d=1n⌊nd⌋≤∑d=1nnd=nH(n)∈O(nlogn)
显然，这比 O(nn) 的平凡估计好上不少. 本例的思路不仅是埃氏筛（Sieve of Eratosthenes）的理论基础，也在杜教筛、快速 Mobius 变换、gcd 卷积[8] 等处出现.
进一步利用此技巧和 p - 级数的估计，我们甚至能在仔细研究 σz(n) 前就得到其前缀和的渐进估计：
Example 3 估计 σz^(n)=∑i=1nσz(i) 的渐进上界.
σz^(n)=∑i=1nσz(i)=∑i=1n∑d∣idz=∑d=1ndz⌊nd⌋≤n∑d=1ndz−1=nP1−z(n)∈{O(nz+1)z>0O(nlogn)z=0O(n)z<0
遗憾的是，对此前缀和做差分并不能得到 σz(n) 的优秀估计.
现在引入一个重要放缩技巧，其在后续估计中屡试不爽.
Proposition 1 ∑d∣nf(d)≤∑i=1nf(⌊ni⌋)
显然，右式比左式多算了 i∤n 的项，因此命题是正确的. 但我们还可以做得更好：
Proposition 2 ∑d∣nf(d)≤∑i=1nf(i)+f(⌊ni⌋)
n 分治. 我们其实已经在 Example 1 估计 σ0(n) 时用过此技巧了.
Example 4 估计 σ1(n) 的渐进上界.
用 Proposition 1： σ1(n)=∑d∣nd≤∑i=1n⌊ni⌋≤nH(n)∈O(nlogn)
可以证明用 Proposition 2 不会得到更优的结果.
我们发现了一个有趣的事实：σ1(n) 和 σ0^(n) 的渐进上界均为 O(nlogn).
Example 5 估计 σz(n) 的渐进上界.
用 Proposition 2 和 p - 级数的性质：
σz(n)=∑d∣ndz≤∑i=1niz+⌊ni⌋z≤{2∑i=1n⌊ni⌋z≤2nz∑i=1ni−z=2nzPz(n)z≥02∑i=1niz=2P−z(n)z<0∈{2nzO(1)z>12nO(logn)z=12nzO(n1−z2)0≤z<12O(n1+z2)−1<z<02O(logn)z=−12O(1)z<−1={O(nz)z>1O(nlogn)z=1O(n1+z2)−1<z<1O(logn)z=−1O(1)z<−1
我们得到了一个相当优秀的渐进上界. 值得关注的是：
当 z=0 时，σ0(n)∈O(n12). 这与 Example 1 的结果一致.
当 z=12 时，σ12(n)∈O(n34)，即 ∑d∣nd∈O(n34). 洛谷 P4980 Polya 定理模板题[9] 的一种比较 trivial 的解法[10] 的时间复杂度证明就来源于此. 我们之后还会在整除分块与杜教筛中见到它.
另外，如果只使用 Proposition 1 ，−1<z<1 部分的渐进上界将只能估计至 O(n). 因此 Proposition 2 是更为优越的.
约数函数更复杂的上限与渐进估计可参考 Wikipedia[7].
2 整除分块
也被称为数论分块. 求 ∑i=1nf(i)g(⌊ni⌋) 我们按 d=⌊ni⌋ 分块求和： ∑dg(d)∑⌊ni⌋=df(i) 可以证明，对一指定的 d，满足 d=⌊ni⌋ 的 i 取遍一连续区间，故若 f 的前缀和能 O(1) 求出，块数量 #{⌊ni⌋}i=1n 即该算法的时间复杂度. 注意到当 i≤n 时，⌊ni⌋ 最多只有 ⌊n⌋ 种取值，而 i≥n 时，1≤⌊ni⌋≤n 表明其也最多只有 ⌊n⌋ 种取值. 因此整除分块的时间复杂度 T1(n)=#{⌊ni⌋}i=1n≤2n∈O(n)
方便起见，后文记 D(n)={⌊ni⌋}i=1n.
2.1 整除分块嵌套
将 Proposition 2 加强，我们有如下通用放缩：
Proposition 3 ∑d∣nf(d)≤∑d∈D(n)f(d)≤∑i=1nf(i)+f(⌊ni⌋)
LHS 成立的关键在于 {d:d∣n}⊂D(n)；而 RHS 的本质就是上述对整除分块块数量上界的估计.
注意到 Proposition 2 是 Example 5 证明的核心，而 Proposition 3 是 Proposition 2 的加强版，故仿造 Example 5 的证明，我们有
Example 6 令 Sz(n)=∑d∈D(n)dz 则前述 Example 5 中 σz(n) 的上界与渐进上界也同样适用于 Sz(n).
现在可以对嵌套整除分块 ∑i=1nf(i)∑j=1⌊ni⌋g(j)h(⌊nij⌋) 的时间复杂度 T2 做出估计了. 对 Example 6 取 z=12，立刻有 T2(n)=∑d∈D(n)T1(d)≤2∑d∈D(n)d=2S12(n)≤4nP12(n)∈O(n34)
我们还可以进一步归纳. 假定 ∀m≥0,∃zm:0≤zm<1,Tm(n)=O(nzm)，我们有
Tm+1(n)=∑d∈D(n)Tm(d)≤C∑d∈D(n)nzm=CSzm(n)∈O(n1+zm2)
因此 zm+1=1+zm2. 边界条件 z0=0，数列递推求得 zm=1−2−m，检验满足条件. 因此 m 重嵌套整除分块的时间复杂度 Tm(n)∈O(n1−2−m)
3 杜教筛
杜教筛可以以低于线性的时间复杂度求解某些数论函数的前缀和. 其思路并不复杂. 设 f 为一数论函数，我们希望快速求得其前缀和 f^(n)=∑i=1nf(i). 考虑数论函数 g 和 h=g∗f， h(n)=∑d∣ng(d)f(nd) 两端做前缀和得 h^(n)=∑i=1nh(i)=∑i=1n∑d∣ig(d)f(id)=∑d=1ng(d)∑i=1⌊nd⌋f(i)=∑d=1ng(d)f^(⌊nd⌋)=g(1)f^(n)+∑d=2ng(d)f^(⌊nd⌋) 因此 f^(n)=1g(1)(h^(n)−∑d=2ng(d)f^(⌊nd⌋))
故若 g、h 的前缀和可 O(1) 算得，根据上式整除分块即可递归地计算出 f 的前缀和.
下面分析算法的复杂度. 注意到 ⌊⌊ni⌋j⌋=⌊nij⌋ 故单轮递归涉及到的自变量均可表示为 d=⌊ni⌋ 的形式. 一个 f^(d) 做整除分块耗时 T1(d)，若采用记忆化递归，由上节分析，算法总时间复杂度为 ∑d∈D(n)T1(d)=T2(n)∈O(n34)
但我们还可以做得更好——考虑先用 O(K) 的时间复杂度线性筛出前 K 个 f(n) 并求前缀和，则递归求解时，d≤K 的 f^(d) 就无需再向下递归了. 为分析此类时间复杂度，对 Proposition 3 做最后一点扩展：
Proposition 4 ∑d∣nd>Kf(d)≤∑d∈D(n)d>Kf(d)≤∑K<i≤nf(i)+∑1≤i≤min{⌊nK⌋,n}f(⌊ni⌋)
特别的，当 K>n 时，有
∑d∣nd>Kf(d)≤∑d∈D(n)d>Kf(d)≤∑1≤i≤⌊nK⌋f(⌊ni⌋)
故用 Proposition 4 ，当 K>n 时，算法在递归部分的时间复杂度降低为
∑d∈D(n)d>KT1(d)=∑1≤i≤⌊nK⌋T1(⌊ni⌋)≤∑1≤i≤⌊nK⌋Cni=Cn∑1≤i≤⌊nK⌋i−12=CnP12(⌊nK⌋)∈nO((nK)12)⊂O(nK−12)
总时间复杂度为 O(K)+O(nK−12)
为最小化时间复杂度，取 K=n23，得到最优时间复杂度 O(n23).
这部分的时间复杂度证明主要参考了文章[11].
References
1. Abuse of notation - wikipedia. (n.d.). https://en.wikipedia.org/wiki/Abuse_of_notation#Function_notation.
2. Graham, R. L., Knuth, D. E., & Patashnik, O. (1994). Concrete mathematics: A foundation for computer science (second, pp. 443–449). Addison-Wesley.
3. Big o notation - wikipedia # family of bachmann–landau notations. (n.d.). https://en.wikipedia.org/wiki/Big_O_notation#Family_of_Bachmann%E2%80%93Landau_notations.
4. 复杂度 - OI wiki. (n.d.). https://oi-wiki.org/basic/complexity/#%E6%B8%90%E8%BF%9B%E7%AC%A6%E5%8F%B7%E7%9A%84%E5%AE%9A%E4%B9%89.
5. Knuth, D. E. (1976). Big omicron and big omega and big theta. SIGACT News, 8(2), 18–24. https://doi.org/10.1145/1008328.1008329
6. Graham, R. L., Knuth, D. E., & Patashnik, O. (1994). Concrete mathematics: A foundation for computer science (second, pp. 47–56). Addison-Wesley.
7. Divisor function - wikipedia # growth_rate. (n.d.). https://en.wikipedia.org/wiki/Divisor_function#Growth_rate.
8. sun123zxy. (2020). sun123zxy’s blog - 原创OI题目 GCD卷积 problem and solution. https://blog.sun123zxy.top/posts/20201206-gcdconv/.
9. P4980 【模板】pólya 定理 - 洛谷 | 计算机科学教育新生态. (n.d.). https://www.luogu.com.cn/problem/P4980.
10. sun123zxy. (2020). sun123zxy’s blog - 等价类计数：Burnside引理 & Polya定理. http://blog.sun123zxy.top/posts/20200321-burnside/#s-4.3.
11. Ander. (2022). 杜教筛. https://zhuanlan.zhihu.com/p/521699400.
