题面
给出N个点,让你画一个最小的包含所有点的圆。
分析
这个问题实际上是找出三个点确定的一个圆,能包含其他所有点。
与三点共圆有密切的关系,先来说三点共圆。
三点共圆有两种情况:三个点共线或不共线,共线则两个最远点构成直径,不共线则圆为三角形外接圆
重点考虑外接圆情形,班经用正弦定理是比较好找的,其实难度是找圆心,也就是用三角形三个顶点坐标来表示外心坐标。
这里代码中的三点共圆圆心计算参考了这篇文章。
其实也有其他的办法如两条线的中垂线交圆心,向量垂直过中点交圆心等,不过代码上用三角形三点坐标直接实现比较容易
下面来谈最小圆覆盖,其实用的是数学归纳法,如果已经覆盖了前 i 个点,就想办法造一个圆让其覆盖 i+1 个点。
下面考虑步骤:
1.已知 1 ~ i-1 个点已经由圆
C
i
−
1
C_{i-1}
Ci−1覆盖,现在欲加入第 i 个点,让它也被圆覆盖
2.如果第 i 个点已经在圆
C
i
−
1
C_{i-1}
Ci−1中,则下一个圆
C
i
=
C
i
−
1
C_i = C_{i-1}
Ci=Ci−1,并且回到第一步继续
3.如果不在,则需要构造新的圆包含 1 ~ i 个点,首先这个新加入的点 i 需要在圆边上(临界情况),这也就确定了圆上的一个点
4.用第一个点和第 i 个点造圆
T
1
T_{1}
T1,这时如果有一个点 j (j<i) 不在
T
1
T_{1}
T1内,就把它取作
C
i
C_i
Ci边上的第二个点。
5.用 i 和 j 两个点造圆
T
2
T_2
T2,同4步,找不在圆
T
2
T_2
T2内的点 k (k<j<i),如果存在,那么k就是圆
C
i
C_i
Ci边上的第三个点。
6.现在 i j k 三个点都在圆
C
i
C_i
Ci 上,且能证明
C
i
C_i
Ci 包含了 1~i 所有点,于是用三点共圆构造
C
i
C_i
Ci,回到1继续迭代。
迭代到
C
n
C_n
Cn即结束
这看似是三层循环,其实期望复杂度 O(n) ,但为了防止数据不合适,在处理前要先打乱给的 n 个点,用 random_shuffle() 可以实现
代码
#include <iostream>
#include<algorithm>
#include<string.h>
#include<string>
#include<math.h>
#include<iomanip>
using namespace std;
struct Point
{
double x, y;
Point() {}
Point(double x,double y):x(x),y(y){}
};//存点结构体
inline double dis(Point &p1, Point &p2)//两个点之间的距离
{
return sqrt((p1.x - p2.x) * (p1.x - p2.x) + (p1.y - p2.y) * (p1.y - p2.y));
}
inline double sqr(double x) { return x * x; }//平方
class min_circle {
private:
Point p[100005];
int n;
public:
Point center;//代表每一次构造出圆的圆心
double r;//代表每一次构造的半径
void init(int n)
{
this->n = n;
for (int i = 0; i < n; i++)cin >> p[i].x >> p[i].y;
random_shuffle(p, p + n);//打乱,避免卡时
}
void circular()//构造圆序列
{
r = 0, center = p[0];
for (int i = 1; i < n; i++)//尝试判定第i个点在不在圆内,如果在则进入下一个点,不在则造一个p[i]在边上的圆
if (dis(center, p[i]) - r > 1e-8)//p[i]在之前造的圆外
{
r=dis(p[0],p[i])/2;//以p[0]和p[i]为直径造一个圆,找出在这个圆外的点p[j],作为第二个在本次构造圆边上的点
center.x = (p[0].x + p[i].x) / 2;
center.y= (p[0].y + p[i].y) / 2;
for (int j = 0; j < i; j++)
{
if (dis(center, p[j]) - r > 1e-8)//发现了这种圆外的j
{
r = dis(p[i], p[j]) / 2;//p[i]和p[j]是两个本轮最终构造圆边上的点,就以这两者为直径构造一个圆,尝试找在其外部的k点
center.x = (p[i].x + p[j].x) / 2;
center.y = (p[i].y + p[j].y) / 2;
for (int k = 0; k < j; k++)
{
if (dis(center, p[k]) - r > 1e-8)//发现这样的k
{
center = build_circle(p[i], p[j], p[k]);//发现p[i],p[j],p[k]都在本轮构造圆边上,圆被确定了
r = dis(center, p[i]);
}
}
}
}
}
}
Point build_circle(Point a, Point b, Point c) {//三点造圆
double a1, a2, b1, b2, c1, c2;
Point ans;
a1 = 2 * (b.x - a.x), b1 = 2 * (b.y - a.y),
c1 = sqr(b.x) - sqr(a.x) + sqr(b.y) - sqr(a.y);
a2 = 2 * (c.x - a.x), b2 = 2 * (c.y - a.y),
c2 = sqr(c.x) - sqr(a.x) + sqr(c.y) - sqr(a.y);
if (fabs(a1-0.0)<1e-8) {
ans.y = c1 / b1;
ans.x = (c2 - ans.y * b2) / a2;
}
else if (fabs(b1-0)<1e-8) {//因为分母上有这两者,这两种情况特判
ans.x = c1 / a1;
ans.y = (c2 - ans.x * a2) / b2;
}
else {
ans.x = (c2 * b1 - c1 * b2) / (a2 * b1 - a1 * b2);
ans.y = (c2 * a1 - c1 * a2) / (b2 * a1 - b1 * a2);
}
return ans;
}
}MC;
int main()
{
ios::sync_with_stdio(false);
int n;
cin >> n;
MC.init(n);//读入 与 随机化,打乱
MC.circular();
cout << fixed << setprecision(10) << MC.r << endl << MC.center.x << " " << MC.center.y;
return 0;
}