\(\mathcal{Description}\)

给定平面上 \(n\) 个点，求最小的能覆盖其中至少 \(m\) 个点的圆半径及一个可能的圆心。

\(n\le500\)，坐标值 \(X\in[0,10^4]\)。

\(\mathcal{Solution}\)

不难想到二分答案 \(r\)，以每个点为圆心，\(r\) 为半径作圆，若 \(r\) 合法则能找到一个被至少 \(m\) 个圆覆盖的点。

但是圆的交极难处理，结合数据范围，考虑通过一些枚举操作来简化问题——钦定最终圆心 \(O\) 在以 \(P_i\) 为圆心的圆内，同时考虑到 \(O\) 的最优性，这一限制可以进一步加强为：圆心 \(O\) 在以 \(P_i\) 为圆心的圆周上，于是可以在圆周上按极角做扫描线。到此，我们得到 \(\mathcal O(n^2\log n\log \frac{X}{\epsilon})\) 的算法，应该能过。

接下来是计算几何算法的常见套路：基于数据随机化降低复杂度。设钦定 \(O\) 在以 \(P_i\) 为圆心的圆周上时，最小的合法半径为 \(r_i\)。假设以 \(i=1..n\) 的顺序依次计算 \(r_i\)，此时我们可以用 \(\mathcal O(n\log n)\) 的时间判断某个 \(r_i\) 是否是前缀最小值。若 \(r_i\) 是前缀最小值，我们还需要 \(\mathcal O(n\log n\log\frac{X}{\epsilon})\) 的时间求出 \(r_i\) 的具体值。由于数据随机，所以期望复杂度为

\[\begin{aligned}
T(n) &= \sum_{i=1}^n\left(\frac{\mathcal O\left(n\log n\log\frac{X}{\epsilon}\right)}{i}+\mathcal O(n\log n)\right)\\
&= \mathcal O\left(n\log^2n\log\frac{X}{\epsilon}\right).
\end{aligned}
\]

相对于期望复杂度的优化，感觉第一步通过适当枚举加强限制，简化问题的思路更为重要。不要一心想着优化复杂度就不敢考虑暴力了。

\(\mathcal{Code}\)

/*+Rainybunny+*/

#include <bits/stdc++.h>

#define rep(i, l, r) for (int i = l, rep##i = r; i <= rep##i; ++i)

#define per(i, r, l) for (int i = r, per##i = l; i >= per##i; --i)

typedef std::pair<double, bool> PDB;

#define fi first

#define se second

const int MAXN = 500;

const double EPS = 1e-7, PI = acos(-1.);

int n, m;

inline double dabs(const double u) { return u < 0. ? -u : u; }

inline int dcmp(const double u) { return dabs(u) < EPS ? 0 : u < 0 ? -1 : 1; }

inline double norm(double u) {

    u = fmod(u, 2. * PI);

    return u < 0. ? u + 2. * PI : u;

}

struct Point {

    double x, y;

    inline Point operator * (const double k) const {

        return { k * x, k * y };

    }

    inline Point operator + (const Point& u) const {

        return { x + u.x, y + u.y };

    }

    inline bool operator == (const Point& u) const {

        return !dcmp(x - u.x) && !dcmp(y - u.y);

    }

    friend inline double dist(const Point& u, const Point& v) {

        return sqrt((u.x - v.x) * (u.x - v.x) + (u.y - v.y) * (u.y - v.y));

    }

} pnt[MAXN + 5];

inline bool check(const Point& O, const double R, Point& resO, double& resR) {

    static std::vector<PDB> evt; evt.clear();

    int cur = 0;

    rep (i, 1, n) if (!(pnt[i] == O)) {

        double d = dist(O, pnt[i]);

        if (d > 2. * R) continue;

        double alp = acos(0.5 * d / R);

        double the = norm(atan2(pnt[i].y - O.y, pnt[i].x - O.x));

        double u = norm(the - alp), v = norm(the + alp);

        if (u > v) ++cur;

        evt.push_back({ u, 1 }), evt.push_back({ v, 0 });

    }

    std::sort(evt.begin(), evt.end());

    if (cur >= m - 1) {

        if (R < resR) resO = O + Point{ 0, resR = R };

        return true;

    }

    for (auto [alp, opt]: evt) {

        if ((cur += opt ? 1 : -1) >= m - 1) {

            if (R < resR) resO = O + Point{ cos(alp), sin(alp) } * (resR = R);

            return true;

        }

    }

    return false;

}

int main() {

    scanf("%d %d", &n, &m);

    rep (i, 1, n) scanf("%lf %lf", &pnt[i].x, &pnt[i].y);

    std::shuffle(pnt + 1, pnt + n + 1, std::mt19937(20120712));

    double ans = 7.5e3; Point ansP = { -1., -1. };

    rep (i, 1, n) {

        if (!check(pnt[i], ans, ansP, ans)) continue;

        double l = 0., r = ans;

        while (l + EPS < r) {

            double mid = 0.5 * (l + r);

            if (check(pnt[i], mid, ansP, ans)) r = mid;

            else l = mid;

        }

    }

    printf("%f\n%f %f\n", ans, ansP.x, ansP.y);

    return 0;

}

巴特西

Solution -「CEOI 2006」「洛谷 P5974」ANTENNA

\(\mathcal{Description}\)

\(\mathcal{Solution}\)

\(\mathcal{Code}\)

最新文章

热门文章