\(\text{powerful number}\)即每個質因子出現次數都\(\geq 2\)的數；一個非常重要的性質是\(\leq n\)的\(\text{powerful number}\)是\(O(\sqrt{n})\)級別的。如果想要求出所有\(\leq n\)的\(\text{powerful number}\)的話，只需要寫一個類似于\(\text{min_25}\)篩第二部分的爆搜即可。大概就是記錄一下當前要找的\(\text{powerful number}\)的最小質數是什么，之后枚舉質數次冪即可。

這個東西比較厲害的一點是可以用來優化求積性函數前綴和的過程。

設\(F(p^q)=p^k\)，且\(F(x)\)是一個積性函數，求\(\sum_{i=1}^nF(i)\)。

考慮找兩個積性函數\(G(x),H(x)\)，滿足\(F(x)=G(x)H(x)\)，同時對于任意質數\(p\)，均滿足\(F(p)=G(p)\)。不難發現，\(F(p)=G(1)H(p)+H(1)G(p)\)，很顯然\(H(1)=G(1)=1\)，于是\(F(p)=H(p)+G(p)\)，又因為\(G(p)=F(p)\)，所以顯然有\(H(p)=0\)。

所以如果\(x\)滿足\(H(x)\neq 0\)，那么\(x\)一定是一個\(\text{powerful number}\)。而\(\displaystyle \sum_{i=1}^nF(i)=\sum_{i\times j\leq n}H(i)G(j)=\sum_{i=1}^nH(i)\sum_{j=1}^{\lfloor \frac{n}{i}\rfloor}G(j)\)，于是我們只需要在\(O(\sqrt{n})\)個\(\text{powerful number}\)處計算即可。

現在的兩個問題是求出\(G(i)\)的前綴和，以及\(H(i)\)的值。

對于這道題，我們令\(G(x)=x^k\)，那么顯然有\(G(p)=F(p)=p^k\)，于是\(G(i)\)的前綴和就是一個自然數冪次方和，可以利用插值等多種多樣的方法求出。

對\(H(i)\)，顯然有\(F(p^q)=\sum_{i=0}^qG(p^i)H(p^{q-i})\)，注意到這實際上是一個多項式卷積的形式，已知\(F,G\)的情況下求\(H\)就是一個多項式求逆問題，一個\(O(q^2)\)的暴力即可。我們在搜\(\text{powerful number}\)的時候實際上在搜質因子組成，于是把各個質因子的函數值乘在一起就好了。

時間復雜度是\(O(\sqrt{n})\)乘上計算\(G(i)\)前綴和的復雜度，對于上面的問題，復雜度就是\(O(k\sqrt{n})\)。

上面那道題的代碼

#include<cmath>
#include<vector>
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#define re register
#define LL long long
const int maxn=3.3e6+5;
const int mod=1e9+7;
std::vector<int> h[300005];
LL n;int k,Sqr,ans,f[maxn],sz[maxn<<1],vis[maxn<<1],p[maxn>>1];
int A[105],B[105],g[105],ifac[105],suf[105],pre[105];
inline int dqm(int x) {return x<0?x+mod:x;}
inline int qm(int x) {return x>=mod?x-mod:x;}
inline int ksm(int a,int b) {
	int S=1;for(;b;b>>=1,a=1ll*a*a%mod)if(b&1)S=1ll*S*a%mod;return S;
}
inline int get(LL nw) {
	int id=nw<=Sqr?nw:Sqr+n/nw;
	if(vis[id]) return sz[id];
	nw%=mod;int tot=0;vis[id]=1;
	pre[0]=nw;for(re int i=1;i<=k;i++)pre[i]=1ll*pre[i-1]*dqm(nw-i)%mod;
	suf[k+1]=1;for(re int i=k;i;--i)suf[i]=1ll*suf[i+1]*dqm(nw-i)%mod;
	for(re int i=1;i<=k;i++) {
		int v=1ll*pre[i-1]*suf[i+1]%mod*g[i]%mod*ifac[i]%mod*ifac[k-i]%mod;
		if((k-i)&1) v=dqm(-v);tot=qm(tot+v);
	}
	return sz[id]=tot;
}
void dfs(int x,int v,LL res) {//表示當前要找的powerful number的最小質數次冪>=第x個質數
	ans=qm(ans+1ll*v*get(res)%mod);
	if(x==p[0]+1)return;if(res/p[x]/p[x]<=0) return;LL R=res/p[x]/p[x];
	for(re int i=x;i<=p[0]&&res/p[i]/p[i]>0;++i,R=res/(p[i]?p[i]:1)/(p[i]?p[i]:1)) 
		for(re int j=2;R>0;++j,R/=p[i]) dfs(i+1,1ll*v*h[i][j]%mod,R);
}
int main() {
	scanf("%lld%d",&n,&k);Sqr=1+sqrt(n);
	for(re int i=2;i<=Sqr;i++) {
		if(!f[i]) p[++p[0]]=i;
		for(re int j=1;j<=p[0]&&p[j]*i<=Sqr;++j) {
			f[p[j]*i]=1;if(i%p[j]==0)break;
		}
	}
	for(re int i=1;i<=p[0];++i) {
		int t=0;LL res=n;while(res>=p[i])res/=p[i],++t;
		A[0]=B[0]=1;A[1]=ksm(p[i],k);
		for(re int j=2;j<=t;++j)A[j]=A[j-1];
		for(re int j=1;j<=t;++j)B[j]=1ll*B[j-1]*A[1]%mod;
		h[i].push_back(1);
		for(re int j=1;j<=t;++j) {
			int nw=A[j];
			for(re int k=0;k<j;++k)nw=dqm(nw-1ll*h[i][k]*B[j-k]%mod);
			h[i].push_back(nw);
		}
	}++k;
	for(re int i=1;i<=k;i++)g[i]=qm(g[i-1]+ksm(i,k-1));ifac[0]=ifac[1]=1;
	for(re int i=2;i<=k;i++)ifac[i]=1ll*(mod-mod/i)*ifac[mod%i]%mod;
	for(re int i=2;i<=k;i++)ifac[i]=1ll*ifac[i-1]*ifac[i]%mod;
	dfs(1,1,n);printf("%d\n",ans);return 0;
}

難點在于構造一個保證\(G(p)=F(p)\)積性函數\(G(i)\)，同時這個函數還能快速求前綴和，所以靈活程度上可能并不如\(\text{min_25}\)篩。