2017 CCPC 哈尔滨 M （ hdu 6242） Geometry Problem (随机化 + 三角形外心)

解题思路：看到1/2这个概率可以很容易想到随机化的算法，每次随机一个三元组，表示从n个点中随机选取三个点，然后判断这三个点的外接圆是不是满足题意的点。可以证明的是选取三个点不能构成这么一个圆的概率是（1-1/8） = 7/8 . 那么选取100次还没找到这个圆的概率已经是1e-6，近似为0.所以算法的复杂度得到证明。

Tips： 比赛的时候一直以为是精度问题导致Wrong Answer，没有注意考虑算法细节，其实是没有考虑小于等于4个点的情况，这种情况下即使几个点共线，选取两个点的中点也是满足题意的。

AC代码：

/*
* @Author: wchhlbt
* @Last Modified time: 2017-11-11
*/

#include <bits/stdc++.h>

#define inf 0x3f3f3f3f
#define pb push_back
#define AA first
#define BB second
#define ONES(x) __builtin_popcount(x)
#define _  << " " <<
using namespace std;

typedef pair<int, int> P;
typedef long long ll ;
int dx[4] = {0,0,1,-1};
int dy[4] = {1,-1,0,0};
const double eps =1e-6;
const int mod = 1000000007;
const double PI = acos(-1.0);
inline int read(){ int num;    scanf("%d",&num);   return num;}
const int maxn = 100007;

//判断x的符号 零、正、负
int sgn(double x)
{
    if(fabs(x) < eps)return 0;
    if(x < 0)return -1;
    else return 1;
}

typedef struct Point
{
    double x,y;
    Point() {} ;
    Point(double _x,double _y):x(_x),y(_y) {};
    //两点相减，返回一个向量; 两向量相减，返回向量的差
    Point operator -(const Point &b) const{
        return Point(x - b.x,y - b.y);
    }
    //两向量相加，返回向量的和
    Point operator +(const Point &b) const{
        return Point(x + b.x,y + b.y);
    }
    bool operator < (const Point &a)const{
        if(a.x==x) return y<a.y;
        return x < a.x;
    }
    bool operator ==(const Point a)const{
        return sgn(x-a.x)==0 && sgn(y-a.y)==0;
    }
    //向量叉积，二维向量叉乘的结果是一个垂直于XY平面的向量，只有z轴分量
    //几何意义：以这两个向量为临边的平行四边形的面积（三角形面积的二倍）
    double operator ^(const Point &b)const{
        return x*b.y - y*b.x;
    }
    //向量点积
    double operator *(const Point &b)const{
        return x*b.x + y*b.y;
    }
    Point operator *(const double a)const{
        return Point(x*a, y*a);
    }
    Point operator /(const double a)const{
        return Point(x/a, y/a);
    }
    void transXY(double B){     //绕原点旋转角度B（弧度值），后x,y的变化
        double tx = x,ty = y;
        x = tx*cos(B) - ty*sin(B);
        y = tx*sin(B) + ty*cos(B);
    }
} Vec;
typedef Point Vec;

//两点间距离
double dist(Point a,Point b){
    return sqrt((a-b)*(a-b));
}
//求向量的长度
double length(Vec A){
    return sqrt(A*A);
}

double dot(Point a,Point b){return a*b;}
double cross(Point a,Point b){return a^b;}

Point CircumCenter(Point a,Point b,Point c){  //三角形的外心
    Point cp;
    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;
    cp.x = a.x + (c1*b2-c2*b1)/d;
    cp.y = a.y + (a1*c2-a2*c1)/d;
    return cp;
}

Point p[maxn];

bool check(int n,Point cp,double r)
{
    int cnt = 0;
    for(int i = 0; i<n; i++){
        if(fabs(dist(cp,p[i])-r) < eps)
            cnt++;
    }
    if(cnt>=(n+1)/2)    return true;
    else    return false;
}

int main()
{
    srand(time(NULL));
    int t = read();
    while(t--){
        int n = read();
        for(int i = 0; i<n; i++){
            double x,y;
            scanf("%lf%lf",&x,&y);
            p[i] = Point(x,y);
        }
        if(n==1){
            printf("%.12f %.12f 0\n",p[0].x,p[0].y);
            continue;
        }
        if(n<=4){
            double x = (p[0].x + p[1].x)/2;
            double y = (p[0].y + p[1].y)/2;
            double r = dist(p[0],p[1])/2;
            printf("%.12f %.12f %.12f\n",x,y,r);
            continue;
        }
        while(1){
            int a = ((rand() << 15 | rand()) % n + n) % n;
            int b = ((rand() << 15 | rand()) % n + n) % n;
            int c = ((rand() << 15 | rand()) % n + n) % n;
            while(b==a) b = ((rand() << 15 | rand()) % n + n) % n;
            while(c==a || c==b) c = ((rand() << 15 | rand()) % n + n) % n;
            Point cp = CircumCenter(p[a],p[b],p[c]);
            if(abs(cp.x)>1e9 || abs(cp.y)>1e9)  continue;
            double r = dist(cp,p[a]);
            if( check(n, cp, r) ){
                //cout << cp.x _ cp.y _ r << endl;
                printf("%.12f %.12f %.12f\n",cp.x,cp.y,r);
                break;
            }
        }
    }
    return 0;
}