P4586 [FJOI2015]最小覆盖双圆问题 （最小圆覆盖+超级计算几何）

最小双圆覆盖问题

问题描述：给定n<=10000个点，找出两个半径相同的圆覆盖掉所有点，求出最小半径。

题解：参考超级大佬

定一思想：确定一个变量，答案确定使，肯定有一条线能将点分割成两部分，使一部分只在R1内，另一部分只在R2内。然后就转变成确定分割线，找两边最小覆盖圆的半径的最大值，就是答案。

问题：怎么确定--二分分割线；分割线不平行于y轴--枚举旋转所有点。

代码：

#include <bits/stdc++.h>
#define ll long long
#define rep(i, a, b) for (ll i = a; i <= b; i++)
using namespace std;
const int N = 1005;
const double an = 1.0 / 100 * acos(-1.0), eps = 1e-9;
const double si = sin(an), co = cos(an);



struct node {
    double x, y;
    node() {}                                     //不知道干嘛
    node(double xx, double yy) : x(xx), y(yy) {}  //使能写成node形式？
    node operator+(const node &b) const { return node(x + b.x, y + b.y); }
    node operator-(const node &b) const { return node(x - b.x, y - b.y); }
    double operator*(const node &b) const { return x * b.y - y * b.x; }
    node operator*(const double &b) const { return node(x * b, y * b); }
    double operator^(const node &b) const { return x * b.x + y * b.y; }
    double len2() { return x * x + y * y; }
    node rot() {
        return node(x * co - y * si, x * si + y * co);
    }                                     //逆时针旋转an度
    node rot90() { return node(-y, x); }  //逆时针旋转90度
    // void shake() { x += reps(), y += reps(); }//扰动
} p[N], st[N], o;
bool cmp(const node &a, const node &b) {
    if (fabs(a.x - b.x) < eps) return a.y - b.y < eps;
    return a.x - b.x < eps;
}
struct Line {
    node p, v;
    inline Line() {}
    inline Line(node pp, node vv) : p(pp), v(vv) {}
    friend node cross(const Line &a, const Line &b) {
        return a.p + a.v * (b.v * (b.p - a.p) / (b.v * a.v));
    }
};
//三点确定一个最小圆
node circle(node a, node b, node c) {
    return cross(Line((a + b) * 0.5, (b - a).rot90()),
                 Line((a + c) * 0.5, (c - a).rot90()));
}
//返回能包含点l1~l2的最小圆的半径，o为圆心
double query(int l1, int l2) {
    if (l1 > l2) return 0;
    ll top = 0;
    double r = 0;
    o = node(0, 0);
    rep(i, l1, l2) st[++top] = p[i];
    random_shuffle(st + 1, st + 1 + top);//
    rep(i, 1, top) if ((st[i] - o).len2() - r > eps) {
        o = st[i], r = 0;
        rep(j, 1, i - 1) if ((st[j] - o).len2() - r > eps) {
            o = (st[i] + st[j]) * 0.5, r = (st[i] - o).len2();
            rep(k, 1, j - 1) if ((st[k] - o).len2() - r > eps)
                o = circle(st[i], st[j], st[k]),
                r = (st[i] - o).len2();
        }
    }
    return sqrt(r);
}
int n;
double res;
int main() {
    // srand(20030719);
    while (cin >> n) {
        res = 1e18;
        if (n == 0) break;
        rep(i, 1, n) cin >> p[i].x >> p[i].y;  //, p[i].shake();
        rep(d, 0, 100) {
            rep(i, 1, n) p[i] = p[i].rot();
            sort(p + 1, p + 1 + n, cmp);
            ll l = 1, r = n;
            while (l <= r) {
                ll mid = (l + r) >> 1;
                double r1 = query(1, mid), r2 = query(mid + 1, n);
                double ans = max(r1, r2);
                if (r1 + r2 - ans - res > eps) break;
                res = min(res, ans);
                if (r1 - r2 < eps)
                    l = mid + 1;
                else
                    r = mid - 1;
            }
        }
        // cout << ">>>";
        printf("%.2lf\n", res);
    }
    return 0;
}