最小圆覆盖

在正式介绍该算法前先来说几句废话（其实有联系滴^-^）

三角关系的探索：

边和角的关系：

正弦定理：

  盗图说明：

有：

证明：

余弦定理：

逐渐进入正题了

：

围绕三点做最小的圆，使得三个点全部在圆的内部或者圆上（最小圆覆盖）。这里存在两种情况：1. 三点全部都在圆上；2. 三个点如果有2个点在圆上，另一个点在圆的内部，那么那两个点一定是直径的两个端点。即：寻找对应三角形内最长的线段--->寻找最长的边--->由正弦定理可知，如果最大的角是钝角，即是第2种情况；如果最大的角是直角，那两个点一定是直径的两个端点，另一个点就在圆上，属于第一种情况；如果最大的角是锐角，那么3点均在圆上，没有边是直径，也属于第1种情况。

已知3个点求它们围成的三角形的面积时除了用叉积，还可以直接用行列式的方法：

对于：

假设圆心在三角形中（及时不设在三角形内部也能由向量的计算推出同样的效果。）那么，行列式计算的结果

  即是3个叉积，加起来就是对应的整个多边形的面积，所以最后要除以2。写成代码：

double area(){
    return fabs((p[0].x*p[1].y+p[2].x*p[0].y+p[1].x*p[2].y-p[2].x*p[1].y-p[1].x*p[0].y-p[0].x*p[2].y)/2.0);
}

find mincircle :三角形两条中垂线：

两式相减处理可以分别求出x和y，

设

那么就有

Y的求解类似上诉过程:

写成程序：

#include <iostream>
#include <cmath>
using namespace std;
struct point{
    double x,y;
};
point p0,p1,p2;
double area(point p0,point p1,point p2){
    return fabs((p0.x*p1.y+p2.x*p0.y+p1.x*p2.y-p2.x*p1.y-p1.x*p0.y-p0.x*p2.y)/2.0);
}
double dis(point p1,point p2){
    return sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y));
}
struct circle{
    point p;
    double r;
};
circle get(point p0,point p1,point p2){
    circle cr;
    double s=area(p0,p1,p2);
    double a=dis(p0,p1), b=dis(p1,p2),c=dis(p0,p2);
    cr.r=a*b*c/4/s;
    double c1=(p0.x*p0.x-p1.x*p1.x+p0.y*p0.y-p1.y*p1.y)/2;
    double c2=(p2.x*p2.x-p1.x*p1.x+p2.y*p2.y-p1.y*p1.y)/2;
    double x21=p2.x-p1.x,  x01=p0.x-p1.x;
    double y21=p2.y-p1.y,  y01=p0.y-p1.y;
    cr.p.x=(c1*y21-c2*y01)/(x01*y21-x21*y01);
    cr.p.y=(c1*x21-c2*x01)/(y01*x21-y21*x01);
    return cr;
}
int main()
{
    p0.x=-3;  p0.y=0;
    p1.x=3;   p1.y=0;
    p2.x=0;   p2.y=4;
    circle ans=get(p0,p1,p2);
    cout<<ans.p.x<<" "<<ans.p.y<<" "<<ans.r<<endl;
    return 0;
}
/*
0 3
0 0
4 0
2 1.5 2.5
-3 0
3 0
0 4
-0 0.875 3.125
Process returned 0 (0x0)   execution time : 0.500 s
Press any key to continue.
*/

寻找多个点的最小覆盖圆可以用上面讨论过的两种情况迭代来写。

献上几道小菜：

hdu 3932 Groundhog Build Home

​​http://acm.hdu.edu.cn/showproblem.php?pid=3932​​

#include <cstdio>
#include <cmath>
#include <iostream>
#include <algorithm>
using namespace std;
const double eps=1e-7;
int x,y,n;
struct point{
    double x,y;
}p[1005];
double area(point p0,point p1,point p2){
    return fabs((p0.x*p1.y+p2.x*p0.y+p1.x*p2.y-p2.x*p1.y-p1.x*p0.y-p0.x*p2.y)/2.0);
}
double dis(point p1,point p2){
    return sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y));
}
point get(point p0,point p1,point p2){
    point p;
    double c1=(p0.x*p0.x-p1.x*p1.x+p0.y*p0.y-p1.y*p1.y)/2;
    double c2=(p2.x*p2.x-p1.x*p1.x+p2.y*p2.y-p1.y*p1.y)/2;
    double x21=p2.x-p1.x,  x01=p0.x-p1.x;
    double y21=p2.y-p1.y,  y01=p0.y-p1.y;
    p.x=(c1*y21-c2*y01)/(x01*y21-x21*y01);
    p.y=(c1*x21-c2*x01)/(y01*x21-y21*x01);
    return p;
}
void min_cover(point p[],point &c,double &r){
    random_shuffle(p,p+n);
    c=p[0];
    r=0;
    for(int i=1;i<n;i++)
        if(dis(c,p[i])+eps>r){
            c=p[i];
            r=0;
            for(int j=0;j<i;j++){
                if(dis(c,p[j])+eps>r){
                    c.x=(p[i].x+p[j].x)/2;
                    c.y=(p[i].y+p[j].y)/2;
                    r=dis(c,p[j]);;
                    for(int k=0;k<j;k++){
                        if(dis(c,p[k])+eps>r){
                            c=get(p[i],p[j],p[k]);
                            r=dis(c,p[k]);
                        }
                    }
                }
            }
        }
}
int main(){
    //freopen("cin.txt","r",stdin);
    while(cin>>x>>y>>n){
        for(int i=0;i<n;i++){
            scanf("%lf%lf",&p[i].x,&p[i].y);
        }
        point cen;
        double r;
        min_cover(p,cen,r);
        printf("(%.1lf,%.1lf).\n%.1lf\n",cen.x, cen.y, r);
    }
    return 0;
}

这种求解外心的方法在其他博客看见的，我没看懂它的原理（||-_-）：

point get(point a,point b,point c){
    double a1 = b.x - a.x, b1 = b.y - a.y, c1 = (a1*a1 + b1*b1)/2;
    double a2 = c.x - a.x, b2 = c.y - a.y, c2 = (a2*a2 + b2*b2)/2;
    double d = a1 * b2 - a2 * b1;
    point p;
    p.x=a.x + (c1*b2 - c2*b1)/d;
    p.y=a.y + (a1*c2 - a2*c1)/d;
    return p;
    //return point(a.x + (c1*b2 - c2*b1)/d,a.y + (a1*c2 - a2*c1)/d);
}

hdu 3007 Buried memory

​​http://acm.hdu.edu.cn/showproblem.php?pid=3007​​

#include <cstdio>
#include <cmath>
#include <iostream>
#include <algorithm>
using namespace std;
const double eps=1e-7;
int n;
struct point{
    double x,y;
}p[505];
double area(point p0,point p1,point p2){
    return fabs((p0.x*p1.y+p2.x*p0.y+p1.x*p2.y-p2.x*p1.y-p1.x*p0.y-p0.x*p2.y)/2.0);
}
double dis(point p1,point p2){
    return sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y));
}
point get(point p0,point p1,point p2){
    point p;
    double c1=(p0.x*p0.x-p1.x*p1.x+p0.y*p0.y-p1.y*p1.y)/2;
    double c2=(p2.x*p2.x-p1.x*p1.x+p2.y*p2.y-p1.y*p1.y)/2;
    double x21=p2.x-p1.x,  x01=p0.x-p1.x;
    double y21=p2.y-p1.y,  y01=p0.y-p1.y;
    p.x=(c1*y21-c2*y01)/(x01*y21-x21*y01);
    p.y=(c1*x21-c2*x01)/(y01*x21-y21*x01);
    return p;
}
void min_cover(point p[],point &c,double &r){
    random_shuffle(p,p+n);
    c=p[0];
    r=0;
    for(int i=1;i<n;i++)
        if(dis(c,p[i])+eps>r){
            c=p[i];
            r=0;
            for(int j=0;j<i;j++){
                if(dis(c,p[j])+eps>r){
                    c.x=(p[i].x+p[j].x)/2;
                    c.y=(p[i].y+p[j].y)/2;
                    r=dis(c,p[j]);;
                    for(int k=0;k<j;k++){
                        if(dis(c,p[k])+eps>r){
                            c=get(p[i],p[j],p[k]);
                            r=dis(c,p[k]);
                        }
                    }
                }
            }
        }
}
int main(){
    //freopen("cin.txt","r",stdin);
    while(cin>>n&&n){
        for(int i=0;i<n;i++){
            scanf("%lf%lf",&p[i].x,&p[i].y);
        }
        point cen;
        double r;
        min_cover(p,cen,r);
        printf("%.2lf %.2lf %.2lf\n",cen.x, cen.y, r);
    }
    return 0;
}

如果没有点的随机化，那么结果不会出错，但是会影响效率：

zoj 1450 Minimal Circle

​​http://acm.zju.edu.cn/onlinejudge/?problemId=450​​

和上一题几乎一样，把点的个数上限改成105即可。