[学习笔记]最小圆覆盖

随机增量算法

图片来自：我

1.random_shuffle

2.枚举增量：点i

圆 C;
for(i=1 to n)
{
    if(P[i] 不在 C 内)
    {
        C = {P[i], 0};
        for(j=1 to i-1)
        {
            if(P[j] 不在 C 内)
            {
                C = {0.5*(P[i]+P[j]), 0.5*dist(P[i], P[j])};
                for(k=1 to j-1)
                {
                    if(P[k] 不在 C 内)
                    {
                        C = 外接圆(P[i], P[j], P[k]);
                    }
                }
            }
        }
    }
}

求外接圆：

 

求出A和B,B和C的中点p1,p2,再把两个向量B-A,C-B旋转90度（顺逆无所谓）

I=p1+t*v1

(p1+t*v1-p2)×v2=0（×表示叉乘）

所以：t=(p2-p1)×v2/(v2×v1)（叉乘乘法分配率）

I=p1+t*v1

注意回来更新半径r的值

 

复杂度：

random_shuffle，所以点的分布比较均匀

每个循环有n/3概率满足是外接圆3点之一

 

#include<bits/stdc++.h>
#define il inline
#define reg register int
#define numb (ch^'0')
using namespace std;
typedef long long ll;
il void rd(int &x){
    char ch;bool fl=false;
    while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true);
    for(x=numb;isdigit(ch=getchar());x=x*10+numb);
    (fl==true)&&(x=-x);
}
namespace Miracle{
const int N=100000+5;
int n;
struct po{
    double x,y;
    po(){x=0;y=0;}
    po(double xx,double yy){
        x=xx;y=yy;
    }
    po friend operator +(po a,po b){
        return po(a.x+b.x,a.y+b.y);
    }
    po friend operator /(po a,double c){
        return po(a.x/c,a.y/c);
    }
    po friend operator -(po a,po b){
        return po(a.x-b.x,a.y-b.y);
    }
    double friend operator *(po a,po b){
        return a.x*b.y-a.y*b.x;
    }
    double len(){
        return x*x+y*y;
    }
}p[N];
po rev(po a){
    swap(a.x,a.y);a.y=-a.y;return a;
}
po circle(const po &A,const po &B,const po &C){
    po lmd=(A+B)/2,rmd=(B+C)/2;
    po vel=rev(B-A),ver=rev(C-B);
    cout<<" mom "<<ver*vel<<endl;
    double t=ver*(rmd-lmd)/(ver*vel);
    vel.x*=t;vel.y*=t;
    return lmd+vel;
}
int main(){
    srand(20011023);
    int n;
    rd(n);
    for(reg i=1;i<=n;++i){
        scanf("%lf%lf",&p[i].x,&p[i].y);
    }
    random_shuffle(p+1,p+n+1);
    po o;
    double r=0.0;
    for(reg i=1;i<=n;++i){
        if((p[i]-o).len()>r){
            o=p[i];r=0;
            for(reg j=1;j<i;++j){
                if((p[j]-o).len()>r){
                    o=(p[j]+p[i])/2;r=(p[j]-o).len();
                    for(reg k=1;k<j;++k){
                        if((p[k]-o).len()>r){
                            o=circle(p[i],p[j],p[k]);
                            r=(p[i]-o).len();
                        }
                    }
                }
            }
        }
//        cout<<" i "<<i<<" : "<<r<<" x "<<o.x<<" y "<<o.y<<endl;
    }
    printf("%.10lf\n%.10lf %.10lf",sqrt(r),o.x,o.y);
    return 0;
}

}
signed main(){
    Miracle::main();
    return 0;
}

/*
   Author: *Miracle*
   Date: 2019/4/30 20:32:29
*/

 单点增量的思路，随机保证“暴力”复杂度。