最小圆覆盖

最小圆覆盖

算法目的：在线性时间复杂度内求出覆盖n个点的最小圆
算法步骤(注意起始最小圆的圆心为第一个点，半径为0，同时1 < i <= n,1 < j < i,1 < k < j)：

1.首先将所有点随机排列；
2.按顺序依次将点加入已经建好的圆中，每添加一个点就进入步骤3；
3.如果点i在当前最小圆的外部，那么说明点i一定在前i个点构成的最小圆的边界上(因为要保证最小圆)，进入步骤4处理；如果点i在当前最小圆的内部，则返回步骤2；
4.此时已经确定点i一定在前i个点所构成的最小圆的边界上，则直接将点i作为圆心，半径为0，重新将前i-1个点加入到这个圆中(类似上面的步骤，不过此时求的是包含点i的前j个点的最小覆盖圆)，每加一个点就进入步骤5；
5.如果当前j点在最小圆外部，则说明点j一定在前j个点(包含点i)构成的最小圆的边界上(与上类似)，进入步骤6；否则不需要更新，返回4；
6.此时已经确定点i，j在前j个点(包括点i)构成的圆的边界上，直接将圆心设为点i和j的中点，直径为两点距离，重新把前j-1个点加入到这个圆中，每加入一个进入7；
7.如果当前点k在最小圆的外部，则点k在前k个点(包括i，j)构成的最小圆的边界，此时i,j,k都在次圆边界，直接3点共圆求出圆心与半径即可；否则无需更新，返回6；

Code

#include<bits/stdc++.h>//O(n)，最小圆覆盖

#define ll long long
#define met(a, x) memset(a,x,sizeof(a))
#define inf 0x3f3f3f3f
#define ull unsigned long long
using namespace std;
const double eps = 1e-12;
struct node {
    double x, y;
} s[500005];
node o;//圆心坐标
double ri;//半径

double dis(node a, node b) {
    return sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y));
}

void getr(node p1, node p2, node p3) {//三个点求三角形圆心坐标和半径
    double a, b, c, d, e, f;
    a = p2.y - p1.y;
    b = p3.y - p1.y;
    c = p2.x - p1.x;
    d = p3.x - p1.x;
    f = p3.x * p3.x + p3.y * p3.y - p1.x * p1.x - p1.y * p1.y;
    e = p2.x * p2.x + p2.y * p2.y - p1.x * p1.x - p1.y * p1.y;
    o.x = (a * f - b * e) / (2 * a * d - 2 * b * c);
    o.y = (d * e - c * f) / (2 * a * d - 2 * b * c);
    ri = dis(o, p1);
}

int main() {
    ios::sync_with_stdio(false);
    cin.tie(0);
    int n;
    cin >> n;
    for (int i = 1; i <= n; i++) {
        cin >> s[i].x >> s[i].y;
    }
    random_shuffle(s + 1, s + n + 1);
    o = s[1];
    ri = 0;
    for (int i = 2; i <= n; i++) {
        if (dis(s[i], o) > ri + eps) {
            o = s[i];
            ri = 0;//第一个点为圆心
            for (int j = 1; j < i; j++) {
                if (dis(o, s[j]) > ri + eps) {
                    o.x = (s[i].x + s[j].x) / 2;
                    o.y = (s[i].y + s[j].y) / 2;
                    ri = dis(o, s[j]);//第一个点和第二个点中点为圆心，距离为直径
                    for (int k = 1; k < j; k++) {
                        if (dis(o, s[k]) > ri + eps) {
                            getr(s[i], s[j], s[k]);//三点定圆
                        }
                    }
                }
            }
        }
    }
    cout << fixed << setprecision(10) << ri << endl;
    cout << fixed << setprecision(10) << o.x << ' ' << fixed << setprecision(10) << o.y << endl;
    return 0;
}