最小圆覆盖问题

//c为圆心，r为半径

const double eps=1e-8;
struct Point {
    double x,y;
} p[505];

double dis(const Point &a,const Point &b) {
    return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
}
Point circumcenter(const Point &a,const Point &b,const Point &c) {
    //返回三角形的外心
    Point ret;
    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;
    ret.x=a.x+(c1*b2-c2*b1)/d;
    ret.y=a.y+(a1*c2-a2*c1)/d;
    return ret;
}



void min_cover_circle(Point *p,int n,Point &c,double &r) { // p为点的集合；c为圆心，r为半径
    random_shuffle(p,p+n); 
    c=p[0];
    r=0;
    for(int i=1; i<n; i++) {
        if(dis(p[i],c)>r+eps) { 
            c=p[i];
            r=0;
            for(int j=0; j<i; j++){
                if(dis(p[j],c)>r+eps) { 
                    c.x=(p[i].x+p[j].x)/2;
                    c.y=(p[i].y+p[j].y)/2;
                    r=dis(p[j],c);
                    for(int k=0; k<j; k++) {
                        if(dis(p[k],c)>r+eps) { 
                        c=circumcenter(p[i],p[j],p[k]);
                        r=dis(p[i],c);
                        }
                    }
                }
            }
        }
    }
}