【题目链接】
题意分析
\[\sum_{i=1}^n\sum_{j=1}^n(ijgcd(i,j))=\sum_{i=1}^n\sum_{j=1}^n\sum_{d=1}^n(dij)[gcd(i,j)=d]=\sum_{d=1}^nd^3\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}(ij)[gcd(i,j)=1]
\]
由于
\[\sum_{i=1}^n\sum_{j=1}^n(ij)[gcd(i,j)=1]=\sum_{i=1}^n\sum_{j=1}^n\sum_{d|i,d|j}μ(d)(ij)=\sum_{d=1}^nμ(d)d^2\sum_{i=1}^{\frac{n}{d}}\sum_{j=1}^{\frac{n}{d}}(ij)=\sum_{d=1}^nμ(d)d^2(1+2+...+\lfloor\frac{n}{d}\rfloor)^2
\]
所以原式转化为
\[\sum_{d=1}^nd^3\sum_{p=1}^{\lfloor\frac{n}{d}\rfloor}μ(p)p^2(1+2+...+\lfloor\frac{n}{pd}\rfloor)^2
\]
令
\[T=pd
\]
那么再一次转化为
\[\sum_{T=1}^n(1+2+...+\lfloor\frac{n}{T}\rfloor)^2\sum_{d|T}d^3(\frac{T}{d})^2μ(\frac{T}{d})=\sum_{T=1}^n(1+2+...+\lfloor\frac{n}{T}\rfloor)^2T^2\sum_{d|T}dμ(\frac{T}{d})
\]
由于
\[Id(u)=u\ \ \ \ \ e(u)=[u=1]\ \ \ \ \ I(u)=1
\]
同时
\[μ*I=e\ \ \ φ*I=Id
\]
\[μ*I*φ=e*φ ⇒ μ*Id=φ
\]
所以原式再一次转化为
\[\sum_{T=1}^n(1+2+...+\lfloor\frac{n}{T}\rfloor)^2T^2φ(T)
\]
左半部分使用等差数列求和公式就可以了 右半部分由于是积性函数 可以使用杜教筛
由于
\[\sum_{d|n}d^2φ(d)(\frac{n}{d})^2=n^2\sum_{d|n}φ(d)=n^3
\]
所以
\[S(n)=\sum_{i=1}^ni^2μ(i)
\]
\[S(n)=\sum_{i=1}^nn^3-\sum_{d=2}^nd^2S(\lfloor\frac{n}{d}\rfloor)
\]
好了 到此结束
CODE:
#include<bits/stdc++.h>
#define M 10000080
using namespace std;
int tot;
long long mod,n,ans;
long long inv2,inv4,inv6;
bool mark[M];
int prime[M],phi[M];
long long sum[M];
map<long long,long long> vis;
long long qpow(long long x,long long y)
{long long tmp=1;for(;y;y>>=1,x=x*x%mod) if(y&1) tmp=tmp*x%mod;return tmp;}
void pre()
{
phi[1]=1;
for(int i=2;i<=10000000;++i)
{
if(!mark[i]) {prime[++tot]=i;phi[i]=i-1;}
for(int j=1;j<=tot&&prime[j]*i<=10000000;++j)
{
mark[prime[j]*i]=1;
if(i%prime[j]==0)
{
phi[prime[j]*i]=prime[j]*phi[i];
break;
}
else phi[prime[j]*i]=(prime[j]-1)*phi[i];
}
}
for(int i=1;i<=10000000;++i) sum[i]=(sum[i-1]+(long long)i*(long long)i%mod*phi[i]%mod+mod)%mod;
inv2=qpow(2,mod-2);
inv4=qpow(4,mod-2);
inv6=qpow(6,mod-2);
}
long long GetSum(long long x)
{
x%=mod;
long long tmp=(1+x)*x%mod*inv2%mod;
return tmp*tmp%mod;
}
long long GetSum2(long long x)
{x%=mod;return x*(x+1)%mod*(2*x%mod+1)%mod*inv6%mod;}
long long GetSum3(long long x)
{x%=mod;return ((x*x%mod)*((x+1)*(x+1)%mod)%mod)*inv4%mod;}
long long Get_Sum(long long x)
{
if(x<=10000000) return sum[x];
if(vis.count(x)) return vis[x];
long long tmp=GetSum3(x);
for(long long l=2,r=0;l<=x;l=r+1)
{
r=x/(x/l);
tmp=(tmp-(GetSum2(r)-GetSum2(l-1)+mod)%mod*Get_Sum(x/l)%mod+mod)%mod;
}
return vis[x]=tmp;
}
int main()
{
scanf("%lld%lld",&mod,&n);pre();
for(long long l=1,r=0;l<=n;l=r+1)
{
r=n/(n/l);
ans=(ans+(GetSum(n/l)*((Get_Sum(r)-Get_Sum(l-1)+mod)%mod))%mod+mod)%mod;
}
printf("%lld\n",ans);
return 0;
}