题意：对于给出的n个询问，每次求有多少个数对(x,y)，满足a≤x≤b，c≤y≤d，
且gcd(x,y) = k，gcd(x,y)函数为x和y的最大公约数。
1≤n≤50000，1≤a≤b≤50000，1≤c≤d≤50000，1≤k≤50000
思路：第一题反演……
利用容斥原理将原询问拆成4个，问题就转化为：
1<=i<=trunc(a div k)，1<=j<=trunc(b div k)，gcd(i,j)=1的(i,j)数对个数
令f(i)表示满足gcd(x,y)=i时(x,y)的对数，F(i)表示满足i|gcd(x,y)的(x,y)的对数
显然F(i)=trunc(n div i)*trunc(m div i)
f(1)=sigma u(d)*n div d*m div d (1<=d<=n)
观察可得n div d*m div d只有2根号n个取值，且每个取值对应的u(i)是连续的一段
然后就可以记录u的前缀和来优化
From http://m.blog.csdn.net/article/details?id=50590197
1 //uses sysutils; 2 const max=50000; 3 var mu,sum,prime:array[0..max]of longint; 4 flag:array[0..max]of longint; 5 a,b,c,d,k,i,j,t,m,cas,v:longint; 6 tmp:double; 7 8 procedure swap(var x,y:longint); 9 var t:int64; 10 begin 11 t:=x; x:=y; y:=t; 12 end; 13 14 function clac(n,m:longint):int64; 15 var i,t1,t2,pos:longint; 16 x,y:int64; 17 begin 18 if n>m then swap(n,m); 19 clac:=0; i:=1; 20 while i<=n do 21 begin 22 x:=n div i; y:=m div i; 23 t1:=n div x; 24 t2:=m div y; 25 if t1<t2 then pos:=t1 26 else pos:=t2; 27 clac:=clac+x*y*(sum[pos]-sum[i-1]); 28 i:=pos+1; 29 end; 30 end; 31 32 begin 33 assign(input,'bzoj2301.in'); reset(input); 34 assign(output,'bzoj2301.out'); rewrite(output); 35 // tmp:=now; 36 read(cas); 37 mu[1]:=1; 38 for i:=2 to max do 39 begin 40 if flag[i]=0 then 41 begin 42 inc(m); prime[m]:=i; 43 mu[i]:=-1; 44 end; 45 j:=1; 46 while (j<=m)and(prime[j]*i<=max) do 47 begin 48 t:=prime[j]*i; 49 flag[t]:=1; 50 if i mod prime[j]=0 then 51 begin 52 mu[t]:=0; 53 break; 54 end; 55 mu[t]:=-mu[i]; 56 inc(j); 57 end; 58 end; 59 for i:=1 to max do sum[i]:=sum[i-1]+mu[i]; 60 for v:=1 to cas do 61 begin 62 read(a,b,c,d,k); 63 dec(a); dec(c); 64 a:=a div k; b:=b div k; c:=c div k; d:=d div k; 65 writeln(clac(b,d)-clac(b,c)-clac(a,d)+clac(a,c)); 66 end; 67 //writeln((now-tmp)*86400:0:2); 68 close(input); 69 close(output); 70 end.
UPD（2018.10.22）：C++
1 #include<cstdio> 2 #include<cstring> 3 #include<string> 4 #include<cmath> 5 #include<iostream> 6 #include<algorithm> 7 #include<map> 8 #include<set> 9 #include<queue> 10 #include<vector> 11 using namespace std; 12 typedef long long ll; 13 typedef unsigned int uint; 14 typedef unsigned long long ull; 15 typedef pair<int,int> PII; 16 typedef vector<int> VI; 17 #define fi first 18 #define se second 19 #define MP make_pair 20 #define N 51000 21 #define M 410000 22 #define eps 1e-8 23 #define pi acos(-1) 24 #define oo 1e9 25 26 int mu[N+10],s[N+10],prime[N+10],flag[N+10]; 27 28 ll calc(int n,int m) 29 { 30 if(n>m) swap(n,m); 31 ll ans=0; 32 int i=1; 33 while(i<=n) 34 { 35 ll x=n/i; 36 ll y=m/i; 37 int t1=n/x; 38 int t2=m/y; 39 int pos=min(t1,t2); 40 ans+=x*y*(s[pos]-s[i-1]); 41 i=pos+1; 42 } 43 return ans; 44 } 45 46 int main() 47 { 48 int cas; 49 scanf("%d",&cas); 50 mu[1]=1; 51 int m=0; 52 for(int i=2;i<=N;i++) 53 { 54 if(!flag[i]) 55 { 56 prime[++m]=i; 57 mu[i]=-1; 58 } 59 for(int j=1;j<=m;j++) 60 { 61 int t=prime[j]*i; 62 if(t>N) break; 63 flag[t]=1; 64 if(i%prime[j]==0) 65 { 66 mu[t]=0; 67 break; 68 } 69 mu[t]=-mu[i]; 70 } 71 } 72 for(int i=1;i<=N;i++) s[i]=s[i-1]+mu[i]; 73 while(cas--) 74 { 75 int a,b,c,d,k; 76 scanf("%d%d%d%d%d",&a,&b,&c,&d,&k); 77 a--; c--; 78 a/=k; b/=k; c/=k; d/=k; 79 ll ans=calc(b,d)-calc(b,c)-calc(a,d)+calc(a,c); 80 printf("%lld ",ans); 81 } 82 return 0; 83 }