P3312 [SDOI2014] 数表- 数论

题解

对于一个询问 \((n,m,a)\)，答案为：

\[\begin{aligned} &\sum_{i=1}^n \sum_{j=1}^m \sigma_1(\gcd(i,j))[\sigma_1(\gcd(i,j))\le a]\\ &=\sum_{i=1}^n \sum_{j=1}^m [\sigma_1(\gcd(i,j))\le a]\sum_{d\mid i,d\mid j} d\\ &=\sum_{d} \sigma_1(d)\cdot [\sigma_1(d)\le a]\cdot \sum_{i=1}^{n/d} \sum_{j=1}^{m/d} [i\perp j]\\ &=\sum_{d} \sigma_1(d)\cdot [\sigma_1(d)\le a]\cdot \sum_{i=1}^{n/d} \sum_{j=1}^{m/d} \epsilon(\gcd(i,j))\\ &=\sum_{d} \sigma_1(d)\cdot [\sigma_1(d)\le a]\cdot \sum_{i=1}^{n/d} \sum_{j=1}^{m/d} \sum_{e\mid i,e\mid j} \mu(e)\\ &=\sum_{d} \sigma_1(d)\cdot [\sigma_1(d)\le a]\cdot \sum_{e} \mu(e) \lfloor \frac{n}{de}\rfloor \lfloor \frac{m}{de}\rfloor\\ &=\sum_T \lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \sum_{de=T} \sigma_1(d)\cdot [\sigma_1(d)\le a]\cdot \mu(e) \end{aligned} \]

设 \(f(n)=\begin{cases} \sigma_1(n) & \sigma_1(n)\le a \\ 0\end{cases}\)，\(g(n)=f*\mu\)。

那么答案为：

\[\sum_{T} \lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor g(T) \]

记 \(N=\max_{\text{all queries}}(n,m)\)。

我们将所有询问离线并按 \(a\) 从小到大排序。用类似扫描线的思想，我们每次会将 \(f\) 的某几个位置的值从 \(0\) 改成 \(\sigma_0\)，总共会有 \(\mathcal{O}(N)\) 次更改。而每次更改 \(f(i)\) 会对 \(g\) 中是 \(i\) 的倍数的位置产生贡献。根据调和级数的结论，我们一共会更改 \(\mathcal{O}(N\log N)\) 次 \(g\)。而整除分块要求我们快速求出 \(g\) 的前缀和，因此我们要用树状数组维护 \(g\)。

复杂度 \(\mathcal{O}(N+T(\sqrt{n}\log n+\log^2 n))\)。

没想到的地方是枚举 \(T=de\)。

代码

#include <iostream>
#include <cstring>
#include <algorithm>
#include <vector>
using namespace std;
#define For(Ti,Ta,Tb) for(int Ti=(Ta);Ti<=(Tb);++Ti)
#define Dec(Ti,Ta,Tb) for(int Ti=(Ta);Ti>=(Tb);--Ti)
#define Debug(...) fprintf(stderr,__VA_ARGS__)
typedef long long ll;
const int Q=2e4+5,N=1e5+5;
const ll Mod=1LL<<31;
struct Fenwick{
	int lim;ll t[N]{};
	void Add(int p,ll k){for(;p<=lim;p+=p&-p) (t[p]+=k)%=Mod;}
	ll Query(int p){ll res=0;for(;p;p-=p&-p) (res+=t[p])%=Mod; return res;}
};
int prcnt,prime[N],notp[N],mu[N];
void Sieve(int mx){
	notp[1]=mu[1]=1;
	For(i,2,mx){
		if(!notp[i]) prime[++prcnt]=i,mu[i]=-1;
		for(int j=1;j<=prcnt&&i*prime[j]<=mx;++j){
			notp[i*prime[j]]=1;
			if(i%prime[j]==0){mu[i*prime[j]]=0;break;}
			mu[i*prime[j]]=-mu[i];
		}
	}
}
int q;ll s0[N],ans[N];
struct Query{int i,n,m,a;}qry[Q];
struct Command{ll a;int d;}cmd[N];
int main(){
	ios::sync_with_stdio(false),cin.tie(nullptr);
	cin>>q;
	int mx=0;
	For(i,1,q){
		cin>>qry[i].n>>qry[i].m>>qry[i].a;
		qry[i].i=i;
		mx=max({mx,qry[i].n,qry[i].m});
	}
	sort(qry+1,qry+q+1,[](const Query &q1,const Query &q2){return q1.a<q2.a;});
	Sieve(mx);
	For(i,1,mx) For(j,1,mx/i) s0[i*j]+=i;
	For(i,1,mx) cmd[i]={s0[i],i};
	sort(cmd+1,cmd+mx+1,[](const Command &c1,const Command &c2){return c1.a<c2.a;});
	Fenwick bit{mx};
	for(int i=1,j=1;i<=q;++i){
		for(;j<=mx&&cmd[j].a<=qry[i].a;++j){
			for(int k=1;cmd[j].d*k<=mx;++k) bit.Add(cmd[j].d*k,cmd[j].a*mu[k]);
		}
		int n=qry[i].n,m=qry[i].m;
		for(int l=1,r;l<=min(n,m);l=r+1){
			r=min(n/(n/l),m/(m/l));
			(ans[qry[i].i]+=(bit.Query(r)-bit.Query(l-1)+Mod)*(n/l)*(m/l))%=Mod;
		}
	}
	For(i,1,q) cout<<ans[i]<<'\n';
	return 0;
}