0%

洛谷 P1742 最小圆覆盖 (随机增量)

题目链接: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)$

代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
#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\n%.10lf %.10lf\n", ans.r, ans.o.x, ans.o.y);
return 0;
}

参考

欢迎关注我的其它发布渠道