Problem. 假定 (f(n)) 为某个积性函数且可以快速筛出前缀和，我们令 (g(n)) 为另一个积性函数，满足：
[g(p^k)=f(p^k)+1 ]现在请求出 (g(n)) 的前缀和。
sol.
对于 (f) 是完全积性函数的情形：
我们首先可以组合意义一发得到：
[g(n)=sum_{dmid n}[gcd(d,frac nd)=1]f(d)
]
然后大力推柿子：
[egin{aligned}
sum_{n=1}^Ng(n)=&sum_{n=1}^Nsum_{dmid n}[gcd(d,frac nd)=1]f(d)\
=&sum_{d=1}^Nf(d)sum_{k=1}^{N/d}[gcd(d,k)=1]\
=&sum_{d=1}^Nf(d)sum_{k=1}^{N/d}sum_{tmid d,tmid k}mu(t)
end{aligned}
]
最后可以得到
[sum_{n=1}^Ng(n)=sum_{t=1}^sqrt{N}mu(t)sum_{d=1}^{N/t^2}f(dt)lfloorfrac{N/t^2}{d}
floor
]
然后使用完全积性函数的性质：
[=sum_{t=1}^sqrt{N}mu(t)f(t)sum_{d=1}^{N/t^2}f(d)lfloorfrac{N/t^2}{d}
floor
]
复杂度为：
[int_1^Nsqrt{frac{N}{x^2}} mathrm dx=mathcal O(sqrt Nlog N)
]
例子：
令 (g(n)) 为积性函数，且满足：
[g(p^k)=p^k+1 ]求 (g) 的前缀和，(nleq 10^{12})
sol.
这就是 (f(n)=n) 的情况，大力爆算即可。
#include<cstdio>
#include<cmath>
typedef long long ll;
const int N = 1E+7;
const ll mod = 998244353;
ll n, ans, Mu[N + 5];
int cnt, prime[N + 5], mu[N + 5];
bool nprime[N + 5];
inline ll S(ll n) { n %= mod; return n * (n + 1) / 2 % mod; }
inline ll calc(ll n, ll k) {
n /= k * k; ll res = 0;
for(ll l = 1, r; l <= n; l = r + 1) {
r = n / (n / l);
(res += (S(r) - S(l - 1)) * (n / l)) %= mod;
} return res;
}
void pre(ll N) {
mu[1] = 1;
for(ll i = 2; i <= N; ++i) {
if(!nprime[i]) prime[++cnt] = i, mu[i] = -1;
for(int j = 1; j <= cnt && i * prime[j] <= N; ++j) {
nprime[i * prime[j]] = 1;
if(i % prime[j]) mu[i * prime[j]] = -mu[i];
else { mu[i * prime[j]] = 0; break; }
}
}
for(int i = 1; i <= N; ++i)
Mu[i] = (Mu[i - 1] + mu[i] * i) % mod;
}
int main() {
scanf("%lld", &n), pre(sqrt(n) + 5);
for(ll l = 1, r; l * l <= n; l = r + 1) {
r = sqrt(n / (n / (l * l)));
(ans += (Mu[r] - Mu[l - 1]) * calc(n, l)) %= mod;
} printf("%d", (ans + mod) % mod);
}
在不是完全积性函数的情况下：
我们定义一个新的函数 (h(n)) 满足：
[h(n)=sum_{dmid n}f(d)
]
然后假定我们有另一个函数 (s(n)) 满足：
[s(n)otimes h(n)=g(n)
]
于是注意到 (s(p)=0)，这个可以大力 Powerful Number 筛：
[egin{aligned}
sum_{n=1}^Ng(n)=&sum_{i=1}^ns(i)sum_{j=1}^{n/i}h(j)\
=&sum_{i=1}^ns(i)sum_{j=1}^{n/i}sum_{dmid j}f(d)\
=&sum_{i=1}^ns(i)sum_{d=1}^{N/i}f(d)lfloorfrac{N/i}{d}
floor\
end{aligned}
]
于是可以快速进行计算，复杂度（绝大多数时候）是 (mathcal O(n^{2/3}))。
一个简单的例子：
令 (g(n)) 为积性函数，且满足：
[g(p^k)=p^k-p^{k-1}+1 ]求 (g) 的前缀和，(nleq 10^{14})
sol.
这就是上面的 (f=varphi) 的特殊情况
这时候，通过上面的柿子可以得到：
[h(n)=nRightarrow sum_{i=1}^ns(i)sum_{j=1}^{n/i}j
]
现在来考虑 (s(p^k)) 的形式：
[g(p^k)=sum_{i=0}^k s(p^i)p^{k-i}=pg(p^{k-1})+s(p^k)
]
那么
[s(p^k)=g(p^k)-pg(p^{k-1})=1-p
]
于是大力计算即可，复杂度 (mathcal O(sqrt n))。
#include<cstdio>
#include<cmath>
typedef long long ll;
const int maxn = 1.1E+7 + 5;
const ll mod = 1E+9 + 7;
ll ans, n;
int prime[maxn], cnt;
bool nprime[maxn];
void pre(int N) {
for(int i = 2; i <= N; ++i) {
if(!nprime[i]) prime[++cnt] = i;
for(int j = 1; j <= cnt && i * prime[j] <= N; ++j) {
nprime[i * prime[j]] = 1;
if(i % prime[j] == 0) break;
}
}
}
inline ll S1(ll n) { n %= mod; return n * (n + 1) / 2 % mod; }
void DFS(int p, ll cur, ll H) {
ll nowp = 1ll * prime[p] * prime[p];
if(nowp > n / cur) return (ans += H * S1(n / cur)) %= mod, void();
DFS(p + 1, cur, H);
for( ; nowp <= n / cur; nowp *= prime[p]) {
DFS(p + 1, cur * nowp, (1 - prime[p]) * H % mod);
if(nowp > n / prime[p]) break;
}
}
int main() {
scanf("%lld", &n), pre(sqrt(n) + 500);
DFS(1, 1, 1), printf("%lld", (ans + mod) % mod);
}
令 (g(n)) 为积性函数，且满足：
[g(p^k)=p^{2k}-p^{2k-1}+1 ]求 (g) 的前缀和，(nleq 10^{10})
sol.
这是上面的 (f=id_1cdot varphi) 的特殊情况，先来考虑 (f) 的前缀和，容易得到：
[fotimes id_1=id_2
]
这个可以杜教筛 (mathcal O(n^{2/3})) 计算。
再来考虑 Powerful Number 筛的部分，我们有：
[h(p^k)=1+sum_{i=1}^kp^i(p^i-p^{i-1})=p^2h(p^{k-1})-p+1
]
现在来考虑 (s(p^k)) 的形式：
[egin{aligned}
g(p^k)&=s(p^k)+sum_{i=1}^k s(p^{k-i})h(p^i)\
&=s(p^k)+sum_{i=1}^ks(p^{k-i})(p^2h(p^{i-1})-p+1)\
&=s(p^k)+p^2g(p^{k-1})-(p-1)sum_{i=0}^{k-1}s(p^{i})\
&=s(p^k)+p^2g(p^{k-1})-(p-1)S(p^{k-1})
end{aligned}
]
于是得到
[s(p^k)=(p-1)S(p^{k-1})+1-p^2
]
于是大力计算即可，复杂度 (mathcal O(n^{2/3}))。
#include<cstdio>
#include<cmath>
#include<unordered_map>
typedef long long ll;
const int maxn = 5E+6 + 5;
const ll mod = 1E+9 + 7;
const ll inv6 = (mod + 1) / 6;
ll ans, N, n, phi[maxn];
int prime[maxn], cnt;
bool nprime[maxn];
std::unordered_map<ll, ll> M;
void pre(int N) {
phi[1] = 1;
for(int i = 2; i <= N; ++i) {
if(!nprime[i]) prime[++cnt] = i, phi[i] = i - 1;
for(int j = 1; j <= cnt && i * prime[j] <= N; ++j) {
nprime[i * prime[j]] = 1;
if(i % prime[j]) phi[i * prime[j]] = phi[i] * (prime[j] - 1);
else { phi[i * prime[j]] = phi[i] * prime[j]; break; }
}
} for(int i = 1; i <= N; ++i)
phi[i] = (phi[i - 1] + i * phi[i]) % mod;
}
inline ll S1(ll n) { n %= mod; return n * (n + 1) / 2; }
inline ll S2(ll n) { n %= mod; return n * (n + 1) % mod * (n * 2 + 1) % mod * inv6 % mod; }
inline ll Phi(ll n) {
if(n <= N) return phi[n];
else if(M.count(n)) return M[n];
ll res = S2(n);
for(ll l = 2, r; l <= n; l = r + 1) {
r = n / (n / l);
(res -= (S1(r) - S1(l - 1)) % mod * Phi(n / l)) %= mod;
} return M[n] = res;
}
inline ll calc(ll n) {
ll res = 0;
for(ll l = 1, r; l <= n; l = r + 1) {
r = n / (n / l);
(res += (Phi(r) - Phi(l - 1)) * (n / l)) %= mod;
} return res;
}
void DFS(int p, ll cur, ll H) {
ll nowp = 1ll * prime[p] * prime[p];
if(nowp > n / cur) return (ans += H * calc(n / cur)) %= mod, void();
DFS(p + 1, cur, H); ll s, S = 1;
for( ; nowp <= n / cur; nowp *= prime[p]) {
s = (S * (prime[p] - 1) + 1 - prime[p] * prime[p]) % mod;
(S += s) %= mod;
DFS(p + 1, cur * nowp, s * H % mod);
if(nowp > n / prime[p]) break;
}
}
int main() {
scanf("%lld", &n), pre(N = pow(n, 0.666) + 500);
DFS(1, 1, 1), printf("%lld", (ans + mod) % mod);
}