BZOJ 3529 【SDOI2014】 数表

题目链接：数表

　　我们一起来膜PoPoQQQ大爷的题解吧Orz

　　首先我们来考虑没有\(a\)的限制该怎么做。显然交换\(n\)，\(m\)答案不变，所以后面默认\(n \le m\)。

　　我们定义两个函数：

\[f(x)=\sum_{d|x}d\]

\[g(x)=\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=x]\]

　　那么显然有：\[ans=\sum_{i=1}^nf(i)g(i)\]

　　\(g\)函数我们可以考虑化简一下：

\begin{aligned}
g(x)&=\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=x] \\
&=\sum_{i=1}^{\lfloor \frac{n}{x} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{x} \rfloor}\sum_{d|i,d|j}\mu(d) \\
&=\sum_{d=1}^{\lfloor \frac{n}{x} \rfloor}\mu(d)\lfloor \frac{n}{dx} \rfloor\lfloor \frac{m}{dx} \rfloor \\
&=\sum_{x|d}\mu(\frac{d}{x})\lfloor \frac{n}{d} \rfloor\lfloor \frac{m}{d} \rfloor
\end{aligned}

　　于是可以得到：

\begin{aligned}
ans=&\sum_{i=1}^nf(i)\sum_{i|d}\mu(\frac{d}{i})\lfloor \frac{n}{d} \rfloor\lfloor \frac{m}{d} \rfloor \\
=&\sum_{d=1}^n\lfloor \frac{n}{d} \rfloor\lfloor \frac{m}{d} \rfloor\sum_{i|d}f(i)\mu(\frac{d}{i})
\end{aligned}

　　我们令\(h(x)=\sum_{i|x}f(i)\mu(\frac{x}{i})\)，那么我们只要有\(h(x)\)的前缀和，就可以在\(O(\sqrt{n})\)的时间内求出\(ans\)了。

　　然后我们来考虑\(a\)的限制。显然，只有\(f(x)\le a\)的\(f(x)\)才会对答案有影响。所以，我们把所有的询问按\(a\)从小到大排好序，并且把所有的\(f(x)\)从小到大依次加进来，使用树状数组维护前缀和，依次处理每个询问，就可以在\(O(Q\sqrt{n}\log{n}+n\log^2 n)\)的时间内解决所有询问了。

　　下面贴代码：

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#define File(s) freopen(s".in","r",stdin),freopen(s".out","w",stdout)
#define maxn 100010

using namespace std;
typedef long long llg;

struct data{
int n,m,a,b;
bool operator < (const data &h)const{return a<h.a;}
}s[maxn];
int n,m,nm,a[maxn],T,c[maxn],ans[maxn];
int mu[maxn],f[maxn],pr[maxn],lp;
const int mod=2147483647;
bool vis[maxn];

int getint(){
int w=0;bool q=0;
char c=getchar();
while((c>'9'||c<'0')&&c!='-') c=getchar();
if(c=='-') c=getchar(),q=1;
while(c>='0'&&c<='9') w=w*10+c-'0',c=getchar();
return q?-w:w;
}

bool cmp(int x,int y){return f[x]<f[y];}
void add(int x,int y){if(y) while(x<=nm) c[x]+=y,x+=x&(-x);}
int sum(int x){
int t=0;
while(x) t+=c[x],x-=x&(-x);
return t;
}

int main(){
File("a");
mu[1]=f[1]=a[1]=1; T=getint();
for(int i=1;s[i].b=i,i<=T;i++){
s[i].n=getint(),s[i].m=getint();
if(s[i].n>s[i].m) swap(s[i].n,s[i].m);
s[i].a=getint(); nm=max(nm,s[i].n);
}
for(int i=2;a[i]=i,i<=nm;i++){
if(!vis[i]) pr[++lp]=i,mu[i]=-1,f[i]=i+1;
for(int j=1;pr[j]*i<=nm;j++){
vis[pr[j]*i]=1;
if(i%pr[j]) mu[pr[j]*i]=-mu[i],f[pr[j]*i]=f[i]*f[pr[j]];
else{ f[pr[j]*i]=f[i]+(f[i]-f[i/pr[j]])*pr[j]; break;}
}
}
sort(s+1,s+T+1); sort(a+1,a+nm+1,cmp);
for(int i=1,now=1,b,la=0;i<=T;i++){
n=s[i].n; m=s[i].m; b=s[i].b;
while(now<=nm && f[a[now]]<=s[i].a){
for(int j=a[now];j<=nm;j+=a[now])
add(j,f[a[now]]*mu[j/a[now]]);
now++;
}
for(int j=1,nt,na;j<=n;j=nt+1){
nt=min(n/(n/j),m/(m/j)); na=sum(nt);
ans[b]+=(na-la)*(n/j)*(m/j); la=na;
}
}
for(int i=1;i<=T;i++) printf("%d\n",ans[i]&mod);
return 0;
}