Description
有一张$n imes m$的数表，其第$i$行第$j$列 $(1≤i≤n,1≤j≤m)$ 的数值为能同时整除$i$和$j$的所有自然数之和。
现在给定$a$，计算数表中不大于$a$的数之和。
Input
输入包含多组数据。
输入的第一行一个整数$Q$表示测试点内的数据组数，接下来Q行，每行三个整数$n,m,a(a≤109)$描述一组数据。
Output
题解：
我数学太水了！！又是一道推公式的题：
egin{aligned}
ans&=sum^n_{i=1}sum^m_{j=1}sum_{dmid gdc(i,j)}d\
&=sum_d^nsum^{lfloorfrac{n}{d}
floor}_{i=1}sum^{lfloorfrac{m}{d}
floor}_{j=1}f(d)[gcd(i,j)=1] &(f(d)=sum_{kmid d}k)\
&=sum_d^n f(d)sum^{lfloorfrac{n}{d}
floor}_{i=1}sum^{lfloorfrac{m}{d}
floor}_{j=1}sum_{kmid gcd(i,j)}mu(k)\
&=sum_k^nsum_d^n f(d)mu(k)lfloorfrac{n}{kd}
floorlfloorfrac{m}{kd}
floor\
&设 T=kd \
ans&=sum_T^nsum_{dmid T}f(d)mu(frac{T}{d})lfloorfrac{n}{T}
floorlfloorfrac{m}{T}
floor\
end{aligned}
前面同样是数论分块，后面依旧是线性筛。
考虑f怎么求
有个约数和定理
若$n=sum_{i=1}^k p_i^{a_i}$，$p_i$为$n$的质因数，那么$n$的约数和 $f(n)$满足$f(n)=prod_{i=1}^ksum_{j=0}^{a_i}p_i^j$
$f$的话首先是个积性函数，我们在筛$mu$的时候想顺便把这个也筛出来
考虑$f(d)$的值，如果说d是质数的话答案显然是$d+1$，下面讨论$d$为合数的情况
设$d=i*p$，其中$p$为质数
$p mid i$，那么$p$和$i$互质，所以$f(d)=f(p)*f(i)$
$p∣i$，设$i=t*p^x$ ，那么根据约数和定理，我们可以得
$$f(i*p) = f(t)f(p^{x+1}) = f(t)sumlimits_{i=0}^{x+1}p^i$$
然后我们把p0（也就是1）拿出来，得到
$$f(i * p) = f(t) + f(t)*sumlimits_{i=1}^{x+1}p^i = f(t) + f(t)*f(p^x)*p$$
然后$i = t * p^x$，所以$f(t) * f(p^x) = f(i)$
所以最后就是$f(d) = f(t) + f(i) * p$
然后就可以筛出$f(d)$啦
剩下的东西
现在加上a的限制，其实就是离线处理
我们先将所有的询问按照a的大小排序，然后从小到大处理
因为分块的时候我们要用到的是g(x)的前缀和，所以用一个树状数组来处理
将f(x)排个序，枚举的时候只枚举到f(x)<a，然后枚举另一个约数求出g，丢到树状数组里面去
求答案的时候直接查询就好了
CODE：
1 #include<iostream> 2 #include<cstdio> 3 #include<algorithm> 4 using namespace std; 5 6 #define lowbit(x) (x&-x) 7 int T,cnt=0,maxn,ans[20005],tmp[100005]; 8 int f[100005],pri[100005],c[100005],mu[100005]; 9 bool vis[100005]; 10 struct Que{ 11 int n,m,a,id; 12 bool operator<(const Que &b)const{ 13 return a<b.a; 14 } 15 }q[20005]; 16 17 bool comp(int a,int b){return f[a]<f[b];} 18 19 void init(){ 20 mu[1]=f[1]=1; 21 for(int i=2;i<=maxn;i++){ 22 if(!vis[i]){ 23 f[i]=i+1,mu[i]=-1; 24 pri[++cnt]=i; 25 } 26 for(int j=1;j<=cnt&&i*pri[j]<=maxn;j++){ 27 vis[i*pri[j]]=true; 28 if(i%pri[j]==0){ 29 mu[i*pri[j]]=0; 30 int t=i; 31 while(t%pri[j]==0)t/=pri[j]; 32 f[i*pri[j]]=f[t]+pri[j]*f[i]; 33 break; 34 }else{ 35 f[i*pri[j]]=f[i]*f[pri[j]]; 36 mu[i*pri[j]]=-mu[i]; 37 } 38 } 39 } 40 } 41 42 void add(int x,int y){ 43 while(x<=maxn){ 44 c[x]+=y; 45 x+=lowbit(x); 46 } 47 } 48 49 int sum(int x){ 50 int ans=0; 51 while(x>=1){ 52 ans+=c[x]; 53 x-=lowbit(x); 54 } 55 return ans; 56 } 57 58 int main(){ 59 scanf("%d",&T); 60 for(int i=1;i<=T;i++){ 61 scanf("%d%d%d",&q[i].n,&q[i].m,&q[i].a); 62 q[i].id=i; 63 if(q[i].n>q[i].m)swap(q[i].n,q[i].m); 64 maxn=max(maxn,q[i].n); 65 } 66 init(); 67 sort(q+1,q+T+1); 68 for(int i=1;i<=maxn;i++)tmp[i]=i; 69 sort(tmp+1,tmp+maxn+1,comp); 70 for(int i=1,now=1;i<=T;i++){ 71 int n=q[i].n,m=q[i].m,a=q[i].a; 72 while(now<=maxn&&f[tmp[now]]<=a){ 73 for(int k=1;k*tmp[now]<=maxn;k++) 74 add(k*tmp[now],f[tmp[now]]*mu[k]); 75 now++; 76 } 77 for(int j=1,pos;j<=n;j=pos+1){ 78 pos=min(n/(n/j),m/(m/j)); 79 ans[q[i].id]+=(n/j)*(m/j)*(sum(pos)-sum(j-1)); 80 } 81 } 82 for(int i=1;i<=T;i++) 83 printf("%d ",ans[i]&0x7fffffff); 84 }