HDU-6242 随机rand()

Problem Description
Alice 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

题目保证有解，若n>=3，则解一定是不共线3点的外接圆
随机选择3点，判断是否符合条件

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=5e4+7;
const int maxp=1e5+7;
const double inf=1e100;
const double eps=1e-8;
const double pi =acos(-1.0);
int sgn(double x){
	if(fabs(x)<eps)return 0;
	return x>0?1:-1;
}
//double hypot(double x,double y){return sqrt(x*x+y*y);}
struct Point{
	double x,y;
	Point(double x=0,double y=0):x(x),y(y){}
	bool operator ==(const Point& b){return !sgn(x-b.x)&&!sgn(y-b.y);}
	bool operator  <(const Point& b)const{return !sgn(x-b.x)?y<b.y:x<b.x;}
	double operator^(const Point& b){return x*b.y-y*b.x;}//叉积
	double operator*(const Point& b){return x*b.x+y*b.y;}//点积
	Point operator +(const Point &b){return Point(x+b.x,y+b.y);}
	Point operator -(const Point& b){return Point(x-b.x,y-b.y);}
	Point operator *(const double &k){return Point(x*k,y*k);}
	Point operator /(const double &k){return Point(x/k,y/k);}
	double len(){return hypot(x,y);}
	double len2(){return x*x+y*y;}
	double rad(Point a,Point b){Point p=*this;return fabs(atan2(fabs((a-p)^(b-p)),(a-p)*(b-p)));}//pa和pb夹角
	double angle(){return atan2(y,x);}//倾斜角
	double angle(Point B){//夹角
		Point A=(*this);
		return acos(A*B/A.len()/B.len());
	}
	double distance(Point b){return hypot(x-b.x,y-b.y);}
	Point trunc(double r){//化为长度为r的向量
		double l=len();
		if(!sgn(l))return *this;
		r/=l;return Point(x*r,y*r);
	}
	Point rotleft(){return Point(-y,x);}//逆时针旋转90 度
	Point rotright(){return Point(y,-x);}//顺时针旋转90 度
	Point rotate(double rad){//rad为弧度 逆时针旋转
		Point A=(*this);
		double c=cos(rad),s=sin(rad);
		return Point(A.x*c-A.y*s,A.x*s+A.y*c);
	}
	Point rotate(Point p,double angle){//绕p逆时针旋转后点
		Point v=((*this)-p).rotate(angle);
		return Point(p.x+v.x,p.y+v.y);
	}
	void input(){scanf("%lf%lf",&x,&y);}
	void show(){printf("(%.2lf %.2lf)\n",x,y);}
};
struct Line{
	Point s,e;
	Line(){};
	Line(Point s,Point e):s(s),e(e){};
	void adjust(){swap(s,e);}
	int relation(Point p){//点和直线关系1在左侧 2在右侧 3在直线上
		int c=sgn((p-s)^(e-s));
		if(c<0)return 1;
		if(c>0)return 2;
		return 3;
	}
	double len(){return s.distance(e);}
	Point crosspoint(Line v){//求两直线的交点//要保证两直线不平行或重合
		double a1=(v.e-v.s)^(s-v.s);
		double a2=(v.e-v.s)^(e-v.s);
		return Point((s.x*a2-e.x*a1)/(a2-a1),(s.y*a2-e.y*a1)/(a2-a1));
	}
	void input(){s.input(),e.input();}
	void show(){s.show(),e.show();}
};
struct circle{
	Point p;//圆心
	double r;//半径
	circle(){}
	circle(Point p,double r):p(p),r(r){}
	circle(Point a,Point b,Point c){//三角形的外接圆//利用两条边的中垂线得到圆心
		Line u=Line((a+b)/2,((a+b)/2)+((b-a).rotleft()));
		Line v=Line((b+c)/2,((b+c)/2)+((c-b).rotleft()));
		p=u.crosspoint(v),r=p.distance(a);
	}
	Point point(double a){return Point(p.x+cos(a)*r,p.y+sin(a)*r);}//根据圆心角算圆上的点
	void input(){p.input();scanf("%lf",&r);}
	void show(){printf("%lf %lf %lf\n",p.x,p.y,r);}
	double area(){return pi*r*r;}
	double circumference(){return 2*pi*r;}//周长
	int relation(Point b){//点和圆的关系//0 圆外//1 圆上//2 圆内
		double dst=b.distance(p);
		if(sgn(dst-r)<0)return 2;
		else if(!sgn(dst-r))return 1;
		return 0;
	}
};
Point P[maxp];
int n;
void solve(){
	if(n==1){
		printf("%lf %lf 1\n",P[0].x+1,P[0].y);
		return;
	}
	if(n<5){
		printf("%lf %lf %lf\n",(P[0].x+P[1].x)/2.0,(P[0].y+P[1].y)/2.0,P[0].distance(P[1])/2.0);
		return;
	}
	while(1){
		int i=rand()%n,j=rand()%n,k=rand()%n;
		if(i==j||i==k||j==k)continue;
		if(!sgn((P[j]-P[i])^(P[k]-P[i])))continue;
		circle O(P[i],P[j],P[k]);
		if(fabs(O.p.x)>1e9||fabs(O.p.y)>1e9||fabs(O.r)>1e9)continue;
		int cnt=0;
		for(int i=0;i<n;++i)cnt+=O.relation(P[i])==1;
		if(cnt<<1>=n){
			O.show();return;
		}
	}
}
int main(){
	int t;
	scanf("%d",&t);
	while(t--){
		srand(time(NULL));
		scanf("%d",&n);
		for(int i=0;i<n;++i)P[i].input();
		solve();
	}
	return 0;
}