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;
}