RT，第一篇题解思路

#include<bits/stdc++.h>
#define int long long
using namespace std;
const long long INF=1LL<<31;
int N=8000000;
int cnt,n;
long long p[8001000],inv2,inv6,ans,zhi[8001000],mod;
bool he[8001000];
map<long long,long long>M;
long long S(long long x){x%=mod;return x*(x+1)%mod*inv2%mod;}
long long Sump(long long x){x%=mod;return x*(x+1)%mod*(x+x+1)%mod*inv6%mod;}
void xxs()
{
	he[1]=p[1]=1;
	for(int i=2;i<=N;i++)
	{
		if(he[i]==0)
		{
			p[i]=(i-1)%mod;
			zhi[++cnt]=i;
		}
		for(int j=1;j<=cnt&&i*zhi[j]<=N;j++)
		{
			he[i*zhi[j]]=true;
			if(i%zhi[j]==0)
			{
				p[i*zhi[j]]=1LL*p[i]*zhi[j]%mod;
				break;
			}
			else
			{
				p[i*zhi[j]]=1LL*p[i]*(zhi[j]-1)%mod;
			}
		}
	}
	for(int i=1;i<=N;i++)p[i]=(p[i-1]+1ll*p[i]*i%mod*i%mod)%mod;
}
long long SF(long long x)
{
    if(x<=N)
		return p[x];
    if(M.find(x)!=M.end())
		return M[x];
    long long ret=S(x);
	ret=ret*ret%mod;
    for(long long i=2,r;i<=x;i=r+1)
    {
        r=x/(x/i);
        long long tt=(Sump(r)-Sump(r-1))%mod;
        ret-=SF(x/i)*tt%mod;
        ret%=mod;
    }
    return M[x]=(ret+mod)%mod;
}
long long quick_pow(long long a,long long b)
{
	long long res=1;
	while(b>0)
	{
		if(b&1)
		{
			res*=a;
			res%=mod;
		}
		a*=a;
		a%=mod;
		b>>=1;
	}
	return res;
}
signed main()
{
	scanf("%lld%lld",&mod,&n);
	inv2=quick_pow(2,mod-2);
	inv6=quick_pow(6,mod-2);
	xxs();
	for(long long i=1,r;i<=n;i=r+1)
    {
        r=n/(n/i);
        long long tt=S(n/i);
		tt=tt*tt%mod;
        long long gg=(SF(r)-SF(i-1))%mod;
        ans+=gg*tt%mod;
        ans%=mod;
    }
    printf("%lld\n",(ans+mod)%mod);
	return 0;
}