求助(应该漏取模了罢QAQ)
查看原帖
求助(应该漏取模了罢QAQ)
335136
LordLaffey楼主2021/9/23 14:31

RT

#include<iostream>
#include<cstdio>
#include<cmath>
#define ll long long
using namespace std;
const ll mod=1e9+7;
const int SIZE=1e6+10;
const ll inv3=333333336;
ll prime[SIZE],S[SIZE],Sp[SIZE],idf[SIZE],ids[SIZE],g[SIZE],gp[SIZE];
//idf 一次项位置   ids 二次项位置 
//g,gp分别为G函数的一次前缀与二次前缀 
bool flag[SIZE];
ll sq,tot,cnt;
ll n,w[SIZE];

ll give(ll i,ll j){
	
	return (w[i]/prime[j]<=sq)?idf[w[i]/prime[j]]:ids[n/(w[i]/prime[j])];
	
}

void next(int x){
	
	if(g[x]<0) g[x]+=mod;
	if(gp[x]<0) gp[x]+=mod;
	
}

void shai(int x){
	
	flag[1]=true;
	for(int i=2;i<=x;i++)
	{
		if(!flag[i])
			prime[++tot]=i,S[tot]=(S[tot-1]+i)%mod,Sp[tot]=(Sp[tot-1]+i*i)%mod;
			
//		int j=1;
//		while(i*prime[j]<=n&&i%prime[j]!=0)
//			flag[i*prime[j]]=true,j++;
		for(int j=1;j<=tot&&i*prime[j]<=x;++j){
			flag[i*prime[j]]=true;
			if (i%prime[j]==0) break;
		}
	}
	
}//step1 质数筛 

ll s(ll x,ll y){
	
	if(prime[y]>=x)	return 0;
	ll k=(x<=sq)?idf[x]:ids[n/x];
//	ll res=(g[k]-h[k]-S[y-1]+y-1)%mod;res=(res+mod)%mod;
	ll res=(gp[k]-g[k]+mod-(Sp[y]-S[y])+mod)%mod;
//	if(y==1) res+=2;
	for(int i=y+1;i<=tot&&prime[i]*prime[i]<=x;i++)
	{
		ll pe=prime[i];
		for(int e=1;pe<=x;e++,pe=pe*prime[i])
		{
			ll a=pe%mod;
			res=(res+a*(a-1)%mod*(s(x/pe,i)+(e!=1)))%mod;
		}

	}
	
	return res%mod;
	
}

int main(){

	scanf("%lld",&n);
	sq=sqrt(n);
	shai(sq);
	
	ll j;
	for(ll i=1;i<=n;i=j+1)
	{
		j=n/(n/i);w[++cnt]=n/i;
		g[cnt]=w[cnt]%mod;
		gp[cnt]=g[cnt]*(g[cnt]+1)/2%mod*(2*g[cnt]+1)%mod*inv3%mod-1;
		g[cnt]=g[cnt]*(g[cnt]+1)/2%mod-1;
		if(w[cnt]<=sq)	idf[w[cnt]]=cnt;
		else	ids[n/(n/i)]=cnt;
//		if(!g[cnt])	g[cnt]+=mod;g[cnt]/=2;
//		if(!gp[cnt])	gp[cnt]+=mod;
	}
	
	for(int j=1;j<=tot;j++)
	{
		for(int i=1;i<=cnt&&prime[j]*prime[j]<=w[i];i++)
		{
			ll k=give(i,j)%mod;
//			g[i]=(g[i]-prime[j]*(g[k]-S[j-1])%mod)%mod; g[i]=(g[i]+mod)%mod;
//			h[i]=(h[i]-h[k]+j-1)%mod; h[i]=(h[i]+mod)%mod;
			g[i]-=prime[j]*(g[k]-S[j-1]+mod)%mod;
//			gp[i]-=((prime[j]*prime[j])%mod)*(gp[k]-Sp[j-1]+mod)%mod;
			gp[i]-=((prime[j]*prime[j])%mod)*(gp[k]-Sp[j-1]+mod)%mod;
			g[i]%=mod;
			gp[i]%=mod;
			next(i);
		}
	}
	cout<<(s(n,0)+1)%mod;
//	printf("%lld",s(n,0)+1);
	
	return 0;
	
	
}

实在找不出来qwq

2021/9/23 14:31
加载中...