min25筛

习题部分:

#include<bits/stdc++.h> 
using namespace std;
const int N = 1e6 + 10,mod = 1e9 + 7;
typedef long long ll;
ll n,sq;
ll v[N],prime[N],sp1[N],sp2[N],cnt;
ll g1[N],g2[N],w[N],id1[N],id2[N],tot;
ll qkp(ll a,ll b)
{
	ll ans = 1;
	for(;b;b>>=1)
	{
		if(b&1) ans = ans * a % mod;
		a = a * a % mod;
	}
	return ans;
}
ll inv2 = qkp(2,mod-2),inv6 = qkp(6,mod-2);
void init()
{
	for(int i=2;i<=sq;++i)
	{
		if(!v[i])
		{
			v[i] = prime[++cnt] = i;
			sp1[cnt] = (sp1[cnt-1] + i) % mod;
			sp2[cnt] = (sp2[cnt-1] + 1ll*i*i %mod) % mod;
		}
		for(int j=1;prime[j]*i<N;++j)
		{
			v[prime[j]*i] = prime[j];
			if(i%prime[j]==0) break;
		}
	}
}
ll g(ll k,ll x)
{
	return (k / (k / x));
}
ll S(ll i,ll j)
{
	if(prime[j] >= i) return 0;
	
	ll p = i <= sq ? id1[i] : id2[n/i];
	//g2存的是x^2,g1存的是x,而要求的是x^2-x,所以公式如下 
	ll ans = ( (g2[p] - g1[p] + mod) % mod - (sp2[j] - sp1[j] + mod) % mod + mod)  % mod;
	
	for(int k=j+1;prime[k]*prime[k]<=i && k<=cnt;k++)
	{
		ll pe = prime[k];
		for(int e=1;pe<=i;++e,pe = pe*prime[k])
		{
			ll x = pe % mod;
			ans = (ans + x*(x-1) % mod * (S(i/pe,k) + (e>1)) % mod) % mod;
		}
	}
	return ans;	
}
int main()
{
	scanf("%lld",&n);
	sq = sqrt(n);
	init();
	//相当于是把积性函数拆开为两个,x和x^2,进行分别运算 
	for(ll l=1,r;l<=n;l=r+1)//预处理g[x,0]
	{
		r = min(n,g(n,l));
		w[++tot] = n / l % mod;
		
		g1[tot] = ( w[tot] * (w[tot]+1) % mod * inv2 % mod -1 + mod) % mod;//对x求一个前缀和,因为1不大于第0个质数,所以要减去 
		//对x^2求一个前缀和,因为1不大于第0个质数,所以要减去 
		g2[tot] = ( w[tot] * (w[tot]+1) % mod * (2*w[tot]+1) % mod * inv6 % mod -1 + mod) % mod;
		
		w[tot] = n / l;
		if(w[tot]<=sq) id1[ w[tot] ] = tot;
		else id2[ n/w[tot] ] = tot;
	}
	//sp1和sp2存的是前面的所有质数的函数值的前缀和 
	for(int j=1;j<=cnt;++j) // g(n,j)
	{
		for(int i=1;i<=tot && prime[j]*prime[j]<=w[i];++i)
		{
			ll tmp = w[i] / prime[j];
			ll p = tmp <= sq ? id1[tmp] : id2[ n/tmp ];
			
			g1[i] = (g1[i] - prime[j] * (g1[p] - sp1[j-1] + mod) % mod + mod) % mod;
			g2[i] = (g2[i] - prime[j] * prime[j] % mod * (g2[p] - sp2[j-1] + mod) % mod + mod) % mod;
		}
	}
	
	printf("%lld\n",S(n,0) + 1);
	
	return 0;
}