http://www.lydsy.com/JudgeOnline/problem.php?id=2820
思路:
ans
=Σp(p为质数) Σgcd(x,y)==p (1<=x<=n,1<=y<=m)
=Σp(p为质数) Σgcd(x,y)==1 (1<=x<=n/p,1<=y<=m/p)
=Σp(p为质数) Σu(d) (d|gcd(x,y)) (1<=x<=n/p,1<=y<=m/p)
=Σp(p为质数) Σu(d) (d|x&&d|y) (1<=x<=n/p,1<=y<=m/p)
=Σp(p为质数) Σd=1..n/p u(d)*floor(n/pd)*floor(m/pd)
k=pd
=Σp(p为质数) Σd=1..n/p u(k/p)*floor(n/k)*floor(m/k)
=Σk=1..n F(k)*floor(n/k)*floor(m/k)
可以预处理F数组
1 #include<algorithm> 2 #include<cstdio> 3 #include<cmath> 4 #include<cstring> 5 #include<iostream> 6 #define N 10000000 7 #define ll long long 8 int mark[N+10],p[N+10],mul[N+10],n,m; 9 ll f[N+10]; 10 int read(){ 11 char ch=getchar();int t=0,f=1; 12 while (ch<'0'||ch>'9'){if (ch=='-') f=-1;ch=getchar();} 13 while ('0'<=ch&&ch<='9'){t=t*10+ch-'0';ch=getchar();} 14 return t*f; 15 } 16 void init(){ 17 mul[1]=1; 18 for (int i=2;i<=N;i++){ 19 if (!mark[i]){ 20 p[++p[0]]=i; 21 mul[i]=-1; 22 } 23 for (int j=1;j<=p[0]&&i*p[j]<=N;j++){ 24 mark[p[j]*i]=1; 25 if (i%p[j]) mul[p[j]*i]=-mul[i]; 26 else{ 27 mul[p[j]*i]=0; 28 break; 29 } 30 } 31 } 32 for (int i=1;i<=p[0];i++){ 33 int pri=p[i]; 34 for (int j=1;j*pri<=N;j++) 35 f[pri*j]+=mul[j]; 36 } 37 f[0]=0; 38 for (int i=1;i<=N;i++) 39 f[i]+=f[i-1]; 40 } 41 int main(){ 42 init(); 43 int T=read(); 44 while (T--){ 45 ll ans=0; 46 n=read(),m=read(); 47 if (n>m) std::swap(n,m); 48 for (int i=1,j;i<=n;i=j+1){ 49 j=std::min(n/(n/i),m/(m/i)); 50 ans+=(ll)(n/i)*(m/i)*(f[j]-f[i-1]); 51 } 52 printf("%lld ",ans); 53 } 54 }