$ans=sum_{i=1}^nsum_{j=1}^nsigma(gcd(i,j))$
枚举gcd为d的所有数得到
$ans=sum_{d<=n}sigma(d)*g(d)$
$g(d)$表示所有(i,j)=d的二元组的数量。
那么可以反演得到$g(i)=sum_{i mid d}mu(lfloor d/i floor )*lfloor n/d floor * lfloor m/d floor$
然后代入然后xjb变换可得
$ans=sum_{d<=n}lfloor n/d floor * lfloor m/d floor sum_{i mid d}mu( lfloor d/i floor ) * sigma(i) $
然后我们要求出$sum_{i mid d}mu(lfloor d/i floor ) *sigma(i) $的前缀和就可以$sqrt n$的时间内解决了
那么我们可以用每个数去暴力更新倍数即可,但是它是一个积性函数,是可以在$Theta(n)$的时间内筛出来的。
但是有A的条件,我们可以去维护前缀和用树状数组,暴力更新倍数即可。
#include <cstdio> #include <cstring> #include <iostream> #include <algorithm> using namespace std; #define F(i,j,k) for (int i=j;i<=k;++i) #define D(i,j,k) for (int i=j;i>=k;--i) #define md 2147483647 #define inf 0x3f3f3f3f #define maxn 100005 struct query{int n,m,k,id,ans;}a[maxn]; struct Bit_Tree{ int x[maxn]; void add(int i,int f) {for (;i<maxn;i+=i&(-i))x[i]+=f;} int gs(int i) { int ret=0; for (;i;i-=i&(-i)) ret+=x[i]; return ret; } }BT; int sigma[maxn],pr[maxn],top,mu[maxn],min_fac_a[maxn],min_fac_sum[maxn],rk[maxn]; void init() { sigma[1]=1;mu[1]=1;rk[1]=1; F(i,2,maxn-1) { rk[i]=i; if (!sigma[i]) { pr[++top]=i; min_fac_a[i]=i; sigma[i]=min_fac_sum[i]=i+1; mu[i]=-1; } F(j,1,top) { if (pr[j]*i>=maxn) break; if (i%pr[j]==0) { sigma[pr[j]*i]=sigma[i]/min_fac_sum[i]* (min_fac_sum[pr[j]*i]=min_fac_sum[i]+min_fac_a[i]*pr[j]); min_fac_a[pr[j]*i]=min_fac_a[i]*pr[j]; mu[pr[j]*i]=0; break; } sigma[pr[j]*i]=sigma[pr[j]]*sigma[i]; min_fac_a[pr[j]*i]=pr[j]; min_fac_sum[pr[j]*i]=pr[j]+1; mu[pr[j]*i]=-mu[i]; } } } int t; bool cmp(query x,query y) {return x.k<y.k;} bool cmp2(query x,query y) {return x.id<y.id;} bool cmp3(int x,int y) {return sigma[x]<sigma[y];} void add(int i) { F(j,1,inf) { if (i*j>=maxn) break; BT.add(i*j,sigma[i]*mu[j]); } } int solve(int n,int m) { int ret=0; for (int i=1,last=0;i<=n;i=last+1) { last=min(n/(n/i),m/(m/i)); ret+=(BT.gs(last)-BT.gs(i-1))*(n/i)*(m/i); } return ret&md; } int main() { init(); sort(rk+1,rk+maxn,cmp3); scanf("%d",&t); F(i,1,t) { scanf("%d%d%d",&a[i].n,&a[i].m,&a[i].k); if (a[i].n>a[i].m) swap(a[i].n,a[i].m); a[i].id=i; } sort(a+1,a+t+1,cmp); int now=0; F(i,1,t) { while (sigma[rk[now+1]]<=a[i].k) add(rk[++now]); a[i].ans=solve(a[i].n,a[i].m); } sort(a+1,a+t+1,cmp2); F(i,1,t) printf("%d ",a[i].ans); }