由于之前做了Luogu P2257 YY的GCD，这里的做法就十分套路了。
建议先看上面一题的推导，这里的话就略去一些共性的地方了。
还是和之前一样设：
[f(d)=sum_{i=1}^a sum_{j=1}^b[gcd(i,j)=d]
]
[F(n)=sum_{n|d} f(d)=lfloorfrac{a}{n}
floorlfloorfrac{b}{n}
floor
]
还是莫比乌斯反演定理推出：
[f(n)=sum_{n|d}mu(lfloorfrac{d}{n}
floor)F(d)
]
这时我们发现，不像上面一题那么繁琐还要对(f(n))求和，这里(ans=f(d))
所以可以直接开始推导答案：
[ans=sum_{d|k}mu(lfloorfrac{k}{d}
floor)F(k)
]
还是考虑换个东西枚举，我们设(lfloorfrac{k}{d} floor=t)，枚举(t)则有：
[ans=sum_{t=1}^{min(lfloorfrac{a}{d}
floor,lfloorfrac{b}{d}
floor)}mu(t)lfloorfrac{a}{tcdot d}
floorlfloorfrac{b}{tcdot d}
floor
]
这个式子已经变成(O(n))的了，还是注意到(lfloorfrac{a}{tcdot d} floorlfloorfrac{b}{tcdot d} floor)可以除法分块，然后只需要对莫比乌斯函数做一个前缀和即可单次(O(sqrt n))。
CODE
#include<cstdio>
#include<cctype>
#define RI register int
using namespace std;
const int P=50005;
int t,n,m,d,prime[P+5],cnt,mu[P+5],sum[P+5]; long long ans; bool vis[P+5];
class FileInputOutput
{
private:
#define tc() (A==B&&(B=(A=Fin)+fread(Fin,1,S,stdin),A==B)?EOF:*A++)
#define pc(ch) (Ftop<S?Fout[Ftop++]=ch:(fwrite(Fout,1,S,stdout),Fout[(Ftop=0)++]=ch))
#define S 1<<21
char Fin[S],Fout[S],*A,*B; int Ftop,pt[25];
public:
FileInputOutput() { A=B=Fin; Ftop=0; }
inline void read(int &x)
{
x=0; char ch; int flag=1; while (!isdigit(ch=tc())) flag=ch^'-'?1:-1;
while (x=(x<<3)+(x<<1)+(ch&15),isdigit(ch=tc())); x*=flag;
}
inline void write(long long x)
{
if (!x) return (void)(pc(48),pc('
')); RI ptop=0;
while (x) pt[++ptop]=x%10,x/=10; while (ptop) pc(pt[ptop--]+48); pc('
');
}
inline void Fend(void)
{
fwrite(Fout,1,Ftop,stdout);
}
#undef tc
#undef pc
#undef S
}F;
#define Pi prime[j]
inline void Euler(void)
{
vis[1]=mu[1]=1; RI i,j; for (i=2;i<=P;++i)
{
if (!vis[i]) prime[++cnt]=i,mu[i]=-1;
for (j=1;j<=cnt&&i*Pi<=P;++j)
{
vis[i*Pi]=1; if (i%Pi) mu[i*Pi]=-mu[i]; else break;
}
}
for (i=1;i<=P;++i) sum[i]=sum[i-1]+mu[i];
}
#undef Pi
inline int min(int a,int b)
{
return a<b?a:b;
}
int main()
{
//freopen("CODE.in","r",stdin); freopen("CODE.out","w",stdout);
for (Euler(),F.read(t);t;--t)
{
F.read(n); F.read(m); F.read(d); ans=0; int lim=min(n/d,m/d);
for (RI l=1,r;l<=lim;l=r+1)
{
r=min(n/(n/l),m/(m/l)); ans+=1LL*(n/(l*d))*(m/(l*d))*(sum[r]-sum[l-1]);
}
F.write(ans);
}
return F.Fend(),0;
}