Description
Solution1
一题标准的莫比乌斯反演,
变式:
ans=∑i=1n∑j|igcd(i,ji)=∑i=1n∑j=1⌊ni⌋gcd(i,j)
设gcd(i,j)=d的数有 f(d) 个, g(d) 表示 d|gcd(i,j) 的个数,
g(d)=∑i=1⌊nd⌋f(id)
,
根据反演:
f(d)=∑i=1⌊nd⌋g(id)μ(i)
那么
ans=∑d=1n√f(d)∗d
(因为gcd(i,j)不可能大于 n−−√ )
展开:
ans=∑d=1n√d∗∑i=1⌊nd⌋g(id)μ(i)
用 T=id
ans=∑T=1ng(T)∑d|Td∗μ(Td)
因为 g(d) 表示 d|gcd(i,j) 的个数,所有可以表示成: gcd(d∗x,d∗y) ,
枚举其中的x,
g(T)=∑x=1⌊nT2⌋⌊nT2∗x⌋
所以当 T>n−−√ 时, g(T)=0 无意义,第一个枚举的界限可以降低到 n−−√
ans=∑T=1n√g(T)∑d|Td∗μ(Td)
g(T) 用分块,后面那一堆暴力算,
复杂度:约 O(n34) (我不会算,CZY用微积分证明了(%%%))
Solution2
ans=∑T=1n√φ(T)∑x=1⌊nT2⌋⌊nT2∗x⌋
因为 φ(n)=∑d|nnd∗μ(d) ;( 详情点这里)
所以我们发现跟上面的式子一模一样。
Code
#include<cstdio>
#include<cmath>
#define fo(i,a,b) for(int i=a;i<=b;i++)
using namespace std;
typedef long long LL;
const int N=400000;
int pr[N/3],mu[N];
int b[N/3];
bool prz[N];
LL g[N],n,m,ans;
void dofirst(int n)
{
mu[1]=1;
fo(i,2,m)
{
if(!prz[i])pr[++pr[0]]=i,mu[i]=-1;
for(int j=1;j<=pr[0] && i*pr[j]<=m;j++)
{
prz[i*pr[j]]=1;
if(!(i%pr[j]))break;
mu[i*pr[j]]=-mu[i];
}
}
}
void zh(int q,int e,int T)
{
if(q>pr[0]){g[T]+=(LL)e*mu[T/e];return;}
fo(i,0,b[q])zh(q+1,e,T),e*=pr[q];
}
LL fk(LL n)
{
LL ans=0,t=1,r;
while(t<=n)
{
r=n/(n/t);
ans+=n/t*(r-t+1);
t=r+1;
}
return ans;
}
int main()
{
freopen("gcd.in","r",stdin);
freopen("gcd.out","w",stdout);
scanf("%lld",&n);m=sqrt(n);
dofirst(m);
fo(I,1,m)
{
pr[0]=0;int i=I;
fo(j,2,sqrt(i))if(!(i%j))
{
pr[++pr[0]]=j;b[pr[0]]=0;
while(!(i%j))i/=j,b[pr[0]]++;
}
if(i-1)pr[++pr[0]]=i,b[pr[0]]=1;
zh(1,1,I);
}
fo(T,1,m)ans+=g[T]*fk(n/T/T);
printf("%lld\n",ans);
return 0;
}