PROGRAM Fractions_polynomiales;
{
Epreuve informatique de l'Ecole Polytechnique           93.157
--------------------------------------------------------------

1°	Soit Fn l'ensemble des applications f(x) de R dans R qui peuvent
s'exprimer sous la forme du quotient de deux polynômes P et Q à
coefficients réels de degré maximum n¯N :


                   P(x)  a0+a1.x+a2.x²+...+an.x^n
             f(x)=-----=--------------------------
                   Q(x)  b0+b1.x+b2.x²+...+bn.x^n

	Ecrire un programme Pascal qui calcule les coefficients de
P(x) et Q(x) à partir de la connaissance des valeurs yi prises par la
fonction f pour k valeurs différentes xi :

         Comment doit-on choisir la valeur de l'entier k ?

Application numérique : n=2  k=5  f(1)= 1 f(2)=2  f(3)=2  f(4)=1  f(5)=0


	
2°	Soit une application g(t) de R+ dans R possédant la propriété (P)
suivante :

(P)                     g(t)=g(1/t)   pour   t>0

Montrer qu'il existe alors une application f(x) dite "transformée" de g
telle que :

                     g(t)=f(t+1/t)  pour t>0

En supposant que l'on dispose d'un moyen de calculer g(t) pour toute valeur
positive de t, écrire un programme Pascal qui trace le graphe de f.



Application numérique :   g(t)=sin(þLog(t))/Log(t)



3°      En supposant que l'on dispose d'un moyen de calculer g possédant
la propriété (P), écrire un programme Pascal qui donne une expression
analytique approchée de sa transformée f grâce à une fraction f~ qui
prend les mêmes valeurs que f pour k valeurs différentes de x.

               f~(xi)=f(xi)   pour   1<=i<=k


==> Applications numériques :

                                                          x+1/x
    g1(x)=x/(x²+1)   g2(x)=x/(x+1)²     g3=[(x+1)²/(x²+1)]

    g4(x)=xsin(1/x)+(sinx)/x

    g5(x)=Arctg(þx)+Arctg(þ/x)    g6(x)=ch(þLog(x))


4°  Dans quels cas peut-on espérer obtenir une représentation exacte de f
par une fraction appartenant à Fn ?
Sinon comment choisir les valeurs xi pour optimiser l'approximation obtenue ?




+------------------------------------------------------------+
¦              Imprimer tous les résultats                   ¦
¦     en indiquant chaque fois à quoi ils correspondent      ¦
+------------------------------------------------------------+


     }


uses crt, modubase;

Const
     maxn=10;
     zoom=4;
     pas=0.1
     nv='x';    { nom de la variable  }


Type
    Points = record
              x,y:real;
             end;
    Vecteurs = record
              x,y:real;
             end;
    Polynome = record
              a:array[0..maxn] of real;
             end;
    Fraction = record
              P,Q:Polynome;
             end;
    Liste_de_points = array[0..maxn] of points;
    Matrice = array[1..maxn,1..maxn] of real;


VAR
   t0,t         : real;
   Ch           : char;
   n            : integer;
   casG         : integer;
   x1,x2,y1,y2  : real;
   Lambda       : real;


Function cosh(t:real):real;  {cosinus hyperbolique}
begin
     cosh:=(exp(t)+exp(-t))/2;
end;

Procedure Zero_P(Var P:polynome);
Var
   i:integer;
begin
     With P do
     begin
          for i:=0 to maxn do
                 a[i]:=0
     end;
end;

Function degre(P:polynome):integer;
Var
   i:integer;
Const
   eps=1E-10;
begin
     degre:=0;
     for i:=1 to maxn do
         if abs(P.a[i])>eps then
            degre:=i;
end;


Function Evalue(P:polynome;t:real):real;
Var
   u,z:real;
   i:integer;
begin
     z:=0;
     u:=1;
     for i:=0 to degre(p) do
         begin
              z:=z+P.a[i]*u;
              u:=u*t;
         end;
     Evalue:=z;
end;


Function g(t:real):real;
{          g0(t)=sin(Lambda*Ln(t))/Ln(t))
           g1(t)=t/(sqr(t)+1)
           g2(t)=x/(x+1)
           g3(t)=(x2+1)/(x+1)2
           g4(t)=xsin(1/x)+(sinx)/x
           gþ(t)=Arctg(þx)+Arctg(þ/x).
}
begin
     case casg of
     0:   if (t<=0) or (t=1) then
             g:=0
          else
             begin
                  Deplace(0,-0.5);
                  Ecris('t=');
                  Ecrisreel(t);
                  Ecris('  ==>   sin(Lambda*Ln(t))/Ln(t)=');
                  Ecrisreel(sin(Lambda*Ln(t))/Ln(t));
                 {WriteLN('sin(Lambda*Ln(t))/Ln(t)');
                 WriteLN('t=',t:9:5,' Lambda*Ln(t)=',Lambda*Ln(t));
                 WriteLN('sin(Lambda*Ln(t))=',sin(Lambda*Ln(t)));
                 WriteLN('sin(Lambda*Ln(t))/Ln(t)=',sin(Lambda*Ln(t))/Ln(t));
                  }
             g:=sin(Lambda*Ln(t))/Ln(t);
             end;
     1:      g:=t/(sqr(t)+1);
     2:      g:=t/sqr(t+1);
     3:      g:=(sqr(t)+1)/sqr(t+1);
     4:      g:=t*sin(1/t)+sin(t)/t;
     5:      g:=arctan(lambda*t)+arctan(lambda/t);
     6:      g:=cosh(lambda*Ln(t));


     end;  {case}
end;

Function Resoud(t:real):real;

{ résoudre   x+1/x=t            x*x-t*x+1=0

si on pose x=exp(u), alors t=2*ch(u)
d'où   u=arg ch(t/2)
et x=exp(arg ch(t/2))}


begin
     resoud:=(t+sqrt(sqr(t)-4))/2;
end;


Function Montre_P(P:polynome;nv:string):string;
Var
   s,v:string;
   i:integer;
begin
     With P do
     begin
          s:='';
          for i:=0 to maxn do
              if a[i]<>0 then
                 begin
                      if (a[i]<0) then
                         s:=s+'-'
                      else
                          if s<>'' then
                             s:=s+'+';

                      str(abs(a[i]):7:4,v);
                      s:=s+v;
                      if i>0 then
                         begin
                              s:=s+'*'+nv;
                              if i>1 then
                                 s:=s+chr(48+i);
                         end;

                 end;
          Montre_P:=s;
     end;   {with P}
end;


Function Montre_F(F:Fraction;nv:string):string;
begin
     With F do
     begin
          Montre_F:='('+Montre_P(P,nv)+')/('+Montre_P(Q,nv)+')';
     end;
end;

Procedure Genere(Var M:Points;F:Fraction;x:real);
begin
     M.x:=x;
     if Evalue(F.Q,x)<>0 then
        M.y:=Evalue(F.P,x)/Evalue(F.Q,x);
end;



Procedure Joindre(M1,M2:points);
begin
     with M1 do Deplace(x,y);
     with M2 do Trace  (x,y);
end ;


Procedure TracePoint(M:Points);
begin
     With M do Croix(x,y);
end;

Procedure Affiche_mat(mat:matrice;n:integer);
var
   i,j:integer;
begin
     WriteLN;
     for i:=1 to n do
         begin
              for j:=1 to n do
                  write(' ',mat[i,j]:10:4);

              WriteLN;
         end;


end;


Procedure Echange2(var x,y:real);
var
   r:real;
begi
     r:=x;
     x:=y;
     y:=r;
end;

Procedure Echange(i1,i2:integer;var mat:matrice;var vec:vecteur;n:integer);
var
   j:integer;
begin
     for j:=1 to n do
         Echange2(mat[i1,j],mat[i2,j]);
     Echange2(vec[i1],vec[i2]);
end;

Procedure Combine(i1,i2:integer;r:real;var mat:matrice;var vec:vecteur;n:integer);
var
   j:integer;
begin
     for j:=1 to n do
         mat[i1,j]:=mat[i1,j]+r*mat[i2,j];
     vec[i1]:=vec[i1]+r*vec[i2];
end;


Procedure Resoud_systeme(mat:matrice;vec:vecteur;n:integer;var sol:vecteur);
                     { Pivot de Gauss    }
var
   i,i1,j,k,i_piv:integer;
   pivot,r:real;
begin
 {    ModeTexte;}

     for i:=1 to n do
         begin
              i_piv:=i;             ;
              for i1:=i to n do
                  if abs(mat[i1,i])>abs(mat[i_piv,i]) then
                          i_piv:=i1;
              Echange(i,i_piv,mat,vec,n);
              if mat[i,i]=0 then
                 begin
                      Ecris('Pas de pivot, donc b0=0');
                      Pause;
                 end
              else
                  for i1:=1 to n d
                      if i1<>i then
                        combine(i1,i,-mat[i1,i]/mat[i,i],mat,vec,n);
{              Affiche_mat(mat,n);Pause;}
         end;
         if mat[i,i]=0 then
              for i:=1 to n do
                  sol[i]:=vec[i]/mat[i,i];
{     ModeGraphique;                      }
end;



Procedure TraceCourbe(F:fraction;nom:string;c:integer)
Const
     n=1000;

Var
   x:real;
   pas : Real;
   M1,M2: points;
begin
     x:=x1;
     pas:=(x2-x1)/n;
     Couleur(vert);
     Deplace(x1,1);
     Ecris(nom);

     Couleur(c);
     Genere(M2,F,x);
     repeat
           x:=x+pas;
           M1:=M2;
           Genere(M2,F,x);
           Joindre(M1,M2);
     until x>x2;
end;

Procedure TraceG(num:integer;nom:string);
Const
     n=1000;

Var
   x:real;
   pas : Real;
   M1,M2: points;

begin

     casg:=num;
     pas:=(x2)/n;
     x:=pas;
     Couleur(Jaune);
     Deplace(x1,1);
     Ecris(nom);

     Couleur(Jaune);
     M2.x:=x;
     M2.y:=g(x);
     repeat
           x:=x+pas;
           M1:=M2;
           M2.x:=x;
           M2.y:=g(x);
           Joindre(M1,M2);
     until x>x2;

     x:=2;
     pas:=(x2)/n;
     Couleur(Brillant);
     Deplace(x1,2);
     Ecris('Courbe f déduite');
     M2.x:=x;
     M2.y:=g(resoud(x));
     repeat
           x:=x+pas;
           M1:=M2;
           M2.x:=x;
           M2.y:=g(resoud(x));
           Joindre(M1,M2);
     until x>x2;
end;

Procedure Init(hx1,hx2,hy1,hy2:real);
Var
   O:Points;
begin
     x1:=hx1;
     x2:=hx2;
     y1:=hy1;
     y2:=hy2;
     ModeGraphique;
      Fenetre(x1,x2,y1,y2);
      Couleur(Bleu);
      Croix(0,0);
      Ecris('o');
      Deplace(0,0);
      Trace(0,1);
      Ecris('y')
      Deplace(0,0);
      Trace(1,0);
      Ecris('x');
end;


Procedure Trace_f_tilde(MI:Liste_de_points;n:integer);
Var
   i,j:integer;
   mat,c:matrice;
   vec,sol:vecteur;
   F:fraction;
begin
            Couleur(Brillant);


     {      construction de la matrice
         A * (a0 a1 a2 ... an b1 b2 ... bn)= (y0 y1 ...  y2n+1)
         Ai1=1  Ai2=xi  Ai3=xi2 .....Ain+2=-yi*xi  Ai2n+1=-yn*xin
          en supposant b0=1
          }

         for i:=0 to 2*n do
          with MI[i] do
           begin
              Croix(x,y);
              for j:=0 to n do
                  mat[1+i,1+j]:=puiss(x,j);
              for j:=1 to n do
                  mat[1+i,1+n+j]:=-y*puiss(x,j);
              vec[1+i]:=y;
           end;
          With F do
                  begin
                       Zero_P(P);
                       Zero_P(Q);
                       Resoud_systeme(mat,vec,2*n+1,sol);
          { Inv_Matrice(mat,c,2*n+1);Mult_Mat_vec(c,vec,sol,2*n+1,2*n+1);}
                       for j:=0 to n do
                           begin
                                P.a[j]:=sol[1+j];
                                if j>0 then
                                   Q.a[j]:=sol[n+1+j]
                                els
                                    Q.a[0]:=1;   {b0=1}
                           end;
                  end;    { with F}
           Deplace(x1,y2*0.6);Ecris('F(x)='+Montre_P(F.P,'x'));
           Deplace(x1,y2*0.5);Ecris('   / '+Montre_P(F.Q,'x'));
           TraceCourbe(F,'f',Brillant);
           Pause;
end;

Procedure Presentation;
begin;
      ModeTexte;
      Efface;
      WriteLN(' ======= Polynomes inverse-symétriques ========');
      WriteLN('                                               ');
      WriteLN('                                               ');
      WriteLN('                                               ');
      WriteLN('   (1) Calcul de P(x)/Q(x)                     ');
      WriteLN('   (2) Tracé de f connaissant g(t)=sin(Lambda*Ln(t))/Ln(t)  ');
      WriteLN('   (3) Calcul de f~                            ');
      WriteLN('                                               ');
      WriteLN('   Tapez votre choix                           ');
      WriteLN('                                               ');
      WriteLN('                                               ');
end;

Procedure Question1;
Const
     M0:Points=(x:1;y:1);
     M1:Points=(x:2;y:2);
     M2:Points=(x:3;y:2);
     M3:Points=(x:4;y:1);
     M4:Points=(x:5;y:0);
Var
   n:integer;    {degré de P et Q}
   MI:Liste_de_points;
   i:integer;

begin
     Efface;
     Writeln('+---------------------- Question n°1 -------------------------------+');
     Writeln('¦                                                                   ¦');
     Writeln('¦   Calcul de P(x)/Q(x)                                             ¦');
     Writeln('¦                                                                   ¦');
     Writeln('¦      P(xi)-yi*Q(xi)=0                                             ¦');
     Writeln('¦    (a0-yi*b0)+(a1-yi*b1)xi+(a2-yi*b2)xi2+....+(an-yi*bn)xin=0     ¦');
     Writeln('¦                                                                   ¦');
     Writeln('¦ <==>  a0+(a1-yi*b1)xi+(a2-yi*b2)xi2+....+(an-yn*bn)xin=b0*yi      ¦');
     Writeln('¦                                                                   ¦');
     Writeln('¦     A * (a0 a1 a2 ... an b1 b2 ... bn)= (b0y0 b0y1 ...  b0y2n+1)  ¦');
     Writeln('¦                                                                   ¦');
     Writeln('¦     Ai1=1  Ai2=xi  Ai3=xi2 .....Ain+2=-yi*xi  Ai2n+1=-yn*xin      ¦');
     Writeln('¦                                                                   ¦')
     Writeln('¦                                                                   ¦');
     Writeln('¦                                                                   ¦');
     Writeln('+-------------------------------------------------------------------+');
     Pause;
     Randomize;
     MI[0]:=M0;
     MI[1]:=M1;
     MI[2]:=M2;
     MI[3]:=M3;
     MI[4]:=M4;
 {    for i:=0 to 4 do
         begin
              MI[i].y:=2*random;
         end;}


     repeat;
            Efface;
            Init(-zoom,zoom*2,-zoom*2,zoom*2);

            n:=2;
            Couleur(Brillant);
           Trace_f_Tilde(MI,n);
           i:=trunc(random*2*n);
           Mi[i].y:=Mi[i].y+random-0.5;
     until false;
end;

Procedure Question2;
Var
   n:integer;    {degré de P et Q}
   MI:Liste_de_points;
   i:integer;

begin
     Efface;
     Writeln('+---------------------- Question n°2 -------------------------------+');
     Writeln('¦                                                                   ¦');
     Writeln('¦            g(t)=sin(Lambda*Log(t))/Log(t)                         ¦');
     Writeln('¦                                                                   ¦');
     Writeln('+-------------------------------------------------------------------+');
     Write('Lambda=');
     ReadLN(Lambda);
     Init(-1,10,-1.2,1.2);
     TraceG(0,'g(t)=sin(Lambda*Log(t))/Log(t)');

     n:=3;
     for i:=0 to 2*n+1 do
      with MI[i] do
          begin
               x:=2+i;
               y:=g(resoud(x));
          end;
      Couleur(Brillant);
      Trace_f_Tilde(MI,n);
      Pause;
end





Procedure Question3;
Var
   t1,t2:real;
   MI:Liste_de_points;
   i:integer;
   n:integer;


begin

  repeat;
     ModeTexte;
     Efface;
     Writeln('+---------------------- Question n°3 -------------------------------+');
     Writeln('¦                                                                   ¦');
     Writeln('¦      (1)  g1(t)=t/(sqr(t)+1)                                      ¦');
     Writeln('¦      (2)  g2(t)=t/sqr(t+1)                                        ¦');
     Writeln('¦      (3)  g3(t)=(sqr(t)+1)/sqr(t+1)                               ¦');
     Writeln('¦      (4)  g4(t)=t*sin(1/t)+(sin(t))/t                             ¦');
     Writeln('¦      (5)  g5(t)=Arctg(Lambda*t)+Arctg(lambda/t)                   ¦');
     Writeln('¦      (6)  g6(t)=ch(Lambda*Log(t))                                 ¦');
     Writeln('¦                                                                   ¦');
     Writeln('+-------------------------------------------------------------------+');
     Writeln;
     Writeln('Tapez votre choix');
          Ch:= ReadKey;
          Case ch of
          '1':  begin
                     n:=1;
                     Init(-1,10,-3,3);
                     TraceG(1,'g1(t)=t/(sqr(t)+1)');

                end;


          '2':  begin
                     n:=1;
                     Init(-1,10,-3,3);
                     TraceG(2,'g2(t)=t/sqr(t+1)');
                end;

          '3' : begin
                     n:=1;
                     Init(-1,10,-3,3);
                     TraceG(3,'  g3(t)=(sqr(t)+1)/sqr(t+1) ');
                end;

          '4' : begin
                     n:=4;
                     Init(-1,10,-3,3);
                     TraceG(4,' g4(x)=t*sin(1/t)+(sin(t))/t ');
                end;

          '5' : begin
                     n:=4;
                     Write('Lambda=');
                     ReadLN(Lambda);
                     Init(-1,10,-3,3);
                     TraceG(5,'gþ(x)=Arctg(þt)+Arctg(þ/t)');
                end;
          '6' : begin
                     n:=4;
                     Write('Lambda=');
                     ReadLN(Lambda);
                     Init(-0.1,10,-1,10);
                     TraceG(6,'g6(t)=ch(Lambda*Log(t))');
                end;
      end ;  {esac}
     for i:=0 to 2*n do
      with MI[i] do
          begin
               x:=2+i;
               y:=g(resoud(x));
          end;
      Couleur(Brillant);
      Trace_f_Tilde(MI,n);
      pause;
  until Ch='0';

end;



Procedure Question4;
     Var
   t1,t2:real;
   h,k:real;

begin
     Efface;
     Writeln('+---------------------- Question n°4 -------------------------------+');
     Writeln('¦                                                                   ¦');
     Writeln('¦     Polynomes aleatoires                                          ¦');
     Writeln('¦         x=P1(t)/(1+t2)     y= P2(t)/(1+t2)                        ¦');
     Writeln('¦                                                                   ¦');
     Writeln('+-------------------------------------------------------------------+');
     Pause;
     Randomize;

     repeat

     until false;

end;

{ Bloc principal }
begin
  Initgraphique;
  ModeTexte;
  Repeat
     Presentation;
     Ch:= ReadKey;
     case ch of
          '1':  Question1;
          '2':  Question2;
          '3' : Question3;
          '4' : Question4;
     end;
     Pause;
  until false;