数论函数求和（下）

发表于 更新于 分类于 ， 阅读次数：

本文共约 10700 字。

Key Words：杜教筛，狄利克雷双曲线法， 贝尔级数，Powerful Number，拟合法，标准积性函数求和，扩展埃氏筛（Min_25筛，洲阁筛），多元积性函数

杜教筛

在接下来的五小节中，我们将注意力集中到一类核心问题上：积性函数前缀和。为了便于讨论，对数论函数 ，记 。

  • 基本和组（又称“块筛”）

    关于 的“基本和组”定义为 的一系列 的值，即我们对 整除分块时所需的值。

例： 关于 的基本和组包含 。

  • 整除定理

结论一： 经历若干次整除，其结果仍然在 内。处理了 关于 的基本和组，就能得到关于 的基本和组。该结论保证了一些递归算法的正确性。

结论二：整除分块的段的边界恰是集合 。注意到 即可证明。

杜教筛利用迪利克雷卷积的性质，求数论函数的基本和组。这一算法最初由杜瑜皓（xydyh） 引入竞赛界，故名。

乘积式

数论函数 满足 ，给出 关于 的基本和组，求 的基本和组。

得到了可整除分块的形式，可以发现，狄利克雷卷积和整除有着自然的关系。只需 的基本和组即可完成整除分块的计算。

若使用整除分块对各个 计算，考虑最大的 个 ，也就是让 ，可得复杂度为 。

若 都是积性函数，则 也是积性函数，可以线性筛。预处理 以内的前缀和，这样就只需要对 即 进行计算。总复杂度为 取 可得最优复杂度 。

除式

数论函数 满足 ，给出 关于 的基本和组，求 的基本和组。

也得到了可整除分块的形式，略有不同的是需要递归计算。复杂度仍是 。

狄利克雷双曲线法*

可以认为是杜教筛的另一种形态。下面推导乘积式，除式大同小异，略去。

给出 求 时：

即是双曲线 下方（含线上）整点对应的函数值的和。

如图，对 的部分（左）、 的部分（中）分别求和，再减去公共部分（右），可得 这容易 计算。该公式中，对称性体现得更加明显，能直接看到基本和组的形式。

例如，当我们计算 的前缀和的时候，有经典结论 用该式计算会比直接乘除分块常数更优。

一些具体函数的基本和组

杜教筛只能起到桥梁的作用，即在不同函数的基本和组之间“转化”，我们必须有最初的“原料”。

对于两个熟悉的完全积性函数， 的前缀和是显然的， 一类的前缀和是自然数幂和，是经典问题（当 较小时有几个常用的小公式）。这两者是我们最初的原料，依据卷积关系，可以用杜教筛进一步合成出其他种类繁多的函数。

对于我们的三个老朋友，有

以上三类常见函数都能被完全积性函数简单表出，也就能杜教筛。下面给出 Luogu4213 的代码。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
#include <algorithm>
#include <cstdio>
#include <map>

using ll = long long;
using uint = unsigned;
const int T = 8e6, MaxT = T+5; 

int p[MaxT/8], totP;
ll mu[MaxT], phi[MaxT];
void sieve(int n) {
    mu[1] = phi[1] = 1;
    for (int i=2; i<=n; i++){
        if (!phi[i]) {
            p[++totP] = i;
            phi[i] = i-1;
            mu[i] = -1;
        }
        for (int j=1; j<=totP; j++){
            int t = i*p[j];
            if (t>n) break;
            if (i%p[j] == 0) {
                phi[t] = phi[i]*p[j];
                break;
            }
            mu[t] = -mu[i];
            phi[t] = phi[i]*(p[j]-1);
        }
    }
    for (int i=2; i<=n; i++){
        mu[i] += mu[i-1];
        phi[i] += phi[i-1];
    }
}

std::map<int, ll> mu2, phi2; 

ll Smu(uint n) {
    if (n<=T) return mu[n];
    if (mu2.find(n) != mu2.end())
        return mu2[n];
    ll ret = 1;
    uint l = 2;
    while(l<=n) {
      	int r = n/(n/l);
      	ret -= Smu(n/l)*(r-l+1);
      	l = r+1;
    }
    return mu2[n] = ret;
}
ll Sphi(uint n) {
    if (n<=T) return phi[n];
    if (phi2.find(n) != phi2.end())
        return phi2[n];
    ll ret = 1ll*n*(n+1)/2;
    uint l = 2;
    while(l<=n) {
      	int r = n/(n/l);
      	ret -= Sphi(n/l)*(r-l+1);
      	l = r+1;
    }
    return phi2[n] = ret;
}

int main() {
    int Q;
    scanf("%d", &Q);
    sieve(T);
    while(Q--) {
        uint N;
        scanf("%u", &N);
        printf("%lld %lld\n", Sphi(N), Smu(N));
    }
    return 0;
}
  • 引理：当 是完全积性函数时，有 。

证明：

接着就能处理下面的函数：

贝尔级数

高维点表示与狄利克雷卷积

回顾第五小节的内容，记 为从小到大第 个质数，根据唯一分解定理，我们可以将任意自然数 写成 这样就能将 对应到无穷序列 ，其中每个数对应一个素数的幂次。

两个自然数相乘时，对应的幂次序列相加，即 因此，如果对每个素数建立一个维度，狄利克雷卷积（乘法卷积）可以视为无穷维加法卷积。也就是说，对于数论函数 ，构造无穷元生成函数 狄利克雷卷积对应该生成函数的乘法。

积性函数的生成函数

对于积性函数 ，它的生成函数可以因式分解：

其中质数 对应的因式为 ，记为 。

在第一节中我们提到一组结论：两个积性函数的狄利克雷卷积仍是积性函数；积性函数的逆仍是积性函数。当时的证明是归纳法，现在，生成函数可以提供一个更简单直接的证明。

积性函数 生成函数能因式分解。考虑积性函数 ，和 ，它们的生成函数是

显然， 也是能因式分解的，故 是积性函数。

对于乘法逆 的情况，证明是类似的，只需要注意到 即可。（，故逆多项式一定存在）

贝尔级数

回顾上一小节的讨论，杜教筛的使用依赖于卷积关系，如何系统地发现积性函数之间的卷积关系呢？

  • 贝尔级数

    对于积性函数 ，定义其关于质数 的贝尔级数为:

    也就是在积性函数的“无穷维生成函数”中抽出一个因式。

  • 定理

    两个积性函数相卷积，其贝尔级数相乘。

贝尔级数把狄利克雷卷积转化为了幂级数卷积。这样，我们可以用幂级数的代数手段来研究积性函数之间的卷积关系。

常见积性函数的贝尔级数

积性函数
贝尔级数

根据这些贝尔级数，可以看出一些数论函数的卷积关系：

  • 将 翻译回数论函数，可以自然地得出 的定义式。

  • 等价于 。

  • 等价于 。

  • 等价于 。

再来看两个比较陌生的积性函数：

  • 的不同质因子数
  • 刘维尔函数：的质因子次数和

写出它们的贝尔级数 可以发现 ，。如果单纯从组合意义出发，很难发现这些卷积关系。

值得留意的是， 甚至是完全积性函数。事实上，我们可以批量制造除 之外的完全积性函数。

  • 完全积性函数的贝尔级数

    对于完全积性函数 ，其贝尔级数

证明：

任取 就能创造相应的完全积性函数。

  • 点积对贝尔级数的影响

    点积完全积性函数 ，对贝尔级数的影响相当于将 代换成 。

证明：写出 的贝尔级数： 恰是 的贝尔级数将 代换成 的结果。

由此，可以便捷地获得下列函数的贝尔级数：

  • 例题 4.2.12. （证明题）

    求证

Solution：

记 ，考虑差分，记 ，则 故 。它的贝尔级数是 翻译为数论函数，知 不难发现 是完全平方数，则 是完全平方数

Powerful Number

方便起见，下文将 powerful number 简称为 PN。

  • Powerful Number

    PN 是所有质因子次数大于等于 的数。

  • 定理

    以内 PN 的个数为 。

该定理是 PN 的力量来源。

证明： 对于一个 PN ，记 此时有 ，这构成一个 到数对 的单射。

计数满足 的数对 的数目，这不会超过 PN 的数目。枚举 考虑满足条件的 的个数，得 ，积分得到 ，证毕。

拟合法

积性函数 对所有质数 满足 ，则称 在素数处互相拟合。现已知 的基本和组，求 的基本和组。

构造函数 满足 (狄利克雷除法)。观察到 ，由于 ，可以得到 。又因为积性，所以 有值的地方都是 PN。

根据杜教筛结论，由 快进到

由于有值的 只有 个，可以直接搜索 PN 进行求和（在搜索的同时得出其积性分解）。为了求出 函数值，需要预处理 以内的质数幂的取值，可以根据 的函数值，在每个质数处分别（暴力）计算幂级数逆，复杂度不超过 。

【例题 4.2.14】给出了搜索 PN 的具体代码。

至此，类似杜教筛，求 的基本和组可以做到 。

由于 PN 的优秀性质，还存在更优的解法。

  • 定理： 关于 的基本和组只有 个本质不同的值。

证明： 在 中，若 ，则 ，这部分只有 个有值的 ，前缀和也就只会变化 次。若 ，显然只有 个。

这样，我们求 的复杂度可以降为 ，经过合理的前缀预处理，总复杂度降为 。

综上，可以总结出拟合法的一般过程：

求解 的基本和组时，构造 满足

  • （即“拟合”）
  • 的基本和组易求。

先求出 的基本和组，然后用 Powerful Number 将 的基本和组修正为 的基本和组。

  • 的前缀和

在第三节中，我们借助 求出了 的前缀和。这个结论需要一定的观察力，下面我们展示如何用拟合法机械地完成这一工作。

的贝尔级数是 ，用 去拟合，其贝尔级数是 。那么，令

很容易看出来 是平方数。套上杜教筛公式 是平方数

这个结果和第三节中得到的相同。

  • 例题 4.2.14. 【模板】Min_25筛

    积性函数 满足 ，求 。

    。

Solution： 注意到 ，构造 ，则 ，且 的基本和组容易杜教筛求出。然后用 PN 计算一个 即可。总复杂度 。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
#include <cstdio>
#include <cmath>
#include <map>

using ll = long long;
const int 
    MaxT = 4.65e6,
    mod = 1000000007,
    inv2 = (mod+1)/2,
    inv6 = (mod+1)/6;

int p[MaxT/8], totP, Hpk[MaxT/8][36], G[MaxT];
int log(int p, ll n) {
    int c = 0;
    while(n) { n/=p; c++; }
    return c;
}
void div(int *f, int *g, int *h, int n) {
    h[0] = 1;
    for (int i=1; i<n; i++) {
        int buf = 0;
        for (int k=1; k<=i; k++)
            buf = (buf + 1ll*g[k]*h[i-k]) %mod;
        h[i] = (f[i]-buf+mod) %mod;
    }
}
void sieve(int T, ll N) {
    G[1]=1;
    for (int i=2; i<=T; i++) {
        if (!G[i]) {
            p[++totP] = i;
            G[i] = i-1;
        }
        for (int j=1; j<=totP; j++) {
            int t = i*p[j];
            if (t>T) break;
            if (i%p[j]==0) {
                G[t] = G[i]*p[j];
                break;
            }
            G[t] = G[i]*(p[j]-1);
        }
    }
    for (int i=1; i<=T; i++)
        G[i] = (G[i-1]+1ll*i*G[i]) %mod;
    for (int i=1; i<=totP; i++) {
        int m = log(p[i], N);
        static int f[36], g[36];
        f[0] = g[0] = 1;
        int pk = 1;
        for (int k=1; k<m; k++) {
            g[k] = 1ll*pk*pk%mod*p[i]%mod*(p[i]-1)%mod;
            pk = 1ll*pk*p[i] %mod;
            f[k] = 1ll*pk*(pk-1) %mod;
        }
        div(f, g, Hpk[i], m);
    }
}

int T;
std::map<ll, int> G2;
int SG(ll n) {
    if (n<=T) return G[n];
    if (G2.find(n) != G2.end())
        return G2[n];
    int ret = 0; 
    ll l = 2;
    while(l<=n) {
        ll r = n/(n/l);
        int sumId = ((r-l+1)%mod)*((l+r)%mod) %mod;
        ret = (ret + 1ll*SG(n/l)*sumId) %mod;
        l = r+1;
    }
    int n0 = n%mod;
    ret = (1ll*n0*(n0+1)%mod*(n0*2+1)%mod*inv6-1ll*ret*inv2) %mod;
    ret = (ret+mod)%mod;
    return G2[n] = (ret+mod) %mod;
}

int ans;
void dfs(ll N, ll d, int h, int t) {
    ans = (ans + 1ll*h*SG(N/d)) %mod;
    for (int i=t; i<=totP; i++) {
        if (d > N/(1ll*p[i]*p[i])) return ;
        int k = 2;
        ll newd = d*p[i]*p[i];
        while (newd<=N) {
            dfs(N, newd, 1ll*h*Hpk[i][k]%mod, i+1);
            newd *= p[i];
            k++;
        }
    } 
}
int main() {
    ll N;
    scanf("%lld", &N);
    T = pow(N, 2./3.);
    sieve(T, N);
    dfs(N, 1, 1, 1);
    printf("%d\n", ans);
    return 0;
}

  • 例题 4.2.15. 「单走一个」

    给定常数 ，有积性函数 ，求 。

    ，。

Solution： 构造函数 ，显然有 。令 。

容易通过插值 求出 的块筛。单个 可以整除分块 计算。

对于单个

对于每个 PN 整除分块计算对应的 ，复杂度为 。

扩展埃氏筛

埃氏筛的思想：从小到大枚举质数，筛去它们的倍数，最终就能得到所有质数。利用该思想，可以开辟解决积性函数求和问题的另一条道路。

埃氏筛的集合递推

约定 指从小到大第 个质数（编号从 开始）， 为 的最小质因子，特殊地，。

记 为删除 的倍数之后 以内剩余的数的集合。即 注意， 会被删除，但 不会被删除。

记 ，即第 轮中删除的数。按照定义有

可以理解为“不含 的质因数”且“至少含一个”，于是有

不难验证，一个“不超过 、不含小于等于 的质因数”的数，乘 之后是“不超过 、不含小于等于 的质因数，至少含一个”，反之也成立，故有双射。

则有集合的递推式

我们已经完成了主要的工作。目前的 会删除 ，这和埃氏筛有一点小差异。尝试去除这一差异，令 熟知埃氏筛只需用 以内的质数来筛除，就能正确地得到 以内所有的质数。因此，当 时，“筛除 的倍数”实际上不会筛掉任何数，即 对于 ，尝试由 的递推式推导 的递推式 综上，我们成功地将埃氏筛刻画为集合递推：

质数 次方和

给出 ，求 以内质数的 次方和。

首先将 内质数筛出来，记个数为 ， 就是 以内素数的集合。

记 即为答案。

根据集合递推式，可以写出 的递推式： 边界：，是自然数幂和。

递推中会用到 ，容易预处理。或者发现 ，可以用 来代替 ，这样就不需要 了。

注意到递推中用到 ，我们需要维护整个基本和组，即对每个 维护 共 个。边界值容易 求出。

质数个数 ，基本和组位置数 ，总状态数为 。若直接用递推式计算，复杂度为 。

使用滚动数组，当 时 ，这部分可以直接不处理。统计需要处理的转移次数：

这给出了一个高效且易于实现的算法。

标准积性函数求和

  • 积性函数 满足 为关于 的多项式（项数较少）， 可以 计算，求 或 关于 的基本和组。这两类问题统称为“标准积性函数求和”。

显然，我们可以用“质数 次方和”处理质数的贡献，接下来处理合数。

在埃氏筛中，我们逐渐将最小质因子为 的数筛去，最终留下质数。接下来我们将这一过程倒过来，逐步恢复最小质因子为 的合数。

  • 洲阁筛

记 为 以内最小素因子大于等于 的数的集合，特别地，其中包括 。注意，其中也包括大于等于 的素数。

记 为 以内最小素因子恰为 的数的集合，则 枚举 的次数，可得 则有 的递推式

接下来将质数从 中去除，令 当 时，“恢复除质数外，最小素因子为 的数”实际上不恢复任何数，即 当 时，有递推式 边界：。

记 则 其中用到的 涉及整除位置，容易用“质数 次方和”预处理。

滚动数组并忽略不变的转移，总复杂度 。

该方法由任之洲在他的集训队论文中提出，故名。

  • Min_25 筛

在只求单个和时适用，可以视为洲阁筛的变种。

枚举最小质因子 ，再枚举最小质因子的次数 ，可得变种递推式 可以发现 时 只能是 ，故

记 ，有 答案是 。需要预处理质数处的和。

不记忆化暴力搜索即可。边界是当 时 ，是质数处的和。

实践表明其在 时效率良好（快于洲阁筛，但理论复杂度没有好的上界）。

这一搜索做法由 Min_25 在他的博客中提出，故名为 Min_25 筛。

多元积性函数

在最后一节，我们将积性函数的概念推广到多元的情况，并借助多元的观点更直观地解决一些问题。

基本概念

  • 二元积性函数 若函数 满足 ，则称为二元积性函数。

  • 积性分解 将 分解为 ，，对于二元积性函数

  • 二元狄利克雷卷积

    两个二元数论函数 的狄利克雷卷积为一个新的二元数论函数，记作 。其满足

类似一元的情况，有如下性质成立：

  • 积性函数的卷积是积性函数。

  • 积性函数的逆是积性函数。

  • 二元贝尔级数

    在素数 处观测积性函数 ，定义贝尔级数为二元形式幂级数

    类似一元的情况，两个函数卷积，对应的贝尔级数相乘。

积性函数矩阵和问题

  • 给出一元积性函数 ，给出 ，求

令 ，可以发现 是二元积性函数。问题转化为二元积性函数的求和。

为了便于书写，对于二元数论函数 ，记 ，即二维前缀和。

  • 低阶拟合法

模仿 Powerful Number 的思路，构造一个数论函数 使得 在素数幂次较低时相等，然后做除法 ， 的非零元素有希望较少。

计算 的矩阵和时，推导类似杜教筛

如果 稀疏，且 的矩阵和易求，我们就得到了一个不错的算法。

不难想到构造 ，这样 ，只要预处理 的基本和组，就能 查询需要的 。

此外， 在低阶处相等：可以发现 。

为了更直观地说明这一点，构造贝尔级数。

由于 的第一行和第一列系数相等，根据幂级数知识， 的第一行和第一列除了 之外，其余系数都是 。

换句话说，对于任意 ，。

  • 密度的分析

    和 Powerful Number 类似，我们很有希望相信 是稀疏的。

    记 ，枚举量是

    其中 ，可以证明 是积性函数（具体过程略）。观察质数幂处的取值

    其中用到了 。

    • 引理：对于数论函数 ，找出最小的 满足 。记 为 的最大值，其余 ，则 。

      其证明超出了本文的范围。

    根据引理，，总枚举量为 。

至此，单组询问可以 解决。

  • 多组询问的一般情况

Luogu3327 的做法利用了 的特殊性（ 只有对角线有值），对于一般的积性函数

设定阈值 ，在询问时，对于 的 ，暴力枚举并计算贡献，复杂度 。

对于 的 ，它对询问 的贡献是 ，分为 个本质不同的矩形。将所有询问离线，用 - 的二维偏序处理这部分的贡献。

二维偏序的总点数不超过

其中 可以按密度估计为 。

总复杂度 ，取 可得复杂度为 。

  • 更快的单组询问：进一步拟合

对于 ， 还不够快，究其根本， 的密度仍然太大。尝试进一步拟合，以降低 的密度。

构造 ，其余为零。令 ，则 ，求答案时

其中

可以整除分块或枚举 计算。（需要用一些筛法预处理 和 的基本和组）

  • 密度的分析

    根据贝尔级数/幂级数知识，对于 我们有 ，不妨设其余都非零。

    记 ，可以证明 是积性函数且 。

    根据之前的引理，，可得枚举量为 。

    其实可以证明枚举量是 。为了更精确的分析，放弃放缩为双曲区域 ，直接考虑二元积性函数 在 中的和。

    （双曲区域中的整点数目比矩形区域要大 倍，这导致数量级估计过大）

    直接分析 较为困难，考虑将其拆成两部分。令 满足 ，其余为 。令 ，易知 。这样做， 占据了主要项，但结构很简单， 的结构复杂但占据了次要项。我们将 拆成了两个更容易分析的函数。

    • 密度的分析

      先考察 的矩阵和， 时总是可以写作 ，枚举 进行计数

    接下来，由

    其中 为满足 的非零 个数。

    注意到 ，根据引理，。 可以根据密度估计为 。

    综上， 的枚举量是 。

  • 对每个非零 计算 的开销

    复杂度是

    利用 ，由对称性，只考虑 。

综上，复杂度可以做到 。可前往 Uoj#885 进行测试。

可以进一步优化到 ，具体方法超出了本文的范围。

习题

见闻*