$\text{Powerful Number}$ 筛
昨天多校出了有关内容的题目,正好学一手。
定义:一个正整数 \(n=\prod_{p_i\in prime} p_i^{c_i}\) 是 \(\text{Powerful Number}\) (以下简称 \(\text{PN}\) )当且仅当 \(\forall i,c_i>1\)。
结论:小于等于正整数 \(N\) 的 \(\text{PN}\) 只有 \(O(\sqrt{N})\) 个。
证明:不会。
于是我们对于一个 \(n\),可以爆搜出所有小于等于它的 \(\text{PN}\)。
那么这东西有什么用呢?可以在有些时候爆杀 \(\text{min25}\)!
考虑对于一个积性函数 \(f(n)\) 去计算它的前缀和。我们可以构造一个积性函数 \(g(n)\) 使得 \(g(p)=f(p)\) 且 \(g\) 的前缀和比较好算。令 \(f=g*h\),显然 \(h(n)\) 也是积性函数。由于 \(f(p)=g(p)h(1)+g(1)h(p)=g(p)\),所以 \(h(p)=0\),也就是说 \(h(n)\) 只有在 \(n\) 是 \(\text{PN}\) 时才非零!由于有
所以我们只要对于 \(i\in \text{PN}\) 都计算出 \(h_i\) 就可以了。又 \(h(n)\) 为积性函数,所以仅需快速求出 \(h(p^c)\)。
我们推一推式子:
即可 \(O(\log n)\) 计算单个 \(h(p^c)\)。这样做时间复杂度为 \(O(\sqrt{n}\log n)\),实际上完全跑不满。当然,如果能推出 \(h(p^c)\) 关于 \(p,c\) 的公式,时间复杂度就正好是 \(O(\sqrt{n})\)。
例题:【模板】Min_25筛
题意:定义 \(f(n)=\begin{cases}1&n=1\\ p^k(p^k-1)&n=p^k \\ f(a)f(b) &n=ab,a\perp b\end{cases}\),求 \(\sum\limits_{i=1}^n f(i)\)。数据范围:$ n\leq 10^{10}$。
可以发现 \(f(p)=p(p-1)=id(p)\varphi(p)\),于是我们可以构造 \(g(n)=id(n)\varphi(n)\)。$ g $ 的前缀和可以杜教筛求出来。现在主要是要求 \(h(p^c)\)。这里提供一种方法。由于有 \(\dfrac{g(p^c)}{g(p^{c-1})}=\dfrac{p^c(p^c-p^{c-1})}{p^{c-1}(p^{c-1}-p^{c-2})}=p^2\),所以有:
有一种特殊情况就是 \(\dfrac{g(p)}{g(1)}=p(p-1) \ne p^2\),因此要把多出的贡献 \(g(p)h(p^{c-1})-p^2g(1)h(p^{c-1})=-ph(p^{c-1})\) 减回去,这样就可以递推计算 \(h(p^c)\)。当然,你可能也可以得到公式 \(h(p^c)=(p-1)(c-1)p^c\)。于是总时间复杂度为 \(O(n^{2/3})\),瓶颈为杜教筛。
代码:
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
inline ll read(){
ll f=1, r=0; char c=getchar();
while (!isdigit(c)) f^=c=='-', c=getchar();
while (isdigit(c)) r=(r<<1)+(r<<3)+(c&15), c=getchar();
return f?r:-r;
}
const int N=5e6+7, mod=1e9+7, inv2=(mod+1)>>1, inv6=(mod+1)/6;
inline void inc(int &x, int y) {x+=y-mod, x+=x>>31&mod;}
inline void dec(int &x, int y) {x-=y, x+=x>>31&mod;}
inline int qpow(int a, int b) {
int res=1;
for (; b; b>>=1, a=(ll)a*a%mod)
if (b&1) res=(ll)res*a%mod;
return res;
}
int s_p, pr[N], sum1[N], phi[N];
bool v[N];
inline void init(int n) {
sum1[1]=phi[1]=1;
for (int i=2; i<=n; i++) {
if (!v[i]) pr[++s_p]=i, phi[i]=i-1;
for (int j=1; j<=s_p && i*pr[j]<=n; j++) {
v[i*pr[j]]=true;
if (i%pr[j]==0) {phi[i*pr[j]]=phi[i]*pr[j]; break;}
phi[i*pr[j]]=phi[i]*phi[pr[j]];
}
sum1[i]=(sum1[i-1]+(ll)phi[i]*i)%mod;
}
}
int sum2[N];
ll n;
inline int& Sum(ll n) {return n<N?sum1[n]:sum2[::n/n];}
inline int S2(ll n) {return n%=mod, n*(n+1)%mod*(2*n+1)%mod*inv6%mod;}
inline int S(ll l, ll r) {return l%=mod, r%=mod, (l+r)*(r+mod-l+1)%mod*inv2%mod;}
inline void S(ll n) {
int &v=Sum(n); if (v) return; v=S2(n);
for (ll l=2, r; l<=n; l=r+1) {
ll d=n/l; r=n/d, S(d);
dec(v, (ll)Sum(d)*S(l, r)%mod);
}
}
void dfs(int i, ll x, int h, int &ans) {
ll pw=(ll)pr[i]*pr[i];
if (i==s_p+1 || x>n/pw) return inc(ans, (ll)h*Sum(n/x)%mod);
dfs(i+1, x, h, ans); int num=(ll)pr[i]*(pr[i]-1)%mod, lst=0;
while (pw*x<=n) {
int a=pw%mod, now=(ll)a*(a+mod-1)%mod;
int nxt=(now+mod-(ll)pr[i]%mod*pr[i]%mod*num%mod+(ll)lst*pr[i])%mod;
dfs(i+1, x*pw, (ll)h*nxt%mod, ans), num=now, pw*=pr[i], lst=nxt;
}
}
int main() {
#ifdef LOCAL
freopen("1.in", "r", stdin);
freopen("1.out", "w", stdout);
#endif
n=read(), init(min(N-1ll, n)), S(n);
int ans=0; dfs(1, 1, 1, ans), printf("%d\n", ans);
return 0;
}