vii.[SDOI2014]数表
仍然是线性筛筛各种东西。我们引出一个东西\(\sigma(x)=\sum\limits_{d|n}d\)，也就是\(x\)的约数和。这个东西明显是积性函数。设\(x=\prod\limits_{i=1}^n(P_i)^{a_i}\)，则\(\sigma(x)=\prod\limits_{i=1}^n(\sum\limits_{j=0}^{a_i}(P_i)^j)\)。
如果我们忽略\(a\)的限制，则我们要求的就是\(\sum\limits_{i=1}^n\sum\limits_{j=1}^m\sigma(\gcd(i,j))\)。
开始推式子。
\(\begin{aligned}\sum_{i=1}^n\sum_{j=1}^m\sigma(\gcd(i,j)) & =\sum_{d=1}^{\min(n,m)}\sum_{i=1}^n\sum_{j=1}^m\sigma(d)[\gcd(i,j)=d]\\& =\sum_{d=1}^{\min(n,m)}\sigma(d)\sum_{i=1}^{n/d}\sum_{j=1}^{m/d}[\gcd(i,j)=1]\\& =\sum_{d=1}^{\min(n,m)}\sigma(d)\sum_{i=1}^{n/d}\sum_{j=1}^{m/d}\sum_{x|i,x|j}\mu(x)\\& =\sum_{d=1}^{\min(n,m)}\sigma(d)\sum_{x=1}^{\min(n/d,m/d)}\mu(x)\left\lfloor\dfrac{n}{dx}\right\rfloor\left\lfloor\dfrac{m}{dx}\right\rfloor\\& =\sum_{T=1}^{\min(n,m)}\left\lfloor\dfrac{n}{T}\right\rfloor\left\lfloor\dfrac{m}{T}\right\rfloor\sum_{x|T}\sigma(x)\mu(\dfrac{T}{x})\end{aligned}\)
假使我们忽略\(a\)的限制，到现在我们已经做完了。但是现在有\(a\)的限制，怎么办呢？
设\(g=\mu*\sigma\)。
将询问按照\(a\)排序，同时用树状数组动态地对于不同的\(a\)更新\(g\)函数。比如说，假设上一个\(a\)是\(5\)，这一个\(a\)是\(6\)，那对于所有\(\sigma(x)=6\)的\(x\)，枚举在\(10^5\)范围内所有\(x\)的倍数\(T\)，更新\(g(T)\)。
而我们最后统计时，就用树状数组维护\(g\)函数的区间和。
哦，对了，关于取模，直接自然溢出即可。
关于线性筛筛各种东西，我推荐一篇blog。里面就有线性筛筛\(\sigma(x)\)的讲解。
代码：
#include<bits/stdc++.h>
using namespace std;
const int N=100000;
int T,pri[100100],mu[100100],sigma[100100],low[100100],sum[100100],t[100100],res[100100];;
pair<int,int>p[100100];
void sieve(){
mu[1]=sigma[1]=1;
for(int i=2;i<=N;i++){
if(!pri[i])low[i]=pri[++pri[0]]=i,mu[i]=-1,sum[i]=sigma[i]=i+1;
for(int j=1;j<=pri[0]&&i*pri[j]<=N;j++){
pri[i*pri[j]]=true;
if(!(i%pri[j])){
mu[i*pri[j]]=0;
low[i*pri[j]]=low[i]*pri[j];
sum[i*pri[j]]=sum[i]+low[i*pri[j]];
sigma[i*pri[j]]=sigma[i]/sum[i]*sum[i*pri[j]];
break;
}
mu[i*pri[j]]=-mu[i];
low[i*pri[j]]=pri[j];
sum[i*pri[j]]=pri[j]+1;
sigma[i*pri[j]]=sigma[i]*sigma[pri[j]];
}
}
for(int i=1;i<=N;i++)p[i]=make_pair(sigma[i],i);
sort(p+1,p+N+1);
}
void add(int x,int y){
while(x<=N)t[x]+=y,x+=x&-x;
}
int ask(int x){
int res=0;
while(x)res+=t[x],x-=x&-x;
return res;
}
struct query{
int n,m,a,id;
friend bool operator <(const query &x,const query &y){
return x.a<y.a;
}
}q[100100];
int solve(int n,int m){
if(n>m)swap(n,m);
int ans=0;
for(int l=1,r;l<=n;l=r+1)r=min(n/(n/l),m/(m/l)),ans+=(n/l)*(m/l)*(ask(r)-ask(l-1));
return ans;
}
int main(){
scanf("%d",&T),sieve();
for(int i=1;i<=T;i++)scanf("%d%d%d",&q[i].n,&q[i].m,&q[i].a),q[i].id=i;
sort(q+1,q+T+1);
for(int i=1,j=1;i<=T;i++){
while(j<=N&&p[j].first<=q[i].a){
for(int k=p[j].second;k<=N;k+=p[j].second)add(k,sigma[p[j].second]*mu[k/p[j].second]);
j++;
}
res[q[i].id]=solve(q[i].n,q[i].m);
}
for(int i=1;i<=T;i++)printf("%d\n",res[i]&(~(1<<31)));
return 0;
}