Description
求
∑i=1n∑d|igcd(d,id)
对于100%的数据, n≤1011 。
Analysis
方法一
这是一个不需要莫比乌斯反演的方法。
∑i=1n∑d|igcd(d,id)
=∑i∗j≤ngcd(i,j)
=∑i∗j≤n∑d|gcd(i,j)φ(d)
这里用到了一个 φ 的性质, ∑d|nφ(d)=n ,具体证明可以看 我的博客。
=∑i∗j≤n∑d|i⋂d|jφ(d)
调换主体,得
=∑d=1n√φ(d)∑i∗j≤⌊nd2⌋1
=∑d=1n√φ(d)∑i=1⌊nd2⌋⌊nd2i⌋
随着 d 的增长,
方法二
这次我们用莫比乌斯反演。
∑i=1n∑d|igcd(d,id)
=∑i=1n∑j=1⌊ni⌋gcd(i,j)
设 f(d) 表示上式中满足 gcd(i,j)=d 的 i,j 有多少对;
g(d) 表示上式中满足 d|gcd(i,j) 的 i,j 有多少对。
首先有
g(d)=∑i=1⌊nd⌋f(d∗i)
反演
f(d)=∑i=1⌊nd⌋g(d∗i)∗μ(i)
而我们要求的
Ans=∑d=1nd∗f(d)
=∑d=1nd∑i=1⌊nd⌋g(d∗i)∗μ(i)
设 T=d∗i ,再调换主体
=∑T=1ng(T)∑i|TTi∗μ(i)
但其实 T 只用枚举到
那么 g(T) 怎么求呢?
若 d|gcd(i,j) ,则可以设 dx,dy 分别为任意的 i,j ,因为必有 d|gcd(dx,dy)
所以我们枚举 x ,因为
g(T)=∑x=1⌊nT2⌋⌊nT2x⌋
所以
Ans=∑T=1n√∑x=1⌊nT2⌋⌊nT2x⌋∗∑i|TTi∗μ(i)
这个东西直接做,预处理一些东西,就可以AC了。
但我们的追求更远。
上面这个式子有没有一点眼熟?
我们把方法一最后得出的结论搬下来
Ans=∑d=1n√φ(d)∑i=1⌊nd2⌋⌊nd2i⌋
,
两个式子比较一下,许多东西非常相似,所以我们可以提出这样一个猜想
φ(n)=∑d|nnd∗μ(d)
证明如下:
设 f(n)=φ(n),g(n)=n
根据上面说的那个欧拉函数的性质,得
g(n)=n=∑d|nφ(d)=∑d|nf(d)
反演
f(n)=φ(n)=∑d|ng(nd)∗μ(d)=∑d|nnd∗μ(d)
得证!
殊途同归!
Code
#include<cstdio>
#include<cmath>
#include<algorithm>
#define fo(i,a,b) for(int i=a;i<=b;i++)
using namespace std;
typedef long long ll;
const int N=316230;
int phi[N],pri[N];
bool bz[N];
void pre(int n)
{
phi[1]=1;
fo(i,2,n)
{
if(!bz[i])
{
phi[i]=i-1;
pri[++pri[0]]=i;
}
fo(j,1,pri[0])
{
int t=i*pri[j];
if(t>n) break;
bz[t]=1;
if(i%pri[j]==0)
{
phi[t]=phi[i]*pri[j];
break;
}
phi[t]=phi[i]*phi[pri[j]];
}
}
}
int main()
{
freopen("gcd.in","r",stdin);
freopen("gcd.out","w",stdout);
ll n,ans=0;
scanf("%lld",&n);
int m=int(sqrt(n));
pre(m);
fo(p,1,m)
{
ll i,l=n/p/p,t=0;
for(ll j=1;j<=l;j=i+1)
{
i=l/(l/j);
t+=(i-j+1)*(l/j);
}
ans+=phi[p]*t;
}
printf("%lld",ans);
fclose(stdin);fclose(stdout);
return 0;
}