题目链接:P1742 最小圆覆盖
题意
给出 N 个点,求最小的包含所有点的圆。
思路
随机增量
最小圆覆盖一般有两种做法:随机增量和模拟退火。随机增量的精确度更高,这里介绍随机增量的做法。
先将所有点随机打乱。
令前 (i - 1) 个点的最小覆盖圆为圆 (O),加入第 (i) 个点。
如果第 (i) 个点在圆 (O) 内或圆 (O) 上,则前 (i) 个点的最小覆盖圆还是圆 (O)。
否则新得到的最小覆盖圆肯定经过第 (i) 个点。然后确定前 (i − 1) 个点中还有哪两个点在最小覆盖圆上。
以第 (i) 个点为圆心,半径为 (0),重复以上过程依次加入点 (P_j)。(圆心为 (frac{P_i+P_j}{2}),半径为 (frac{|P_iP_j|}{2}))
固定两个点之后再重复以上步骤找第三个点。(因为需要三个点来确定一个圆)
遍历完所有点之后,所得到的圆就是最小覆盖圆。
时间复杂度为 (O(n))
代码
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef double db;
const db eps = 1e-8;
const db pi = acos(-1.0);
const ll maxn = 1e5 + 10;
inline int dcmp(db x) {
if(fabs(x) < eps) return 0;
return x > 0? 1: -1;
}
class Point {
public:
double x, y;
Point(double x = 0, double y = 0) : x(x), y(y) {}
void input() {
scanf("%lf%lf", &x, &y);
}
bool operator<(const Point &a) const {
return (!dcmp(x - a.x))? dcmp(y - a.y) < 0: x < a.x;
}
bool operator==(const Point &a) const {
return dcmp(x - a.x) == 0 && dcmp(y - a.y) == 0;
}
db dis2(const Point a) {
return pow(x - a.x, 2) + pow(y - a.y, 2);
}
db dis(const Point a) {
return sqrt(dis2(a));
}
db dis2() {
return x * x + y * y;
}
db dis() {
return sqrt(dis2());
}
Point operator+(const Point a) {
return Point(x + a.x, y + a.y);
}
Point operator-(const Point a) {
return Point(x - a.x, y - a.y);
}
Point operator*(double p) {
return Point(x * p, y * p);
}
Point operator/(double p) {
return Point(x / p, y / p);
}
db dot(const Point a) {
return x * a.x + y * a.y;
}
db cross(const Point a) {
return x * a.y - y * a.x;
}
db ang(Point a) {
return acos((a.dis() * dis()) / dot(a));
}
};
typedef Point Vector;
class Circle {
public:
Point o;
db r;
Circle() {}
Circle(Point o, db r):o(o), r(r){}
// 三点定圆
Circle(Point A, Point B, Point C) {
double a1 = B.x - A.x, b1 = B.y - A.y, c1 = (a1 * a1 + b1 * b1) / 2;
double a2 = C.x - A.x, b2 = C.y - A.y, c2 = (a2 * a2 + b2 * b2) / 2;
double d = a1 * b2 - a2 * b1;
o.x = A.x + (c1 * b2 - c2 * b1) / d;
o.y = A.y + (a1 * c2 - a2 * c1) / d;
r = o.dis(A);
}
Point point(db a) {
return Point(o.x + cos(a) * r, o.y + sin(a) * r);
}
// 点在圆内
bool PinC(Point p) {
db d = p.dis(o);
return dcmp(d - r) < 0;
}
// 点在圆外
bool PoutC(Point p) {
db d = p.dis(o);
return dcmp(d - r) > 0;
}
};
vector<Point> p;
// 最小圆覆盖
Circle min_circle(vector<Point> p) {
int sz = p.size();
random_shuffle(p.begin(), p.end());
Circle ans(p[0], 0.0);
for(int i = 0; i < sz; ++i) {
if(ans.PoutC(p[i])) {
ans = Circle(p[i], 0);
for(int j = 0; j < i; ++j) {
if(ans.PoutC(p[j])) {
ans = Circle((p[i] + p[j]) / 2.0, p[i].dis(p[j]) / 2.0);
for(int k = 0; k < j; ++k) {
if(ans.PoutC(p[k])) {
ans = Circle(p[i], p[j], p[k]);
}
}
}
}
}
}
return ans;
}
int main() {
int n;
scanf("%d", &n);
for(int i = 0; i < n; ++i) {
Point tmp;
tmp.input();
p.push_back(tmp);
}
Circle ans = min_circle(p);
printf("%.10lf
%.10lf %.10lf
", ans.r, ans.o.x, ans.o.y);
return 0;
}