题解：P3517 [POI 2011] WYK-Plot

· · 题解

::::warning 这篇题解使用了洛谷的 2025 年 7 月渲染更新功能！ ::::

题意

给定平面上 n 个点 (x_i,y_i)，你需要将其按照编号划分成不超过 m 段，对于划分成 [l,r] 的段的代价如下：

  • 先选定一个点 (x_0,y_0)；
  • 该段的代价即为 \max_{i\in[l,r]}\operatorname{dis}\{(x_0,y_0),(x_i,y_i)\}，其中 \operatorname{dis}\{(x_0,y_0),(x_1,y_1)\} 表示平面上点 (x_0,y_0) 到 (x_1,y_1) 的直线距离。

求最小的划分总代价，并且输出每一段你选择的 (x_0,y_0)。精度要求 10^{-6}。

题解

首先，对于确定的一段 [l,r]，最优的 (x_0,y_0) 的选择就是最小圆覆盖问题，即用一个半径最小的圆覆盖 (x_l,y_l),(x_{l+1},y_{l+1}),\cdots,(x_r,y_r) 这些点。

::::info[最小圆覆盖的做法] 假设要求 (x_1,y_1),(x_2,y_2),\cdots,(x_k,y_k) 这些点的最小圆覆盖。将所有点随机打乱，从小到大枚举 i，假设已经维护出 1\sim i-1 的最小圆覆盖 C，如果 i 在 C 里那么最小圆覆盖不变；否则 i 一定在 1\sim i 的最小圆覆盖 C' 上（否则可以调整得到更优的圆覆盖）。此时我们把新的 1\sim i 的最小圆覆盖 C' 设成圆心为 (x_i,y_i)、半径为 0 的圆，然后枚举 j\in[1\sim i-1]，如果 j 在 C' 里那么枚举下一个 j；否则 i,j 一定一起在最小圆覆盖上。此时我们就把 C' 设为以连接 (x_i,y_i) 和 (x_j,y_j) 的线段为直径的圆，接着枚举 k\in[1,j-1]，如果 k 在 C' 里那么枚举下一个 k，否则 i,j,k 三个点就能确定一个圆了，把 C' 设成过这三个点的圆作为 C'。

时间复杂度期望 O(n)，大概有一个 9 的常数。 :::: 考虑二分答案，二分最大半径 R，每次贪心地选择尽可能长的段划分，最后判断是否划分了 \le m 即可。现在考虑怎么知道从 i 开始最远能够划到哪里。

由于最小圆覆盖依赖点集随机化，因此我们不能从小到大枚举右端点并且动态维护最小圆覆盖。因此考虑再二分最远的右端点 j，判断 [i,j] 的最小圆覆盖是否 \le R。可是这样的复杂度非常高，原因在于可能最后的右端点离左端点非常近，但是我们却从很远的地方开始二分。

那么我们考虑先进行一次倍增，即枚举 j=i+2^0,i+2^1,\cdots i+2^k 时 [i,j] 的最小圆覆盖的半径是否 \le R；假设我们发现当 j=i+2^k 时，[i,i+2^{k-1}] 可以划分但是 [i,i+2^k] 就无法划分，这就意味着 j 只可能在 [i+2^{k-1},i+2^k-1] 的范围里，也就是我们只需要在这个范围里二分即可。

这么做的时间复杂度可以这么分析：

所以总时间复杂度就是 O(n\log n\log V)。

Code

#include <bits/stdc++.h>
bool MemoryST; using namespace std;
#define ll long long
#define mk make_pair
#define open(x) freopen(#x".in", "r", stdin), freopen(#x".out", "w", stdout)
#define lowbit(x) ((x) & (-(x)))
#define lson l, mid, rt << 1
#define rson mid + 1, r, rt << 1 | 1
#define BCNT __builtin_popcount
#define cost_time (1e3 * clock() / CLOCKS_PER_SEC) << "ms"
#define cost_space (abs(&MemoryST - &MemoryED) / 1024.0 / 1024.0) << "MB"
const int inf = 0x3f3f3f3f; 
const ll linf = 1e18; 
mt19937 rnd(random_device{}());
template<typename T> void chkmax(T& x, T y) { x = max(x, y); }
template<typename T> void chkmin(T& x, T y) { x = min(x, y); }
template<typename T> T abs(T x) { return (x < 0) ? -x : x; }
struct Point { double x, y; Point ( double x0 = 0, double y0 = 0) { x = x0, y = y0; }};
struct Circle { double R; Point O; Circle(double R0 = 0, Point O0 = Point()) { R = R0, O = O0; }};
const double eps = 1e-8;
bool MemoryED; int main() { 
    int n, m; cin >> n >> m; vector<Point> a(n);
    for (auto &[x, y] : a) cin >> x >> y;
    auto sqr = [&](double x) -> double { return x * x; };
    auto cmp = [&](double x, double y) -> bool { return x - y < eps && x - y > -eps; };
    auto dis = [&](Point A, Point B) -> double { return sqrt((A.x - B.x) * (A.x - B.x) + (A.y - B.y) * (A.y - B.y)); };
    auto getcir = [&](Point A, Point B, Point C) -> Circle { // 求经过三个点的圆的圆心 & 半径
        double a1, a2, b1, b2, c1, c2, R; Point O;
        a1 = (B.x - A.x) * 2.0, b1 = (B.y - A.y) * 2.0;
        c1 = sqr(B.x) - sqr(A.x) + sqr(B.y) - sqr(A.y);
        a2 = (C.x - A.x) * 2.0, b2 = (C.y - A.y) * 2.0;
        c2 = sqr(C.x) - sqr(A.x) + sqr(C.y) - sqr(A.y);
        if (cmp(a1, 0)) O.y = c1 / b1, O.x = (c2 - O.y * b2) / a2;
        else if (cmp(b1, 0)) O.x = c1 / a1, O.y = (c2 - O.x * a2) / b2;
        else O = Point((c2 * b1 - c1 * b2) / (a2 * b1 - a1 * b2), (c2 * a1 - c1 * a2) / (b2 * a1 - b1 * a2));
        return Circle(dis(O, A), O);
    }; 
    auto solve = [&](vector<Point> p) -> Circle { // 求最小圆覆盖
        shuffle(p.begin(), p.end(), rnd); Point O = p[0]; double R = 0; int len = (int)p.size();
        for (int i = 0; i < len; i ++) if (dis(O, p[i]) - R > eps) {
            O = Point((p[0].x + p[i].x) / 2, (p[0].y + p[i].y) / 2), R = dis(p[0], p[i]) / 2;
            for (int j = 1; j < i; j ++) if (dis(O, p[j]) - R > eps) {
                O = Point((p[j].x + p[i].x) / 2, (p[j].y + p[i].y) / 2), R = dis(p[i], p[j]) / 2;
                for (int k = 0; k < j; k ++) if (dis(O, p[k]) - R > eps) {
                    Circle nw = getcir(p[i], p[j], p[k]);
                    O = nw.O, R = nw.R;
                }
            }
        } 
        return Circle(R, O);
    }; vector<Point> Stand;
    double l = 0, r = 1e7, Ans = 1e7; for (int _ = 0; _ <= 50 && r - l > eps; _ ++) {
        double mid = (l + r) / 2;
        auto check = [&](double R) -> pair<bool, vector<Point> > {
            // cout << "[check] " << R << '\n';
            int cnt = 0; vector<Point> ans;
            for (int i = 0, j; i < n && cnt <= m; cnt ++) {
                Point stand = a[i];
                for (j = 2; i + j <= n; j <<= 1) {
                    vector<Point> now; for (int k = i; k < i + j; k ++) now.push_back(a[k]);
                    Circle c = solve(now);
                    if (c.R - R > eps) break;
                    stand = c.O;
                } int nxt = min(n - 1, i + (j >> 1) - 1);
                // cout << "expected segment: " << i + (j >> 1) << ' ' << i + j - 1 << '\n';
                for (int l = min(n - 1, i + (j >> 1)), r = min(n - 1, i + j - 1); l <= r; ) {
                    int mid = (l + r) >> 1; vector<Point> now;
                    for (int k = i; k <= mid; k ++)
                        now.push_back(a[k]);
                    Circle c = solve(now);
                    if (c.R - R > eps) r = mid - 1;
                    else stand = c.O, nxt = mid, l = mid + 1;
                } 
                // cout << "Part: " << i << ' ' << nxt << '\n';
                ans.push_back(stand), i = nxt + 1;
            } // cout << R << ' ' << cnt << '\n';
            return mk(cnt <= m, ans);
        }; auto [ok, tmp] = check(mid);
        if (ok) Ans = mid, Stand = tmp, r = mid;
        else l = mid;
    } cout << fixed << setprecision(8);
    cout << Ans << '\n' << (int)Stand.size() << '\n';
    for (auto [x, y] : Stand) cout << x << ' ' << y << '\n';
    return 0;
}