powerful number integrability function prefix and

\(\text{powerful number}\) is the number of times each prime factor occurs \(\geq 2\); a very important property is \(\leq n\) at the \(text{powerful number}\) level \(O(\sqrt{n})\.If you want to find all \(\leq n\) of \(\text{powerful number}\), you only need to write something similar to \(\text{min_25}\) Screen the second part of the explosive search.It is probably to record what the minimum prime number of the \(\text{powerful number}\) you are looking for, and then enumerate the prime powers.

One of the great things about this is that it can be used to optimize the process of prefix sum for quadratic functions.

Set \(F(p^q)=p^k\) and \(F(x)\) is a product function, requiring \(\sum_{i=1}^nF(i)\).

Consider finding two product functions \(G(x),H(x)\), satisfying \(F(x)=G(x)H(x)\), and satisfying \(F(p)=G(p)\ for any prime number.Not difficult to find, \(F(p)=G(1)H(p)+H(1)G(p)\), obviously\(H(1)=G(1)=1\), so \(F(p)=H(p)+G(p)\), and since \(G(p)=F(p)\, there is obviously \(H(p)=0\).

So if \(x\) satisfies \(H(x)\neq 0\), then \(x\) must be a \(\text{powerful number}\).And \(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)\ so we just need to calculate at \(O(\sqrt{n})\\) number(powerfultext{n}).

The two problems now are finding the prefix sum of \(G(i)\) and the value of \(H(i)\).

For this question, let's make \(G(x)=x^k\), so obviously there is \(G(p)=F(p)=p^k\), so the prefix of (G(i)\ is a power of natural numbers, which can be solved by various methods such as interpolation.

For \(H(i)\), there is obviously \(F(p^q)=\sum_{i=0}^qG(p^i)H(p^{q-i})\), notice that this is actually a polynomial convolution, and that \(H\) is a polynomial inversion problem when \(F,G\) is known, and that a \(O(q^2)\) violence is sufficient.When we search for \(\text{powerful number}\), we are actually looking for the composition of the mass factors, so multiply the function values of each mass factor together.

Time complexity is the complexity of \(O(\sqrt{n})\) multiplied by the calculation \(G(i)\) prefix and, for the above problem, the complexity is \(O(k\sqrt{n})\.

Code for the question above

#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) {//Indicates the smallest prime power of the powerful number to be found >=the Xth prime number
	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;
}

The difficulty is to construct a guaranteed \(G(p)=F(p)\) product function \(G(i)\), which can also quickly sum prefixes, so it may not be as flexible as (\text{min_min)25}\) sieve.

Posted on Wed, 17 Jun 2020 21:09:49 -0400 by datona