2016 Multi-University Training Contest 1 1011 tetrahedron 四面体内切圆

tetrahedron

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

Total Submission(s): 408    Accepted Submission(s): 157


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

 
给四点求 四面体内切圆 
公式：

#include<stdio.h>
#include<string.h>
#include<algorithm>
#include<iostream>
#include<stdlib.h>
#include<math.h>
const double eps=1e-8;
#define zero(x) (((x) > 0 ? (x) : (-x)) < eps)

struct point//点也是向量
{
double x,y,z;
point(double x=0,double y=0,double z=0):x(x),y(y),z(z){}

};
int isline(point p,point s1,point s2) ;

//重载
point operator + (point a,point b)
{
return point(a.x+b.x,a.y+b.y,a.z+b.z);
}

point operator - (point a,point b)
{
return point(a.x-b.x,a.y-b.y,a.z-b.z);
}

point operator * (point a,double b)
{
return point(a.x*b,a.y*b,a.z*b);
}

point operator / (point a,double b)
{
return point(a.x/b,a.y/b,a.z/b);
}
//×乘
point operator * (point a,point b)
{
return point(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 pdp(point a,point b)
{
return a.x*b.x+a.y*b.y+a.z*b.z;
}

//向量长度 AKA 模长
double plen(point a)
{
return sqrt(pdp(a,a));
}
//两点距离
double pplen(point a,point b)
{
return sqrt(fabs(pdp(a-b,a-b)));
}
//
//int main()
//{
//	point a = point(0,0,0);
//	point b = point(0,0,1);
//	point c = point(0,1,0);
//	point d = point(1,0,0);
//	point e = point(1,1,1);
//
//	 printf("%.2f\n",pplen(a,e));
//}

//平面法向量 Plane normal vector   两边叉乘
point  pnv(point a,point b,point c)
{
return (a-b)*(a-c);
}
//int main()
//{
//	point a = point(0,0,0);
//	point b = point(1,0,0);
//	point c = point(0,1,0);
//	point d = point(0,0,1);
//	point e = point(1,1,1);
//	point f = pnv(a,b,c);
//	 printf("%.2f %.2f %.2f\n",f.x,f.y,f.z);
//}

//点到平面的距离  p surface len（法向量）乘点（到平面上任意一点的距离）除以（法向量模长）
double psurlen(point p,point s1,point s2,point s3)
{
point f = pnv(s1,s2,s3);

return fabs((pdp(f,p-s1))/plen(f));

}
//int main()
//{
//	point a = point(0,0,0);
//	point b = point(1,0,0);
//	point c = point(0,1,0);
//	point d = point(0,0,1);
//	point e = point(1,1,1);
//
//	 printf("%.2f\n",psurlen(a,b,c,d)*3);// = 根号3 (b,c,d)这个平面截正方体对角线三分之一
//
//}

// 四点共面  到平面的距离==零
int issur(point p,point s1,point s2,point s3)
{
if(isline(p,s1,s2)||isline(s3,s1,s2)
||isline(p,s1,s3)||isline(p,s3,s2))
return 1;//任意三点共线
return zero(psurlen(p,s1,s2,s3));
}
//int main()
//{
//	point a = point(0,0,0);
//	point b = point(1,0,0);
//	point c = point(0,1,0);
//	point d = point(0,0,1);
//	point e = point(1,1,1);
//	point f = point(1,1,0);
//
//	 printf("%d\n",issur(a,b,c,f));
//
//}

//三点共线 任意两点连线 叉乘为0
int isline(point p,point s1,point s2)
{
return zero(plen((s1-s2)*(s2-p)));
}
// int main()
//{
//	point a = point(0,0,0);
//	point b = point(1,0,0);
//	point c = point(0,1,0);
//	point d = point(2,2,2);
//	point e = point(1,1,1);
//	printf("%.2f\n",pllen(e,d,a));
//	 printf("%d\n",isline(e,d,a));
//
//}

//点到直线距离  p line len（法向量）乘点（到线上任意一点的距离）除以（法向量模长）
double pllen(point p,point s1,point s2)
{
if(isline(p,s1,s2)) //共线为零
return 0;
point f = s1-s2;
if(zero((pdp(f,p-s1))))//s1为投影点的话
return plen(p-s1);
//	printf("%.2f\n",pdp(f,p-s1));
return fabs((plen(f*(p-s1)))/plen(f)); //直线p-s1 叉乘直线s1-s2 除以 直线s1-s2模长

}
// int main()
//{
//	point a = point(0,0,0);
//	point b = point(1,0,0);
//	point c = point(0,1,0);
//	point d = point(0,0,1);
//	point e = point(1,1,1);
//	point f = point(2,2,2);
//
//	 printf("%.2f\n",pllen(d,b,c));
//
//}

//三点面积  square of surface
double S(point p1,point p2,point p3)
{
if(isline(p1,p2,p3))
return 0;
//	 printf("%.2f %.2f\n",pplen(p2,p3),pllen(p1,p2,p3)*2);
return fabs(pplen(p2,p3)*pllen(p1,p2,p3)/2);
}
// int main()
//{
//	point a = point(0,0,0);
//	point b = point(1,0,0);
//	point c = point(0,1,0);
//	point d = point(0,0,1);
//	point e = point(1,1,1);
//	point f = point(2,2,2);
//
//	 printf("%.2f\n",S(d,b,c)*2);
//
//}

//四面体体积 V4
double v4(point p1,point p2,point p3,point p4)
{
if(issur(p1,p2,p3,p4))
return 0;//共面为零
return S(p1,p2,p3)*psurlen(p4,p1,p2,p3)/3;//底乘高/3
return fabs(pdp((p2-p1)*(p3-p1),(p4-p1))/6);//混合积几何意义
}
// int main()
//{
//	point a = point(0,0,0);
//	point b = point(1,0,0);
//	point c = point(0,1,0);
//	point d = point(0,0,1);
//	point e = point(1,1,1);
//	point f = point(2,2,2);
//	printf("%f\n",pdp((b-a)*(c-a),(d-a))/6);
//	 printf("%f\n",v4(d,b,c,a));
//
//}

int main()
{
point a,b,c,d;
while(~scanf("%lf%lf%lf",&a.x,&a.y,&a.z))
{
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);
//	printf("%d\n",issur(a,b,c,d));

if(issur(a,b,c,d))
{
printf("O O O O\n");
continue;
}
double s1 = S(d,b,c);
double s2 = S(a,d,c);
double s3 = S(a,b,d);
double s4 = S(a,b,c);
double sum = s1+s2+s3+s4;
double r = 3*v4(a,b,c,d)/(sum);
double xx = (s1*a.x+s2*b.x+s3*c.x+s4*d.x)/sum;
double yy = (s1*a.y+s2*b.y+s3*c.y+s4*d.y)/sum;
double zz = (s1*a.z+s2*b.z+s3*c.z+s4*d.z)/sum;
printf("%.4f %.4f %.4f %.4f\n",xx,yy,zz,r);
}
}