tel/TPmatlab/integration/TP3g/calcul_poly_gauss_symb.m

   1 function P=calcul_poly_gauss_symb(c,n)\r
   2 \r
   3 %       calcul_poly_gauss_symb : calcul des polynômes de Gauss en symbolique (récurrence)\r
   4 %\r
   5 % *********************************************************\r
   6 %\r
   7 % P=calcul_poly_gauss_symb(c,n) renvoie les polynômes de Gauss P_0,P_1,...,P_n\r
   8 %   pour les quatre méthodes génériques (Legendre, Tchebytchev, Hermite \r
   9 %   et Laguerre; calcul par récurrence.\r
  10 %\r
  11 %       variables d'entrées : \r
  12 %    *  n : entier naturel.\r
  13 %     * c définit la méthode : \r
  14 %            c=1 : Gauss-Legendre,\r
  15 %            c=2 : Gauss-Tchebytchev,\r
  16 %            c=3 : Gauss-Hermite,\r
  17 %            c=4 : Gauss-Laguerre.\r
  18 %       variables de sortie\r
  19 %     * P : tableau de taille (n+1)*(n+1) tel que \r
  20 %           pour tout i dans {1,...,n+1}, L(i,1:i) contient\r
  21 %           le polynôme L_{i-1} (coefficient dominant en tête)\r
  22 %           à coefficients symboliques\r
  23 % \r
  24 % \r
  25 % ************ Fonctions auxiliaires utilisées ************\r
  26 %\r
  27 %       aucune\r
  28 %\r
  29 % *********************************************************\r
  30 %\r
  31 \r
  32 \r
  33 % Contrôles d'entrée\r
  34 % nombre d'arguments\r
  35 if nargin~=2\r
  36    error('nombre d''arguments de la fonction incorrect');\r
  37 end\r
  38 % autres tests éventuels \r
  39 if (fix(n)~=n) | (n<0)\r
  40    error('l''entier n doit être un entier naturel');\r
  41 end\r
  42 if ~((c==1) | (c==2) | (c==3) | (c==4))\r
  43    error('c doit être un entier égal à 1, 2, 3 ou 4');\r
  44 end\r
  45 \r
  46 % Corps d'algorithme.\r
  47 P=zeros(n+1,n+1);\r
  48 P=sym(P);\r
  49 if (c~=2)\r
  50    P(1,1)=1;\r
  51    Lb=[];\r
  52    La=[1];\r
  53    for k=0:n-1      \r
  54       sk=sym(k);\r
  55       switch(c)\r
  56          case 1\r
  57             auxi0=1/(sk+1);\r
  58             auxi1=(2*sk+1)*auxi0;\r
  59             auxi2=-sk*auxi0;\r
  60             Lc=[auxi1*La 0]+[0 0 auxi2*Lb];\r
  61          case 3\r
  62             Lc=[2*La 0]-[0 0 2*sk*Lb];\r
  63          case 4\r
  64             auxi0=1/(sk+1);\r
  65             auxi1=(2*sk+1)*auxi0;\r
  66             auxi2=-sk*auxi0;\r
  67             Lc=[-auxi0*La 0]+[0 auxi1*La]+[0 0 auxi2*Lb];             \r
  68       end\r
  69       Lb=La;\r
  70       La=Lc;\r
  71       P(k+2,1:k+2)=Lc;    \r
  72    end\r
  73 else   \r
  74    P(1,1)=1;\r
  75    if (n>=1)\r
  76       P(2,1:2)=[1 0];\r
  77       Lb=1;\r
  78       La=[1 0];\r
  79       for k=1:n-1      \r
  80          Lc=[2*La 0]-[0 0 Lb];\r
  81          Lb=La;\r
  82          La=Lc;\r
  83          P(k+2,1:k+2)=Lc;    \r
  84       end\r
  85    end\r
  86 end\r