【JZOJ4872】太阳神（莫比乌斯反演）（数论分块）

题意：

给定一个$n$，求有多少对$(a,b)$满足$a\le n,b\le n,lcm(a,b)>n$，$n\le 1e10$

---

题解：

首先转化一下，计算$a\le n,b\le n,lcm(a,b)\le n$的对数。

$$\begin{array}{rl}Ans=& \sum _{a=1}^n\sum _{b=1}^n[lcm(a,b)\le n]\\ =& \sum _{a=1}\sum _{b=1}[\frac{a\cdot b}{gcd(a,b)}\le n]\\ =& \sum _{d=1}\sum _{a=1}\sum _{b=1}[gcd(a,b)=d][\frac{a\cdot b}{d}\le n]\\ =& \sum _{d=1}\sum _{a=1}\sum _{b=1}[gcd(a,b)=1][abd\le n]\\ =& \sum _{k=1}\mu (k)\sum _{d=1}\sum _{a=1}\sum _{b=1}[abd\le \frac{n}{k^2}]\end{array}$$

所以我们需要计算的就是$abc\le m$的有序三元组$(a,b,c)$个数。

令$a<b<c$，相等的情况可以特殊处理。枚举$a\in [1,\sqrt[3]{m}],b\in [a+1,\sqrt{n/a}]$，$c$可以直接统计。

复杂度可以通过积分来估计。

对于$m$进行一次这样的计算的复杂度为：

$$\begin{array}{rl}T(m)=& {\int }_1^{\sqrt[3]{m}}\sqrt{\frac{m}{x}}dx\\ =& \sqrt{m}\cdot {\int }_1^{\sqrt[3]{m}}{x}^{-\frac{1}{2}}dx\\ =& m^{\frac{1}{2}}\cdot O(m^{\frac{1}{6}}-o(1))\\ =& O(m^{\frac{2}{3}})\end{array}$$

那么总复杂度为：

$$\begin{array}{rl}Q(n)=& {\int }_1^{\sqrt{n}}(\frac{n}{x^2}{)}^{\frac{2}{3}}dx\\ =& n^{\frac{2}{3}}{\int }_1^{\sqrt{n}}{x}^{-\frac{4}{3}}dx\\ =& n^{\frac{2}{3}}\cdot O(n^{-\frac{1}{6}}+o(1))\\ =& O(n^{\frac{2}{3}})\end{array}$$

---

代码：

#include<bits/stdc++.h>
#define ll long long
#define re register
#define cs const

using std::cerr;
using std::cout;

cs int mod=1e9+7;

cs int N=1e5+7;
int p[N],pc;
bool mark[N];
int mu[N];

inline void linear_sieves(int lim){
	mu[1]=1;
	for(int re i=2;i<=lim;++i){
		if(!mark[i])p[++pc]=i,mu[i]=-1;
		for(int re j=1;i*p[j]<=lim;++j){
			mark[i*p[j]]=true;
			if(i%p[j])mu[i*p[j]]=-mu[i];
			else break;
		}
	}
} 

ll last,res;

inline ll solve(ll n){
	if(n==last)return res;last=n,res=0;
	for(int a=1;(ll)a*a*a<=n;++a){
		res+=3ll*(n/a/a-a)+1,res%=mod;
		for(int re b=a+1,lb=sqrt(n/a);b<=lb;++b)
		res+=6ll*(n/a/b-b)+3,res%=mod;
	}
	return res;
}

signed main(){
#ifdef zxyoi
	freopen("sun.in","r",stdin);
#endif
	ll n;scanf("%lld",&n);
	int lim=sqrt(n);linear_sieves(lim);
	ll ans=0;last=-1;
	for(int re i=1;i<=lim;++i)if(mu[i])ans+=mu[i]*solve(n/i/i),ans%=mod;
	n%=mod;
	ans=((n*n-ans)%mod+mod)%mod;
	cout<<ans<<"\n";
	return 0;
}