【Luogu3768】簡單的數學題（莫比烏斯反演，杜教篩）

【Luogu3768】簡單的數學題（莫比烏斯反演，杜教篩）

題面

洛谷
\[求\sum_{i=1}^n\sum_{j=1}^nijgcd(i,j)\]
$ n<=10^9$

題解

很明顯的把\(gcd\)提出來
\[\sum_{d=1}^nd\sum_{i=1}^n\sum_{j=1}^nij[gcd(i,j)==d]\]
習慣性的提出來
\[\sum_{d=1}^nd^3\sum_{i=1}^{n/d}\sum_{j=1}^{n/d}ij[gcd(i,j)==1]\]
後面這玩意很明顯的來一發莫比烏斯反演
\[\sum_{d=1}^nd^3\sum_{i=1}^{n/d}\mu(i)i^2(1+2+...[\frac{n}{id}])^2\]
寫起來好麻煩呀
我就設\(sum(x)=1+2+3+...x\)
令\(T=id\)
提出來！

\[\sum_{T=1}^nsum(\frac{n}{T})^2\sum_{d|T}d^3\frac{T}{d}^2\mu(\frac{T}{d})\]

有些\(d\)能夠約掉
\[\sum_{T=1}^nsum(\frac{n}{T})^2T^2\sum_{d|T}d\mu(\frac{T}{d})\]

如今若是把後面給篩出來
能夠\(O(\sqrt n)\)求啦
如今，問題來了
\[T^2\sum_{d|T}d\mu(\frac{T}{d})\]怎麼算？？

考慮一個式子：
\[(id*\mu)(i)=\varphi(i)\]
也就是說，\(\mu\)和\(id(x)=x\)的狄利克雷卷積等於\(\varphi(i)\)
太神奇啦！！！

因此說，
\[T^2\sum_{d|T}d\mu(\frac{T}{d})=T^2\varphi(T)\]

令\[f(i)=i^2\varphi(i)\]
\[S(n)=\sum_{i=1}^nf(i)\]

杜教篩套路的式子拿出來
\[g(1)S(n)=\sum_{i=1}^n(g*f)(i)-\sum_{i=2}^ng(i)S(\frac{n}{i})\]
仍是發現有\(\varphi(i)\)的項
想到\[\sum_{d|i}\varphi(d)=i\]
因此令\(g(x)=x^2\)
因此
\[S(n)=\sum_{i=1}^n(g*f)(i)-\sum_{i=2}^ng(i)S(\frac{n}{i})\]

\[(g*f)(i)=\sum_{d|i}f(d)g(\frac{i}{d})=\sum_{d|i}d^2\varphi(d)\frac{i}{d}^2\]
\[=i^2\sum_{d|i}\varphi(d)=i^3\]

因此
\[S(n)=\sum_{i=1}^ni^3-\sum_{i=2}^ni^2S(\frac{n}{i})\]
根據小學奧數的經驗：
\(1^3+2^3+....n^3=(1+2+....n)^2=sum(n)^2\)

因此如今有：
\[ans=\sum_{T=1}^nsum(\frac{n}{T})^2\ T^2\sum_{d|T}d\mu(\frac{T}{d})\]
前面能夠數論分塊
後面用杜教篩能夠再非線性時間裏面求出前綴和
這道題目就搞定啦

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<set>
#include<map>
#include<vector>
#include<queue>
using namespace std;
int MAX=8000000;
#define MAXN 8000000
#define ll long long
inline ll read()
{
    ll x=0,t=1;char ch=getchar();
    while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
    if(ch=='-')t=-1,ch=getchar();
    while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
    return x*t;
}
ll MOD,n,inv6,inv2;
int pri[MAXN],tot;
ll phi[MAXN+10];
bool zs[MAXN+10];
map<ll,ll> M;
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;
}
void pre()
{
    zs[1]=true;phi[1]=1;
    for(int i=2;i<=MAX;++i)
    {
        if(!zs[i])pri[++tot]=i,phi[i]=i-1;
        for(int j=1;j<=tot&&i*pri[j]<=MAX;++j)
        {
            zs[i*pri[j]]=true;
            if(i%pri[j])phi[i*pri[j]]=1ll*phi[i]*phi[pri[j]]%MOD;
            else{phi[i*pri[j]]=1ll*phi[i]*pri[j]%MOD;break;}
        }
    }
    for(int i=1;i<=MAX;++i)phi[i]=(phi[i-1]+1ll*phi[i]*i%MOD*i%MOD)%MOD;
}
ll Sum(ll x){x%=MOD;return x*(x+1)%MOD*inv2%MOD;}
ll Sump(ll x){x%=MOD;return x*(x+1)%MOD*(x+x+1)%MOD*inv6%MOD;}
ll SF(ll x)
{
    if(x<=MAX)return phi[x];
    if(M[x])return M[x];
    ll ret=Sum(x);ret=ret*ret%MOD;
    for(ll i=2,j;i<=x;i=j+1)
    {
        j=x/(x/i);
        ll tt=(Sump(j)-Sump(i-1))%MOD;
        ret-=SF(x/i)*tt%MOD;
        ret%=MOD;
    }
    return M[x]=(ret+MOD)%MOD;
}
int main()
{
    MOD=read();n=read();
    MAX=min(1ll*MAX,n);
    inv2=fpow(2,MOD-2);
    inv6=fpow(6,MOD-2);
    pre();
    ll ans=0;
    for(ll i=1,j;i<=n;i=j+1)
    {
        j=n/(n/i);
        ll tt=Sum(n/i);tt=tt*tt%MOD;
        ll gg=(SF(j)-SF(i-1))%MOD;
        ans+=gg*tt%MOD;
        ans%=MOD;
    }
    printf("%lld\n",(ans+MOD)%MOD);
    return 0;
}