传送门
可以先不管(a)的限制
显然问题可以转化，对于一个数(i*j)，能同时整除(i)并且整除(j)的数(x)，也就是满足(x|i&&y|j)，则(x|gcd(i,j))，也就是(x)为(gcd(i,j))的约数，设(g(n))为(n)的约数和，那么
[ans=sum_{d=1}^{min(n,m)}sum_{i=1}^{n}sum_{j=1}^{m}[gcd(i,j)==d]g(d)\
ans=sum_{d=1}^{min(n,m)}g(d)sum_{i=1}^{n}sum_{j=1}^{m}[gcd(i,j)==d]\
]
然后设
[f(d)=sum_{i=1}^{n}sum_{j=1}^{m}[gcd(i,j)==d]\
F(d)=sum_{d|a}f(a)=lfloor frac{n}{d}
floorlfloor frac{m}{d}
floorsum_{i=1}^{n}sum_{j=1}^{m}[d|gcd(i,j)]\
f(x)=sum_{x|d}F(d)mu(frac{d}{x})=sum_{x|d}lfloor frac{n}{d}
floorlfloor frac{m}{d}
floormu(frac{d}{x})\
]
然后(ans)可以转化了
[ans=sum_{d=1}^{min(n,m)}g(d)f(d)\
ans=sum_{d=1}^{min(n,m)}g(d)sum_{d|y}lfloor frac{n}{y}
floorlfloor frac{m}{y}
floormu(frac{y}{d})\
ans=sum_{d=1}^{min(n,m)}lfloor frac{n}{d}
floorlfloor frac{m}{d}
floor sum_{d|y}g(d)mu(frac{y}{d})\
]
然后约数和可以线筛的时候处理出来
这个时候应该拾起那个没有考虑的(a)了
考虑离线，当(g(d)<=a)的时候会产生贡献，排序一下
然而数论分块是需要前缀和的，有什么快速维护前缀和并且支持单点修改的数据结构呢
树状数组是个不错的选择
代码：
#include<cstdio>
#include<iostream>
#include<algorithm>
using namespace std;
void read(int &x) {
char ch; bool ok;
for(ok=0,ch=getchar(); !isdigit(ch); ch=getchar()) if(ch=='-') ok=1;
for(x=0; isdigit(ch); x=x*10+ch-'0',ch=getchar()); if(ok) x=-x;
}
#define rg register
#define lowbit(i) (i&(-i))
const int maxn=1e5+10,mod=(1<<31)-1;
int T,ans[maxn],mu[maxn],pri[maxn],tot,mx,s[maxn],f[maxn];bool vis[maxn];
struct o{int n,m,a,id;}a[maxn];
struct oo{int x,y;}g[maxn];
bool cmp(oo a,oo b){return a.x<b.x;}
bool Cmp(o a,o b){return a.a<b.a;}
void prepare(int x)
{
mu[1]=1,s[1]=1;
for(rg int i=2;i<=x;i++)
{
if(!vis[i])pri[++tot]=i,mu[i]=-1,s[i]=i+1;
for(rg int j=1;j<=tot&&pri[j]*i<=x;j++)
{
vis[pri[j]*i]=1;
if(!(i%pri[j])){s[i*pri[j]]=s[i]*s[pri[j]]-s[i/pri[j]]*pri[j];break;}
else s[i*pri[j]]=s[i]*s[pri[j]],mu[i*pri[j]]=-mu[i];
}
}
for(rg int i=1;i<=x;i++)g[i].x=s[i],g[i].y=i;
sort(g+1,g+x+1,cmp);
}
void add(int x,int y){for(rg int i=x;i<=mx;i+=lowbit(i))f[i]+=y;}
int get(int x){int ans=0;for(rg int i=x;i;i-=lowbit(i))ans+=f[i];return ans;}
int solve(int n,int m)
{
if(n>m)swap(n,m);int ans=0;
for(rg int i=1,j;i<=n;i=j+1)
{
j=min(n/(n/i),m/(m/i));
int t=(n/i)*(m/i)*(get(j)-get(i-1));
ans+=t;
}
return ans;
}
int main()
{
read(T);
for(rg int i=1;i<=T;i++)read(a[i].n),read(a[i].m),read(a[i].a),a[i].id=i,mx=max(mx,max(a[i].n,a[i].m));
prepare(mx);sort(a+1,a+T+1,Cmp);int now=1;
for(rg int i=1;i<=T;i++)
{
while(now<=mx&&g[now].x<=a[i].a)
{
for(rg int j=g[now].y;j<=mx;j+=g[now].y)
add(j,g[now].x*mu[j/g[now].y]);
now++;
}
ans[a[i].id]=solve(a[i].n,a[i].m);
}
for(rg int i=1;i<=T;i++)printf("%d
",ans[i]&((1<<31)-1));
}