【SDOI2014】数表【莫比乌斯反演】【树状数组】

传送门

传送门

题意：

$T$组询问给定$n,m,a$,求

$$\sum _{i=1}^N\sum _{j=1}^M[(\sum _{d\in N^{\ast }}[d\mid i][d\mid j]d)\le a]\sum _{d\in N^{\ast }}[d\mid i][d\mid j]d$$

$1\le n,m\le 1e5,1\le a\le 1e9,1\le T\le 2e4$

对$2^{31}$取模

反演神题

假设没有前面的限制

$$\sum _{i=1}^N\sum _{j=1}^M\sum _{d\in N^{\ast }}[d\mid i][d\mid j]$$

记$F(x)$表示约数和

$$\sum _{i=1}^N\sum _{j=1}^{M}F(gcd(i,j))$$

套路性地枚举$gcd$

$$\sum _d\sum _{i=1}^N\sum _{j=1}^M[gcd(i,j)=d]F(d)$$

套路性地提到前面

$$\sum _{d}F(d)\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=d]$$

后面是个套路性的反演

设

$$f(d)=\sum _{i=1}^N\sum _{j=1}^M[gcd(i,j)=d]$$

$$g(d)=\sum _{d\mid n}f(n)=\lfloor \frac{N}{d}\rfloor \lfloor \frac{M}{d}\rfloor$$

$$f(d)=\sum _{d\mid n}g(n)\mu (\frac{n}{d})=\sum _{d\mid n}\lfloor \frac{N}{n}\rfloor \lfloor \frac{M}{n}\rfloor \mu (\frac{n}{d})$$

代回去

$$\sum _{d}F(d)\sum _{d\mid n}\lfloor \frac{N}{n}\rfloor \lfloor \frac{M}{n}\rfloor \mu (\frac{n}{d})$$

$$\sum _{n=1}^{min(N,M)}\lfloor \frac{N}{n}\rfloor \lfloor \frac{M}{n}\rfloor \sum _{d\mid n}F(d)\mu (\frac{n}{d})$$

左边是个套路性的整除分块 不管

右边并没有很好的性质，只能筛出来暴力求……

等下，原题是有限制$F(d)\le a$的

设

$$G(n)=\sum _{d\mid n}F(d)\mu (\frac{n}{d})$$

我们发现只有不超过$a$的$F(d)$才对答案有贡献

而$T,N$都很小

所以可以把$F(d)$和$a$排个序，拿个指针指一下

跑到一个$a$把前面的$F(d)$枚举倍数统计更新$G$,再拿整除分块算答案

需要维护单点加，区间询问，用个树状数组即可

复杂度大概$O(NlogN+TlogT\sqrt{T})$

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cctype>
#include <algorithm>
#define MAXN 100005
using namespace std;
const int N=100000;
int np[MAXN],pl[MAXN],cnt;
int mu[MAXN],f[MAXN],mp[MAXN],pt[MAXN];
inline int qpow(int a,int p)
{
	int ans=1;
	while (p)
	{
		if (p&1) ans*=a;
		a*=a;p>>=1;
	}
	return ans;
}
void init()
{
	np[1]=mu[1]=f[1]=1;
	for (int i=2;i<=N;++i)
	{
		if (!np[i])
		{
			pl[++cnt]=i;
			mu[i]=-1;
			mp[i]=i;pt[i]=1;
		}
		int x;
		for (int j=1;(x=i*pl[j])<=N;++j)
		{
			np[x]=1;mp[x]=pl[j];
			if (i%pl[j]==0)
			{
				mu[x]=0;pt[x]=pt[i]+1;
				break;
			}
			mu[x]=-mu[i];pt[x]=1;
		}
	}
	for (int i=2;i<=N;i++)
	{
		int t=qpow(mp[i],pt[i]);
		if (i==t) f[i]=f[i/mp[i]]+t;
		else f[i]=f[t]*f[i/t];
	}
}
int p[MAXN];
inline bool cmp(const int& a,const int& b){return f[a]<f[b];}
struct query{int n,m,a,id;}q[MAXN];
inline bool operator <(const query& x,const query& y){return x.a<y.a;}
struct BIT
{
	int s[MAXN];
	inline int lowbit(const int& x){return x&-x;}
	inline void modify(int x,const int& v){for (;x<=N;s[x]+=v,x+=lowbit(x));}
	inline int query(int x){int ans=0;for (;x;ans+=s[x],x-=lowbit(x));return ans;}
}g;
int res[MAXN];
//int calc(int x,int y,int a)
//{
//	int ans=0;
//	for (int i=1;i<=x&&i<=y;i++)
//		if (x%i==0&&y%i==0)
//			ans+=i;
//	if (ans<=a) return ans;
//	else return 0;
//}
int main()
{
	init();
	for (int i=1;i<=N;i++) p[i]=i;
	sort(p+1,p+N+1,cmp);
	int pos=0;
	int T;
	scanf("%d",&T);
	for (int i=1;i<=T;i++) 
	{
		scanf("%d%d%d",&q[i].n,&q[i].m,&q[i].a);q[i].id=i;
//		int t=0;
//		for (int x=1;x<=q[i].n;x++)
//			for (int y=1;y<=q[i].m;y++)
//				t+=calc(x,y,q[i].a);
//		printf("%d\n",t&(~(1<<31)));
	}
	sort(q+1,q+T+1);
	for (int i=1;i<=T;i++)
	{
		while (pos<N&&f[p[pos+1]]<=q[i].a)
		{
			++pos;
			for (int d=1;d*p[pos]<=N;++d) g.modify(d*p[pos],f[p[pos]]*mu[d]);
		}
		if (q[i].n>q[i].m) swap(q[i].n,q[i].m);
		int ans=0;
		for (int l=1,r;l<=q[i].n;l=r+1)
		{
			r=min(q[i].n/(q[i].n/l),q[i].m/(q[i].m/l));
			ans+=(q[i].n/l)*(q[i].m/l)*(g.query(r)-g.query(l-1));
		}
		res[q[i].id]=ans;
	}
	for (int i=1;i<=T;i++) printf("%d\n",res[i]&(~(1<<31)));
	return 0;
}