51nod 1238 最小公倍数之和 V3

我们首先计算

f(n)===i=1nlcm(i,n)d|ni=1ndnidd[gcd(i,nd)=1]nd|ndφ(d)+e(d)2

于是
==i=1nf(i)12d=1ndφ(d)+e(d)i=1ndid12d=1nd2φ(d)+e(d)i=1ndi

我们只需要杜教筛求出 id2φ 的前缀和就可以了。
(id2φ)×id2=id2(φ×1)=id3
因此
=i=1ni2φ(i)i=1ni3i=2ni2d=1nid2φ(d)

求出 f <script id="MathJax-Element-984" type="math/tex">f</script>的前缀和之后不难得到答案。

#include<cstdio>
#include<cmath>
#include<algorithm>
using namespace std;
#define LL long long
const int p=1000000007,maxn=5000010,mod=1000007,inv2=500000004,inv6=166666668;
int phi[maxn],prm[maxn],sum[maxn],tot,m;
LL n;
/*void pause()
{
    int x;
    x=1;
}*/
int inc(int x,int y)
{
    x+=y;
    return x>=p?x-p:x;
}
int dec(int x,int y)
{
    x-=y;
    return x<0?x+p:x;
}
struct hashset
{
    int fir[mod+10],ne[mod*5],ans[mod*5],tot;
    LL val[mod*5];
    int find(LL n)
    {
        for (int i=fir[n%mod];i;i=ne[i])
            if (val[i]==n) return ans[i];
        return -1;
    }
    void ins(LL n,int x)
    {
        if (tot>mod*5-5) return;
        int y=n%mod;
        ne[++tot]=fir[y];
        fir[y]=tot;
        val[tot]=n;
        ans[tot]=x;
    }
}h;
int sum1(LL n)
{
    n%=p;
    return n*inc(n,1)%p*inv2%p;
}
int sum2(LL n)
{
    n%=p;
    return n*inc(n,1)%p*inc(inc(n,n),1)%p*inv6%p;
}
int sum3(LL n)
{
    int x=sum1(n);
    return (LL)x*x%p;
}
int calc(LL n)
{
    if (n<=m) return sum[n];
    int ret=h.find(n),ans,x,y;
    if (ret!=-1) return ret;
    ans=sum3(n);
    for (LL i=2,j;i<=n;i=j+1)
    {
        //if (n==10000000000LL&&j==3333333333LL) pause();
        j=n/(n/i);
        x=dec(sum2(j),sum2(i-1));
        y=calc(n/i);
        ans=dec(ans,(LL)x*y%p);
        //if (x<0||y<0||ans<0) printf("no:%lld,%lld,%lld\n",n,i,j);
    }
    h.ins(n,ans);
    return ans;
}
int main()
{
    int ans=0,x,y;
    scanf("%lld",&n);
    m=pow(n,2.0/3);
    //m=100;
    phi[1]=sum[1]=1;
    for (int i=2;i<=m;i++)
    {
        if (!phi[i])
        {
            prm[++tot]=i;
            phi[i]=i-1;
        }
        for (int j=1;j<=tot&&(LL)i*prm[j]<=m;j++)
            if (i%prm[j]) phi[i*prm[j]]=phi[i]*(prm[j]-1);
            else
            {
                phi[i*prm[j]]=phi[i]*prm[j];
                break;
            }
        sum[i]=inc(sum[i-1],(LL)phi[i]*i%p*i%p);
    }
    for (LL i=1,j;i<=n;i=j+1)
    {
        j=n/(n/i);
        x=inc(dec(calc(j),calc(i-1)),i==1);
        y=sum1(n/i);
        ans=inc(ans,(LL)x*y%p);
    }
    printf("%d\n",dec(ans,sum1(n)));
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值