{ Центр вписанной в тетраэдр сферы }
type vector = array[1..3] of real;
{-----------------------------------------------------------------}
procedure vpro(a,b: vector; var c: vector);
begin
c[1]:=a[2]*b[3] - a[3]*b[2];
c[2]:=a[3]*b[1] - a[1]*b[3];
c[3]:=a[1]*b[2] - a[2]*b[1];
end;
{-----------------------------------------------------------------}
function spro(a,b: vector): real;
begin spro:=a[1]*b[1] + a[2]*b[2] + a[3]*b[3]; end;
{-----------------------------------------------------------------}
procedure norm(var a: vector);
{ Нормирует вектор a }
var n: real;
begin n:=sqrt(spro(a,a));
a[1]:=a[1]/n; a[2]:=a[2]/n; a[3]:=a[3]/n; end;
{-----------------------------------------------------------------}
procedure vadd(a,b: vector; t: real; var c: vector);
{ c = a + t*b }
begin
c[1]:=a[1]+t*b[1]; c[2]:=a[2]+t*b[2]; c[3]:=a[3]+t*b[3]; end;
{-----------------------------------------------------------------}
procedure nbiss(a,b, c,d: vector; var n: vector);
{ Нормаль биссектрисы двугранного угла, проходящего через ребро AB }
var n1,n2, ab,bc,bd: vector;
begin
vadd(b,a,-1, ab); vadd(c,b,-1, bc); vadd(d,b,-1, bd);
vpro(ab,bc, n1); norm(n1);
vpro(ab,bd, n2); norm(n2);
vadd(n1,n2,1, n); end;
{-----------------------------------------------------------------}
function dist(a,b,c, x: vector): real;
{ Расстояние от точки x до плоскости, содержащей точки a,b,c }
var n, ab,ac,ax: vector;
begin
vadd(b,a,-1, ab); vadd(c,a,-1, ac); vadd(x,a,-1, ax);
vpro(ab,ac, n); norm(n);
dist:=spro(ax,n); end;
{-----------------------------------------------------------------}
{-----------------------------------------------------------------}
function det(a1,a2,a3: vector): real;
var x: vector;
begin vpro(a2,a3, x); det:=spro(a1,x); end;
{-----------------------------------------------------------------}
procedure kramer(a1,a2,a3, f: vector; var x: vector);
{ Метод Крамера: a1,a2,a3 - строки системы, f - правые части }
var d: real; b1,b2,b3: vector;
begin
d:=det(a1,a2,a3);
b1:=a1; b2:=a2; b3:=a3;
b1[1]:=f[1]; b2[1]:=f[2]; b3[1]:=f[3];
x[1]:=det(b1,b2,b3)/d;
b1:=a1; b2:=a2; b3:=a3;
b1[2]:=f[1]; b2[2]:=f[2]; b3[2]:=f[3];
x[2]:=det(b1,b2,b3)/d;
b1:=a1; b2:=a2; b3:=a3;
b1[3]:=f[1]; b2[3]:=f[2]; b3[3]:=f[3];
x[3]:=det(b1,b2,b3)/d;
end;
{-----------------------------------------------------------------}
{-----------------------------------------------------------------}
procedure vrand(var x: vector);
begin x[1]:=random; x[2]:=random; x[3]:=random; end;
{-----------------------------------------------------------------}
{=================================================================}
var a,b,c,d, x, a1,a2,a3, f: vector;
{=================================================================}
begin
vrand(a); vrand(b); vrand(c); vrand(d);
nbiss(a,b,c,d, a1); nbiss(b,c,d,a, a2); nbiss(c,d,a,b, a3);
f[1]:=spro(a1,a); f[2]:=spro(a2,b); f[3]:=spro(a3,c);
kramer(a1,a2,a3, f, x);
writeln('=====================================');
writeln(dist(a,b,c, x):21:15, dist(a,b,c, d):21:15);
writeln(dist(a,b,d, x):21:15, dist(a,b,d, c):21:15);
writeln(dist(a,c,d, x):21:15, dist(a,c,d, b):21:15);
writeln(dist(b,c,d, x):21:15, dist(b,c,d, a):21:15);
readln;
end.