2016 Multi-University Training Contest 1 1011 hdu 5733 四面体内切球

链接：戳这里

tetrahedron

Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 65536/65536 K (Java/Others)

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

思路：

s1,s2,s3,s4表示四面体的四个面的面积，V表示四面体的体积

R=3V/(s1+s2+s3+s4)

球心计算公式：

代码：

/*
S1 -- S4 是四面体的四个面积,内切球的半径为高,这样的4个体积加起来就是总体积
内切球半径的求法 r = 3V / (S1+S2+S3+S4)
算四面体内切球的球心
http://www.docin.com/p-504197705.html?qq-pf-to=pcqq.c2c
*/
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<string>
#include<vector>
#include <ctime>
#include<queue>
#include<set>
#include<map>
#include<stack>
#include<iomanip>
#include<cmath>
#define mst(ss,b) memset((ss),(b),sizeof(ss))
#define maxn 0x3f3f3f3f
#define MAX 1000100
///#pragma comment(linker, "/STACK:102400000,102400000")
typedef long long ll;
typedef unsigned long long ull;
#define INF (1ll<<60)-1
using namespace std;
const double eps=1e-8;
#define zero(x) (((x)>0?(x):-(x))<eps)
struct point3{
    double x,y,z;
    point3(double x=0,double y=0,double z=0):x(x),y(y),z(z){}
};
typedef point3 vec3;
vec3 operator + (vec3 a,vec3 b){
    return vec3(a.x+b.x,a.y+b.y,a.z+b.z);
}
vec3 operator - (vec3 a,vec3 b){
    return vec3(a.x-b.x,a.y-b.y,a.z-b.z);
}
vec3 operator * (vec3 a,double p){
    return vec3(a.x*p,a.y*p,a.z*p);
}
vec3 operator / (vec3 a,double p){
    return vec3(a.x/p,a.y/p,a.z/p);
}
point3 operator * (point3 a,point3 b){
    return point3(a.y*b.z-a.z*b.y,a.z*b.x-a.x*b.z,a.x*b.y-a.y*b.x);
}
struct line3{
    point3 a,b;
};
struct plane3{
    point3 a,b,c;
};
point3 xmult(point3 u,point3 v){
    point3 ret;
    ret.x=u.y*v.z-v.y*u.z;
    ret.y=u.z*v.x-u.x*v.z;
    ret.z=u.x*v.y-u.y*v.x;
    return ret;
}
double dmult(point3 u,point3 v){
    return u.x*v.x+u.y*v.y+u.z*v.z;
}
point3 subt(point3 u,point3 v){
    point3 ret;
    ret.x=u.x-v.x;
    ret.y=u.y-v.y;
    ret.z=u.z-v.z;
    return ret;
}
/*****四面体体积*********/
double volume(point3 a, point3 b, point3 c, point3 d) {
    return fabs(dmult( (b - a) * (c - a) , (d - a) ) ) / 6.0;
}
/*****平面法向量*********/
point3 pvec(point3 s1,point3 s2,point3 s3){
    return xmult(subt(s1,s2),subt(s2,s3));
}
/*****两点距离***********/
double distance(point3 p1,point3 p2){
    return sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y)+(p1.z-p2.z)*(p1.z-p2.z));
}
/*****向量长度***********/
double vlen(point3 p){
    return sqrt(p.x*p.x+p.y*p.y+p.z*p.z);
}
/*****点到平面距离***************/
double ptoplane(point3 p,point3 s1,point3 s2,point3 s3){
    return fabs(dmult(pvec(s1,s2,s3),subt(p,s1)))/vlen(pvec(s1 ,s2,s3));
}
/*****判四点共面*********/
int dots_onplane(point3 a,point3 b,point3 c,point3 d){
    return zero(dmult(pvec(a,b,c),subt(d,a)));
}
/*****判三点共线*********/
int dots_inline(point3 p1,point3 p2,point3 p3){
    return vlen(xmult(subt(p1,p2),subt(p2,p3)))<eps;
}
/*****点到平面的投影*****/
point3 shade_ptoplane(point3 p, point3 a, point3 b, point3 c) {
    point3 nor = (b - a) * (c - a);
    point3 nor0 = ( nor * dmult(nor, p - a) ) / vlen(nor) / vlen(nor);
    return (p - nor0);
}
double Dot(point3 a,point3 b){
    return a.x*b.x+a.y*b.y+a.z*b.z;
}
double Length(point3 a){
    return sqrt(Dot(a,a));
}
point3 Cross(point3 a,point3 b){
    return point3(a.y*b.z-a.z*b.y,a.z*b.x-a.x*b.z,a.x*b.y-a.y*b.x);
}
double area2(point3 a,point3 b,point3 c){
    return Length(Cross(b-a,c-a))*0.5;
}
point3 A,B,C,D,O;
int main(){
    while(scanf("%lf%lf%lf",&A.x,&A.y,&A.z)!=EOF){
        scanf("%lf%lf%lf",&B.x,&B.y,&B.z);
        scanf("%lf%lf%lf",&C.x,&C.y,&C.z);
        scanf("%lf%lf%lf",&D.x,&D.y,&D.z);
        if(dots_onplane(A,B,C,D) || dots_inline(A,B,C)
           || dots_inline(A,C,D) || dots_inline(A,B,D) || dots_inline(B,C,D)){
                printf("O O O O\n");
                continue;
        }
        double s1=area2(B,C,D);
        double s2=area2(A,C,D);
        double s3=area2(A,B,D);
        double s4=area2(A,B,C);
        double R=3.0*volume(A,B,C,D)/(s1+s2+s3+s4);
        double x1=(s1*A.x+s2*B.x+s3*C.x+s4*D.x)/(s1+s2+s3+s4);
        double y1=(s1*A.y+s2*B.y+s3*C.y+s4*D.y)/(s1+s2+s3+s4);
        double z1=(s1*A.z+s2*B.z+s3*C.z+s4*D.z)/(s1+s2+s3+s4);
        printf("%.4f %.4f %.4f %.4f\n",x1,y1,z1,R);
    }
    return 0;
}