杜教筛是用来求一类积性函数的前缀和,利用数论分块的思想来降低复杂度
假设我们现在要求 (S(n) = sum_{i = 1}^n f(i)) ,(f(i)) 为积性函数, (n leqslant 10^{12}) ,假设有另一个积性函数 (g)。我们来求它们狄利克雷卷积的前缀和
[egin{align}
&sum_{i = 1}^n (g * f) = sum_{i = 1}^n sum_{d mid i} g(d) f(frac{i}{d})
\
=& sum_{d = 1}^n g(d) sum_{d|i} f(frac{i}{d})
\
=& sum_{d = 1}^n g(d) sum_{i = 1}^{frac{n}{d}}f(i)
\
=& sum_{d = 1}^n g(d) S(frac{n}{d})
end{align}
]
容斥处理一下,得到
[g(1)S(n) = sum_{d = 1}^n g(d)S(frac{n}{d}) - sum_{d = 2}^n g(d)S(frac{n}{d})
]
至于 (g) 的选择,需要看具体情况
原函数 (f(n)) | 选取的函数 (g(n)) | 递归式 |
---|---|---|
(mu(n)) | (1) | (S(n) = 1 - sum_{d = 2}^n S(frac{n}{d})) |
(varphi(n)) | (1) | (S(n) = frac{n (n + 1)}{2} - sum_{d = 2}^n S(frac{n}{d})) |
(ncdot varphi(n)) | (S(n)=sum_{i=1}^{n}i^2-sum_{d=2}^{n}dcdot S(lfloorfrac{n}{d} floor)) |
注意询问数量较多时可以采用 unordered_map
做记忆化来加速。
卡常什么的太讨厌了
#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 7.5e6+5;
map<signed,int> mp,mm;
bool isNotPrime[N + 5];
int mu[N + 5], phi[N + 5];
signed primes[N + 5], cnt;
inline void euler() {
isNotPrime[0] = isNotPrime[1] = true;
mu[1] = 1;
phi[1] = 1;
for (int i = 2; i <= N; i++) {
if (!isNotPrime[i]) {
primes[++cnt] = i;
mu[i] = -1;
phi[i] = i - 1;
}
for (int j = 1; j <= cnt; j++) {
int t = i * primes[j];
if (t > N) break;
isNotPrime[t] = true;
if (i % primes[j] == 0) {
mu[t] = 0;
phi[t] = phi[i] * primes[j];
break;
} else {
mu[t] = -mu[i];
phi[t] = phi[i] * (primes[j] - 1);
}
}
}
}
int Mu(int n) {
if(n<=N) return mu[n];
if(mm[n]) return mm[n];
register signed l=2,r;
int ans=1;
while(l<=n) {
r=n/(n/l);
ans-=(r-l+1)*Mu(n/l);
l=r+1;
}
mm[n]=ans;
return ans;
}
int Phi(int n) {
if(n<=N) return phi[n];
if(mp[n]) return mp[n];
register signed l=2,r;
int ans=n*(n+1)/2;
while(l<=n) {
r=n/(n/l);
ans-=(r-l+1)*Phi(n/l);
l=r+1;
}
mp[n]=ans;
return ans;
}
signed main() {
ios::sync_with_stdio(false);
euler();
for(int i=1;i<=N;i++) mu[i]+=mu[i-1], phi[i]+=phi[i-1];
int n,t;
cin>>n;
while(n--) {
cin>>t;
cout<<Phi(t)<<" "<<Mu(t)<<endl;
}
}