Description
Ans=∑i=ab∑j=1ilcm(i,j)i
Solution 【莫比乌斯反演|【51NOD 1227】平均最小公倍数】题目要我们求的就是这个嘛:
Ans=∑i=1n∑j=1ijgcd(i,j)
转化一下
Ans=∑i=1n(φ(i)?i/2)??ni?
( 先枚举两个互质的数,再算它们的倍数)(前面有phi的括号是互质的数的和)
用杜教筛处理 ∑xi=1φ(i)?i (下面的n和上面的不一样)
∑d=1n∑i|d(φ(i)?i)?di
∑d=1nd?∑i|dφ(i)
∑d=1nd2
设 Sd=d?(d+1)?(2?d+1)/6
回到刚才的式子,转化一下(i*j即为之前的d)
∑i=1ni?∑j=1?ni?φ(j)?j=Sn
∑i=1nφ(i)?i=Sn?∑i=2ni?∑j=1?ni?φ(j)?j
回到Ans的式子,用杜教筛+分块即可
复杂度: O(n23)
Code
#include
#include
#define fo(i,a,b) for(int i=a;
i<=b;
i++)
using namespace std;
typedef long long LL;
const int N=3876656,mo=1000000007,M=200397;
const LL eni=500000004,sni=166666668;
int n,n1;
bool prz[N+10];
int pr[N/2];
int phi[N+10];
int Hx[M][2];
int HX(int q)
{
int i=q%M;
while(Hx[i][0]&&Hx[i][0]!=q)i=(i+1)%M;
return i;
}
LL SM(LL s,LL t){return (s+t)*(t-s+1)%mo*eni%mo;
}
LL Gphi(int q)
{
if(q<=N)return phi[q];
int t=HX(q);
if(Hx[t][0])return Hx[t][1];
Hx[t][0]=q;
LL ans=0;
for(int i=2,nx;
i<=q;
i=nx+1)
{
nx=q/(q/i);
ans=(ans+SM(i,nx)*Gphi(q/i)%mo)%mo;
}
q%=mo;
return Hx[t][1]=(LL)q*(q+1)%mo*(2*q+1)%mo*sni%mo-ans;
}
LL Gans(int n)
{
int ans=n%mo;
for(int i=2,nx;
i<=n;
i=nx+1)
{
nx=n/(n/i);
ans=(ans+(LL)(n/i)*(Gphi(nx)-Gphi(i-1))%mo*eni%mo)%mo;
}
return ans;
}
int main()
{
phi[1]=1;
fo(i,2,N)
{
if(!prz[i])pr[++pr[0]]=i,phi[i]=i-1;
fo(j,1,pr[0])
{
int t=pr[j]*i;
if(t>N)break;
prz[t]=1;
phi[t]=phi[i]*pr[j];
if(i%pr[j]==0)break;
phi[t]=phi[i]*(pr[j]-1);
}
}
fo(i,2,N)phi[i]=((LL)phi[i]*((LL)i)%mo+phi[i-1])%mo;
scanf("%d%d",&n1,&n);
printf("%lld\n",(Gans(n)-Gans(n1-1)+mo)%mo);
return 0;
}
推荐阅读
- 数论|Codeforces 235E Number Challenge 莫比乌斯反演+数论
- 杜教筛|hdu 6706 huntian oy 杜教筛
- 数论|HDU 6706 huntian oy(杜教筛)
- 刷题小结|51nod 1227 平均最小公倍数
- 筛法|[JZOJ5134][SDOI省队集训2017]三元组
- Code|CodeForces 839 D.Winter is here(莫比乌斯反演+组合数学)
- 莫比乌斯反演|Codeforces 235 E Number Challenge(莫比乌斯反演)
- 杜教筛|[学习笔记] 杜教筛 (51nod 1244+1227 +1237 +1238+1239) - 数论
- [51Nod 1238] 最小公倍数之和 (恶心杜教筛)