$\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}\) 时才非零!由于有

\[\begin{aligned} \sum_{i=1}^n f_i&=\sum_{i=1}^n \sum_{d|i} h_d g_{d/i}\\ &=\sum_{d=1}^n h_d\sum_{i=1}^{n/d} g_i\\ \end{aligned} \]

所以我们只要对于 \(i\in \text{PN}\) 都计算出 \(h_i\) 就可以了。又 \(h(n)\) 为积性函数,所以仅需快速求出 \(h(p^c)\)

我们推一推式子:

\[\begin{aligned} f(p^c)&=\sum_{i=0}^c g(p^i)h(p^{c-i})\\ &=g(1)h(p^c)+\sum_{i=1}^c g(p^i)h(p^{c-i})\\ g(1)h(p^c)&=f(p^c)-\sum_{i=1}^c g(p^i)h(p^{c-i}) \end{aligned} \]

即可 \(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\),所以有:

\[\begin{aligned} h(p^c)&=f(p^c)-\sum_{i=1}^c g(p^i)h(p^{c-i})\\ &=f(p^c)-p^2\sum_{i=1}^c g(p^{i-1})h(p^{c-i})\\ &=f(p^c)-p^2 f(p^{c-1}) \end{aligned} \]

有一种特殊情况就是 \(\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;
}

推荐习题:简单的函数HDU7186,Project Euler 639。

posted @ 2022-08-03 22:13  b1ts  阅读(99)  评论(3编辑  收藏  举报