Program ep95_301;

{
     Ecole Polytechnique  ...........................   95.301





1/ Ecrire un programme Pascal qui calcule le volume d'un tétraèdre dont
trois arêtes concourantes ont respectivement pour mesures  x, y, z tandis
que les trois arêtes restantes forment un triangle équilatéral de côté d.
Soit f(x,y,z) ce volume.

2/ Quelles propriétés géométriques peut-on prévoir a priori pour la surface
d'équation  f(x,y,z)=Cte   ?
Illustrer l'allure de cette surface par une famille de coupes obtenues
en assignant à la variable z une série de valeurs fixées.
Que peut-on dire par exemple des coupes obtenues dans un diagramme (x,y)
lorsque la valeur assignée à z devient très petite ? très grande ?


[FIGURE]

3/ On se place maintenant à l'intérieur d'un plan Euclidien.
Un point M se déplace à l'intérieur d'un triangle équilatéral ABC de côté d.
On considère les distances respectives du point M aux sommets du triangle :

             MA=x      MB=y     MC=z

Montrer qu'il existe entre les quatre quantités x,y,z et d une relation du
type  d=g(x,y,z) où g est une fonction que l'on explicitera.
Ecrire un programme Pascal qui calcule la valeur de d connaissant celles
de x, y et z.
(Application numérique proposée : x=73, y=88, z=95)


4/ Utiliser les résultats précédents pour rechercher systématiquement les
cas de figure où x, y et z et d peuvent tous les quatre s'exprimer exactement
par des nombres entiers.


:--------------------------------------------------------------:
:    Imprimer tous les résultats                               :
:     en indiquant chaque fois à quoi ils correspondent        :
:--------------------------------------------------------------:
                            -=-=-=-
 }

{ Indications de solution :

Question 1


Soit H la projection orthogonale de M dans le triangle équilatéral  ABC
On cherche u, v tels que   AH=uAB+vBC

Pythagore généralisé donne
AH.AB = AM.AB = (AM²+AB²-BM²)/2=(x²+d²-y²)/2
AH.AC = AM.AC = (AM²+AC²-CM²)/2=(x²+d²-z²)/2

Or
AH.AB = (uAB+vBC)AB=d²u+d²v/2
AH.AC = (uAB+vBC)AC=d²u/2+d²v

En combinant :
u=(2d²+2x²+2y²-3z²)/6d²
v=(d²+x²-2y²+z²)/3d²

Donc :
AH²=u²AB²+v²AC²+2uvAB.AC=d²u²+d²v²+d²uv=d²(u²+v²+uv)

MH²=MA²-AH²=x²-d²(u²+v²+uv)=
      = [36d²x²-(2d²+2x²+2y²-3z²)²-4(d²+x²-2y²+z²)²-2(2d²+2x²+2y²-3z²)(d²+x²-
2y²+z²)]/36d²
Le volume est égal à MH.Aire(ABC)/3=


V= d. Racine(d²(x²+y²+z²)+(x²y²+y²z²+z²x²) -(x^4+y^4+z^4) -d^4)/2

----------------------------------------------------------------------
}


uses
    modubase;
var
   a,b,c,d:real;
   x:real;
     Car         : Char ;
     Question    : Char;

type
    Myreal=real;
    points=record
       x,y,z:Myreal;
    end;

var
   PA,PB,PC,PM:Points;

procedure definir(VAR p:points;x,y,z:MyReal);
begin
     p.x:=x;
     p.y:=y;
     p.z:=z;
end;

function distance(A,B:points):MyReal;
begin
     distance:=sqrt(sqr(A.x-B.x)+sqr(A.y-B.y)+sqr(A.z-B.z));
end;

function volume6(AB,AC,AD,BC,CD,BD:Myreal):Myreal;
{calcul du volume d'un tétraèdre ABCD en fonction des 6 arêtes}
var
   rad:MyReal;
begin
     ab:=AB*AB; ac:=AC*AC; ad:=AD*AD; bc:=BC*BC; cd:=CD*CD;bd:=BD*BD;
     rad:=ab*cd*(ac+ad+bc+bd-ab-cd)
         +ac*bd*(ab+ad+bc+cd-ac-bd)
         +ad*bc*(ab+ac+bd+cd-ad-bc)
         -bc*cd*bd-cd*ad*ac-ab*bd*ad-ab*bc*ac;
     if rad<0 then
        Volume6:=-1     { la construction du tétraèdre est impossible }
     else
        Volume6:=sqrt(rad)/12;
end;

function volume4(x,y,z,d:Myreal):Myreal;
begin
     volume4:=volume6(x,y,z,d,d,d);
end;

function volume(x,y,z,d:Myreal):Myreal;
var
    x2,y2,z2,d2:MyReal;
   sigma,tau,q:myreal;
begin
     x2:=x*x;
     y2:=y*y;
     z2:=z*z;
     d2:=d*d;
     sigma:=x2+y2+z2;
     tau:=x2*y2+y2*z2+z2*x2;

     q:=3*tau-sigma*sigma+d2*(sigma-d2);

     if q<0 then
       begin
           { WriteLN('Construction du tétraèdre impossible');}
            volume:=-1;  { valeur impossible }
       end
     else
         volume:=d*sqrt(q)/12;
end;

procedure verification(a,b,c:integer);
var
   q1,q2:Longint;
   a2,b2,c2:Longint;
   det:Longint;
   x:Longint;
begin
     a2:=a*a;
     b2:=b*b;
     c2:=c*c;
     q1:=a2*a2+b2*b2+c2*c2-(a2*b2+b2*c2+c2*a2);;
     q2:=a2+b2+c2;
     det:=q2*q2-4*q1;
     x:=round(sqrt(det));
     Writeln('q1=',q1, '   q2=',q2,' det=',det,' x²=',x*x);
     pause;
end;

function solution (a,b,c:integer):real;
var
   a2,b2,c2:Longint;
   q1,q2:real;
   d1:real;
begin
     a2:=a*a;
     b2:=b*b;
     c2:=c*c;
     q1:=a2*a2+b2*b2+c2*c2-(a2*b2+b2*c2+c2*a2);
     q2:=a2+b2+c2;
     d1:=q2*q2-4*q1;


     if d1>=0 then
        solution:=sqrt((q2+sqrt(d1))/2)
     else
         solution:=-0.5;
end;

function g(x,y,z:MyReal):real;
var
   x2,y2,z2:MyReal;
   q1,q2:real;
   d1:real;

begin
     x2:=x*x;
     y2:=y*y;
     z2:=z*z;
     q1:=x2*x2+y2*y2+z2*z2-(x2*y2+y2*z2+z2*x2);
     q2:=x2+y2+z2;
     d1:=q2*q2-4*q1;

     if d1<0 then
       begin
            WriteLN('Pas de solution');
            g:=0;
       end
     else
        g:=sqrt((q2+sqrt(d1))/2);
end;

function g1(x,y,z:MyReal):real;
var
   x2,y2,z2:MyReal;
   q1,q2:real;
   d1:real;
begin
     x2:=x*x;
     y2:=y*y;
     z2:=z*z;
     q1:=x2*x2+y2*y2+z2*z2-(x2*y2+y2*z2+z2*x2);
     q2:=x2+y2+z2;
     d1:=q2*q2-4*q1;


     if d1<0 then
       begin
            WriteLN('Pas de solution');
            g1:=0;
       end
     else
        g1:=sqrt((q2-sqrt(d1))/2);
end;

procedure VariationVolume(x,y,z:Myreal);
var
   xmin,xmax,ymin,ymax:Myreal;
   d,v:MyReal;
   pas:MyReal;
begin
           Efface;
           xmin:=-0.1;
           xmax:=3;
           ymin:=-0.1;
           ymax:=1;
           Fenetre(xmin,xmax,ymin,ymax);
           Couleur(Vert);
           x_axe(0,0,1);
           y_axe(0,0,1);
           Deplace(0,ymax*0.9);
           Ecris('Volume du Tétraèdre V=f(x,y,z,d) en fonction de d');
           Deplace(0,ymax*0.50);Trace(x,ymax*0.50);Ecris('x=');EcrisReel(x);
           Deplace(0,ymax*0.40);Trace(y,ymax*0.40);Ecris('y=');EcrisReel(y);
           Deplace(0,ymax*0.30);Trace(z,ymax*0.30);Ecris('z=');EcrisReel(z);
           Couleur(Brillant);
           pas:=0.01;
           d:=pas;
           Deplace(d,Volume(x,y,z,d));
           repeat
                 d:=d+pas;
                 v:=Volume(x,y,z,d);
                 if v<0 then
                     Deplace(d,0)
                 else
                     Trace(d,v);
           until d>xmax;
end;

procedure Presentation;
begin;
      Efface;
      WriteLN(' ====== EP95.301 Volume du Tétraèdre  =====');
      WriteLN('                                  ');
      WriteLN('                                                  ');
      WriteLN(' (1)  Calcul V=f(x,y,z,d)     ');
      WriteLN(' (2)  Famille de coupes (x,y) pour V et z constant ');
      WriteLN(' (3)  Relation d=g(x,y,z)   ');
      WriteLN(' (4)  recherche des solutions entières en (x,y,z,d) ');
      WriteLN('                              ');
      Write('   Tapez votre Choix  :  ');
end;

{ Tracé 3d (c) Nicolas Lejeune 1995 ! ;^) }
Const
  COEF = 0.5 ;
Procedure Trace3d (x, y, z : Real) ;
Begin
  Trace (x + z*COEF, y + z*COEF) ;
End ;

Procedure Deplace3d (x, y, z : Real) ;
Begin
  Deplace (x + z*COEF, y + z*COEF) ;
End ;

procedure Question1;
    procedure test(x,y,z,d:MyReal);
    begin
          WriteLN('f(',x:4:1,',',y:4:1,',',z:4:1,',',d:4:1,')=',Volume4(x,y,z,d),'=',volume(x,y,z,d));
    end;
begin
     Efface;
     Writeln('+---------------------- Question n°1 -------------------------------+');
     Writeln('¦                                                                   ¦');
     Writeln('¦      Calcul V=f(x,y,z,d)                                          ¦');
     Writeln('+-------------------------------------------------------------------+');

     Test(1,1,1,1);
     Test(1,1,1,1.5);
     Test(2,2,3,3.4);
     Test(2,2,2,1);
     Test(2,3,4,2);
     Test(2,3,4,3);
     Test(2,3,4,2.5);
     Test(1,2,3,2.5);
     Pause;

     ModeGraphique;
     VariationVolume(1,2,2);

     PAuse;
     repeat
           VariationVolume(random,random,random);
           Pause;
     until false;

end;


procedure Question2;
var
    x,y,z,z0:MyReal;
    theta:MyReal;
    r:Myreal;
    procedure Calcule(theta:real);
    begin
         Definir(PM,PC.x+r*cos(theta),PC.y+r*sin(theta),z);
         x:=distance(PM,PA);
         y:=distance(PM,PB);

    end;
begin
     Efface;
     Writeln('+---------------------- Question n°2 -------------------------------+');
     Writeln('¦                                                                   ¦');
     Writeln('¦      Famille de coupes (x,y) pour V et z constant                 ¦');
     Writeln('+-------------------------------------------------------------------+');
     d:=10;
     Definir(PA,d/2,0,0);
     Definir(PB,-d/2,0,0);
     Definir(PC,0,d*sqrt(3)/2,0);
     ModeGraphique;
     Isofenetre(-1,40,-1);
     X_axe(0,0,10);
     Y_axe(0,0,10);
     z0:=2;
     z:=z0;
     repeat
           Definir(PM,0,0,z0);
           z:=z+0.2;
           r:=sqrt(sqr(z)-sqr(z0));
           couleur(1+random(10));
           theta:=2*pi;
           Calcule(theta);
           Deplace3d(x,y,z);
           repeat
                 theta:=theta-0.01;
                 Calcule(theta);
                 Trace3d(x,y,z);
           until theta<0;
     until z>100;
     pause;
end;

procedure Question3;
var
   x,y,z:MyReal;
begin
      Efface;
     Writeln('+---------------------- Question n°3 -------------------------------+');
     Writeln('¦                                                                   ¦');
     Writeln('¦      Relation d=g(x,y,z)                                          ¦');
     Writeln('+-------------------------------------------------------------------+');
     WriteLN;
     repeat
           WriteLN('Entrer x,y,z=');
           Readln(x,y,z);
           WriteLN('Le côté du triangle équilatéral est =',g(x,y,z));
           WriteLN(' ( Solution extérieure non demandée :',g1(x,y,z),')');

     until false;
end;

procedure Question4;
var
   a,b,c,s:integer;
begin
      Efface;
     Writeln('+---------------------- Question n°4 -------------------------------+');
     Writeln('¦                                                                   ¦');
     Writeln('¦      recherche des solutions entières en (x,y,z,d)                ¦');
     Writeln('+-------------------------------------------------------------------+');
     WriteLN;
     verification(57,65,73);
     verification(73,88,95);
     verification(43,147,152);
     verification(69,115,118);
     x:=solution(57,73,65);
     a:=57;
     b:=73;
     c:=65;
     s:=190;
     repeat
           inc(s);
           write('.',s);
           a:=0;
           repeat
                 Inc(a);
                 b:=a;
                 repeat
                       Inc(b);
                        c:=s-a-b;
                        x:=solution(a,b,c);
                        if abs(x-round(x))<1E-06 then
                          if (a+b-x)>0.5 then
                           if x>c then
                             begin
                                  WriteLN;
                                  writeLN('a=',a:8,' b=',b:8,' c=',c:8,' -> =',solution(a,b,c):15:9);
                             end;
                 until c<b;
           until 3*a>s;
     until false;
     readln;
end;



begin
  Initgraphique;
  Randomize;
  while true do
  begin
     Modetexte;
     Presentation;
     write('Question choisie ? ');
     readln(question);
     case question of
          '1':  Question1;
          '2':  Question2;
          '3' : Question3;
          '4' : Question4;
     end;
     Pause;
  end;
end.