【SPOJ】DIVCNTK min_25筛

时间：2021-10-18 00:41:34

题目大意

　　给你 \(n,k\)，求

\[S_k(n)=\sum_{i=1}^n\sigma_0(i^k)
\]

　　对 \(2^{64}\) 取模。

题解

　　一个min_25筛模板题。

　　令 \(f(n)=\sigma_0(n^k)\)，那么 \(S_k(n)=\sum_{i=1}^nf(i)\)，而且

\[\begin{cases}
f(1)&=1\\
f(p)&=k+1\\
f(p^c)&=kc+1
\end{cases}
\]

　　直接上min_25筛就好了。

　　时间复杂度：\(O(\frac{n^\frac{3}{4}}{\log n})\)

　　某些实现可能会有一点微小的细节。

代码

#include<cstdio>

#include<cstring>

#include<algorithm>

#include<cstdlib>

#include<ctime>

#include<utility>

#include<functional>

#include<cmath>

#include<assert.h>

using namespace std;

typedef long long ll;

typedef unsigned long long ull;

typedef pair<int,int> pii;

typedef pair<ll,ll> pll;

void open(const char *s){

#ifndef ONLINE_JUDGE

	char str[100];sprintf(str,"%s.in",s);freopen(str,"r",stdin);sprintf(str,"%s.out",s);freopen(str,"w",stdout);

#endif

}

int rd(){int s=0,c,b=0;while(((c=getchar())<'0'||c>'9')&&c!='-');if(c=='-'){c=getchar();b=1;}do{s=s*10+c-'0';}while((c=getchar())>='0'&&c<='9');return b?-s:s;}

void put(int x){if(!x){putchar('0');return;}static int c[20];int t=0;while(x){c[++t]=x%10;x/=10;}while(t)putchar(c[t--]+'0');}

int upmin(int &a,int b){if(b<a){a=b;return 1;}return 0;}

int upmax(int &a,int b){if(b>a){a=b;return 1;}return 0;}

const int M=100010;

int pri[M],cnt,b[M];

ll f1[M],f2[M];

void init()

{

	for(int i=2;i<=100000;i++)

	{

		if(!b[i])

			pri[++cnt]=i;

		for(int j=1;j<=cnt&&i*pri[j]<=100000;j++)

		{

			b[i*pri[j]]=1;

			if(i%pri[j]==0)

				break;

		}

	}

	pri[cnt+1]=100001;

}

ll n,k;

int m;

ll dfs(ll x,int y)

{

	if(x<=1||x<pri[y])

		return 0;

	if(pri[y]>m)

		return (x<=m?f1[x]:f2[n/x])-f1[m];

	ll s=(x<=m?f1[x]:f2[n/x])-f1[pri[y]-1];

	for(int i=y;i<=cnt&&(ll)pri[i]*pri[i]<=x;i++)

	{

		ll x1=x/pri[i];

		for(int j=1;x1>=pri[i];j++,x1/=pri[i])

			s+=dfs(x1,i+1)*(j*k+1)+((j+1)*k+1);

	}

	return s;

}

void solve()

{

	scanf("%lld%lld",&n,&k);

	m=sqrt(n)+0.5;

	int mx=n/(m+1);

	for(int i=2;i<=m;i++)

		f1[i]=i-1;

	for(int i=1;i<=mx;i++)

		f2[i]=n/i-1;

	for(int i=1;i<=cnt&&pri[i]<=m;i++)

	{

		ll x=f1[pri[i]-1];

		int n1=min((ll)mx/pri[i],n/pri[i]/pri[i]);

		int n2=min((ll)mx,n/pri[i]/pri[i]);

		for(int j=1;j<=n1;j++)

			f2[j]-=f2[j*pri[i]]-x;

		for(int j=n1+1;j<=n2;j++)

			f2[j]-=f1[n/j/pri[i]]-x;

		for(int j=m;j>=(ll)pri[i]*pri[i];j--)

			f1[j]-=f1[j/pri[i]]-x;

	}

	for(int i=2;i<=m;i++)

		f1[i]*=k+1;

	for(int i=1;i<=mx;i++)

		f2[i]*=k+1;

	ll ans=dfs(n,1);

	ans++;

	printf("%llu\n",ans);

}

int main()

{

	open("divcntk");

	int t;

	scanf("%d",&t);

	init();

	while(t--)

		solve();

	return 0;

}

相关文章