【BZOJ3529】数表（莫比乌斯反演，树状数组）

题解

首先不管\(A\)的范围的限制

要求的东西是

其中\(\sigma(x)\)表示\(x\)的约数之和

约数之和是一个积性函数，可以线性筛

具体的做法请参考皮皮亮的Blog

根据常见的套路

把\(gcd\)给提出来

后面那个东西不用说了吧。。

我原来已经推过很多遍啦

可以参考原来推过的式子

所以我直接写啦

还是一样的，设\(T=id\)

然后提出来

如果没有\(A\)的限制

这就可以直接筛后面的那个狄利克雷卷积

然后\(O(\sqrt n)\)的回答每一组询问

现在有了\(A\)的限制，貌似筛出来也没有什么用了

考虑到\(n,m<=10^5\)应该可以乱搞了

所以干脆就先筛出\(\sigma\)和\(\mu\)

然后，这个限制还是很不好搞？？？

范围这么小

那我就暴力呀

离线询问后按照\(A\)排序

然后把\(\sigma\)排序

依次插到相对应的位置

因为要求前缀，就搞一个树状数组来算前缀和就行啦

最后是取膜，为啥要用\(long\ long\)

\(int\)自然溢出就行啦

暴力出奇迹

#include<iostream>

#include<cstdio>

#include<cstdlib>

#include<cstring>

#include<cmath>

#include<algorithm>

#include<set>

#include<map>

#include<vector>

#include<queue>

using namespace std;

#define MAX 100000

inline int read()

{

	int x=0,t=1;char ch=getchar();

	while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();

	if(ch=='-')t=-1,ch=getchar();

	while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();

	return x*t;

}

bool zs[MAX+1];

int pri[MAX],tot;

int mu[MAX+1],sig[MAX+1],sumd[MAX+1],powd[MAX+1];

int T,N,ans[MAX];

struct Ask{int n,m,A,id;}q[MAX];

int p[MAX+1],c[MAX+1];

bool operator<(Ask a,Ask b){return a.A<b.A;}

void pre(int N)

{

	zs[1]=true;sig[1]=mu[1]=p[1]=1;

	for(int i=2;i<=N;++i)

	{

		p[i]=i;

		if(!zs[i])pri[++tot]=i,mu[i]=-1,sig[i]=i+1,sumd[i]=i+1,powd[i]=i;

		for(int j=1;j<=tot&&i*pri[j]<=N;++j)

		{

			zs[i*pri[j]]=true;

			if(i%pri[j])

			{

				mu[i*pri[j]]=-mu[i];

				sig[i*pri[j]]=sig[i]*sig[pri[j]];

				sumd[i*pri[j]]=1+pri[j];

				powd[i*pri[j]]=pri[j];

			}

			else

			{

				powd[i*pri[j]]=powd[i]*pri[j];

				sumd[i*pri[j]]=sumd[i]+powd[i*pri[j]];

				sig[i*pri[j]]=sig[i]/sumd[i]*sumd[i*pri[j]];

				break;

			}

		}

	}

}

bool cmp(int a,int b){return sig[a]<sig[b];}

inline int lowbit(int x){return x&(-x);}

void Add(int x,int w){while(x<=N)c[x]+=w,x+=lowbit(x);}

int Query(int x){int ret=0;while(x)ret+=c[x],x-=lowbit(x);return ret;}

int Solve(int x,int y)

{

	int now=0,lst=0,ret=0;

	for(int i=1,j;i<=x;i=j+1)

	{

		j=min(x/(x/i),y/(y/i));

		now=Query(j);

		ret+=(x/i)*(y/i)*(now-lst);

		lst=now;

	}

	return ret;

}

int main()

{

	T=read();

	for(int i=1;i<=T;++i)

	{

		q[i].id=i,q[i].n=read(),q[i].m=read(),q[i].A=read();

		if(q[i].n>q[i].m)swap(q[i].n,q[i].m);

		N=max(N,q[i].n);

	}

	pre(N);

	sort(&q[1],&q[T+1]);sort(&p[1],&p[N+1],cmp);

	for(int i=1,j=1;i<=T;++i)

	{

		for(;j<=N&&sig[p[j]]<=q[i].A;++j)

			for(int k=p[j];k<=N;k+=p[j])

				if(mu[k/p[j]])Add(k,sig[p[j]]*mu[k/p[j]]);

		ans[q[i].id]=Solve(q[i].n,q[i].m);

	}

	for(int i=1;i<=T;++i)

	{

		if(ans[i]<0)ans[i]+=2147483647,ans[i]++;

		printf("%d\n",ans[i]);

	}

	return 0;

}