GDSOI 2016 T1 互补约数

Description

求∑i=1n∑d|igcd(d,id)

n<=10^11

Solution

首先，我们发现gcd中的两个东西是所有乘积不超过n的数对，即

Ans=∑i∑j,i∗j<=ngcd(i,j)

然后Ans=∑i=1n∑j=1⌊ni⌋gcd(i,j)

咦？

长得那么像莫比乌斯的经典形式。只不过j的上界是会变化的。

那么我们可以尝试一下。

设F(d)表示上式中gcd(i,j)=d的个数，G(d)表示d|gcd(i,j)的个数。

根据反演，F(d)=∑i=1⌊nd⌋G(id)μ(i)

那么Ans=∑d=1nF(d)∗d

展开Ans=∑d=1nd∑i=1⌊nd⌋G(id)μ(i)

考虑枚举T=i∗d

Ans=∑T=1nG(T)∑i|TTiμ(i)

似乎复杂度一点都没有优化？

我们考虑一下如何求G

如果两个数的gcd=T，那么这两个数可以表示成XT,YT.

那么我们可以枚举X,G(T)=∑x=1⌊nT2⌋⌊nT2x⌋

这个东西的复杂度大概是O(N12)

由于n的下面有个平方，于是这个东西的复杂度肯定没问题。

于是发现原来的式子中，T只需要枚举到n√就好了。

那么我们可以考虑预处理出∑i|TTiμ(i)设为a(T)

那么原式就变成了优美的Ans=∑T=1n√G(T)a(T)

时间复杂度大约为O(N34)

实际远远不到，O(跑得过)O(∩_∩)O~

Code

#include<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define fo(i,a,b) for(int i=a;i<=b;i++)
#define N 1000005
using namespace std;
typedef long long ll;
ll a
,p
,n,ans;
int mu
,Mx;
bool bz
;
int main() {
freopen("gcd.in","r",stdin);
freopen("gcd.out","w",stdout);
scanf("%lld",&n);mu[1]=1;Mx=sqrt(n);
fo(i,2,Mx) {
if (!bz[i]) p[++p[0]]=i,mu[i]=-1;
fo(j,1,p[0]) {
ll k=i*p[j];if (k>Mx) break;
bz[k]=1;if (!(i%p[j])) break;
mu[k]=-mu[i];
}
}
fo(i,1,Mx) fo(j,1,Mx/i) a[i*j]+=mu[i]*j;
fo(i,1,Mx) {
ll m=n/(ll)i/(ll)i,sum=0;
for(ll l=1,r=0;l<=m;l=r+1)
r=m/(m/l),sum=sum+(r-l+1)*(m/l);
ans=ans+a[i]*sum;
}
printf("%lld",ans);
}