P3768 简单的数学题 洛谷 （杜教筛）

#include <bits/stdc++.h>
#define ll long long
#define int long long 
using namespace std;
template <class cl>
void read(cl &x)
{
    x = 0;
    int f = 0;
    char ch;
    ch = getchar();
    while (!isdigit(ch))
    {
        f = f | (ch == '-'), ch = getchar();
    }
    while (isdigit(ch))
    {
        x = (x << 1) + (x << 3) + (ch ^ 48);
        ch = getchar();
    }
    x = f ? -x : x;
    return;
}
template <class cl>
void put(cl x)
{
    if (x < 0)
        putchar('-'), x = -x;
    if (x > 9)
        put(x / 10);
    putchar(x % 10 + '0');
    return;
}
int mod;
int kpow(int a,int b)
{
   a=a%mod;
   int ans=1;
   while(b)
   {
       if(b&1)
       {
          ans=ans%mod*(a%mod)%mod;
       }
       b=b>>1;
       a=(a%mod)*(a%mod)%mod;
   }
   return ans;
}
int cnt=0;
// 害，搞了半天是把平方和和立方和用混了。
const int N=5e6+10;
const int nn=5e6;
int prime[N];
int phi[N];
bitset<N>g;
void fip()
{
    g[1]=g[0]=1;phi[1]=1;
    for(int i=2;i<=nn;i++)
    {
        if(g[i]==0)
        {
            prime[++cnt]=i;
            phi[i]=i-1;
        }
        for(int j=1;j<=cnt&&prime[j]<=nn/i;j++)
        {
            g[i*prime[j]]=1;
            if(i%prime[j]==0)
            {
                phi[i*prime[j]]=prime[j]*phi[i];
                break;
            }
            phi[i*prime[j]]=phi[i]*(prime[j]-1);
        }
    }
    for(int i=1;i<=nn;i++)
    {
        phi[i]=(phi[i-1]+i*i%mod*phi[i]%mod)%mod;
    }
}
int inv6;
unordered_map<int,int>mp;
int sum3(int x)//立方和公式
{
    x=x%mod;
    int tt=(x*(x+1)/2)%mod;
    return tt*tt%mod;
}
int sum2(int x)//平方和公式
{
    x=x%mod;
    return x*(x+1)%mod*(2*x+1)%mod*inv6%mod;
}
int du(int x)
{
   if(x<=nn) return phi[x];
   if(mp[x]) return mp[x];
   int res=sum3(x);
   int ree=0;
   for(int l=2,r;l<=x;l=r+1)
   {
       r=x/(x/l);
       ree=(ree+(sum2(r)-sum2(l-1)+mod)%mod*du(x/l)%mod)%mod;
   }
   res=(res-ree+mod)%mod;
   return mp[x]=res;
}
void to()
{
    int n;
    read(mod);
    read(n);
    inv6=kpow(6,mod-2);
    fip();
    int ans=0;
    for(int l=1,r;l<=n;l=r+1)
    {
        r=n/(n/l);
        ans=(ans+(du(r)-du(l-1)+mod)%mod*sum3(n/l)%mod)%mod;
    }
    put(ans);
    puts("");return ;
}
signed main()
{
    to();
    return 0;
}