BZOJ3564: [SHOI2014]信号增幅仪（最小圆覆盖）

传送门

题意：
求一个固定半长轴方向和半长轴半短轴比值的最小椭圆覆盖平面上所有点。

题解：
最小圆覆盖。
首先将坐标旋转，再把所有横坐标除以比值得到最小圆覆盖的经典模板。

最小圆覆盖的做法：
首先，在已经确定最小圆的点集中加入一个在最小圆外面的点，那么这个点一定在新点集的最小圆上。

算法：
初始化圆形为$p_1$，半径为$0$
1.从左到右依次加入点，若在当前圆内则跳过，否则问题转化为在一个已经最小圆确定的点集中加入另外一个在圆外的点，可以肯定这个点在最小圆上。
2.确定$1～i$的最小圆上的一个点之后，依次枚举$1～i-1$的点，分别确定$\{1～j\}\cup i$的最小圆。对于一个点$j$，如果在$\{1～j-1\}\cup i$的最小圆中则跳过，否则确定了$\{1～j\}\cup i$的第二个点。
3.同样，要确定$\{1～j\}\cup i$的第三个点，那么要确定$\{1～j-1\}\cup j \cup i$的最小圆，同2对于每个j枚举k即可。

时间复杂度看起来像$n^3$,实际上随机点顺序之后是$O(n)$的。

#include<bits/stdc++.h>
using namespace std;
inline int read(){
    char ch=getchar();int i=0,f=1;
    while(!isdigit(ch)){if(ch=='-')f=-1;ch=getchar();}
    while(isdigit(ch)){i=(i<<1)+(i<<3)+ch-'0';ch=getchar();}
    return i*f;
}
const double PI=acos(-1.0);
const int Maxn=5e5+50;
int n;
double R;
struct point{
    double x,y;
    point(double x=0,double y=0):x(x),y(y){}
    friend inline point operator +(const point &a,const point &b){return point(a.x+b.x,a.y+b.y);}
    friend inline point operator -(const point &a,const point &b){return point(a.x-b.x,a.y-b.y);}
    friend inline double operator *(const point &a,const point &b){return a.x*b.y-a.y*b.x;}
    friend inline double dis(const point &a,const point &b){return (a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y);}
    friend inline point mid(const point &a,const point &b){return point((a.x+b.x)/2.0,(a.y+b.y)/2.0);}
    inline point normalvector(){return point(-y,x);}
}p[Maxn],C;
struct line{
    point a,b;
    line(point a,point b):a(a),b(b){}
    friend inline point inter(const line &a,const line &b){
        double k1=(a.b-b.b)*(a.a-b.b);
        double k2=(a.a-b.a)*(a.b-b.a);
        #define t k1/(k1+k2)
        return point(b.b.x+(b.a.x-b.b.x)*t,b.b.y+(b.a.y-b.b.y)*t);
        #undef t
    }
};
inline void calc(point &C,double &R,const point &a,const point &b,const point &c){
    C=inter(line(mid(a,b),mid(a,b)+(b-a).normalvector()),line(mid(b,c),mid(b,c)+(c-b).normalvector()));
    R=dis(C,a);
}
int main(){
    srand(19491001);n=read();
    for(int i=1;i<=n;i++){
        p[i].x=read(),p[i].y=read();
    }
    double theta=360.0-read();(theta/=180.0)*=PI;double P=read();
    for(int i=1;i<=n;i++){
        double x=p[i].x,y=p[i].y;
        p[i].x=x*cos(theta)-y*sin(theta);
        p[i].y=x*sin(theta)+y*cos(theta);
        p[i].x/=P;
    }
    random_shuffle(p+1,p+n+1);
    C=p[1],R=0;
    for(int i=2;i<=n;i++){
        if(dis(p[i],C)>R){
            C=p[i],R=0;
            for(int j=1;j<i;j++){
                if(dis(p[j],C)>R){
                    R=dis(p[j],p[i])/4.0,C=mid(p[j],p[i]);
                    for(int k=1;k<j;k++){
                        if(dis(p[k],C)>R)
                            calc(C,R,p[i],p[j],p[k]);
                    }
                }
            }
        }
    }
    printf("%.3f",sqrt(R));
}