【bzoj1336/1337/2823】[Balkan2002]Alien最小圆覆盖 随机增量法

题目描述

给出N个点，让你画一个最小的包含所有点的圆。

输入

先给出点的个数N,2<=N<=100000，再给出坐标Xi,Yi.(-10000.0<=xi,yi<=10000.0)

输出

输出圆的半径，及圆心的坐标

样例输入

6
8.0 9.0
4.0 7.5
1.0 2.0
5.1 8.7
9.0 2.0
4.5 1.0

样例输出

5.00
5.00 5.00

---

题解

随机增量法求最小圆覆盖裸题

求法：设初始圆为某空圆，先枚举第一个点，如果不在当前圆内，则令当前圆为这一个点的最小圆覆盖并枚举第二个点，如果不在则变为这两个点的最小圆覆盖并枚举第三个点，如果不在则变为这三个点的最小圆覆盖。

看上去是三重循环，但是实际上时间复杂度为期望$O(n)$的，证明参见 http://blog.csdn.net/lthyxy/article/details/6661250。

需要先将点随机排序以防止被刻意卡掉。

另外求三个点的公共圆时可以直接套用坐标公式，参见代码。

#include <cmath>
#include <cstdio>
#include <algorithm>
#define N 100010
using namespace std;
const double eps = 1e-15;
double x[N] , y[N];
int id[N];
inline double squ(double x)
{
	return x * x;
}
int main()
{
	srand(20011011);
	int n , i , j , k;
	double px = 0 , py = 0 , r = 0 , x1 , x2 , x3 , y1 , y2 , y3;
	scanf("%d" , &n);
	for(i = 1 ; i <= n ; i ++ ) scanf("%lf%lf" , &x[i] , &y[i]) , id[i] = i;
	random_shuffle(id + 1 , id + n + 1);
	for(i = 1 ; i <= n ; i ++ )
	{
		if(squ(px - x[id[i]]) + squ(py - y[id[i]]) > r + eps)
		{
			px = x[id[i]] , py = y[id[i]] , r = 0;
			for(j = 1 ; j < i ; j ++ )
			{
				if(squ(px - x[id[j]]) + squ(py - y[id[j]]) > r + eps)
				{
					px = (x[id[i]] + x[id[j]]) / 2 , py = (y[id[i]] + y[id[j]]) / 2 , r = (squ(x[id[i]] - x[id[j]]) + squ(y[id[i]] - y[id[j]])) / 4;
					for(k = 1 ; k < j ; k ++ )
					{
						if(squ(px - x[id[k]]) + squ(py - y[id[k]]) > r + eps)
						{
							x1 = x[id[i]] , x2 = x[id[j]] , x3 = x[id[k]];
							y1 = y[id[i]] , y2 = y[id[j]] , y3 = y[id[k]];
							px = (x1 * x1 * (y2 - y3) + x2 * x2 * (y3 - y1) + x3 * x3 * (y1 - y2) - (y1 - y2) * (y2 - y3) * (y3 - y1)) / (x1 * (y2 - y3) + x2 * (y3 - y1) + x3 * (y1 - y2)) / 2;
							py = ((x1 - x2) * (x2 - x3) * (x3 - x1) - y1 * y1 * (x2 - x3) - y2 * y2 * (x3 - x1) - y3 * y3 * (x1 - x2)) / (x1 * (y2 - y3) + x2 * (y3 - y1) + x3 * (y1 - y2)) / 2;
							r = squ(px - x1) + squ(py - y1);
						}
					}
				}
			}
		}
	}
	printf("%.15lf
%.15lf %.15lf
" , sqrt(r) , px , py);
	return 0;
}