luoguP3768 简单的数学题

简单的数学题

题意

求 \((\sum_{i=1}^n\sum_{j=1}^n ijgcd(i,j))~mod~p\)

解法

很容易可以得出式子：

$ \sum_{d=1}^{n} \varphi(d) d^2 ~\sum_{i=1}^{n} i^3$

考虑如何求 \(\sum_{d=1}^{n} \varphi(d) d^2\)

以下内容来自这里

令\(f(i)= i^2 \varphi(i)\)

\(S(n)=\sum_{i=1}^nf(i)\)

杜教筛套路的式子拿出来

\(S(n)=\sum_{i=1}^n(g*f)(i)-\sum_{i=2}^ng(i)S(\frac{n}{i})g(1)\)

还是发现有 \(\varphi(i)\) 的项

想到 \(\sum_{d|i}\varphi(d)=id\)

所以令 \(g(x)=x^2\)

因为 \(S(n)=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})\)

代码

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn=5e6+15;
const int lim=5e6;
ll mod;
ll poww(ll a,int b){
    ll ans=1;
    while(b){
        if(b&1)ans=ans*a%mod;
        a=a*a%mod;
        b>>=1;
    }
    return ans;
}

bool nop[maxn];
vector<int> p;
int phi[maxn];
ll s[maxn],inv_6;

void _init(int n) {
    inv_6=poww(6,mod-2);
    nop[1]=phi[1]=1;
    
    for(int i=2;i<=n;i++) {
        if(!nop[i]) phi[i]=i-1,p.push_back(i);
        for(int j=0;i*p[j]<=n&&j<p.size();j++) {
            nop[i*p[j]]=1;
            if(i%p[j]==0) {phi[i*p[j]]=phi[i]*p[j];break;}
            phi[i*p[j]]=phi[i]*(p[j]-1);
        }
    }
    for(int i=1;i<=n;i++) s[i]=(s[i-1]+1ll*i*1ll*i%mod*phi[i]%mod)%mod;
}

ll s2(ll n) {
    n%=mod;
    return (n*(n+1)%mod*((2*n%mod+1)%mod)%mod*inv_6)%mod;
}

ll s3(ll n) {
    n%=mod;
    ll k=n*(n+1)/2%mod;
    return k*k%mod;
}

map<ll,ll> du;
ll _du(ll n) {
    if(n<=lim) return s[n];
    if(du.count(n)) return du[n];
    ll ans=s3(n);
    for(ll l=2,r;l<=n;l=r+1) {
        r=n/(n/l);
        ans=(ans- _du(n/l)*((s2(r)-s2(l-1)+mod)%mod)%mod +mod)%mod;
    }
    return du[n]=ans;
}

ll n;
int main() {
    scanf("%lld%lld",&mod,&n);	
    _init(lim);
    ll ans=0;
    for(ll l=1,r;l<=n;l=r+1) {
        r=n/(n/l);
        ans=(ans+ s3(n/l)*((_du(r)-_du(l-1)+mod)%mod)%mod )%mod;
    }
    printf("%lld\n",ans);
    return 0;
}