学习笔记 - 最小圆覆盖
最小圆覆盖问题是这样的：
给出 \(N\) 个点，让你画一个最小的包含所有点的圆。
首先，答案一定是一个由这些点里面的三个点所确定的圆，或者由其中的两个点作为直径的圆。
所以我们很容易得到一个 \(O(n^4)\) 的算法 - 枚举三个在圆上的点并判断有没有包含。显然只要存在一个包含了所有点的圆，就是最小圆了。（也就是说不需要判断是不是最小，
下面是一种在期望情况下为 \(O(n)\) 的算法 - 随机增量法。
这个算法基于一个这样的事实：对于任意的 \(N\) 个点，其中的一个点如果没有出现在其它的点的最小覆盖圆上，那么这个点一定出现在了所有 \(N\) 个点的最小覆盖圆的圆周上。
对于第 \(i\) 个点 \(P_i\)，我们假设已经得到了一个包含前 \(i - 1\) 个点的圆 \(C\)。
如果 \(P_i\) 已经在 \(C\) 中了，那么不需要考虑，直接跳过就可以了。
否则，\(P_i\) 一定在包含前 \(i\) 个点的圆的圆周上。我们将 \(C\) 重置为以 \(P_i\) 为圆心，\(0\) 为半径的圆，然后枚举 \(j(1 \leq j < i)\)。
此时圆 \(C\) 应该表示的是包含前 \(j-1\) 个点和 \(P_i\) 的圆。
如果点 \(P_j\) 在圆 \(C\) 中，那么依然不需要考虑，继续枚举下一个 \(j\) 就可以了。
否则点 \(P_j\) 一定在包含前 \(j\) 个点和 \(P_i\) 的圆的圆周上。将圆 \(C\) 重置为以 \(P_iP_j\) 为直径的圆。
然后枚举 \(k(1\leq k < j)\)，此时圆 \(C\) 表示的应该是包含前 \(k-1\) 个点和点 \(P_j, P_i\) 的圆。
如果 \(P_k\) 在 \(C\) 中，那么有也不要考虑。否则点 \(P_k\) 一定在包含前 \(k\)点和点 \(P_j, P_i\) 的圆的圆周上，直接将圆 \(C\) 重置为由点 \(P_i, P_j, P_k\) 确定的圆就可以了。
上面的的算法的复杂度看起来是 \(O(n^3)\) 的，实际上，我们可以证明，在期望情况下，它的复杂度是 \(O(n)\) 的。
首先，对于一个长度为 \(n\) 的顺序随机的点序列 \(P\)，假设它们的最小覆盖圆为 \(C\)，则对于 \(\forall P_i\) 出现在 \(C\) 的圆周上的概率均等，又一共会有 \(3\) 个点出现在圆周上，那么每一个的出现在圆周上的概率为 \(\frac 3n\)。
然后，对于第 \(i\) 个点，若它没有出现在由前 \(i-1\) 个点的最小覆盖圆中，那么它一定在前 \(i\) 个点的最小覆盖圆的圆周上，如上文所述，这个概率为 \(\frac 3i\)。因此，需要枚举 \(j\) 的情况只有 \(\frac 3i\)。同理，对于 \(j\) 来说，需要枚举 \(k\) 的情况也只有 \(\frac 3j\) 的概率。
因此在期望情况下，总的时间 \(T(n) = \sum\limits_{i=1}^n \frac 3i \sum\limits_{j=1}^{i=1} \frac 3j \cdot j = \sum\limits_{i=1}^n \frac 3i \cdot 3i = 9n\)。于是总期望时间复杂度为 \(O(n)\)。
模板题 - bzoj1336 [Balkan2002]Alien最小圆覆盖
直接用板子就可以了。
#include<bits/stdc++.h>
#define fec(i, x, y) (int i = head[x], y = g[i].to; i; i = g[i].ne, y = g[i].to)
#define dbg(...) fprintf(stderr, __VA_ARGS__)
#define File(x) freopen(#x".in", "r", stdin), freopen(#x".out", "w", stdout)
#define fi first
#define se second
#define pb push_back
template<typename A, typename B> inline char smax(A &a, const B &b) {return a < b ? a = b , 1 : 0;}
template<typename A, typename B> inline char smin(A &a, const B &b) {return b < a ? a = b , 1 : 0;}
typedef long long ll; typedef unsigned long long ull; typedef std::pair<int, int> pii;
template<typename I>
inline void read(I &x) {
int f = 0, c;
while (!isdigit(c = getchar())) c == '-' ? f = 1 : 0;
x = c & 15;
while (isdigit(c = getchar())) x = (x << 1) + (x << 3) + (c & 15);
f ? x = -x : 0;
}
const int N = 100000 + 7;
const double eps = 1e-10;
int n, m;
inline int dcmp(const double &x) { return fabs(x) < eps ? 0 : (x < 0 ? -1 : 1); }
struct Point {
double x, y;
inline Point(const double &x = 0, const double &y = 0) : x(x), y(y) {}
} a[N];
inline double dist(const Point &a, const Point &b) { return sqrt((a.x - b.x) * (a.x - b.x) + (a.y - b.y) * (a.y - b.y)); }
inline std::pair<Point, double> get_C(const Point &p1, const Point &p2, const Point &p3) {
double a1 = 2 * (p1.x - p2.x), b1 = 2 * (p1.y - p2.y), c1 = (p1.x * p1.x + p1.y * p1.y) - (p2.x * p2.x + p2.y * p2.y);
double a2 = 2 * (p2.x - p3.x), b2 = 2 * (p2.y - p3.y), c2 = (p2.x * p2.x + p2.y * p2.y) - (p3.x * p3.x + p3.y * p3.y);
double x = (b1 * c2 - b2 * c1) / (a2 * b1 - a1 * b2), y = (a1 * c2 - a2 * c1) / (a1 * b2 - a2 * b1);
double r = sqrt((p1.x - x) * (p1.x - x) + (p1.y - y) * (p1.y - y));
return std::make_pair(Point(x, y), r);
}
inline void work() {
// std::mt19937 rnd(time(0) + (ull)new char);
// std::shuffle(a + 1, a + n + 1, rnd);
srand(time(0) + (ull)new char);
std::random_shuffle(a + 1, a + n + 1);
Point O = a[1];
double r = 0;
for (int i = 2; i <= n; ++i) if (dcmp(dist(a[i], O) - r) > 0) {
O = a[i], r = 0;
for (int j = 1; j < i; ++j) if (dcmp(dist(a[j], O) - r) > 0) {
O = Point((a[i].x + a[j].x) / 2, (a[i].y + a[j].y) / 2), r = dist(a[i], a[j]) / 2;
for (int k = 1; k < j; ++k) if (dcmp(dist(a[k], O) - r) > 0) {
const std::pair<Point, double> &tmp = get_C(a[i], a[j], a[k]);
O = tmp.fi, r = tmp.se;
}
}
}
printf("%.2lf\n", r);
printf("%.2lf %.2lf\n", O.x, O.y);
}
inline void init() {
scanf("%d", &n);
for (int i = 1; i <= n; ++i) scanf("%lf%lf", &a[i].x, &a[i].y);
}
int main() {
#ifdef hzhkk
freopen("hkk.in", "r", stdin);
#endif
init();
work();
fclose(stdin), fclose(stdout);
return 0;
}
别的一些题目
[https://lydsy.com/JudgeOnline/problem.php?id=2280](BZOJ2280 [Poi2011]Plot) 题解
