最小圆覆盖学习笔记
前置知识：平面几何初步。
最小圆覆盖
最小圆覆盖是什么，顾名思义：给出N个点，让你画一个最小的包含所有点的圆，要求出圆的半径，及圆心的坐标。
三点定圆大家都知道吧，唯一要求是三点不共线。
求最小圆，就要枚举这些点，去确定圆。
对于可（du）敬（liu）的出题人，他可不想让你这么简单的做完这个题。
这时，我们就要用到随机增量法。
用这个东西：random_shuffle(a+1,a+n+1);把出题人善意给你的数据顺序打乱。
（你要问为什么？一会儿会解决吧）
下面我们可以开始求了。三点定圆，当然要三层循环，枚举三个点。
对于第一个点，我们让其为圆心。对于第二点我们和第一个点连线,取其中点为圆心.连线为直径,对于第三个点，我们进行三点定圆.
来了第四个点,我们判断其是否在圆中,如果在圆中,我们继续判断,如果不在圆中，我们取一二个点,和第四个点进行三点定圆.
因为第四个点不在前三个点的圆中.所以第四个点定的新圆一定比前三个定的圆大.
显然，第三个点的时间复杂度是o(n)的，表面上时间复杂的是O(n^3)
的,实际上在第一次执行完最内层循环.已经确定了一个圆。再执行外边的两层循环,就会涉及到概率.
均摊时间复杂度O(n),具体我也不会证.
cq(i,2,n){
if(dis(a[i],o)>ri+eps){
o=a[i];ri=0;
cq(j,1,i-1){
if(dis(o,a[j])>ri+eps){
o.x=(a[i].x+a[j].x)/2;
o.y=(a[i].y+a[j].y)/2;
ri=dis(o,a[j]);
cq(k,1,j-1){
if(dis(o,a[k])>ri+eps){
tt(a[i],a[j],a[k]);
}
}
}
}
}
}
三点定圆
前提是三点不共线。
圆的一般式方程大家都知道吧：(x−x_0)2+(y−y_0)2=r^2
我们设三个点为\((x_1,y_1)\),\((x_2,y_2)\),\((x_3,y_3)\)带入到方程里
公式（1）（2）相减，（1）（3）相减之后经过化简可以得到：
\((x_1−x_2)x_0+(y_1−y_2)y_0=\frac{(x^2_1−x^2_2)−(y^2_2−y^2_1)}{2}\)
\((x_1−x_3)x_0+(y_1−y_3)y_0=\frac{(x^2_1−x^2_3)−(y^2_3−y^2_1)}{2}\)
上面的式子太复杂，换一个好看一点的。
设
\(a=x_1−x_2\)
\(b=y_1−y_2\)
\(c=x_1−x_3\)
\(d=y_1−y_3\)
\(e=\frac{(x^2_1−x^2_2)−(y^2_2−y^2_1)}{2}\)
\(f=\frac{(x^2_1−x^2_3)−(y^2_3−y^2_1)}{2}\)
那么 解得
de−bf
\(x_0=\frac{de−bf}{bc-ad}\)
\(y_0=\frac{af−ce}{bc-ad}\)
有了 x_0 和 y_0 的值后，带入(1) 式就可以得到 r的值。至此，三点确定圆的问题就解决了。
没了
代码如下
void tt(Point p1,Point p2,Point p3){
double a,b,c,d,e,f;
a=p2.y-p1.y;
b=p3.y-p1.y;
c=p2.x-p1.x;
d=p3.x-p1.x;
f=p3.x*p3.x+p3.y*p3.y-p1.x*p1.x-p1.y*p1.y;
e=p2.x*p2.x+p2.y*p2.y-p1.x*p1.x-p1.y*p1.y;
o.x=(a*f-b*e)/(2*a*d-2*b*c);
o.y=(d*e-c*f)/(2*a*d-2*b*c);
ri=dis(o,p1);
}
完整代码
#include<iostream>
#include<cstdio>
#include<cmath>
#include<algorithm>
#define cq(i,s,n) for(int i=s;i<=n;i++)
using namespace std;
const double eps=1e-12;
struct Point{
double x,y;
}a[500005];
Point o;
int n;
double ri;
double dis(Point a,Point b){
return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
}
void tt(Point p1,Point p2,Point p3){
double a,b,c,d,e,f;
a=p2.y-p1.y;
b=p3.y-p1.y;
c=p2.x-p1.x;
d=p3.x-p1.x;
f=p3.x*p3.x+p3.y*p3.y-p1.x*p1.x-p1.y*p1.y;
e=p2.x*p2.x+p2.y*p2.y-p1.x*p1.x-p1.y*p1.y;
o.x=(a*f-b*e)/(2*a*d-2*b*c);
o.y=(d*e-c*f)/(2*a*d-2*b*c);
ri=dis(o,p1);
}
int main(){
scanf("%d",&n);
cq(i,1,n){
scanf("%lf%lf",&a[i].x,&a[i].y);
}
random_shuffle(a+1,a+n+1);
o=a[1];ri=0;
cq(i,2,n){
if(dis(a[i],o)>ri+eps){
o=a[i];ri=0;
cq(j,1,i-1){
if(dis(o,a[j])>ri+eps){
o.x=(a[i].x+a[j].x)/2;
o.y=(a[i].y+a[j].y)/2;
ri=dis(o,a[j]);
cq(k,1,j-1){
if(dis(o,a[k])>ri+eps){
tt(a[i],a[j],a[k]);
}
}
}
}
}
}
printf("%.10lf\n%.10lf %.10lf",ri,o.x,o.y);
return 0;
}
\({\Huge\color{Salmon}{习题}}\)（话说可以算成双倍经验）
P2533 [AHOI2012]信号塔
P1742 最小圆覆盖
