Geometry Problem HDU - 6242 （随机化+计算几何）

lice is interesting in computation geometry problem recently. She found a interesting problem and solved it easily. Now she will give this problem to you :

You are given N distinct points (Xi,Yi) on the two-dimensional plane. Your task is to find a point P and a real number R , such that for at least ⌈N2⌉ given points, their distance to point P is equal to R

.

Input

The first line is the number of test cases.

For each test case, the first line contains one positive number N(1≤N≤105) .

The following N lines describe the points. Each line contains two real numbers Xi and Yi (0≤|Xi|,|Yi|≤103) indicating one give point. It's guaranteed that N

points are distinct.

Output

For each test case, output a single line with three real numbers XP,YP,R , where (XP,YP) is the coordinate of required point P . Three real numbers you output should satisfy 0≤|XP|,|YP|,R≤109 .

It is guaranteed that there exists at least one solution satisfying all conditions. And if there are different solutions, print any one of them. The judge will regard two point's distance as R if it is within an absolute error of 10−3 of R

.

Sample Input

1
7
1 1
1 0
1 -1
0 1
-1 1
0 -1
-1 0

Sample Output

0 0 1

思路：

对于这道题因为一半的点都在一个圆上，我们可以想到随机化，我们展开组合数公式，可以算出每一次不成功的概率为1-1/8 那么有公式1-（1-p)^n 可以得出随机100次差不多就够了，所以每次随机出3个点，特判n<=4的情况，三点确定一圆，可以通过三角形外心算出

ac代码：

#include<bits/stdc++.h>
#define LL long long
#define INF 0x3f3f3f3f
#define INFLL 0x3f3f3f3f3f3f3f
#define lson rt<<1
#define rson rt<<1|1
using namespace std;
const int maxn = 1e5+50;
int T;
int n;
double x,y,r;
struct node
{
    double x,y;
}N[maxn];
int Rand(int L, int R) {//区间内随机数生成函数
    return (LL)rand() * rand() % (R - L + 1) + L;
}
double dis2(int a)
{
    return sqrt((N[a].x-x)*(N[a].x-x)+(N[a].y-y)*(N[a].y-y));
}
void cal(int a,int b,int c)//求外心
{
     double a1 = N[b].x - N[a].x, b1 = N[b].y - N[a].y, c1 = (a1*a1 + b1*b1)/2;
    double a2 = N[c].x - N[a].x, b2 = N[c].y - N[a].y, c2 = (a2*a2 + b2*b2)/2;
    double d = a1 * b2 - a2 * b1;
    x = N[a].x + (c1*b2 - c2*b1)/d,y = N[a].y + (a1*c2 - a2*c1)/d;
    r  = dis2(a);
}
double dis(int a,int b)
{
    return sqrt((N[a].x-N[b].x)*(N[a].x-N[b].x)+(N[a].y-N[b].y)*(N[a].y-N[b].y));
}

int check(int a,int b,int c)
{
    cal(a,b,c);
    int cnt = 0;
    for(int i = 0;i<n;i++){
        if(fabs(dis2(i)-r)<=1e-6){
            cnt++;
        }
    }
    //cout<<cnt<<endl;
    if(cnt>=(n+1)/2)
        return 1;
    else
        return 0;
}
int main()
{
    scanf("%d",&T);
    while(T--)
    {
        scanf("%d",&n);
        for(int i = 0;i<n;i++){
            scanf("%lf%lf",&N[i].x,&N[i].y);
        }
        if(n==1){
            printf("%lf %lf %lf\n",N[0].x,N[0].y-1,1.0);
        }
        else if(n<=4)
        {
            double ansx = (N[0].x+N[1].x)/2;
            double ansy = (N[0].y+N[1].y)/2;
            double ansr =  dis(0,1);
            printf("%lf %lf %lf\n",ansx,ansy,ansr/2);
        }
        else{
            while(1){
            int a = Rand(0,n-1);
            int b = Rand(0,n-1);
            int c = Rand(0,n-1);
            while(a==b){
                b = Rand(0,n-1);
            }
            while(a==c||b==c){
                c = Rand(0,n-1);
            }
            int falg = check(a,b,c);
            if(fabs(x)>1e9||fabs(y)>1e9||fabs(r)>1e9) continue;
            //cout<<a<<' '<<b<<' '<<c<<endl;
            if(falg){
                printf("%lf %lf %lf\n",x,y,r);
                break;
            }
            }
        }
    }
}