[BZOJ1337/1336/2823] 最小圆覆盖

Description

给出N个点,让你画一个最小的包含所有点的圆。

Input

先给出点的个数N,2<=N<=100000,再给出坐标Xi,Yi.(-10000.0<=xi,yi<=10000.0)

Output

输出圆的半径,及圆心的坐标

Sample Input

6
8.0 9.0
4.0 7.5
1.0 2.0
5.1 8.7
9.0 2.0
4.5 1.0

Sample Output

5.00
5.00 5.00

题目分析

随机增量法太神辣 (其实刷只是为了PK+三倍经验)

以下内容转自commonc神犇的博客


算法目的:最小圆覆盖算法可以在线性时间复杂度内求出覆盖n个点的最小圆

算法步骤:
     ①首先现将所有点随机排列

     ②按顺序把点一个一个的加入(一步一步的求前i个点的最小覆盖圆),每加入一个点就进入③

     ③如果发现当前i号点在当前的最小圆的外面,那么说明点i一定在前i个点的最小覆盖圆边界上,我们转到④来进一步确定这个圆,否则前i个点的最小覆盖圆与前i-1个点的最小覆盖圆是一样的,则不需要更新,返回②

     ④此时已经确认点i一定在前i个点的最小覆盖圆的边界上了,那么我们可以把当前圆的圆心设为第i个点,半径为0,然后重新把前i-1个点加入这个圆中(类似上面的步骤,只不过这次我们提前确定了点i在圆上,目的是一步一步求出包含点i的前j个点的最小覆盖圆),每加入一个点就进入⑤

     ⑤如果发现当前j号点在当前的最小圆的外面,那么说明点j也一定在前j个点(包括i)的最小覆盖圆边界上,我们转到⑥来再进一步确定这个圆,否则前j个点(包括i)的最小覆盖圆与前i-1个点(包括i)的最小覆盖圆是一样的,则不需要更新,返回④

     ⑥此时已经确认点i,j一定在前j个点(包括i)的最小覆盖圆的边界上了,那么我们可以把当前圆的圆心设为第i个点与第j的点连线的中点,半径为到这两个点的距离(就是找一个覆盖这两个点的最小圆),然后重新把前j-1个点加入这个圆中(还是类似上面的步骤,只不过这次我们提前确定了两个点在圆上,目的是求出包含点i,j的前k个点的最小覆盖圆),每加入一个点就进入⑦

     ⑦如果发现当前k号点在当前的最小圆的外面,那么说明点k也一定在前k个点(包括i,j)的最小覆盖圆边界上,我们不需要再进一步确定这个圆了(因为三个点能确定一个圆!),直接求出这三点共圆,否则前k个点(包括i,j)的最小覆盖圆与前k-1个点(包括i,j)的最小覆盖圆是一样的,则不需要更新。

时间复杂度:O(N)
空间复杂度:O(N)


最后 最关键的部分 如何求出三点共圆的圆心坐标?
公式:
假设三点坐标 圆心坐标
则:

#include <cstdio>
#include <cstring>
#include <set>
#include <map>
#include <vector>
#include <cmath>
#include <queue>
#include <algorithm>
using namespace std;
int n;
struct your
{
    double x,y;
}a[500010];
double R;
your O;
double sqr(double x) { return x*x; }
double dis(your x,your y)
{
    return sqrt(sqr(x.x-y.x)+sqr(x.y-y.y));
}
int incir(your x)
{
    if(dis(O,x)<=R) return 1;
    return 0;
}
your slove(your x,your y,your z)
{
    your tmp;
    tmp.x=(sqr(z.x)+sqr(z.y)-sqr(x.x)-sqr(x.y))*(x.y-y.y)-(sqr(y.x)+sqr(y.y)-sqr(x.x)-sqr(x.y))*(x.y-z.y);
    tmp.x/=2;
    tmp.x/=(x.x-y.x)*(x.y-z.y)-(x.y-y.y)*(x.x-z.x);
    tmp.y=(sqr(z.x)+sqr(z.y)-sqr(x.x)-sqr(x.y))*(x.x-y.x)-(sqr(y.x)+sqr(y.y)-sqr(x.x)-sqr(x.y))*(x.x-z.x);
    tmp.y/=2;
    tmp.y/=(x.y-y.y)*(x.x-z.x)-(x.x-y.x)*(x.y-z.y);
    return tmp;
}
int main()
{   
    scanf("%d",&n);
    for(int i=1;i<=n;i++) scanf("%lf%lf",&a[i].x,&a[i].y);
    random_shuffle(a+1,a+n+1);
    for(int i=1;i<=n;i++)
        if(!incir(a[i]))
        {
            O.x=a[i].x,O.y=a[i].y,R=0;
            for(int j=1;j<i;j++)
                if(!incir(a[j]))
                {
                    O.x=(a[i].x+a[j].x)/2;
                    O.y=(a[i].y+a[j].y)/2;
                    R=dis(O,a[i]);
                    for(int k=1;k<j;k++)
                        if(!incir(a[k]))
                        {
                            O=slove(a[i],a[j],a[k]);
                            R=dis(a[i],O);
                        }
                }
        }
    printf("%.10lf\n%.10lf %.10lf",R,O.x,O.y);
    return 0;
}

发表评论

邮箱地址不会被公开。 必填项已用*标注