program ep95_316;

{
   Epreuve informatique de l'Ecole Polytechnique         95.316
   -------------------------------------------------------------



On considère la fonction réelle :

                f(x)=a.x^b.exp(-x/g)

où a, b et g sont des paramètres réels positifs.

1/ Discuter l'allure de cette fonction selon la valeur des paramètres.
Montrer qu'elle admet un extremum. Ecrire un programme Pascal qui trace
son graphe pour les valeurs positives de la variable x, et qui cherche
la position de l'extremum.


2/ On connaît la valeur de f(x) pour trois valeurs particulières, différentes
et strictement positives de x :

	f(x1)=y1    f(x2)=y2     f(x3)=y3

Ecrire un programme Pascal qui calcule les trois paramètres a, b et g connaissant x1,
x2, x3 et y1, y2, y3.
(Application numérique suggérée :  x1=1 x2=2 x1=3 et  y1=1 y2=3 y3=2)


3/ Ecrire un progarmme Pascal qui calcule l'intégrale

      SOMME de 0 à u [ x^b.exp(-x).dx]

         pour tout couple (u,b) de R²+

4/ En considérant b comme un paramètre, tracer le graphe de l'application :

      u ->  g(u,b) = SOMME de 0 à u [ x^b.exp(-x).dx]

Discuter son allure en fonction des valeurs de b .
Montrer que lorsque u tend vers l'infini, g(u,b) tend vers une limite, que
l'on notera  Gamma(b).
Ecrire une programme Pascal qui calcule cette limite.
Remarquer qu'elle prend une valeur entière lorsque le réel b est un nombre
entier. Y a-t-il une explication ?


5/ En considérant maintenant  u comme un paramètre, tracer le graphe de
l'application :

              b ->  g(u,b) = SOMME de 0 à u [ x^b.exp(-x).dx]

En faisant tendre u vers l'infini, retrouver le résultat Gamma(b) obtenu à
la question précédente.
Montrer que ce graphe de l'application b -> g(u,b) est une courbe convexe,
dont on calculera le point extremum E.
Tracer le lieu du point E lorsque le paramètre u varie.


:--------------------------------------------------------------:
:    Imprimer tous les résultats                               :
:     en indiquant chaque fois à quoi ils correspondent        :
:--------------------------------------------------------------:
                            -=-=-=-
 }


uses
    modubase, crt, printer;


const
     maxn=10;
type
  Myint=Longint;
  Myreal = Double;
  Distance1=myreal;
  Distance2=myreal;
  Distance3=myreal;
  Distance4=myreal;
  Distance6=myreal;
  Distance8=myreal;
  Distance9=myreal;
  Distance10=myreal;
  Angle=MyReal;

  Points =record
      x,y:Distance1;
      t:string
  end;
var
   ch:char;
   n:integer;
   alpha,beta,gamma:MyReal;


procedure Comparer(s:string;x,y:Myreal);
va
   ecart:real;
const
    eps=0.01;
begin
     Ecris(s+' ');
     ecart:=abs(x-y)/max(eps,max(abs(x),abs(y)));
     if ecart<eps then
      Ecris(' OK.')
      else
       begin
             Ecris('! écart de ');
             EcrisEntier(round(ecart*100));
             Ecris('%');
       end;
end;

function distance(P1,P2:points):Distance1;
begin
     distance:=sqrt(sqr(p1.x-p2.x)+sqr(p1.y-p2.y));
end;

procedure NewPoint(var p:points;t:string;x,y:Distance1);
begin
    p.t:=t;
    p.x:=x;
    p.y:=y;
end;

procedure MontrePoint(p:points);
begin
    Croix(p.x,p.y);
    Ecris(p.t);
end;

procedure MontreSegment(p1,p2:points);
begin
     Deplace(p1.x,p1.y);
     Trace(p2.x,p2.y);
end;

function f(x:MyReal):MyReal;
begi
     f:=alpha*Puissreal(x,beta)*exp(-x/gamma);
end;

procedure extremum(VAR x,y:Myreal);
{ methode :
La dérivée est f'(x)=alpha/gamma*(beta*x^(beta-1)-x^beta)*exp(-x/gamma))
        Elle s'annule quand x=beta}



var
  dx,a,b,c,f1,f2,f3 : MyReal;

begin
     x:=beta;
     y:=f(x);
end;



procedure Tracer;
var
   x,y,x1,y1,dx,xmax:Myreal;
begin
     xmax:=20;
     ModeGraphique;
     Fenetre(-1,xmax,-1,20);
     Couleur(Vert);
     x_axe(0,0,1);
     Y_axe(0,0,1);
     Couleur(Brillant);
     dx:=0.1;
     beta:=0;
     repeat
     beta:=beta+1;
     Deplace(1,15);Ecris('beta=');EcrisReel(beta);
     x:=0;
     y:=f(x);
     repeat
           x1:=x;
           y1:=y;
           x:=x+dx;
           y:=f(x);
           Deplace(x1,y1);
           Trace(x,y);
     until x>xmax;
     Extremum(x,y);
     Croix(x,y);
     Pause;
     Inc(n);
     until false;
end;

procedure Passage(P1,P2,P3:Points);
var
   x,y,x1,y1,dx:Myreal;
   det:MyReal;
   xmax:MyReal
begin
     xmax:=20;
     Fenetre(-1,xmax,-1,20);
     Couleur(Vert);
     x_axe(0,0,1);
     Y_axe(0,0,1);
     Couleur(Brillant);
     with P1 do begin Croix(x,y);Ecris(t);end;
     with P2 do begin Croix(x,y);Ecris(t);end;
     with P3 do begin Croix(x,y);Ecris(t);end;


    det:= (P1.x-P3.x)*Ln(P2.x)
         +(P2.x-P1.x)*Ln(P3.x)
         +(P3.x-P2.x)*Ln(P1.x)
    WriteLN('det=',det:12:6);Pause;
    beta:=((P1.x-P3.x)*Ln(P2.y)
           +(P2.x-P1.x)*Ln(P3.y)
           +(P3.x-P2.x)*Ln(P1.y))/det;
    gamma:=det/
           (Ln(P1.x/P3.x)*Ln(P2.y)
           +Ln(P2.x/P1.x)*Ln(P3.y)
           +Ln(P3.x/P2.x)*Ln(P1.y));
    WriteLN('beta=',beta:12:6);
    WriteLN('gamma=',gamma:12:6);


    alpha:=P1.y*exp(P1.x/gamma-beta*Ln(P1.x));
    WriteLN('(1) alpha=',alpha:12:6);
    alpha:=P2.y*exp(P2.x/gamma-beta*Ln(P2.x));
    WriteLN('(2) alpha=',alpha:12:6);
    alpha:=P3.y*exp(P3.x/gamma-beta*Ln(P3.x));
    WriteLN('(3) alpha=',alpha:12:6);
     dx:=0.1;
     Deplace(1,15);Ecris('beta=');EcrisReel(beta);
     x:=0;
     y:=f(x);
     repeat
           x1:=x;
           y1:=y;
           x:=x+dx;
           y:=f(x);
           Deplace(x1,y1);
           Trace(x,y);
     until x>xmax;
     Extremum(x,y);
     Croix(x,y);
     Pause;
end;

procedure Integrale(umax:Myreal);
var
   x,y,x1,y1,dx:MyReal;
begin
     alpha:=1;
     gamma:=1;
     beta:=0;
     ModeGraphique;
     Fenetre(-1,umax,-1,20);
     Couleur(Vert);
     x_axe(0,0,1);
     Y_axe(0,0,1);
     Couleur(Brillant);
     Deplace(2,18);Ecris('Courbes  u -> g(u,beta)');
   repeat
     Couleur(round(beta)+1);
     dx:=0.01;
     x:=0;
     y:=0;
     repeat
           x1:=x;
           y1:=y;
           x:=x+dx;
           y:=y+f(x)*dx;
           Deplace(x1,y1);
           Trace(x,y);
     until x>umax;
     Deplace(umax-2,y);Ecris('beta=');EcrisEntier(round(beta));
     Ecris(' ginf=');EcrisReel(y);
      beta:=beta+1;
until beta>10;
Pause;
end;

function Gbeta(b,dx,xmax:MyReal):Myreal;
var
   x,y:Myreal;
begin
     Beta:=b;
     x:=0;
     y:=0;
     repeat
           x:=x+dx;
           y:=y+f(x)*dx;
     until x>xmax;
     Gbeta:=y;

end;

function G(u,b,dx:MyReal):Myreal;
var
   x,y:Myreal;
begin
     Beta:=b;
     x:=0;
     y:=0;
     repeat
           x:=x+dx;
           y:=y+f(x)*dx;
     until x>u;
     G:=y;

end;

procedure TraceGamma;
var
   xmax, umax,x,y,x1,y1,dx:Myreal;
   Xextr,Yextr:MyReal;
   Xextr1,Yextr1:MyReal;
begin
     Alpha:=1;
     GAmma:=1;
     ModeGraphique;
     xmax:=2;
     Fenetre(-0.1,xmax,-0.1,1.1);
     Couleur(Vert);
     x_axe(0,0,1);
     Y_axe(0,0,1);
     Couleur(Brillant);
     Croix(1,1);
     Deplace(0.1,1);Ecris('The Euler function')
     Deplace(0.1,0.9);Ecris('beta -> g(u,beta)');
     umax:=0.2;

  repeat
     Couleur(Rouge);
     dx:=0.01;
     x:=0;
     y:=Gbeta(x,0.01,umax);
     Xextr1:=Xextr;
     Yextr1:=Yextr;
     Yextr:=1;
     repeat
           x1:=x;
           y1:=y;
           x:=x+dx;
           y:=Gbeta(x,0.01,umax);
           Deplace(x1,y1);
           Trace(x,y);
           if y<Yextr then
             begin
                  Xextr:=x;
                  Yextr:=y;
                  DEplace(0.1,0.1);
                  Ecris('Extr=(');EcrisReel(Xextr);
                  Ecris(',');EcrisReel(Yextr);Ecris(')');

             end;
     until x>xmax;
     Couleur(Brillant);
     if umax>0.3 then
      begin
        Deplace(Xextr,Yextr);
        Trace(Xextr1,Yextr1);
      end;
     Deplace(1,y);
     Ecris('u=');EcrisReel(umax);
     Ecris('beta=');EcrisReel(umax);
     Ecris(' ginf=');EcrisReel(y);
     umax:=umax+1;
  until false;
     Pause;

end;




procedure Question_1
begin
     n:=3;
     alpha:=1;
     beta :=2;
     gamma:=1;
    Tracer;
end;



procedure Question_2;
begin
    Efface;
end;

procedure Question_3;
VAr
  P1,P2,P3:Points;
begin
     ModeGraphique;
     RAndomize;
     repeat
           Efface;
           NewPoint(P1,'P1',1,1+Random);
           NewPoint(P2,'P2',2,3+Random);
           NewPoint(P3,'P3',3,3+Random);
   {  alpha:=1;
               beta:=1;
               gamma:=1;
     with P1 do begin x:=1;y:=f(x);end;
     with P2 do begin x:=2;y:=f(x);end;
     with P3 do begin x:=3;y:=f(x);end;}

          Passage(P1,P2,P3);
     until false;
end;


procedure Question_4;
begin
     Efface;
     Integrale(20);
end;

procedure Question_5;
begin
     TraceGamma;
end;

procedure Question_6;
var
   x,y,x1,y1,x2,y2,x3,y3:Myreal;
   xold,yold:MyReal;
   a,b:MyReal;
   u:MyReal;
const
    pas=0.1;
    umax=10;

begin
     alpha:=1;
     gamma:=1;
     ModeGraphique;
     Fenetre(-0.1,2,-0.1,1.1);
     Couleur(Vert);
     x_axe(0,0,1);
     Y_axe(0,0,1);
     Couleur(Rouge);
     x:=0;
     y:=Gbeta(x,0.01,umax);
     repeat
           x1:=x;
           y1:=y;
           x:=x+pas;
           y:=Gbeta(x,0.01,umax);
           Deplace(x1,y1)
           Trace(x,y);
     until x>2;


     Couleur(Brillant);
     Croix(1,1);
     Deplace(0.1,1);Ecris('Lieu des extrema');
    u:=0.1;
    repeat
     Deplace(0.1,0.9);Ecris('u=');EcrisReel(u);
     x1:=0;
     x2:=0.5;
     x3:=4;
     y1:=G(u,x1,pas);
     y2:=G(u,x2,pas);
     y3:=G(u,x3,pas);
     repeat

           a:=((y3-y2)/(x3-x2)-(y2-y1)/(x2-x1))/(x3-x1);
           b:=(y2-y1)/(x2-x1)-a*(x1+x2);


           x:=-b/(2*a);
           y:=G(u,x,pas);
           Couleur(jaune);
           Point(x,y);
           if y1<y3 then
             begin
                if x>x2 then begin x3:=x;y3:=y; end
                        else begin x3:=x2;y3:=y2;x2:=x;y2:=y; end;
             end
           else
             begin
                  if x<x2 then begin x1:=x;y1:=y; end
                        else begin x1:=x2;y1:=y2;x2:=x;y2:=y; end;
             end;

    until (x3-x2)<1.E-8;
    Couleur(Brillant);
    Croix(x,y);
    if u>0.3 then
      begin
         Deplace(xold,yold);
         Trace(x,y);
      end;
     Deplace(0.1,0.1);
     Ecris(' x=');EcrisReel(x);
     Ecris(' y=');EcrisReel(y);

      xold:=x;
      yold:=y;
     u:=u+0.01;
    until false

end;


procedure Presentation;
begin;
      Efface;
      WriteLN(' ====== f(x)=x^n.exp--x/alpha) =========');
      WriteLN('                                        ');
      WriteLN('                                        ');
      WriteLN('   (1) Tracé de la courbe f(akpha,beta,gamma)  ');
      WriteLN('   (2) Lieu des extrema de f(alpha,beta,gamma)  ');
      WriteLN('   (3) Passage par trois points  ');
      WriteLN('   (4) Courbes u    -> g(u,beta) ');
      WriteLN('   (5) Courbes beta -> g(u,beta) ');
      WriteLN('   (6) Lieu des extrema ');

      WriteLN('                              ');
      Write('   Tapez votre Choix  :  ');
end;


begi
  Randomize;
  InitGraphique;
    repeat
           ModeTexte;
           Presentation;
          Read(Ch);
          case ch of
               '1':  Question_1;
               '2':  Question_2;
               '3':  Question_3;
               '4':  Question_4;
               '5':  Question_5;
               '6':  Question_6;
               '0' : Halt
          end;
    until false;
end.
