shit!

#include <cstdio>
#include <iostream>
#include <cmath>
using namespace std;
const int N = 200000;
const int M = N+5;
const int MOD = 1e9+7;
const int inv = 166666668;
#define int long long
int read()
{
	int x=0,f=1;char c;
	while((c=getchar())<‘0‘ || c>‘9‘) {if(c==‘-‘) f=-1;}
	while(c>=‘0‘ && c<=‘9‘) {x=(x<<3)+(x<<1)+(c^48);c=getchar();}
	return x*f;
}
int n,cnt,tot,p[M],vis[M],w[M],sp1[M],sp2[M];
int sqr,id1[M],id2[M],g1[M],g2[M];
void init(int n)
{
	for(int i=2;i<=n;i++)
	{
		if(!vis[i])
		{
			p[++cnt]=i;
			sp1[cnt]=(sp1[cnt-1]+i)%MOD;
			sp2[cnt]=(sp2[cnt-1]+i*i)%MOD;
		}
		for(int j=1;j<=cnt && i*p[j]<=n;j++)
		{
			vis[i*p[j]]=1;
			if(i%p[j]==0) break;
		}
	}
}
int S(int x,int y)
{
	if(x<=p[y]) return 0;
	int k=x<=sqr?id1[x]:id2[n/x];
	int ans=(g2[k]-g1[k]-(sp2[y]-sp1[y]))%MOD;
	for(int i=y+1;i<=cnt && p[i]*p[i]<=x;i++)
	{
		int pr=p[i];
		for(int e=1;pr<=x;e++,pr*=p[i])
		{
			int xx=pr%MOD;
			ans=(ans+xx*(xx-1)%MOD*(S(x/pr,i)+(e!=1)))%MOD;
		}
	}
	return ans;
}
signed main()
{
	init(N);
	n=read();sqr=sqrt(n);
	for(int l=1,r;l<=n;l=r+1)
	{
		r=n/(n/l);
		w[++tot]=n/l;
		g1[tot]=w[tot]%MOD;//要模要模 
		g2[tot]=g1[tot]*(g1[tot]+1)%MOD*(2*g1[tot]+1)%MOD*inv%MOD-1;
		g1[tot]=g1[tot]*(g1[tot]+1)/2%MOD-1;
		if(n/l<=sqr) id1[n/l]=tot;//编号方式是x 
		else id2[n/(n/l)]=tot;//编号方式是n/x
	}
	for(int i=1;i<=cnt;i++)
		for(int j=1;j<=tot && p[i]*p[i]<=w[j];j++)
		{
			int k=w[j]/p[i]<=sqr?id1[w[j]/p[i]]:id2[n/(w[j]/p[i])];
			g1[j]=(g1[j]-p[i]*(g1[k]-sp1[i-1]))%MOD;
			g2[j]=(g2[j]-p[i]*p[i]%MOD*(g2[k]-sp2[i-1]))%MOD;
		}
	printf("%lld\n",(S(n,0)+1+MOD)%MOD);
}

shit!

原文：https://www.cnblogs.com/C202044zxy/p/14355392.html