洛谷 P3768 简单的数学题 （莫比乌斯反演）

题意：求$(\sum_{i=1}^{n}\sum_{j=1}^{n}ijgcd(i,j))mod p$（p为质数，n<=1e10）

很显然，推式子。

$\sum_{i=1}^{n}\sum_{j=1}^{n}ijgcd(i,j)$

=$\sum_{d=1}^{n}\sum_{i=1}^{n}\sum_{j=1}^{n}ijd[gcd(i,j)==d]$

=$\sum_{d=1}^{n}d^3\sum_{i=1}^{\lfloor \frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor \frac{n}{d}\rfloor}ij[gcd(i,j)==d]$

=$\sum_{d=1}^{n}d^3\sum_{i=1}^{\lfloor \frac{n}{d}\rfloor}\mu(i)i^2S({\lfloor \frac{n}{id}\rfloor})^2,S(n)=(n+1)*n/2$

=$\sum_{T=1}^{n}S({\lfloor \frac{n}{T}\rfloor})^2\sum_{d|T}d^3(\frac{T}{d})^2\mu(\frac{T}{d})$

=$\sum_{T=1}^{n}S({\lfloor \frac{n}{T}\rfloor})^2T^2\sum_{d|T}d\mu(\frac{T}{d})$

由$\mu*id=\varphi $可得$\sum_{T=1}^{n}S({\lfloor \frac{n}{T}\rfloor})^2T^2\varphi (T)$

前面整除分块，只需要预处理$T^2\varphi(T)$ 前缀和即可。

由于n有1e10那么大，就需要用到非线性的求前缀和的方法，这里用到杜教筛，见代码。

#include <bits/stdc++.h>
#define ll long long
using namespace std;
const int N=1e6+5;
int pri[N],tot,phi[N],sum[N];
bool p[N];
ll n,MD,ans,inv6;
unordered_map<ll,int> w;
void init() {
    phi[1]=1;
    for(int i=2;i<N;i++) {
        if(!p[i]) phi[i]=i-1,pri[tot++]=i;
        for(int j=0;j<tot&&pri[j]*i<N;j++) {
            p[i*pri[j]]=true;
            if(i%pri[j]==0) {
                phi[i*pri[j]]=phi[i]*pri[j];
                break;
            }
            else phi[i*pri[j]]=phi[i]*phi[pri[j]];
        }
    }
    for(int i=1;i<N;i++) sum[i]=(sum[i-1]+1LL*i*i%MD*phi[i]%MD)%MD;
}
ll pre_3(ll x) {
    x%=MD;
    ll t=x*(x+1)/2%MD;
    return t*t%MD;
}
ll pre_2(ll x) {
    x%=MD;
    return x*(x+1)%MD*(2*x+1)%MD*inv6%MD;
}
int quick_pow(int x,int y) {
    int ans=1;
    while(y) {
        if(y&1) ans=1LL*ans*x%MD;
        y>>=1;
        x=1LL*x*x%MD;
    }
    return ans;
}
int cal(ll x) {
    if(x<N) return sum[x];
    if(w[x]) return w[x];
    ll ans=pre_3(x);
    for(ll l=2,r;l<=x;l=r+1) {
        r=x/(x/l);
        ans=(ans-(pre_2(r)-pre_2(l-1)+MD)%MD*cal(x/l)%MD+MD)%MD;
    }
    return w[x]=ans;
}
int main() {
    scanf("%lld%lld",&MD,&n);
    inv6=quick_pow(6,MD-2),init();
    for(ll l=1,r;l<=n;l=r+1) {
        r=n/(n/l);
        ans=(ans+pre_3(n/l)*(cal(r)-cal(l-1)+MD)%MD)%MD;
    }
    printf("%lld\n",ans);
    return 0;
}