为了好写式子，先不管$a$的限制
设$facs$为因子和，那么有
$ans=sumlimits_{i=1}^nsumlimits_{j=1}^mfacs(gcd(i,j))$
再设$f(k)=sumlimits_{i=1}^nsumlimits_{j=1}^m[gcd(i,j)==k]$
熟悉的东西，再写一遍=。=
$f(k)=sumlimits_{i=1}^nsumlimits_{j=1}^m[gcd(i,j)==k]$
$=sumlimits_{i=1}^{min(leftlfloorfrac{n}{k} ight floor,leftlfloorfrac{m}{k} ight floor)}μ(i)leftlfloorfrac{n}{ik} ight floorleftlfloorfrac{m}{ik} ight floor$
$sumlimits_{i=1}^{min(n,m)}[k|i]μ(frac{i}{k})leftlfloorfrac{n}{i} ight floorleftlfloorfrac{m}{i} ight floor$
那么
$ans=sumlimits_{i=1}^{min(n,m)}facs(i)f(i)$
$=sumlimits_{i=1}^{min(n,m)}facs(i)sum_{d|i}μ(frac{i}{d})leftlfloorfrac{n}{i} ight floorleftlfloorfrac{m}{i} ight floor$
$=sumlimits_{i=1}^{min(n,m)}leftlfloorfrac{n}{i} ight floorleftlfloorfrac{m}{i} ight floorsum_{d|i}facs(d)μ(frac{i}{d})$
预处理后面的那个东西，前面的每次$O(sqrt n)$回答
等等还有$a$的限制
先读进来按$a$排序，然后依次插入后面那个函数值，也就是单点修改+区间查询，树状数组解决
1 #include<cstdio> 2 #include<cstring> 3 #include<algorithm> 4 using namespace std; 5 const int N=100005,M=1e5,Mod=2147483647; 6 struct a 7 { 8 int n,m,a; 9 int idx,ans; 10 }qry[N]; 11 struct b 12 { 13 int f,idx; 14 }sum[N]; 15 int npr[N],pri[N],mul[N],bit[N],T,p,cnt; 16 bool cmp(a x,a y) 17 { 18 return x.a<y.a; 19 } 20 bool com(b x,b y) 21 { 22 return x.f<y.f; 23 } 24 bool cpr(a x,a y) 25 { 26 return x.idx<y.idx; 27 } 28 void Add(int pos,int tsk) 29 { 30 while(pos<=M) 31 bit[pos]+=tsk,pos+=pos&-pos; 32 } 33 int Query(int pos) 34 { 35 int ret=0; 36 while(pos) 37 ret+=bit[pos],pos-=pos&-pos; 38 return ret; 39 } 40 void Prework() 41 { 42 npr[1]=true,mul[1]=1,p=1; 43 for(int i=2;i<=M;i++) 44 { 45 if(!npr[i]) pri[++cnt]=i,mul[i]=-1; 46 for(int j=1;j<=cnt&&i*pri[j]<=M;j++) 47 { 48 npr[i*pri[j]]=true; 49 if(i%pri[j]==0) break; 50 else mul[i*pri[j]]=-mul[i]; 51 } 52 } 53 for(int i=1;i<=M;i++) sum[i].idx=i; 54 for(int i=1;i<=M;i++) 55 for(int j=i;j<=M;j+=i) sum[j].f+=i; 56 } 57 int Solve(int n,int m) 58 { 59 int ret=0; 60 if(n>m) swap(n,m); 61 for(int i=1,j;i<=n;i=j+1) 62 { 63 j=min(n/(n/i),m/(m/i)); 64 ret+=(n/i)*(m/i)*(Query(j)-Query(i-1)); 65 } 66 return ret&Mod; 67 } 68 int main() 69 { 70 Prework(); 71 scanf("%d",&T); 72 for(int i=1;i<=T;i++) 73 scanf("%d%d%d",&qry[i].n,&qry[i].m,&qry[i].a),qry[i].idx=i; 74 sort(qry+1,qry+1+T,cmp); 75 sort(sum+1,sum+1+M,com); 76 for(int i=1;i<=T;i++) 77 { 78 while(p<=M&&sum[p].f<=qry[i].a) 79 { 80 for(int j=sum[p].idx;j<=M;j+=sum[p].idx) 81 Add(j,mul[j/sum[p].idx]*sum[p].f); p++; 82 } 83 qry[i].ans=Solve(qry[i].n,qry[i].m); 84 } 85 sort(qry+1,qry+1+T,cpr); 86 for(int i=1;i<=T;i++) printf("%d ",qry[i].ans); 87 return 0; 88 }