P2257 YY的GCD

博客详细介绍了如何运用莫比乌斯反演和数论分块技巧解决计算所有素数对$(i, j)$在给定区间内出现次数的问题。通过对求和次序的灵活变换,将复杂度降低到$O(n log log n)$,并提供了相应的埃氏筛和预处理算法实现。

Label

灵活变换求和次序的普通莫比乌斯反演

Description

给定T(T=104)T(T=10^4)T(T=104)n,m(1≤n,m≤107)n,m(1\leq n,m\leq 10^7)n,m(1n,m107),求:

∑i=1n∑j=1m[(i,j)∈prime]\sum_{i=1}^{n}\sum_{j=1}^{m}[(i,j)\in prime]i=1nj=1m[(i,j)prime]

Solution

∑i=1n∑j=1m[(i,j)∈prime]\sum_{i=1}^{n}\sum_{j=1}^{m}[(i,j)\in prime]i=1nj=1m[(i,j)prime]

=∑p∈primep≤min(n,m)∑i=1n∑j=1m[(i,j)=p]=\sum_{p\in prime}^{p\le min(n,m)}\sum_{i=1}^{n}\sum_{j=1}^{m}[(i,j)=p]=pprimepmin(n,m)i=1nj=1m[(i,j)=p]

=∑p∈primep≤min(n,m)∑i=1⌊np⌋∑j=1⌊mp⌋[(i,j)=1]=\sum_{p\in prime}^{p\le min(n,m)}\sum_{i=1}^{\lfloor\frac{n}{p}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{p}\rfloor}[(i,j)=1]=pprimepmin(n,m)i=1pnj=1pm[(i,j)=1]

=∑p∈primep≤min(n,m)∑i=1⌊np⌋∑j=1⌊mp⌋∑d∣(i,j)μ(d)=\sum_{p\in prime}^{p\le min(n,m)}\sum_{i=1}^{\lfloor\frac{n}{p}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{p}\rfloor}\sum_{d|(i,j)}\mu(d)=pprimepmin(n,m)i=1pnj=1pmd(i,j)μ(d)

=∑p∈primep≤min(n,m)∑d=1min(⌊np⌋,⌊mp⌋)μ(d)∑i=1⌊np⌋∑j=1⌊mp⌋[d∣(i,j)]=\sum_{p\in prime}^{p\le min(n,m)}\sum_{d=1}^{min(\lfloor\frac{n}{p}\rfloor,\lfloor\frac{m}{p}\rfloor)}\mu(d)\sum_{i=1}^{\lfloor\frac{n}{p}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{p}\rfloor}[d|(i,j)]=pprimepmin(n,m)d=1min(pn,pm)μ(d)i=1pnj=1pm[d(i,j)]

=∑p∈primep≤min(n,m)∑d=1min(⌊np⌋,⌊mp⌋)μ(d)∑i=1⌊npd⌋∑j=1⌊mpd⌋=\sum_{p\in prime}^{p\le min(n,m)}\sum_{d=1}^{min(\lfloor\frac{n}{p}\rfloor,\lfloor\frac{m}{p}\rfloor)}\mu(d)\sum_{i=1}^{\lfloor\frac{n}{pd}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{pd}\rfloor}=pprimepmin(n,m)d=1min(pn,pm)μ(d)i=1pdnj=1pdm

=∑p∈primep≤min(n,m)∑d=1min(⌊np⌋,⌊mp⌋)μ(d)⌊npd⌋⌊mpd⌋(1)=\sum_{p\in prime}^{p\le min(n,m)}\sum_{d=1}^{min(\lfloor\frac{n}{p}\rfloor,\lfloor\frac{m}{p}\rfloor)}\mu(d)\lfloor\frac{n}{pd}\rfloor\lfloor\frac{m}{pd}\rfloor(1)=pprimepmin(n,m)d=1min(pn,pm)μ(d)pdnpdm1

我们按照常规反演处理方法化简到这里时不难发现:后面的这个式子虽然可以数论分块处理,但此题的数据范围来讲,最前面求和符号代表的枚举素数操作的复杂度无法令人接受。此处,我们必须将素数求和号拿到式子后面去。

事实上,如果把ddd的求和号提到最前面,后面的⌊npd⌋⌊mpd⌋\lfloor\frac{n}{pd}\rfloor\lfloor\frac{m}{pd}\rfloorpdnpdm的计算还是无法避开枚举素数。事实上,最好的处理办法是:将原式改写为第一层求和符号枚举pdpdpd的式子。

考虑在当前最简式=∑p∈primep≤min(n,m)∑d=1min(⌊np⌋,⌊mp⌋)μ(d)⌊npd⌋⌊mpd⌋=\sum_{p\in prime}^{p\le min(n,m)}\sum_{d=1}^{min(\lfloor\frac{n}{p}\rfloor,\lfloor\frac{m}{p}\rfloor)}\mu(d)\lfloor\frac{n}{pd}\rfloor\lfloor\frac{m}{pd}\rfloor=pprimepmin(n,m)d=1min(pn,pm)μ(d)pdnpdmp,dp,dp,d的含义:枚举素数的同时,我们又枚举所有ppp的倍数pdpdpd。所以,考虑当pd=ipd=ipd=i时,p,dp,dp,d需满足什么条件才使得其对⌊ni⌋⌊mi⌋\lfloor\frac{n}{i}\rfloor\lfloor\frac{m}{i}\rfloorinim有贡献,我们可以直接将(1)式变换成如下式子:

∑i=1n⌊ni⌋⌊mi⌋∑p∈prime,p∣ip≤iμ(⌊ip⌋)(2)\sum_{i=1}^{n}\lfloor\frac{n}{i}\rfloor\lfloor\frac{m}{i}\rfloor\sum_{p\in prime,p|i}^{p\leq i}\mu(\lfloor\frac{i}{p}\rfloor)(2)i=1ninimpprime,pipiμ(pi)2

(其实上述式子pd(=i)pd(=i)pd(=i)下界应为2,但pd=1pd=1pd=1时对应的项值为0,对答案无影响。事实上,(从某篇文章中看到的)我们通过变换求和次序的好处之一就在于:大大省略了对答案无影响的0项的枚举计算。)

对于∑p∈prime,p∣ip≤iμ(⌊ip⌋)\sum_{p\in prime,p|i}^{p\leq i}\mu(\lfloor\frac{i}{p}\rfloor)pprime,pipiμ(pi),在O(n)O(n)O(n)预处理出μ(1)∼μ(n)\mu(1)\sim \mu(n)μ(1)μ(n)的值后,可利用未经优化的埃氏筛预处理。之后便可对(2)式进行数论分块。

算法时间复杂度:O(nloglogn)O(nloglogn)O(nloglogn)(瓶颈在于埃氏筛)。

Code

#include<cstdio>
#include<iostream>
#define ri register int
#define ll long long
using namespace std;

const int MAXN=1e7+20;
int T,N,M,prime[MAXN],cnt,mu[MAXN];
ll sum[MAXN],ans;
bool notprime[MAXN];

void Eulasieve()
{
	mu[1]=1,notprime[1]=true;
	for(ri i=2;i<=MAXN;++i)
	{
		if(!notprime[i]) prime[++cnt]=i,mu[i]=-1;
		for(ri j=1;j<=cnt&&i*prime[j]<=MAXN;++j)
		{
			notprime[i*prime[j]]=true;
			if(i%prime[j]==0) break;
			else mu[i*prime[j]]=-mu[i];
		}
	}
}

void Estsieve()
{
	sum[1]=0;
	for(ri i=2;i<=MAXN;++i)
		if(!notprime[i])
			for(ri j=1;i*j<=MAXN;++j) sum[i*j]+=(ll)mu[j];
	for(ri i=1;i<=MAXN;++i) sum[i]=sum[i-1]+sum[i];
}

int main()
{
	std::ios::sync_with_stdio(false);
	Eulasieve(),Estsieve();
	cin>>T;
	for(ri op=1;op<=T;++op)
	{
		cin>>N>>M;
		if(N>M) swap(N,M);
		ans=0LL;
		for(ri l=2,r;l<=N;l=r+1)
		{
			r=min(N/(N/l),M/(M/l));
			ans+=(ll)(N/l)*(ll)(M/l)*(sum[r]-sum[l-1]);
		}
		cout<<ans<<'\n';
	}
	return 0;
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值