Description
有一张n×m的数表，其第i行第j列（1<=i<=n，1<=j<=m）的数值为能同时整除i和j的所有自然数之和。给定a，计算数表中不大于a的数之和。
Input
输入包含多组数据。
输入的第一行一个整数Q表示测试点内的数据组数，接下来Q行，每行三个整数n，m，a(|a| < =10^9)描述一组数据。
Output
对每组数据，输出一行一个整数，表示答案模2^31的值。
Sample Input
2
4 4 3
10 10 5
Sample Output
20
148
HINT
1<=n,m<=10^5 , 1<=Q<=2×10^4
Solution
莫比乌斯反演，膜拜PoPoQQQ
先不管a的限制，我们直接求数表中所有数之和(ans')
设(F(i))为i的约数和，(g(i))为在限制n和m范围内gcd为i的对数的个数，那么
[ans'=sum_{i=1}^{min(n,m)}F(i)g(i)
]
(g(i))的反演已经是老套路了，直接得出(g(i)=sum_{i|d}mu(frac{d}{n})lfloor frac{n}{d}
floor lfloor frac{m}{d}
floor)
继续推
[ans'=sum_{i=1}^{min(n,m)}F(i)sum_{i|d}mu(frac{d}{i})lfloor frac{n}{d}
floor lfloor frac{m}{d}
floor=sum_{i=1}^{min(n,m)}F(i)sum_{k=1}^{lfloor frac{min(n,m)}{i}
floor}mu(k)lfloor frac{n}{ik}
floor lfloor frac{m}{ik}
floor
]
设(T=ik)
[ans'=sum_{T=1}^{min(n,m)}lfloor frac{n}{T}
floor lfloor frac{m}{T}
floorsum_{k=1}^{lfloor frac{min(n,m)}{T}
floor}F(k)mu(lfloor frac{T}{k}
floor)
]
前面部分整除分块，我们只要处理出后半部分的前缀和就行了
对于(F(k))直接(O(nlog_2n))暴力枚举约数和它的倍数
对于(mu)这就不用说了吧。。。（看之前的文章）
好的，我们把改装之后的问题解决了，可是原问题呢？
回到真正的(ans)，因为有a的限制，所以那些(F(i))大于a的是不能加贡献的
那么我们把询问按a排序，然后用树状数组维护(F(k)mu(lfloor frac{T}{k}
floor))的前缀和，只有当前询问的a大于扫描到的(F(i))，才把(F(i))能产生的贡献加入树状数组。这样就保证了不改加入的贡献不会被加入
#include<bits/stdc++.h>
#define ll long long
const int MAXT=20000+10,MAXN=100000+10,Mod=0x7fffffff;
int T,cnt,vis[MAXN],prime[MAXN],mu[MAXN],C[MAXN],ans[MAXT],limit;
struct question{
int n,m,a;
int id;
inline bool operator < (const question &A) const{
return a<A.a;
};
};
question Q[MAXT];
struct node{
int s;
int id;
inline bool operator < (const node &A) const{
return s<A.s;
};
};
node F[MAXN];
template<typename T> inline void read(T &x)
{
T data=0,w=1;
char ch=0;
while(ch!='-'&&(ch<'0'||ch>'9'))ch=getchar();
if(ch=='-')w=-1,ch=getchar();
while(ch>='0'&&ch<='9')data=((T)data<<3)+((T)data<<1)+(ch^'0'),ch=getchar();
x=data*w;
}
template<typename T> inline void write(T x,char c=' ')
{
if(x<0)putchar('-'),x=-x;
if(x>9)write(x/10);
putchar(x%10+'0');
if(c!=' ')putchar(c);
}
template<typename T> inline void chkmin(T &x,T y){x=(y<x?y:x);}
template<typename T> inline void chkmax(T &x,T y){x=(y>x?y:x);}
template<typename T> inline T min(T x,T y){return x<y?x:y;}
template<typename T> inline T max(T x,T y){return x>y?x:y;}
inline void init()
{
for(register int i=1;i<=limit;++i)
for(register int j=i;j<=limit;j+=i)F[j].s+=i;
for(register int i=1;i<=limit;++i)F[i].id=i;
std::sort(F+1,F+limit+1);
memset(vis,1,sizeof(vis));
vis[0]=vis[1]=0;
mu[1]=1;
for(register int i=2;i<=limit;++i)
{
if(vis[i])
{
prime[++cnt]=i;
mu[i]=-1;
}
for(register int j=1;j<=cnt&&i*prime[j]<=limit;++j)
{
vis[i*prime[j]]=0;
if(i%prime[j])mu[i*prime[j]]=-mu[i];
else break;
}
}
}
inline int lowbit(int x)
{
return x&(-x);
}
inline int sum(int x)
{
int res=0;
while(x)
{
res+=C[x];
x-=lowbit(x);
}
return res;
}
inline void add(int x,int k)
{
while(x<=limit)
{
C[x]+=k;
x+=lowbit(x);
}
}
inline int solve(int n,int m)
{
int res=0;
for(register int i=1;;)
{
if(i>min(n,m))break;
int j=min(n/(n/i),m/(m/i));
res+=(n/i)*(m/i)*(sum(j)-sum(i-1));
i=j+1;
}
return res&Mod;
}
int main()
{
read(T);
for(register int i=1;i<=T;++i)
{
read(Q[i].n);read(Q[i].m);read(Q[i].a);
Q[i].id=i;
chkmax(limit,min(Q[i].n,Q[i].m));
}
init();
std::sort(Q+1,Q+T+1);
for(register int i=1,j=1;i<=T;++i)
{
while(j<=limit&&F[j].s<=Q[i].a)
{
for(register int p=F[j].id,k=1;p<=limit;p+=F[j].id,++k)add(p,F[j].s*mu[k]);
++j;
}
ans[Q[i].id]=solve(Q[i].n,Q[i].m);
}
for(register int i=1;i<=T;++i)write(ans[i],'
');
return 0;
}