tetrahedron

Problem Description

Given four points ABCD, if ABCD is a tetrahedron（四面体）, calculate the inscribed（内接） sphere（球） of ABCD.

Input

Multiple test cases (test cases ≤100 ). Each test cases contains a line of 12 integers [−1e6,1e6] indicate the coordinates of four vertices of ABCD. Input ends by EOF.

Output

Print the coordinate of the center of the sphere and the radius, rounded to 4 decimal places. If there is no such sphere, output “O O O O”.

Sample Input

0 0 0 2 0 0 0 0 2 0 2 0 0 0 0 2 0 0 3 0 0 4 0 0

Sample Output

0.4226 0.4226 0.4226 0.4226 O O O O

思路：

首先用叉积判断是否四点共面，若四点共面，则不能形成四面体。随后，利用四面体体积公式，海伦公式（求三角形面积），四面体内心公式，求出四面体的内切圆半径和内心。

设四个点是A、B、C、D，那么内切球的圆心可以表示为

||||代表欧几里德的范式，本质是求距离。 注意： 在||(b-a)(c-a)||d中，b-a是求两个点b和a的差，返回值ba点={b.x-a.x, b.y-a.y,b.z-a.z}。同理，c-a是求两个点c和a的差，返回值ca点={c.x-a.x, c.y-a.y,c.z-a.z}。(b-a)(c-a)是两个向量ba点和ca点的叉积，返回值baca点= {(b.y-a.y)(c.z-a.z)-(c.y-a.y)(b.z-a.z), (b.z-a.z)(c.x-a.x)-(b.x-a.x)(c.z-a.z), (b.x-a.x)(c.y-a.y)-(b.y-a.y)(c.x-a.x)} ||(b-a)(c-a)||是求向量baca点的长度，返回值=， = sqrt（((b.y-a.y)(c.z-a.z)-(c.y-a.y)(b.z-a.z))^2 +((b.z-a.z)(c.x-a.x)-(b.x-a.x)(c.z-a.z))^2 +( (b.x-a.x)(c.y-a.y)-(b.y-a.y)(c.x-a.x))^2），结果是一个常数e ||(b-a)(c-a)||d是求向量d的倍数乘积，返回值={e*d.x, e*d.y, e*d.z} 还差的函数：常数乘向量，向量加法。

半径r可以表示为

程序代码：

#include <iostream> #include <cstdio> #include <cstdlib> #include <cmath> #include <map> #include <vector> #include <queue> #include <cstring> #include <string> #include <algorithm> using namespace std; typedef long ll; typedef unsigned long ull; #define MM(a,b) memset(a,b,sizeof(a)); #define inf 0x7f7f7f7f #define FOR(i,n) for(int i=1;i<=n;i++) #define CT continue; #define PF printf #define SC scanf const int mod=1000000007; const int N=1e6+100; int n,m,c[N],pre[N],sum[N],a[N],ans[N]; struct Point{ ll x,y,z; void read() { scanf("%lld%lld%lld",&x,&y,&z); } }p[6]; /*两点相减，在||（b-a）*(c-a)||*d中，b-a是求两个点b和a的差，返回值ba点={b.x-a.x,b.y-a.y,b.z-a.z} Point operator-(Point a,Point b) { return (Point){b.x-a.x,b.y-a.y,b.z-a.z}; }*/ Point operator-(Point a,Point b) { Point m; m.x=b.x-a.x; m.y=b.y-a.y; m.z=b.z-a.z; return m; } /*两点叉积，（b-a）*(c-a)，ba，ca的叉积，{a.y*b.z-b.y*a.z,a.z*b.x-a.x*b.z,a.x*b.y-a.y*b.x} Point cross(Point a,Point b) { return (Point){a.y*b.z-b.y*a.z,a.z*b.x-a.x*b.z,a.x*b.y-a.y*b.x}; }*/ Point cross(Point a,Point b) { Point n; n.x=a.y*b.z-b.y*a.z; n.y=a.z*b.x-a.x*b.z; n.z=a.x*b.y-a.y*b.x; return n; } //点到原点的距离，||（b-a）*(c-a)||是求向量baca点的长度，结果是一个常数，||·||欧几里得的范式，本质是求距离 double dis(Point a) { return sqrt(a.x*a.x+a.y*a.y+a.z*a.z); } //两个向量的点积 double dot(Point a,Point b) { return a.x*b.x+a.y*b.y+a.z*b.z; } //判断是否共面，共面dot(m,c-a)为0 double pointtoface(Point c,Point a,Point b,Point d) { Point m=cross(b-a,d-a); return dot(m,c-a)/dis(m); }//三维几何中点到面的距离利用 向量a*向量b=|a|*|b|cos(ang) int main() { long double s[5]; //输入四个点的坐标 while(~scanf("%lld%lld%lld",&p[1].x,&p[1].y,&p[1].z)) { p[2].read(); p[3].read(); p[4].read(); //判断是否共面，（b-a）*(c-a)与d求点积，点积为0则共面 if(dot(cross(p[2]-p[1],p[3]-p[1]),p[4])==0) {printf("O O O O\n"); CT;} double ts=0; //求面积 //(c-b)*(d-b) s[1]=dis(cross(p[3]-p[2],p[4]-p[2]))/2; //(c-a)*(d-a) s[2]=dis(cross(p[3]-p[1],p[4]-p[1]))/2; //(b-a)*(d-a) s[3]=dis(cross(p[2]-p[1],p[4]-p[1]))/2; //(c-a)*(b-a) s[4]=dis(cross(p[3]-p[1],p[2]-p[1]))/2; printf("%f %f %f %f \n",s[1],s[2],s[3],s[4]); //求面积和 for(int i=1;i<=4;i++) ts+=s[i]; //求高 double h=pointtoface(p[3],p[1],p[2],p[4]); double r=fabs(s[3]/ts*h); //求坐标 double x=(s[1]*p[1].x+s[2]*p[2].x+s[3]*p[3].x+s[4]*p[4].x)/ts; double y=(s[1]*p[1].y+s[2]*p[2].y+s[3]*p[3].y+s[4]*p[4].y)/ts; double z=(s[1]*p[1].z+s[2]*p[2].z+s[3]*p[3].z+s[4]*p[4].z)/ts; printf("%.4f %.4f %.4f %.4f %.4f\n",ts,x,y,z,r); } return 0; }