【洛谷4466】[国家集训队] 和与积（莫比乌斯反演）

点此看题面

• 给定\(n\)，求有多少对\(a,b\)满足\(1\le a<b\le n\)且\((a+b)|ab\)。
• \(n\le2^{31}-1\)

初始推式子

说实话，最近翻了翻数论小蓝本，于是看到这种式子首先就想到令\(d=gcd(a,b)\)得到\(a_0=\frac ad,b_0=\frac bd\)互质。

然后就发现原式相当于是：

\[(a_0+b_0)d|a_0b_0d^2\Leftrightarrow (a_0+b_0)|a_0b_0d \]

由于\(gcd(a_0+b_0,a_0)=gcd(b_0,a_0)=1\)（同理\(gcd(a_0+b_0,b_0)=1\)），所以我们可以直接消去整除符号右边的\(a_0,b_0\)，得到一个非常漂亮的式子：

\[(a_0+b_0)|d \]

考虑我们现在的几个限制：

\[\begin{cases} b=b_0d\le n\Rightarrow d\le\lfloor\frac n{b_0}\rfloor\\ (a_0+b_0)|d\Rightarrow a_0+b_0\le d \end{cases} \]

综合两个条件得到\(a_0+b_0\le\lfloor\frac n{b_0}\rfloor\)。

所以说，我们发现，\(a_0,b_0,d\)都应该是不超过\(\sqrt n\)的。

莫比乌斯反演

根据上面的结论，答案的计算式就应该是：

\[\sum_{a_0=1}^{\sqrt n}\sum_{b_0=a_0+1}^{\sqrt n}\lfloor\frac{\lfloor\frac{n}{b_0}\rfloor}{a_0+b_0}\rfloor[gcd(a_0,b_0)=1] \]

经典莫比乌斯反演，枚举\(gcd\)得到：

\[\sum_{d=1}^{\sqrt n}\mu(d)\sum_{a_0=1}^{\lfloor\frac{\sqrt n}d\rfloor}\sum_{b_0=a_0+1}^{\lfloor\frac{\sqrt n}d\rfloor}\lfloor\frac{\lfloor\frac n{b_0d}\rfloor}{(a_0+b_0)d}\rfloor \]

其中大分母中的\(d\)可以移到分子的分母中：

\[\sum_{d=1}^{\sqrt n}\mu(d)\sum_{a_0=1}^{\lfloor\frac{\sqrt n}d\rfloor}\sum_{b_0=a_0+1}^{\lfloor\frac{\sqrt n}d\rfloor}\lfloor\frac{\lfloor\frac n{b_0d^2}\rfloor}{a_0+b_0}\rfloor \]

我们发现\(a_0\)只出现在分母\(a_0+b_0\)这一项中，因此我们不如不枚举\(a_0\)，而是去枚举分母，得到：

\[\sum_{d=1}^{\sqrt n}\mu(d)\sum_{b_0=2}^{\lfloor\frac{\sqrt n}d\rfloor}\sum_{t=b_0+1}^{2b_0-1}\lfloor\frac{\lfloor\frac n{b_0d^2}\rfloor}{t}\rfloor \]

发现此时的\(\sum_{t=b_0+1}^{2b_0-1}\lfloor\frac{\lfloor\frac n{b_0d^2}\rfloor}{t}\rfloor\)已经成为一个标准的除法分块形式。

那么我们只要以共计\(O(\sqrt nlog\sqrt n)\)的复杂度枚举\(d,b_0\)，然后除法分块\(O(\sqrt{\sqrt n})\)枚举\(t\)即可解决此题。

代码：\(O(n^{\frac34}\log\sqrt n)\)

#include<bits/stdc++.h>
#define Tp template<typename Ty>
#define Ts template<typename Ty,typename... Ar>
#define Reg register
#define RI Reg int
#define Con const
#define CI Con int&
#define I inline
#define W while
#define SN 47000
#define LL long long
using namespace std;
int n,sn;struct LinearSieve
{
	int Pt,P[SN+5],mu[SN+5];I int operator [] (CI x) {return mu[x];}
	I LinearSieve()//线性筛预处理μ
	{
		mu[1]=1;for(RI i=2,j;i<=SN;++i) for(!P[i]&&(mu[P[++Pt]=i]=-1),
			j=1;j<=Pt&&i*P[j]<=SN;++j) if(P[i*P[j]]=1,i%P[j]) mu[i*P[j]]=-mu[i];else break;
	}
}Mu;
I LL Calc(CI L,CI R,CI n)//除法分块
{
	RI l,r;LL t=0;for(l=L;l<=min(n,R);l=r+1) r=min(R,n/(n/l)),t+=1LL*(n/l)*(r-l+1);return t;
}
int main()
{
	RI i,j;LL t=0;for(scanf("%d",&n),sn=sqrt(n),i=1;i<=sn;++i)//枚举d
		for(j=2;i*j<=sn;++j) t+=Mu[i]*Calc(j+1,2*j-1,n/i/i/j);return printf("%lld\n",t),0;//枚举b0，然后除法分块
}