洛谷3678：简单的数学题（画柿子+杜教筛）

最近看了《魔法少女伊莉雅~雪下的誓言》，
这部番除了百合与萌以外还有其他东西，
就和当初看《fate》的我一样。

题面
题意：求

∑i=1n∑j=1nij(i,j)

n≤1010

这就是数论题，自己数学课上半个小时推出来，看Kscla的博客一个下午才懂。

根据我两页草稿纸的结果，原式为

∑i=1ni3∑j=1⌊ni⌋j2φ(j)

对右边进行杜教筛就可以了。
（我今天本来想学洲阁筛的，结果沉迷到这题了）

#include <iostream>#include <fstream>#include <algorithm>#include <cmath>#include <ctime>#include <cstdio>#include <cstdlib>#include <cstring>#include <map>using namespace std;#define mmst(a, b) memset(a, b, sizeof(a))#define mmcp(a, b) memcpy(a, b, sizeof(b))typedef long long LL;const int N=10100000,nn=10000000;LL n,p,ans;int phi[N],prime[N],num;bool b[N];map<LL,LL> mp;LL get1(LL x){    LL y=x+1;    if(!(y&1))    y=y/2;    else    x=x/2;    x%=p;    y%=p;    return x*x%p*y%p*y%p;}LL get2(LL x){    LL y=x+1,z=2*x+1;    bool f2=false,f3=false;    if ( !f2 && !(x%2) ) x/=2,f2=true;    if ( !f3 && !(x%3) ) x/=3,f3=true;    if ( !f2 && !(y%2) ) y/=2,f2=true;    if ( !f3 && !(y%3) ) y/=3,f3=true;    if ( !f2 && !(z%2) ) z/=2,f2=true;    if ( !f3 && !(z%3) ) z/=3,f3=true;    x%=p;    y%=p;    z%=p;    return x*y%p*z%p;}LL work(LL x){    if(x<=nn)    return phi[x];    if(mp.count(x))    return mp[x];    LL res=get1(x),last;    for(LL i=2;i<=x;i=last+1)    {        last=x/(x/i);        res=(res-(get2(last)-get2(i-1)+p)%p*work(x/i)%p+p)%p;    }    mp[x]=res;    return res;}int main(){    cin>>p>>n;    phi[1]=1;    for(int i=2;i<=nn;i++)    {        if(!b[i])        {            phi[i]=i-1;            prime[++num]=i;        }        for(int j=1;j<=num&&prime[j]*i<=nn;j++)        {            b[i*prime[j]]=1;            if(i%prime[j]==0)            {                phi[i*prime[j]]=phi[i]*prime[j];                break;            }            phi[i*prime[j]]=phi[i]*(prime[j]-1);        }    }    for(int i=1;i<=nn;i++)    phi[i]=((LL)phi[i-1]+(LL)phi[i]*i%p*i%p)%p;    LL last;    for(LL i=1;i<=n;i=last+1)    {        last=n/(n/i);        ans=(ans+(get1(last)-get1(i-1)+p)%p*work(n/i)%p)%p;    }    cout<<ans<<endl;    return 0;}