#树状数组，莫比乌斯反演，整除分块#bzoj 3529 JZOJ 3623 洛谷 3312 数表

题目

求$$\sum _{i=1}^n\sum _{j=1}^m{\sigma }_{gcd(i,j)}[{\sigma }_{gcd(i,j)}\le w]$$
多组数据

---

分析

先把后面条件扔掉,默认$n\le m$
$$\sum _{i=1}^n\sum _{j=1}^m{\sigma }_{gcd(i,j)}=\sum _{d=1}^n\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }{\sigma }_d[gcd(i,j)==1]$$
$$=\sum _{d=1}^n\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }\sum _{k|i,k|j}\mu (k)$$
$$=\sum _{d=1}^n\sum _{k=1}^{\lfloor \frac{n}{d}\rfloor }\mu (k)\lfloor \frac{n}{kd}\rfloor \lfloor \frac{m}{kd}\rfloor$$
令$i=kd$，那么
$$=\sum _{i=1}^n\lfloor \frac{n}{i}\rfloor \lfloor \frac{m}{i}\rfloor \sum _{d|i}{\sigma }_d\mu (\frac{i}{d})$$
主要是后面这个部分，由于${\sigma }_d\le w$的限制，考虑离线询问，用树状数组实时维护，用枚举倍数的方式，时间复杂度
$$O(q\sqrt{n}logn+nlo{g}^{2}n)$$

---

代码

#include <cstdio>
#include <cctype>
#include <algorithm>
#define rr register
#define min(a,b) ((a)<(b)?(a):(b))
using namespace std;
const int N=100101,M=20101;
struct rec{int n,m,w,rk;}b[M];
int si[N],mu[N],rk[N],v[N],prime[M>>1],t,cnt,ans[M],c[N];
inline signed iut(){
	rr int ans=0; rr char c=getchar();
	while (!isdigit(c)) c=getchar();
	while (isdigit(c)) ans=(ans<<3)+(ans<<1)+(c^48),c=getchar();
	return ans;
}
bool cmp1(int x,int y){return si[x]<si[y];}
inline void init(int n){
	si[1]=mu[1]=rk[1]=1;
	for (rr int i=2;i<=n;++i){
		rk[i]=i;
		if (!v[i]) v[i]=si[i]=i+1,mu[prime[++cnt]=i]=-1;
		for (rr int j=1;j<=cnt&&prime[j]*i<=n;++j){
			if (!(i%prime[j])){
				v[i*prime[j]]=v[i]*prime[j]+1;
				si[i*prime[j]]=si[i]/v[i]*v[i*prime[j]];
				break;
			}
			mu[i*prime[j]]=-mu[i];
			v[i*prime[j]]=prime[j]+1;
			si[i*prime[j]]=si[i]*v[i*prime[j]];
		}
	}
	sort(rk+1,rk+1+n,cmp1);
}
bool cmp(rec x,rec y){return x.w<y.w;}
inline void add(int x,int y){while (x<=N) c[x]+=y,x+=-x&x;}
inline signed anss(int x){
	rr int ans=0;
	while (x) ans+=c[x],x-=-x&x;
	return ans;
}
inline signed answ(int n,int m){
	rr int ans=0,ls=0,now;
	for (rr int l=1,r;l<=n;l=r+1,ls=now){
		r=min(n/(n/l),m/(m/l)),now=anss(r);
		ans+=(n/l)*(m/l)*(now-ls);
	}
	return ans;
}
void print(int ans){
	if (ans>9) print(ans/10);
	putchar(ans%10+48);
}
signed main(){
	init(N); rr int t=iut();
	for (rr int i=1;i<=t;++i){
		rr int n=iut(),m=iut(),w=iut();
		if (n>m) n^=m,m^=n,n^=m;
		b[i]=(rec){n,m,w,i};
	}
	sort(b+1,b+1+t,cmp);
	for (rr int i=1,j=1;i<=t;++i){
		for (;j<=N&&si[rk[j]]<=b[i].w;++j)
		    for (rr int k=1;k*rk[j]<=N;++k)
		        if (mu[k]) add(k*rk[j],mu[k]*si[rk[j]]);
		ans[b[i].rk]=answ(b[i].n,b[i].m);
	}
	for (rr int i=1;i<=t;++i) print(ans[i]+(ans[i]<0)*(1ll<<31)),putchar(10);
	return 0;
}