[luoguP3768]简单的数学题

传送门

惯用套路

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

然后简单的莫比乌斯反演一下，得到

$$f(d)=d^3\sum_{T=1}^{n/d}\mu(T)T^2\sum_{i=1}^{n/Td}\sum_{j=1}^{n/Td}ij\\ $$

然后知道答案式

$$ans=\sum_{d=1}^{n}f(d)=\sum_{d=1}^{n}d^3\sum_{T=1}^{n/d}\mu(T)T^2\sum_{i=1}^{n/Td}\sum_{j=1}^{n/Td}ij\\ $$

这个时候去算，复杂度是$O(n)$的，无法通过

然后设$k=Td$

$$ans=\sum_{k=1}^{n}\sum_{i=1}^{n/k}\sum_{j=1}^{n/k}ij\sum_{d|k}\mu(\frac{k}{d})(\frac{k}{d})^2d^3\\ ans=\sum_{k=1}^{n}\sum_{i=1}^{n/k}\sum_{j=1}^{n/k}ijk^2\sum_{d|k}\mu(\frac{k}{d})d\\ $$

然后现在只要能够筛出$k^2\sum_{d|k}\mu(\frac{k}{d})d$，就可以$O(\sqrt{n})$解决了

我们设$g(k)=k^2\sum_{d|k}\mu(\frac{k}{d})d$

现在考虑杜教筛，我们有一个很经典的式子

$$\mu*id=\varphi\\ g(k)=k^2\varphi(k)$$

则我们设$S(n)=\sum_{i=1}^{n}g(i)$，套上杜教筛的式子

$$S(n)h(1)=\sum_{i=1}^{n}(h*g)-\sum_{i=2}^{n}S(\lfloor\frac{n}{i}\rfloor)h(i)$$

发现

$$h*g=\sum_{d|n}h(\lfloor\frac{n}{d}\rfloor)g(d)=\sum_{d|n}h(\lfloor\frac{n}{d}\rfloor)d^2\varphi(d)\\ $$

考虑如何把$\varphi(d)$消去

因为

$$\sum_{d|n}\varphi(d)=n$$

则

$$h*g=\sum_{d|n}h(\lfloor\frac{n}{d}\rfloor)d^2n\\ $$

然后我们可以考虑将$h(i)$设为$i^2$

则

$$h*g=n^3\\ S(n)=\sum_{i=1}^{n}i^3-\sum_{i=2}^{n}i^2S(\lfloor\frac{n}{i}\rfloor)$$

这个就已经可以处理了

然后对于这个答案式太复杂了

$$ans=\sum_{k=1}^{n}\sum_{i=1}^{n/k}\sum_{j=1}^{n/k}ijk^2\sum_{d|k}\mu(\frac{k}{d})d\\ $$

我们设$s(x)=\sum_{i=1}^{x}\sum_{j=1}^{x}ij$，这个是可以用等差数列解决的

然后答案式就简洁多了

$$ans=\sum_{k=1}^{n}s(\lfloor\frac{n}{k}\rfloor)g(k)$$

代码：

#include<cstdio>
#include<algorithm>
#include<iostream>
#include<cstring>
#include<tr1/unordered_map>
using namespace std;
void read(int &x)
{
    char ch;bool ok;
    for(ok=0,ch=getchar();!isdigit(ch);ch=getchar())if(ch=='-')ok=1;
    for(x=0;isdigit(ch);x=x*10+ch-'0',ch=getchar());if(ok)x=-x;
}
#define rg register
#define ll long long
const int maxn=3e6+10;
int mod,ans,phi[maxn],pri[maxn],tot,inv2,inv6,f[maxn];
ll n;bool vis[maxn];
tr1::unordered_map<ll,int>mp;
int mul(int x,int y){return 1ll*x*y-1ll*x*y/mod*mod;}
int add(int x,int y){return x+y>=mod?x+y-mod:x+y;}
int del(int x,int y){return x-y<0?x-y+mod:x-y;}
int mi(int a,int b){
    int ans=1;
    while(b){
        if(b&1)ans=mul(ans,a);
        b>>=1,a=mul(a,a);
    }
    return ans;
}
void prepare()
{
    phi[1]=1;
    for(rg int i=2;i<=3e6;i++){
        if(!vis[i])pri[++tot]=i,phi[i]=i-1;
        for(rg int j=1;j<=tot&&pri[j]*i<=3e6;j++){
            vis[pri[j]*i]=1;
            if(!(i%pri[j])){phi[pri[j]*i]=mul(pri[j],phi[i]);break;}
            else phi[pri[j]*i]=mul(phi[pri[j]],phi[i]);
        }
    }
    for(rg int i=1;i<=3e6;i++)f[i]=mul(mul(i,i),phi[i]);
    for(rg int i=1;i<=3e6;i++)f[i]=add(f[i],f[i-1]);
}
int s(int n) {return mul(mul(mul(n*2+1,n+1),n),inv6);}
int square(ll l,ll r){return del(s(r%mod),s((l-1)%mod));}
int sum(ll n){return mul(mul((n+1)%mod,n%mod),inv2);}
int sum1(ll n){return mul(mul(n%mod,(n+1)%mod),inv2);}
int get_phi(ll n){
    if(n<=3e6)return f[n];
    if(mp[n])return mp[n];
    int ans=mul(sum1(n),sum1(n));
    for(rg ll i=2,j;i<=n;i=j+1)j=n/(n/i),ans=del(ans,mul(square(i,j),get_phi(n/i)));
    return mp[n]=ans;
}
int main()
{
    read(mod),scanf("%lld",&n);
    prepare(),inv2=mi(2,mod-2),inv6=mi(6,mod-2);
    for(rg ll i=1,j,t;i<=n;i=j+1){
        j=n/(n/i),t=mul(sum(n/i),sum(n/i));
        ans=add(ans,mul(t,del(get_phi(j),get_phi(i-1))));
    }
    printf("%d\n",ans);
}