http://acm.split.hdu.edu.cn/showproblem.php?pid=5733
题意:给你4个定点,这个4个定点构成一个四边形,求这个四边形的内切圆的圆心和半径
没有的话输出O O O O
思路:混合积可以求四边形的体积
r = v*3/(s1+s2+s3+s4)
s1,s2,s3,s4分别是四个面的面积
然后圆心
x = a.x*s1+b.x*s2+c.x*s3+d.x*s4;
y = a.y*s1+b.y*s2+c.y*s3+d.y*s4;
z = a.z*s1+b.z*s2+c.z*s3+d.z*s4;
s1 是a点不构成的面积
s2 是b点不构成的面积
1 #include <stdio.h> 2 #include <string.h> 3 #include <math.h> 4 #include <iostream> 5 using namespace std; 6 7 struct point3 8 { 9 double x,y,z; 10 }; 11 12 double lenth(point3 u){ 13 return sqrt(u.x*u.x+u.y*u.y+u.z*u.z); 14 } 15 point3 xmult(point3 u,point3 v) 16 { 17 point3 ret; 18 ret.x=u.y*v.z-v.y*u.z; 19 ret.y=u.z*v.x-u.x*v.z; 20 ret.z=u.x*v.y-u.y*v.x; 21 return ret; 22 } 23 24 double dmult(point3 u,point3 v) 25 { 26 return u.x*v.x+u.y*v.y+u.z*v.z; 27 } 28 29 30 point3 operator - ( point3 A , point3 B ) 31 { 32 point3 ret; 33 ret.x = A.x-B.x; 34 ret.y = A.y-B.y; 35 ret.z = A.z-B.z; 36 return ret; 37 } 38 double area( point3 A , point3 B , point3 C ) 39 { 40 return lenth(xmult(B-A,C-A))/2 ; 41 } 42 43 double volume( point3 A , point3 B , point3 C , point3 D ) 44 { 45 // printf("%lf ",xmult(B-A,C-A)); 46 //printf("%lf ",dmult(D-A,xmult(B-A,C-A))); 47 return fabs(dmult(D-A,xmult(B-A,C-A)))/6 ; 48 } 49 50 51 52 int main() 53 { 54 55 point3 a,b,c,d; 56 while(~scanf("%lf%lf%lf%lf%lf%lf%lf%lf%lf%lf%lf%lf",&a.x,&a.y,&a.z,&b.x,&b.y,&b.z,&c.x,&c.y,&c.z,&d.x,&d.y,&d.z)) 57 { 58 double s1 = area(b,c,d); 59 double s2 = area(a,c,d); 60 double s3 = area(a,b,d); 61 double s4 = area(a,b,c); 62 double v = volume(a,b,c,d); 63 if(v<1e-8) 64 { 65 printf("O O O O "); 66 continue; 67 } 68 double s = s1+s2+s3+s4; 69 double r = v*3/s; 70 double x = a.x*s1+b.x*s2+c.x*s3+d.x*s4; 71 double y = a.y*s1+b.y*s2+c.y*s3+d.y*s4; 72 double z = a.z*s1+b.z*s2+c.z*s3+d.z*s4; 73 printf("%.4lf %.4lf %.4lf %.4lf ",x/s,y/s,z/s,r); 74 } 75 return 0; 76 }