51nod 1220 约数之和

最新推荐文章于 2022-03-02 20:24:41 发布

fyc_kabuto

最新推荐文章于 2022-03-02 20:24:41 发布

阅读量341

点赞数 1

分类专栏：莫比乌斯反演

本文链接：https://blog.csdn.net/qq_36808030/article/details/79035343

版权

莫比乌斯反演专栏收录该内容

17 篇文章 0 订阅

订阅专栏

题意：

求：

\sum i n \sum j n d (i j)

$\sum_i^n\sum_j^nd(ij)$

题解：

首先有这个东西： $nm$ 的因子是 $\sum_{p|n}\sum_{q|m}\frac npq[(p,q)=1]$
所以要求

\sum i n \sum j n \sum p | i \sum q | j n p q [(p, q) = 1]

$\sum_i^n\sum_j^n\sum_{p|i}\sum_{q|j}\frac npq[(p,q)=1]$
然后上反演。

\sum d n μ (d) d \sum p ⌊ n d ⌋ s u m (⌊ n p d ⌋) \sum q ⌊ n d ⌋ q ⌊ n d q ⌋

$\sum_d^n\mu(d)d\sum_p^{\lfloor\frac nd\rfloor}sum(\lfloor\frac n{pd}\rfloor)\sum_q^{\lfloor\frac nd\rfloor}q\lfloor\frac n{dq}\rfloor$
挺好推的不写过程了。
然后杜教筛

∑ndμ(d)d $\sum_d^n\mu(d)d$
code：

#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#define LL long long
const LL mod=1000000007;
using namespace std;
LL mu[1000010],prime[1000010],pr=0,n;
bool v[1000010];
LL Sum[100010];
void pre()
{
    mu[1]=1;memset(v,true,sizeof(v));
    for(LL i=2;i<=1000000;i++)
    {
        if(v[i]) prime[++pr]=i,mu[i]=-1;
        for(LL j=1;j<=pr&&i*prime[j]<=1000000;j++)
        {
            v[i*prime[j]]=false;
            if(i%prime[j]==0){mu[i*prime[j]]=0;break;}
            mu[i*prime[j]]=-mu[i];
        }
    }
    for(LL i=1;i<=1000000;i++) mu[i]=(mu[i-1]+mu[i]*i)%mod;
}
LL sum(LL x) {return x*(x+1)%mod*500000004%mod;}
LL solve(LL x)
{
    if(x<=1000000) return mu[x];
    if(Sum[n/x]!=0) return Sum[n/x];
    LL ans=1,j;
    for(LL i=2;i<=x;i=j+1)
    {
        j=x/(x/i);
        ans-=(sum(j)-sum(i-1))%mod*solve(x/i)%mod;
        ans%=mod;
    }
    Sum[n/x]=ans;
    return ans;
}
LL calcp(LL n)
{
    LL ans=0,j;
    for(LL i=1;i<=n;i=j+1)
    {
        j=n/(n/i);
        ans=(ans+(j-i+1)*sum(n/i)%mod)%mod;
    }
    return ans;
}
LL calcq(LL n)
{
    LL ans=0,j;
    for(LL i=1;i<=n;i=j+1)
    {
        j=n/(n/i);
        ans=(ans+(sum(j)-sum(i-1))%mod*(n/i)%mod)%mod;
    }
    return ans;
}
int main()
{
    pre();
    scanf("%lld",&n);
    LL ans=0,j;
    for(LL i=1;i<=n;i=j+1)
    {
        j=n/(n/i);
        ans=(ans+(solve(j)-solve(i-1))%mod*calcp(n/i)%mod*calcq(n/i)%mod)%mod;
    }
    printf("%lld",(ans+mod)%mod);
}