hdu 5733 tetrahedron（2016 Multi-University Training Contest 1计算几何）

题目链接http://acm.hdu.edu.cn/showproblem.php?pid=5733

题意：一道计算几何的题，给出四面体的四个顶点的空间坐标，计算四面体内切球的球心坐标和半径。

一个赤裸裸的公式题

利用海伦公式求出四个面的面积，然后根据公式算出球心(x,y,z)的坐标

之后计算球心到任意一个面的距离，即半径。任取三点求出空间平面方程，然后利用点到平面的距离公式即可算出

最后判断是否存在的情况，其实只需要考虑四个点是否能形成一个四面体即可，特判四点共面的情况，即点到平面的距离为0

#include <iostream>
#include <cstdio>
#include <cmath>
using namespace std;

template<class T> T gcd(T a, T b) {
    return b ? gcd(b, a % b) : a;
}

#define lson rt<<1,l,mid
#define rson rt<<1|1,mid+1,r

typedef long long ll;
const int INF=0x3f3f3f3f;
const int maxn=1e5+10;
int T,n,m;

struct Point{
	double x;
	double y;
	double z;
}p[5];

double getArea(Point p,Point q,Point m){
	double a=sqrt((p.x-q.x)*(p.x-q.x)+(p.y-q.y)*(p.y-q.y)+(p.z-q.z)*(p.z-q.z));
	double b=sqrt((m.x-q.x)*(m.x-q.x)+(m.y-q.y)*(m.y-q.y)+(m.z-q.z)*(m.z-q.z));
	double c=sqrt((p.x-m.x)*(p.x-m.x)+(p.y-m.y)*(p.y-m.y)+(p.z-m.z)*(p.z-m.z));
	double r=(a+b+c)/2;
	return sqrt(r*(r-a)*(r-b)*(r-c));
}

int main(){
#ifndef ONLINE_JUDGE
	freopen("test.in","r",stdin);
	freopen("test.out","w",stdout);
#endif
	while(~scanf("%lf%lf%lf%lf%lf%lf%lf%lf%lf%lf%lf%lf",&p[1].x,&p[1].y,&p[1].z,&p[2].x,&p[2].y,&p[2].z,&p[3].x,&p[3].y,&p[3].z,&p[4].x,&p[4].y,&p[4].z)){
		double a=p[1].y*p[2].z-p[1].y*p[3].z-p[2].y*p[1].z+p[2].y*p[3].z+p[3].y*p[1].z-p[3].y*p[2].z;
		double b=-p[1].x*p[2].z+p[1].x*p[3].z+p[2].x*p[1].z-p[2].x*p[3].z-p[3].x*p[1].z+p[3].x*p[2].z;
		double c=p[1].x*p[2].y-p[1].x*p[3].y-p[2].x*p[1].y+p[2].x*p[3].y+p[3].x*p[1].y-p[3].x*p[2].y;
		double d=-p[1].x*p[2].y*p[3].z+p[1].x*p[3].y*p[2].z+p[2].x*p[1].y*p[3].z-p[2].x*p[3].y*p[1].z-p[3].x*p[1].y*p[2].z+p[3].x*p[2].y*p[1].z;
		double tmp=fabs(a*p[4].x+b*p[4].y+c*p[4].z+d)/sqrt(a*a+b*b+c*c);
		
		if(!(tmp<0||tmp>0)){
			printf("O O O O\n");
			continue;
		}

		double s1=getArea(p[2],p[3],p[4]);
		double s2=getArea(p[1],p[3],p[4]);
		double s3=getArea(p[1],p[2],p[4]);
		double s4=getArea(p[1],p[2],p[3]);
		double x=(s1*p[1].x+s2*p[2].x+s3*p[3].x+s4*p[4].x)/(s1+s2+s3+s4);
		double y=(s1*p[1].y+s2*p[2].y+s3*p[3].y+s4*p[4].y)/(s1+s2+s3+s4);
		double z=(s1*p[1].z+s2*p[2].z+s3*p[3].z+s4*p[4].z)/(s1+s2+s3+s4);

		double r=fabs(a*x+b*y+c*z+d)/sqrt(a*a+b*b+c*c);
		printf("%.4f %.4f %.4f %.4f\n",x,y,z,r);
	}
	return 0;
}