Bzoj2280 [POI2011] WYK-Plot

, , , ,

Description

https://www.luogu.com.cn/problem/P3517

Solution

倍增典题，但是我不会。

一个显然的想法是二分以后，对于每个左端点再二分一下，里面做最小圆覆盖，找到这个左端点的最右的右端点。

但考虑一种所有只能分成长度为 $1$ 的段的情况，此时复杂度是 $O(n^2\log n)$。

然后我就不会了啊。发现倍增可以和二分很好的结合在一起。

考虑先倍增出一个区间长度 $2^k$，表示这个左端点的区间长度至少为 $2^k$，最多为 $2^{k+1}-1$。然后在这里面二分。

由于倍增的总长度是 $O(n)$ 的，而二分的总长度与倍增的总长度同阶。所以复杂度 $O(n\log^2 n)$。然而我被卡精度了，放一下 loj 97 的代码。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
#include<bits/stdc++.h>
#define ll long long
#define double long double
#define ull unsigned long long
#define maxn 100005
#define put() putchar('\n')
#define Tp template<typename T>
#define Ts template<typename T,typename... Ar>
using namespace std;
Tp void read(T &x){
    int f=1;x=0;char c=getchar();
    while (c<'0'||c>'9') {if (c=='-') f=-1;c=getchar();}
    while (c>='0'&&c<='9') {x=x*10+c-'0';c=getchar();}
    x*=f;
}
namespace Debug{
	Tp void _debug(char* f,T t){cerr<<f<<'='<<t<<endl;}
	Ts void _debug(char* f,T x,Ar... y){while(*f!=',') cerr<<*f++;cerr<<'='<<x<<",";_debug(f+1,y...);}
	#define gdb(...) _debug((char*)#__VA_ARGS__,__VA_ARGS__)
}using namespace Debug;
#define fi first
#define se second
#define mk make_pair
const int mod=1e9+7;
int power(int x,int y=mod-2) {
	int sum=1;
	while (y) {
		if (y&1) sum=sum*x%mod;
		x=x*x%mod;y>>=1;
	}
	return sum;
}
int n,m;
mt19937 rnd(181766);
struct node {
	double x,y;
}a[maxn],C,b[maxn],g[maxn];
double r;
const double eps=1e-9;
double dis(node x,node y) {return sqrt((x.x-y.x)*(x.x-y.x)+(x.y-y.y)*(x.y-y.y));}
void getcircle(node x,node y,node z) {
	double a1=2*(z.x-x.x),a2=2*(z.x-y.x),b1=2*(z.y-x.y),b2=2*(z.y-y.y);
	if (abs(a1*b2-a2*b1)<=eps) {
		C.x=(min({x.x,y.x,z.x})+max({x.x,y.x,z.x}))/2;
		C.y=(min({x.y,y.y,z.y})+max({x.y,y.y,z.y}))/2;
		r=max({dis(C,x),dis(C,y),dis(C,z)});
		// assert(r<=dis(C,x));
		// assert(r<=dis(C,y));
		// assert(r<=dis(C,z));
		// assert(0);
		return ;
	}
	double c1=z.x*z.x+z.y*z.y-x.x*x.x-x.y*x.y,c2=z.x*z.x+z.y*z.y-y.x*y.x-y.y*y.y;
	C.y=(a1*c2-a2*c1)/(a1*b2-a2*b1);
	C.x=(b1*c2-b2*c1)/(b1*a2-b2*a1);
	r=dis(C,x);
	// assert((r-dis(C,x))<=1e-5);
	// assert((r-dis(C,y))<=1e-5);
	// assert((r-dis(C,z))<=1e-5);
}
int times=0;
double solve(int L,int R) {
	int i,j,k;
	for (i=L;i<=R;i++) a[i-L+1]=b[i];
	int n=R-L+1;
	shuffle(a+1,a+1+n,rnd);
	C.x=C.y=r=0;
	for (i=1;i<=n;i++) if (dis(a[i],C)>r) {
		C=a[i],r=0;
		for (j=1;j<i;j++) if (dis(a[j],C)>r) {	
			C.x=(a[i].x+a[j].x)/2;
			C.y=(a[i].y+a[j].y)/2;
			r=dis(a[i],C);
			for (k=1;k<j;k++) if (dis(a[k],C)>r) {
				getcircle(a[i],a[j],a[k]);
				// ++times;
			}
		}
	}
	
	return r;
}
int nex[maxn],K;
bool check(double val) {
	int sum=0,i;times=0;
	for (i=1;i<=n;i++) {
		int k;
		for (k=0;i+(1<<k+1)-1<=n;k++) {
			if (solve(i,i+(1<<k+1)-1)<=val) ;
			else break;
		}
		int L=i+(1<<k)-1,R=min(n,i+(1<<k+1))+1;
		while (L+1<R) {
			int mid=L+R>>1;
			if (solve(i,mid)<=val) L=mid;
			else R=mid;
		}
		// gdb(i,L);
		solve(i,L);
		nex[i]=L;g[i]=C;
		sum++;
		// if (sum>m) return 0;
		i=L;
	}
	// gdb(times);
	K=sum;
	return sum<=m;
}
signed main(void){
	int i,j,k;
	read(n);read(m);
	for (i=1;i<=n;i++) scanf("%Lf%Lf",&b[i].x,&b[i].y);
	double l=0,r=2e6,mid;
	// solve(1,4);gdb(C.x,C.y,::r);
	// printf("%d\n",check(9.8));return 0;
	while (l+eps<r) {
		mid=(l+r)/2;
		if (check(mid)) r=mid;
		else l=mid;  
		// gdb(l,r);
	}
	check(r);
	printf("%.12Lf\n",r);
	printf("%d\n",K);
	for (i=1;i<=n;i++) {
		printf("%.12Lf %.12Lf\n",g[i].x,g[i].y);
		i=nex[i];
	}
	return 0;
}
//i=begin && g++ $i.cpp -o $i -std=c++14 && ./$i