深刻认识到自己是彩笔的这一事实。
Powerful Number：所有质因子次数都不为 1 的数。
有一个结论是小于等于 $n$ 的 Powerful Number 个数是 $O(\sqrt n)$ 的。
可以用来求积性函数前缀和。假设我们正在求 $F$ 的前缀和，然后现在我们有一个积性函数 $G$ 满足 $\forall p \in \mathbb{P},F§=G§$，且 $G$ 的前缀和可次线性求。然后我们现在有个函数 $H$ 满足 $F=G*H$（狄利克雷卷积），显然 $H$ 也是积性函数。然后根据杜教筛那套推法，我们有：
$$
\sum_{i=1}^nF(n)=\sum_{i=1}^nH(n)\sum_{j=1}^{\left\lfloor\frac ni\right\rfloor}G(j)\tag 1
$$
我们观察到：
$$
F(p)=G(1)H(p)+G(p)H(1)
$$
又因为 $F§=G§$，则我们有 $H§=0$。然后容易发现 $H$ 只在 Powerful Number 处值非 0，在非 Powerful Number 的地方由于 $H$ 是积性函数则它被分解为了 $H(\frac np)\cdot H§=H(\frac np)\cdot 0=0$。于是我们可以直接 dfs 出所有的 Powerful Number 并代到 $(1)$ 式里面算。现在问题是怎么算 $H$。我们发现：
$$
\begin{aligned}
F(n)&=\sum_{d|n}G(d)\cdot H\left(\frac nd\right)\\
F(n)&=H(n)+\sum_{d|n,d\ne1}G(d)\cdot H\left(\frac nd\right)\\
H(x)&=F(n)-\sum_{d|n,d\ne1}G(d)\cdot H(\frac nd)
\end{aligned}
$$
实际上我们只需要算 $H(p^k)$，毕竟我们有个枚举 Powerful Number 的过程，乘起来就行。
$$
H(p^k)=F(p^k)-\sum_{j=1}^kG(p^j)H(p^{k-j})
$$
好算。就是写着不怎么方便
P5325 【模板】Min_25 筛
$f$ 是积性函数，且 $f(p^k)=p^k(p^k-1)$。求
$$
\sum_{i=1}^nf(i)
$$
对 $10^9+7$ 取模。
$1\le n\le 10^{10}$
我们设 $G(x)=x\varphi(x)$（注意 $G(x)=x(x-1)$ 不是积性函数），然后按照上面一样做就好了。$G(x)$ 的前缀和可以用杜教筛。
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
|
#include <cmath>
#include <cstdio>
typedef long long qe;
const int $L = 5624000;
const int $p = $L;
const int $pk = 35;
const int mod = 1000000007;
const int inv6 = 166666668;
inline int add(int x, int y) { return (x += y) >= mod? x - mod: x; }
inline void Add(int &x, int y) { if ((x += y) >= mod) x -= mod; }
inline int sub(int x, int y) { return (x -= y) < 0? x + mod: x; }
inline void Sub(int &x, int y) { if ((x -= y) < 0) x += mod; }
inline int div2(int x) { return ((x & 1)? x + mod: x) >> 1; }
qe n;
int L, ans;
bool co[$L];
int pr[$p], G[$L], Hp[$L][$pk];
inline void sieve() {
G[1] = 1;
for (int i = 2; i <= L; ++i) {
if (!co[i]) { pr[++pr[0]] = i; G[i] = i - 1; }
for (int j = 1, k; j <= pr[0] && (k = i * pr[j]) <= L; ++j) {
co[k] = 1;
if (i % pr[j]) G[k] = (qe)G[i] * (pr[j] - 1) % mod;
else { G[k] = (qe)G[i] * pr[j] % mod; break; }
}
G[i] = (G[i - 1] + (qe)G[i] * i) % mod;
}
qe tmp[$pk];
for (int i = 1; i <= pr[0]; ++i) {
qe w = (qe)pr[i] * pr[i]; if (w > n) break;
const int del = w % mod;
auto &cur = Hp[i];
cur[0] = tmp[0] = 1;
tmp[1] = (qe)pr[i] * (pr[i] - 1) % mod;
for (int j = 2; w <= n; ++j, w *= pr[i]) {
tmp[j] = (qe)tmp[j - 1] * del % mod;
cur[j] = (qe)w % mod * (w % mod - 1) % mod;
for (int k = 1; k <= j; ++k)
Sub(cur[j], (qe)tmp[k] * cur[j - k] % mod);
}
}
}
inline int sum1(qe x) { x %= mod; return x * (x + 1) / 2 % mod; }
inline int sum2(qe x) { x %= mod; return x * (x + 1) % mod * (x * 2 + 1) % mod * inv6 % mod; }
int S(qe n) {
static int mem[2005];
if (n <= L) return G[n];
const int index = ::n / n;
if (mem[index]) return mem[index];
int ans = sum2(n), lst = 1;
for (qe l = 2, r; l <= n; l = r + 1) {
const qe t = n / l; r = n / t;
const int cur = sum1(r);
Sub(ans, (qe)sub(cur, lst) * S(t) % mod);
lst = cur;
}
return mem[index] = ans;
}
void dfs(qe x, int H, int lst) {
Add(ans, (qe)H * S(n / x) % mod);
const qe lim = n / x;
for (int i = lst + 1; i <= pr[0]; ++i) {
if ((qe)pr[i] * pr[i] > lim) return;
qe cur = (qe)x * pr[i] * pr[i];
for (int j = 2; cur <= n; ++j, cur *= pr[i])
dfs(cur, (qe)H * Hp[i][j] % mod, i);
}
}
int main() {
scanf("%lld", &n);
L = std::max((int)std::pow(n, 0.675) + 5, 100);
sieve();
dfs(1, 1, 0);
printf("%d\n", ans);
}
|
