题目大意
给你 (n)，求
[sum_{a=1}^nsum_{b=1}^nsum_{c=1}^n[frac{1}{a}+frac{1}{b}=frac{1}{c}][gcd(a,b,c)=1]\
]
(nleq {10}^{12})
题解
[egin{align}
&sum_{a=1}^nsum_{b=1}^nsum_{c=1}^n[frac{1}{a}+frac{1}{b}=frac{1}{c}][gcd(a,b,c)=1]\
=&sum_{a=1}^nsum_{b=1}^nsum_{c=1}^n[c(a+b)=ab][gcd(a,b,c)=1]\
=&sum_{a=1}^nsum_{b=1}^n[(a+b)mid ab][gcd(a,b,c)=1]\
end{align}
]
通过打表可以发现，一对数 (a,b(aleq b)) 满足条件的充要条件是 (bleq n) 且 (frac{a}{gcd(a,b)}+frac{b}{gcd(a,b)}=gcd(a,b))
证明：
若 (gcd(a,b)=1)，则 (gcd(a+b,ab)=1)
记 (g=gcd(a,b),a=ga',b=gb')，则
[a+b=g(a'+b')\
ab=g^2a'b'\
frac{ab}{a+b}=frac{g^2a'b'}{g(a'+b')}=frac{ga'b'}{a'+b'}\
]
所以 ((a'+b')mid g)
若 (a'+b' eq g)，则 (gcd(a,b,c)=frac{g}{a'+b'})，所以 (a'+b'=g) 所以 (a=a'(a'+b'),b=b'(a'+b'),c=a'b')
解法一
记 (x=a',y=b')
[egin{align}
ans&=sum_{x=1}^nsum_{y=x}^n[xy+y^2leq n]gcd(x,y)=1\
&=sum_{d=1}^sqrt{n}mu(d)sum_{x=1}sum_{y=x}[xy+y^2leq frac{n}{d^2}]\
&=sum_{d=1}^sqrt{n}mu(d)sum_{y=1}minleft(leftlfloorfrac{leftlfloorfrac{n}{d^2}
ight
floor-y^2}{y}
ight
floor,y
ight)
end{align}
]
时间复杂度：(O(sqrt nlog n))
解法二
(a+b=g^2)
枚举 (g)，那么 (gcd(g,a')=1)
当 (g) 比较小的时候 ((1leq gleq sqrt n))，(a') 有 (varphi(g)) 种取值，可以直接筛
当 (g) 比较大的时候 ((sqrt n<g<sqrt{2n}))，
[sum_{i=1}^frac{n}{g}[gcd(g,i)=1]\
=sum_{dmid g}mu(d)lfloorfrac{n}{gd}
floor
]
可以暴力枚举因子。
时间复杂度：(O(sqrt nlog n))
代码
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cstdlib>
#include<ctime>
#include<utility>
#include<functional>
#include<cmath>
#include<vector>
//using namespace std;
using std::min;
using std::max;
using std::swap;
using std::sort;
using std::reverse;
using std::random_shuffle;
using std::lower_bound;
using std::upper_bound;
using std::unique;
using std::vector;
typedef long long ll;
typedef unsigned long long ull;
typedef double db;
typedef std::pair<int,int> pii;
typedef std::pair<ll,ll> pll;
void open(const char *s){
#ifndef ONLINE_JUDGE
char str[100];sprintf(str,"%s.in",s);freopen(str,"r",stdin);sprintf(str,"%s.out",s);freopen(str,"w",stdout);
#endif
}
void open2(const char *s){
#ifdef DEBUG
char str[100];sprintf(str,"%s.in",s);freopen(str,"r",stdin);sprintf(str,"%s.out",s);freopen(str,"w",stdout);
#endif
}
int rd(){int s=0,c,b=0;while(((c=getchar())<'0'||c>'9')&&c!='-');if(c=='-'){c=getchar();b=1;}do{s=s*10+c-'0';}while((c=getchar())>='0'&&c<='9');return b?-s:s;}
void put(int x){if(!x){putchar('0');return;}static int c[20];int t=0;while(x){c[++t]=x%10;x/=10;}while(t)putchar(c[t--]+'0');}
int upmin(int &a,int b){if(b<a){a=b;return 1;}return 0;}
int upmax(int &a,int b){if(b>a){a=b;return 1;}return 0;}
int gcd(int a,int b)
{
return b?gcd(b,a%b):a;
}
int check(int x,int y)
{
return x*y%(x+y)==0;
}
const int N=1000010;
int b[N];
int pri[N];
int miu[N];
int cnt;
int main()
{
// open("loj6482");
ll n;
scanf("%lld",&n);
miu[1]=1;
for(int i=2;i<=1000000;i++)
{
if(!b[i])
{
pri[++cnt]=i;
miu[i]=-1;
}
for(int j=1;j<=cnt&&i*pri[j]<=1000000;j++)
{
b[i*pri[j]]=1;
if(i%pri[j]==0)
break;
miu[i*pri[j]]=-miu[i];
}
}
ll ans=0;
for(int i=1;(ll)i*i<=n;i++)
if(miu[i])
{
ll s=0;
for(int j=1;;j++)
{
int z=min((n/i/i-(ll)j*j)/j,(ll)j);
if(z<=0)
break;
s+=z;
}
ans+=miu[i]*s;
}
ans=ans*2-1;
printf("%lld
",ans);
// int ans=0;
// for(int i=1;i<=n;i++)
// for(int j=1;j<=n;j++)
// if(gcd(i,j)==1&&i*(i+j)<=n&&j*(i+j)<=n)
// ans++;
// printf("%d
",ans);
// return 0;
// int n=1000;
// for(int i=1;i<=n;i++)
// for(int j=i;j<=n;j++)
// if(gcd(i,j)==1&&check(i*(i+j),j*(i+j))&&gcd(gcd(i*(i+j),j*(i+j)),i*j)==1)
// if(j*(i+j)<=n)
// printf("%d %d %d
",i*(i+j),j*(i+j),i*j);
// return 0;
// for(int i=1;i<=n;i++)
// for(int j=i;j<=n;j++)
// if(i*j%(i+j)==0)
// {
// int z=i*j/(i+j);
// if(gcd(gcd(i,j),z)==1)
// printf("%d %d %d %d %d %d
",i,j,z,gcd(i,j),i/gcd(i,j),j/gcd(i,j));
// }
// return 0;
}