积性函数的线性筛法总结
内容
原理
线性筛质数原理就是找到一个质数来更新范围内的所有以这个质数为最小质因子的合数。
积性函数性质$f(ab)=f(a)f(b)\;\;\;\;(a,b)=1$
那么我们可以想到线筛过程中可以把积性函数筛出来，要讨论一下：
首先是质数，一般积性函数的质数都有特殊的性质，因为它不能被两个数乘出来。。($\mu(p) = -1, \varphi(p) = p – 1$)
然后就是筛的时候的讨论，比如筛$\varphi$的时候
void init() {
phi[1] = 1;
for (int i = 2; i <= N; ++i) {
if (!p[i])
pr[++pr[0]] = i, phi[i] = i - 1;
for (int j = 1; j <= pr[0] && i * pr[j] <= N; ++j) {
p[i * pr[j]] = 1;
if (i % pr[j] == 0) {
phi[i * pr[j]] = phi[i] * pr[j];
break;
}
phi[i * pr[j]] = phi[i] * phi[pr[j]];
}
}
}
当i % pr[j] != 0的时候，这两个数就是互质的，那么可以直接相乘。
当i % pr[j] == 0的时候，每个积性函数都不是相同的，都需要化下式子。比如$\varphi$来说，可以得到$\varphi(ab) = \varphi(a)b\;\;(a \%b=0)$
其实讨论的内容就是知道了$f(p_i^{k_i})$和$p_i$，要计算$f(p_i^{k_i+1})$。稍微的把这两个函数拆开，相减/相除一下再讨论下贡献就ok了。。实在想不出来就打表找规律。。。。
例题
bzoj 1968
题意
求$$\sum_{i = 1}^nd(i)$$
题解
设$x = p_1^{k_1} \times p_2^{k_2} \times \cdots \times p_n^{k_n}$
那么$d(x) = (k_1 + 1)\times (k_2+1)\times \cdots \times(k_n+1)$，是按照组合的思路求出来，下面的筛法也是一样的思路
这里筛的时候还需要维护一个最小质因子指数函数$t$
假如是质数，约数和为2（1和它自己），t函数为1
假如互质，直接相乘，t函数为1（现在这个是最小的质因子）
假如不互质，考虑这一个质因子会做出多少贡献，算下得出d[i * pr[j]] = d[i] / (t[i] + 1) * (t[i] + 2)
这样就可以筛了：
#include <cstdio>
#include <algorithm>
using namespace std;
const int N = 1000000;
bool p[N + 10];
int n, t[N + 10], pr[N + 10];
long long k, ans, d[N + 10];
void init(int n) {
d[1] = 1, t[1] = 0;
for (int i = 2; i <= n; ++i) {
if (!p[i])
pr[++pr[0]] = i, d[i] = 2, t[i] = 1;
for (int j = 1; j <= pr[0] && i * pr[j] <= n; ++j) {
p[i * pr[j]] = 1;
if (i % pr[j] == 0) {
d[i * pr[j]] = d[i] / (t[i] + 1) * (t[i] + 2);
t[i * pr[j]] = t[i] + 1;
break;
}
d[i * pr[j]] = d[i] * d[pr[j]];
t[i * pr[j]] = 1;
}
}
}
main() {
scanf("%d", &n);
init(n);
for (int i = 1; i <= n; ++i)
ans += d[i];
printf("%lld\n", ans);
}
BZOJ2813
题意
求约数平方和$r(x)$
题解
假如x是质数，那么$r(x) = x^2 + 1^2$
假如互质的话，用原来的所有的约数乘上新的质数，那么新的平方和就是$r(x) * p^2$，化个式子得到$r(x * p) = r(x)*(p^2+1)=r(x)*r(p)$，也就证明了这是个积性函数
假如不互质，那么不能直接乘出来构造新的约数，这样会乘出来相同的约数。那么我们可以再维护一个函数$g$表示去掉最小质因子的部分，那么可以先用所有约数乘上$p$构造出还有$p$的约数，再加上$(g(rx))$累加所有不含它的约数，那就是$r(x) = r(g(x)) + r(x) * p^2$。。
#include <cstdio>
#include <algorithm>
using namespace std;
const int N = 10000000, mod = 1000000007;
bool p[N + 10];
int n, T, a, b, c, t[N + 10], pr[N + 10];
long long k, ans1, ans2, d[N + 10], r[N + 10], g[N + 10];
void init(int n) {
r[1] = d[1] = 1;
for (int i = 2; i <= n; ++i) {
if (!p[i])
pr[++pr[0]] = i, d[i] = 2, t[i] = g[i] = 1,
r[i] = (1ll * i * i + 1ll) % mod;
for (int j = 1; j <= pr[0] && i * pr[j] <= n; ++j) {
p[i * pr[j]] = 1;
if (i % pr[j] == 0) {
d[i * pr[j]] = (1ll * d[i] / (t[i] + 1ll) * (t[i] + 2ll)) % mod;
t[i * pr[j]] = t[i] + 1;
g[i * pr[j]] = g[i];
r[i * pr[j]] = (r[g[i]] + 1ll * r[i] * pr[j] % mod * pr[j] % mod) % mod;
break;
}
d[i * pr[j]] = d[i] * d[pr[j]];
t[i * pr[j]] = 1;
g[i * pr[j]] = i % mod;
r[i * pr[j]] = (r[i] + 1ll * r[i] * pr[j] % mod * pr[j] % mod) % mod;
}
}
}
main() {
scanf("%d%d%d%d%d", &T, &n, &a, &b, &c);
init(c);
while (T--) {
long long tmp;
tmp = d[n] + ((n & 1) ? 1 : 0);
ans1 = (ans1 + tmp) % mod;
tmp = r[n] + ((n & 1) ? 4 : 0);
ans2 = (ans2 + tmp) % mod;
n = (1ll * n * a + b) % c + 1;
}
printf("%lld\n%lld\n", ans1, ans2);
}
BZOJ4407
题意
$$\sum_{i = 1}^n\sum_{j = 1}^m(i,j)^k$$
题解
$$
\large
\begin{aligned}
&\;\;\;\;\; \sum_{i = 1}^n\sum_{j = 1}^m (i,j)^k \\
&= \sum_d d^k \sum_{i = 1}^n\sum_{j = 1}^m [(i,j)=b] \\
&= \sum_d d^k \sum_{i = 1}^{n \over d}\sum_{j = 1}^{m \over d}[(i,j)=1] \\
&= \sum_d d^k \sum_{i = 1}^{n \over d}\sum_{j = 1}^{m \over d}\sum_{p|(i,j)} \mu(p) \\
&= \sum_d d^k \sum_p \mu(p) {n\over dp} {m\over dp} \\
&= \sum_d d^k \sum_p \mu(p) {n\over T} {m\over T} \\
&= \sum_{T = 1}^{min(n,m)} {n\over T} {m\over T} \sum_{d|T}d^k \mu({T\over d})
\end{aligned}
$$
关键就在于这个函数
$$f(T)=\sum_{d|T}d^k\mu({T\over d})$$
这个是两个积性函数的狄利克雷卷积，它肯定是一个积性函数
假如是质数，那么$f(i) = i ^k – 1$
假如互质，直接乘
假如不互质，对$f(i∗p)$有影响的$t$一定是与对$f(i)$有影响的$t$中，因为其余的因数都至少包含$p^2$，$\mu=0$。而此时$d$增加了$p$倍，故$f(i∗p)=f(i)∗p^k$。
#include <cstdio>
#include <algorithm>
using namespace std;
const int N = 5000000, mod = 1000000007;
bool p[N + 10];
int T, pr[N + 10], n, m;
long long k, f[N + 10], ans;
long long fast\_pow(long long a, int p)
{
long long ans=1;
for (; p; p >>= 1, a = a * a % mod)
if (p & 1) ans = ans * a % mod;
return ans;
}
void init() {
f[1] = 1;
for (int i = 2; i <= N; ++i) {
if (!p[i])
pr[++pr[0]] = i, f[i] = (fast\_pow(i, k) - 1 + mod) % mod;
for (int j = 1; j <= pr[0] && i * pr[j] <= N; ++j) {
p[i * pr[j]] = 1;
if (i % pr[j] == 0) {
f[i * pr[j]] = f[i] * fast\_pow(pr[j], k) % mod;
break;
}
f[i * pr[j]] = f[i] * f[pr[j]] % mod;
}
}
for (int i = 2; i <= N; ++i) f[i] = (f[i] + f[i - 1]) % mod;
}
main() {
scanf("%d%d", &T, &k);
init();
while (T--) {
scanf("%d%d", &n, &m);
if (n > m) swap(n, m);
ans = 0;
for (int i = 1, las; i <= n; i = las + 1) {
las = min(n / (n / i), m / (m / i));
ans += (f[las] - f[i - 1]) * (n / i) % mod * (m / i) % mod;
ans = (ans % mod + mod) % mod;
}
printf("%lld\n", ans);
}
}
还有一个例题就是bzoj2154