版权声明:本文为博主原创文章,未经博主允许不得转载。 https://blog.csdn.net/hzj1054689699/article/details/83347031
Description
求
i=1∑nj=1∑nlcm(i,j)
n≤1010
Solution
将
lcm拆开
原式化为
i=1∑nj=1∑ngcd(i,j)ij
令
d=gcd(i,j),i=i/d,j=j/d
方便起见以下的除法默认为下取整
=d=1∑ndi=1∑dnj=1∑dnij[gcd(i,j)=1]
考虑对称性,有
=d=1∑nd×⎝⎛2⎝⎛i=1∑dnij=1∑ij[gcd(i,j)=1]⎠⎞−1⎠⎞
(因为1被计算了2次,因此最后要-1)
我们知道,小于等于n且与n互质的数的和是
2φ(n)∗n(考虑一个小于等于n且与n互质的数i,显然n-i也与n互质,这样总共有
2φ(n)组)
1除外,小于等于1且与1互质的数的和是1,因此还要加上一个
21
因此又有
=d=1∑nd×⎝⎛2⎝⎛i=1∑dn2φ(i)∗i2⎠⎞−1+21×2⎠⎞
=d=1∑nd×i=1∑dnφ(i)×i2
设
f(n)=φ(n)⋅n2
现在要做的就是求
f(n)的前缀和
S(n)
考虑杜教筛
我们发现
f∗id2=id3(其中
∗为狄利克雷卷积),即
d∣n∑φ(d)⋅d2⋅(dn)2=n2d∣n∑φ(d)=n3
那么
i=1∑ni3=j=1∑nd∣j∑f(dn)d2
接下来就是杜教筛套路做法,就不赘述了,直接贴上最终式子
令
sk(n)=i=1∑nik
S(n)=s3(n)−d=2∑nd2S(dn)
这就直接分块做,根据数学知识(狗头)
有自然数平方和公式
s2(n)=6n(n+1)(2n+1),s3(n)=(2n(n+1))2
那么最终就是
Ans(n)=d=1∑nd⋅S(dn)
Code
#include <cstdio>
#include <iostream>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <cstdlib>
#define fo(i,a,b) for(int i=a;i<=b;++i)
#define fod(i,a,b) for(int i=a;i>=b;--i)
#define LL long long
#define mo 1000000007
#define N 5000005
using namespace std;
LL n,sm[N],ny6,ny2,h[N/10];
int pr[N/10],phi[N];
bool bz[N];
LL ksm(LL k,LL n)
{
LL s=1;
for(;n;n>>=1,k=k*k%mo) if(n&1) s=s*k%mo;
return s;
}
LL squ(LL x)
{
x%=mo;
return(x*x%mo*(x+1)%mo*(x+1)%mo*ny2%mo*ny2%mo);
}
LL sqr(LL x)
{
x%=mo;
return (x*(x+1)%mo*(2*x+1)%mo*ny6%mo);
}
void prp(int n)
{
sm[1]=phi[1]=1;
fo(i,2,n)
{
if(!bz[i]) pr[++pr[0]]=i,phi[i]=i-1;
for(int j=1;j<=pr[0]&&i*pr[j]<=n;j++)
{
bz[i*pr[j]]=1;
if(i%pr[j]==0)
{
phi[i*pr[j]]=phi[i]*pr[j];
break;
}
phi[i*pr[j]]=phi[i]*(pr[j]-1);
}
sm[i]=(sm[i-1]+(LL)phi[i]*(LL)i%mo*(LL)i)%mo;
}
}
LL get(LL m)
{
if(m<=N-5) return sm[m];
int p=n/m;
if(h[p]) return h[p];
h[p]=squ(m);
LL d=2,d1,q;
while(d<=m)
{
q=m/d,d1=m/q;
h[p]=(h[p]-get(q)*(sqr(d1)-sqr(d-1))%mo+mo)%mo;
d=d1+1;
}
return h[p];
}
int main()
{
cin>>n;
prp(N-5);
ny6=ksm(6,mo-2),ny2=ksm(2,mo-2);
LL d=1,d1,ans=0;
while(d<=n)
{
LL m=n/d,d1=n/m;
ans=(ans+(d+d1)%mo*(d1-d+1)%mo*ny2%mo*((get(m)+mo)%mo))%mo;
d=d1+1;
}
printf("%lld\n",ans);
}