P8058 [BalkanOI2003] Farey 序列 题解

分析

考虑二分答案。

对于当前二分的答案 $x$，设 $cnt$ 表示 Farey 序列中 $\frac{p}{q}\le x$ 的满足条件的数量。对于一组 $(i,j)$，若 $\frac{j}{i}\le x$，则 $j\le \lfloor i\times x\rfloor$。得到暴力式子：

$$cnt=\sum _{i=1}^n\sum _{j=1}^{\lfloor i\times x\rfloor }[gcd(i,j)=1]$$

定义 $f_i=\sum _{j=1}^{\lfloor i\times x\rfloor }[gcd(i,j)=1]$。根据容斥，有：

$$f_i=\lfloor i\times x\rfloor -\sum _{j=1}^{\lfloor i\times x\rfloor }[gcd(i,j)>1]$$

将 $gcd(i,j)$ 的值表示出来：

$$f_i=\lfloor i\times x\rfloor -\sum _{j=1}^{\lfloor i\times x\rfloor }[gcd(i,j)=k\wedge k>1]$$

把 $k$ 提出来：

$$f_i=\lfloor i\times x\rfloor -\sum _{k=2}^i\sum _{j=1}^{\lfloor \lfloor \frac{i}{k}\rfloor \times x\rfloor }[gcd(i,j)=1]$$

不难发现，第二个求和等价于 $f_{\lfloor \frac{i}{k}\rfloor }$，则有：

$$f_i=\lfloor i\times x\rfloor -\sum _{k=1}^i{f}_{\lfloor \frac{i}{k}\rfloor }$$

化简得：

$$f_i=\lfloor i\times x\rfloor -\sum _{k|i\wedge 1<k<i}^{}{f}_k$$

$cnt=\sum _{i=1}^n{f}_i$。枚举 $i$ 的约数，则能够在 $O(n\sqrt{n})$ 的复杂度求出 $cnt$，然后复杂度就是 $O(n\sqrt{n}\log n)$。

对于求 $p,q$ 的值，暴力枚举即可。

代码

#include<bits/stdc++.h>
using namespace std;
//#define int long long
#define re register
#define il inline
#define pii pair<int,int>
#define x first
#define y second
#define gc getchar()
#define rd read()
#define debug() puts("------------")

namespace yzqwq{
	il int read(){
		int x=0,f=1;char ch=gc;
		while(ch<'0'||ch>'9'){if(ch=='-') f=-1;ch=gc;}
		while(ch>='0'&&ch<='9') x=(x<<1)+(x<<3)+(ch^48),ch=gc;
		return x*f;
	}
	il int qmi(int a,int b,int p){
		int ans=1;
		while(b){
			if(b&1) ans=ans*a%p;
			a=a*a%p,b>>=1;
		}
		return ans;
	}
	il auto max(auto a,auto b){return (a>b?a:b);}
	il auto min(auto a,auto b){return (a<b?a:b);}
	il int gcd(int a,int b){
		if(!b) return a;
		return gcd(b,a%b);
	}
	il int lcm(int a,int b){
		return a/gcd(a,b)*b;
	}
	il void exgcd(int a,int b,int &x,int &y){
		if(!b) return x=1,y=0,void(0);
		exgcd(b,a%b,x,y);
		int t=x;
		x=y,y=t-a/b*x;
		return ;
	}
	mt19937 rnd(time(0));
}
using namespace yzqwq;

const int N=4e4+10;
double eps=1e-10;
int n,k,f[N];

il bool check(double x){
	int cnt=0;
	for(re int i=1;i<=n;++i){
		f[i]=i*x;
		for(re int j=2;j*j<=i;++j){
			if(i%j==0){
				f[i]-=f[j];
				if(i/j!=j) f[i]-=f[i/j];
			}
		}
		cnt+=f[i];
	}
	return cnt>=k;
}

il void solve(){
	n=rd,k=rd;
	double l=0.0,r=1.0,ans=0.0;
	while((r-l)>eps){
		double mid=(l+r)/2.0;
		if(check(mid)) ans=mid,r=mid;
		else l=mid;
	}
	int p=0,q=0;
	double c=40001.0;
	for(re int Q=1;Q<=40001;++Q){
		int P=Q*ans;
		if(fabs(double(P*1.0/Q)-ans)<c) c=fabs(double(P*1.0/Q)-ans),p=P,q=Q;
	}
	printf("%lld %lld\n",p,q);
	return ;
}

signed main(){
//	freopen(".in","r",stdin);
//	freopen(".out","w",stdout);
	int t=1;while(t--)
	solve();
	return 0;
}