[SPOJ20174]DIVCNT3 - Counting Divisors (cube):Min_25筛

分析

首先,STO ywy OTZ,ywy TQL%%%!

说一下这道题用min_25筛怎么做。

容易发现,对于所有质数\(p\),都满足\(f(p)=4\),于是我们就可以直接通过\([1,x]\)内的质数的个数\(h(x)\)来求出\(g(x)=\sum_{i=1}^{x}f(i) \times [i \in prime]\)了,即\(g(x)\)可以等价地表示为\(g(x)=4 \times h(x)\)。如何求\(h(x)\)是min_25筛的基本操作就不过多赘述了。而且进一步分析我们可以得出,\(f(p^q)=3 \times q+1,p \in prime\),这个结论在实现min_25筛时也很重要。

代码里的g[x]其实是我上面写的\(h(x)\)

说一句题外话,min_25筛的话,使用哈希表实现更直观,但是如果用id1,id2两个数组实现的话常数会显著减小(哈希表的地方我注释掉了)。

代码

#include <bits/stdc++.h>
#define rin(i,a,b) for(register int i=(a);i<=(b);++i)
#define irin(i,a,b) for(register int i=(a);i>=(b);--i)
#define trav(i,a) for(register int i=head[a];i;i=e[i].nxt)
typedef long long LL;
using std::cin;
using std::cout;
using std::endl;

inline LL read(){
	LL x=0,f=1;char ch=getchar();
	while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=getchar();}
	while(ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
	return x*f;
}

const LL MAXN=1e11+5;
const int MAXR=1e6+5;
const int HMOD=3e6-1;
LL n,num[MAXR],g[MAXR];
int cnt,tot,prm[MAXR],id1[MAXR],id2[MAXR];
bool vis[MAXR];

/*
struct Hash_Map{
	int siz,head[HMOD+5],nxt[MAXR],to[MAXR],sta[MAXR],top;
	LL val[MAXR];
	inline void insert(LL x,int y){
		++siz;
		nxt[siz]=head[x%HMOD];
		val[siz]=x;
		to[siz]=y;
		head[x%HMOD]=siz;
		sta[++top]=x%HMOD;
	}
	inline int operator [] (LL x){
		for(register int i=head[x%HMOD];i;i=nxt[i])
			if(val[i]==x) return to[i];
	}
	void clear(){
		siz=0;
		while(top) head[sta[top--]]=0;
	}
}mp;
*/

void pre_process(int n){
	rin(i,2,n){
		if(!vis[i]) prm[++cnt]=i;
		rin(j,1,cnt){
			if(i*prm[j]>n) break;
			vis[i*prm[j]]=true;
			if(i%prm[j]==0) break;
		}
	}
}

LL solve(LL x,int y){
//	int k=mp[x];
	int k=x<=MAXR-5?id1[x]:id2[n/x];
	LL ret=(g[k]<<2)-((y-1)<<2);
	rin(i,y,cnt){
		if(1ll*prm[i]*prm[i]>x) break;
		LL temp1=prm[i],temp2=1ll*prm[i]*prm[i];
		for(register int j=1;temp2<=x;++j,temp1=temp2,temp2*=prm[i])
			ret+=(3*j+1)*solve(x/temp1,i+1)+3*j+4;
	}
	return ret;
}

int main(){
	int T=read();
	pre_process(500000);
	while(T--){
		tot=0;
//		mp.clear();
		n=read();
		LL nxti;
		for(register LL i=1;i<=n;i=nxti){
			nxti=n/(n/i)+1;
			num[++tot]=n/i;
//			mp.insert(num[tot],tot);
			if(num[tot]<=MAXR-5) id1[num[tot]]=tot;
			else id2[nxti-1]=tot;
			g[tot]=num[tot]-1;
		}
		rin(i,1,cnt){
			if(1ll*prm[i]*prm[i]>n) break;// 这个一定要加,否则会T第一个点!
			rin(j,1,tot){
				if(1ll*prm[i]*prm[i]>num[j]) break;
//				int k=mp[num[j]/prm[i]];
				int k=num[j]/prm[i]<=MAXR-5?id1[num[j]/prm[i]]:id2[n/(num[j]/prm[i])];
				g[j]-=g[k]-(i-1);
			}
		}
		printf("%lld\n",solve(n,1)+1);
	}
	return 0;
}

posted on 2019-01-11 19:01  ErkkiErkko  阅读(265)  评论(0编辑  收藏  举报