JZOJ6383. 【NOIP2019模拟2019.10.07】果实摘取

Description

• 小 D 的家门口有一片果树林，果树上果实成熟了，小 D 想要摘下它们。
  为了便于描述问题，我们假设小 D 的家在二维平面上的 (0, 0) 点，所有坐标范围的绝对值不超过 N 的整点坐标上都种着一棵果树。（(0, 0) 这个点没有果树）
• 小 D 先站在 (0, 0) 处，正对着 (1, 0) 的方向。
• 每次摘果实时，小 D 会逆时针选择他能看到的第 K 棵还未摘取果实的果树，然后向着这个方向走去，在行走的过程中摘下沿路的所有的果树上的果树果实，直到走到果树林的边缘。
• 接下来，小 D 回到 (0, 0) 处，正对着上一次摘果实的果树的方向。
• 小 D 会重复这个过程，直到所有的果实都被摘取，小 D 感兴趣的是，最后一棵被摘下果实的果树是哪一棵？
• 注意小 D 不能看到被任何其他果树遮挡着的果树。
• n,K<=1e5

Solution

• 考虑一共有多少条果树，也就是
  $(sum_{i=1}^{n}phi (i)+1)*8$
• 接下来每一次跳K就是一个经典的约瑟夫问题了。
• 有一种$O(n)$的递推做法
• $f[i]=(f[i-1]+K) mod i ,f[1]=0$
• 意义就是考虑长度为i的约瑟夫问题，钦定从0开始，选择掉了K-1，剩下的是一个i-1的子问题，集合是${k,k+1...n-1,0,1...,k-2}$，这个子问题的编号加上k就是当前的编号了。
• 由于n（此n非读入的n）非常大，发现如果多次加K都没有取模得话，我们可以将这些压在一起做，推一推算一算就变成log的了。
• 那么问题转化成为求斜率第x个的互质的坐标。
• 先可以给x对于一个象限里的个数去一个模。
• 然后二分一个小数表示最大斜率mid，暴力枚举横坐标，计算纵坐标满足在斜率内的互质点的个数。
• 即对于每一个i，求出：
  $sum_{j=1}^{i*mid}[gcd(i,j)=1]$
• 容(fan)斥(yan)一下就可以得到
  $sum_{j=1}^{i*mid} sum_{d|gcd(i,j)} mu[d]$
• 将d提前
  $sum_{d|i}mu[d]*left lfloor left lfloor i*mid ight floor /d ight floor$
• 对于每一i枚举它的因子，总复杂度即为n/1+n/2+n/3+n/4…
• 即$O(n logn )$
• 加上二分$O(n log ^2 n )$
• 很短很好打

#include<cstdio>
#include<cmath>
#include<cstring>
#include<algorithm>
#define maxn 100005
#define db double 
#define E 3e-10
#define min(a,b) ((a<b)?a:b)
#define ll long long 
using namespace std;

int n,m,i,j,tp,xx,yy;
int tot,pri[maxn],bz[maxn],phi[maxn],u[maxn];
ll sum,sum0,num,f,x,d,k;
db l,r,mid;

int gcd(int x,int y){return (x%y==0)?y:gcd(y,x%y);}

void Getphi(){
	u[1]=1;
	for(i=2;i<=n;i++) {
		if (!bz[i]) bz[i]=1,u[i]=-1,phi[i]=i-1,pri[++tot]=i;
		for(j=1;j<=tot&&i*pri[j]<=n;j++){
			bz[i*pri[j]]=1;
			if (i%pri[j]==0){
				phi[i*pri[j]]=phi[i]*pri[j];
				u[i*pri[j]]=0;
				break;
			}  else {
				phi[i*pri[j]]=(pri[j]-1)*phi[i];
				u[i*pri[j]]=-u[i];
			}
		}
	}
}

int main(){
	scanf("%d%d",&n,&m);
	Getphi();
	for(i=1;i<=n;i++) sum+=phi[i];
	sum++,sum0=sum*8;
	f=0,x=1;
	while (x<sum0){
		d=min(sum0-x,(x-f)/(m-1)+((x-f)%(m-1)>0));
		f=(f+d*m)%(x+d),x+=d;
	} 
	tp=f/(sum*2),num=f%(sum*2);
	if (num==0) xx=n,yy=0; else {
		l=0,r=n;
		while (abs(r-l)>E){
			 mid=(l+r)/2;
			 sum=0;
			 for(d=1;d<=n;d++) for(i=d;i<=n;i+=d){
			 	k=i*mid; k=min(k,n);
				sum+=k/d*u[d];
			 }
			 if (sum<num) l=mid; else 
			 if (sum>num) r=mid; else {
			 	xx=n,yy=0;
				for(i=1;i<=n;i++) {
					k=i*mid; k=min(k,n);
					if (1.0*k/i>1.0*yy/xx) 
						xx=i,yy=k;
				}
				break;
			 }
		}
	}
	while (tp--) swap(xx,yy),xx=-xx;
	if (xx&&yy){
		k=min(n/abs(xx),n/abs(yy));
		xx*=k,yy*=k;
	}
	printf("%d %d",xx,yy);
}