LOJ #2159. 「POI2011 R1」Plot

好难写啊！

这题如果保证数据随机，那么可以直接跑一个最小圆覆盖，先二分半径，再贪心覆盖。
但是出题人显然不会这么善良。
于是就可以倍增，\(1,2,4,8,16...\)，这样尝试长度，找到最大可行二进制长度（即最高位）后，再逐位确定。
复杂度\(O(nlog^2(n))\)
但是写完之后又被卡了精度，改随机数种子才可以过。

#include <bits/stdc++.h>
using namespace std;
const int N=100010;
typedef double ld;
const ld EPS=1e-8;
struct point{
    ld x,y;
    ld distance(){
	return sqrt(x*x+y*y);
    }
    point operator *(const ld &z) const{
	return {z*x,z*y};
    }
    point operator /(const ld &z) const{
	return {x/z,y/z};
    }
    ld operator ^(const point &_) const{
	return x*_.x+y*_.y;
    }
    ld operator *(const point &_) const{
	return x*_.y-y*_.x;
    }
    point operator +(const point &_) const{
	return {x+_.x,y+_.y};
    }
    ld sqr() const{
	return x*x+y*y;
    }
    point operator -(const point &_) const{
	return {x-_.x,y-_.y};
    }
    point rotate(const ld &alpha) const{
	ld tc=cos(alpha),ts=sin(alpha);
	return {x*tc-y*ts,x*ts+y*tc};
    }
    point normal() const{
	return {-y,x};
    }
}a[N],c[N];
struct line{
    point x,y;
    point generate(const ld &c) const{
	return x+y*c;
    }
    point cross(const line &u) const{
	ld s1=y*u.y;
	ld s2=u.y*(x-u.x);
	//cerr<<s1<<" "<<s2<<" "<<u.y.x<<" "<<u.y.y<<endl;
	//cerr<<x.x<<" "<<x.y<<" "<<y.x<<" "<<y.y<<endl;
	return generate(s2/s1);
    }
};
struct cir{
    point o;
    ld r;
    bool in(const point &d) const{
	return (o-d).sqr()<=r;
    }
};
bool rec,re;
int n,m;
int num;
point out[N];
cir solve(const point &x,const point &y,const point &z){
    //line l1{(x+y)/2,(x-y).rotate(M_PI_2)};
    //line l2{(y+z)/2,(z-y).rotate(M_PI_2)};
    line l1{(x+y)/2,(x-y).normal()};
    line l2{(y+z)/2,(z-y).normal()};
    if (abs(l1.y*l2.y)<EPS){
	ld len=max(max((x-y).sqr(),(x-z).sqr()),(y-z).sqr());
	if (len==(x-y).sqr()){
	    return {(x+y)/2,sqrt(len)/2};
	}
	if (len==(x-z).sqr()){
	    return {(x+z)/2,sqrt(len)/2};
	}
	return {(y+z)/2,len/2};
    }
    point tmp=l1.cross(l2);
    return {tmp,(tmp-x).sqr()};
}
ld check(int l,int r){
    if (r-l+1==1){
	if (re) out[++num]=a[l];
	return 0;
    }
    for (int i=l; i<=r; ++i) c[i]=a[i];
    random_shuffle(c+l,c+r+1);
    cir tmp={c[l],0};
    for (int i=l+1; i<=r; ++i){
	if (tmp.in(c[i])) continue;
	tmp={c[i],0};
	for (int j=l; j<i; ++j)
	    if (!tmp.in(c[j])){
		tmp={(c[i]+c[j])/2,(c[i]-c[j]).sqr()/4};
		for (int k=l; k<j; ++k)
		    if (!tmp.in(c[k])){
			//cerr<<i<<" "<<j<<" "<<k<<endl;
			//cerr<<c[i].x<<" "<<c[i].y<<endl;
			//cerr<<c[j].x<<" "<<c[j].y<<endl;
			//cerr<<c[k].x<<" "<<c[k].y<<" "<<l<<" "<<r<<endl;
			tmp=solve(c[i],c[j],c[k]);
		    }
	    }
    }
    if (re){
	out[++num]=tmp.o;
    }
    return tmp.r;
}
int getnext(const int &x,const ld &mid){
    int i;
    for (i=2; x+i-1<=n; i<<=1) if (check(x,x+i-1)>mid*mid) break;
    for (int j=(i>>=1); j; j>>=1)
	if (x+i+j-1<=n&&check(x,x+i+j-1)<=mid*mid) i+=j;
    if (rec){
	re=1;
	check(x,x+i-1);
	re=0;
    }
    return x+i;
}
bool maybe(const ld &mid){
    //cerr<<"maybe:"<<mid<<endl;
    int rest=m;
    for (int i=1; i<=n; i=getnext(i,mid))
	if (--rest<0) return 0;
    //cerr<<"SUCC"<<endl;
    return 1;
}
int main(){
    srand(19260817);
    scanf("%d%d",&n,&m);
    for (int i=1; i<=n; ++i) scanf("%lf%lf",&a[i].x,&a[i].y);
    ld ret=-1;
    for (ld l=0,r=3e6,mid=(l+r)/2; l+EPS<r; mid=(l+r)/2)
	if (maybe(mid)) r=ret=mid; else l=mid;
    cout<<fixed<<setprecision(10);
    cout<<ret<<endl;
    rec=1;
    maybe(ret);
    cout<<num<<endl;
    for (int i=1; i<=num; ++i) cout<<out[i].x<<" "<<out[i].y<<'\n';
}