P3768 简单的数学题 解题报告
1.前置知识：
欧拉筛，莫比乌斯反演，狄利克雷卷积，杜教筛
2.莫反+狄利克雷卷积+欧拉筛
看到这道题的瞬间，按照DNA来一个莫反
\[\begin{aligned}
&\sum_{i=1}^{n}\sum_{j=1}^{n}{ij\gcd(i,j)}\\
=&\sum_{k=1}^{n}k\sum_{i=1}^{n}\sum_{j=1}^{n}{ij[\gcd(i,j)=k]}\\
=&\sum_{k=1}^{n}{k^3}\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{k}\rfloor}{ij[\gcd(i,j)=1]}\\
=&\sum_{k=1}^{n}{k^3}\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{k}\rfloor}{ij\sum_{d|\gcd(i,j)}{\mu(d)}}\\
=&\sum_{k=1}^{n}{k^3}\sum_{d=1}^{\lfloor\frac{n}{k}\rfloor}{\mu(d)}\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{k}\rfloor}{ij[d|\gcd(i,j)]}\\
=&\sum_{k=1}^{n}{k^3}\sum_{d=1}^{\lfloor\frac{n}{k}\rfloor}{\mu(d)d^2}\sum_{i=1}^{\lfloor\frac{n}{kd}\rfloor}i\sum_{j=1}^{\lfloor\frac{n}{kd}\rfloor}{j}\\
\end{aligned}
\]
根据莫反惯例，设\(T=kd,sum(x)=\sum_{i=1}^{n}i\)
则原式可化为：
\[\begin{aligned}
&\sum_{T=1}^{n}sum(\lfloor\frac{n}{T}\rfloor)^2\sum_{d|T}{d^2\mu(d)(\frac{T}{d})^3}\\
=&\sum_{T=1}^{n}sum(\lfloor\frac{n}{T}\rfloor)^2T^2\sum_{d|T}{\mu(d)(\frac{T}{d})}\\
\end{aligned}
\]
设
\[F(T)=T^2\sum_{d|T}{\mu(d)(\frac{T}{d})}
\]
根据狄利克雷卷积可知，
\[\begin{aligned}
F(T)&=T^2\sum_{d|T}{\mu(d)(\frac{T}{d})}\\
&=T^2*((id\times\mu)(T))\\
&=T^2*((I\times\varphi\times\mu)(T))\\
&=T^2*((I\times\mu\times\varphi)(T))\\
&=T^2*((I\times\mu\times\varphi)(T))\\
&=T^2*((\epsilon\times\varphi)(T))\\
&=T^2 \varphi(T)\\
\end{aligned}
\]
原式可以化成
\[\sum_{T=1}^{n}sum(\lfloor\frac{n}{T}\rfloor)^2F(T)
\]
F函数可以用欧拉筛线性求出，求出F的前缀和后式子可以在\(O(\sqrt{n})\)的复杂度内求出
现在，这道题可以做到\(O(n)\)求解，能拿60分。
3.杜教筛优化
但是,n的范围达到了1e10,因此我们必须通过比欧拉筛的方式求出F的前缀和，比如说杜教筛
设
\[S(n)=\sum_{i=1}^nf(i)
\]
\[g(n)=n^2
\]
\[h(n)=(f\times g)(n)=\sum_{d|n}{d^2\varphi(d)(\frac{n}{d})^2} =n^2*\sum_{d|n}\varphi(d)=n^3
\]
\[sump(n)=\sum_{i=1}^ni^2
\]
由伟大的杜教筛式子得：
\[S(n)g(1)=\sum_{i=1}^nh(i)-\sum_{i=2}^nS(\lfloor \frac{n}{i}\rfloor)g(i)
\]
代入定义可得
\[S(n)=\sum_{i=1}^ni^3-\sum_{i=2}^ni^2S(\lfloor \frac{n}{i}\rfloor)
\]
由一个小学奥数和两个定积分得：
\[sum(n)=\frac{n(n+1)}{2}
\]
\[sump(n)=\frac{n^2(n+1)^2}{4}
\]
\[\sum_{i=1}^ni^3=sum(n)^2
\]
那么我们现在可以用杜教筛求出S(n)了。
杜教筛的时间复杂度为\(O(n^\frac{2}{3})\),属于O(能过)的范围内。
4.代码
#include<bits/stdc++.h>
using namespace std;
#define LL long long
#define maxn 5000005
bool vis[maxn];
int cnt,prime[maxn],phi[maxn];
LL s1[maxn],n,mod,inv6,inv2;
map <LL,LL> s2;
LL sum(LL x)
{
x%=mod;
return x*(1+x)%mod*inv2%mod;
}
LL sump(LL x)
{
x%=mod;
return x*(1+x)%mod*(1+x+x)%mod*inv6%mod;
}
void oula(LL n)
{
phi[1]=1;
for(int i=2;i<=n;i++)
{
if(!vis[i])
{
prime[++cnt]=i;
phi[i]=i-1;
}
for(int j=1;j<=cnt&&i*prime[j]<=n;j++)
{
vis[i*prime[j]]=1;
if(i%prime[j])
{
phi[i*prime[j]]=phi[i]*(prime[j]-1);
}else
{
phi[i*prime[j]]=phi[i]*prime[j];
}
}
}
for(int i=1;i<=n;i++)
{
s1[i]=s1[i-1]+1ll*i*i%mod*phi[i]%mod;
s1[i]%=mod;
}
}
LL dujiao(LL n)
{
if(n<=maxn-5)return s1[n];
if(s2[n])return s2[n];
LL ret = sum(n);
ret=ret*ret%mod;
for(LL l=2,r=0;l<=n;l=r+1)
{
r=n/(n/l);
LL tt=(sump(r)-sump(l-1))%mod;
ret-=tt*dujiao(n/l)%mod;
ret%=mod;
}
ret=(ret+mod)%mod;
s2[n]=ret;
return ret;
}
LL solve(LL n)
{
LL ret=0;
for(LL l=1,r=0;l<=n;l=r+1)
{
r=(n/(n/l));
LL tt=sum(n/l);tt=tt*tt%mod;
LL gg=dujiao(r)-dujiao(l-1);gg%=mod;
ret+=gg*tt%mod;
ret%=mod;
}
ret=(ret+mod)%mod;
return ret;
}
LL fpow(LL a,LL b)
{
LL s=1;
while(b){if(b&1)s=s*a%mod;a=a*a%mod;b>>=1;}
return s;
}
int main()
{
scanf("%lld%lld",&mod,&n);
inv2=fpow(2,mod-2);
inv6=fpow(6,mod-2);
oula(maxn-5);
printf("%lld",solve(n));
return 0;
}
