$ans=\sum_{i=1}^n\sum_{j=1}^n\sigma(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 \rfloor )*\lfloor n/d \rfloor * \lfloor m/d \rfloor$
然后代入然后xjb变换可得
$ans=\sum_{d<=n}\lfloor n/d \rfloor * \lfloor m/d \rfloor \sum_{i \mid d}\mu( \lfloor d/i \rfloor ) * \sigma(i) $
然后我们要求出$\sum_{i \mid d}\mu(\lfloor d/i \rfloor ) *\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\n",a[i].ans);
}