scilab/modules/differential_equations/help/fr_FR/ode.xml

   1 <?xml version="1.0" encoding="UTF-8"?>
   2 <refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:id="ode" xml:lang="fr">
   3     <refnamediv>
   4         <refname>ode</refname>
   5         <refpurpose>solveur d'équations différentielles ordinaires</refpurpose>
   6     </refnamediv>
   7     <refsynopsisdiv>
   8         <title>Séquence d'appel</title>
   9         <synopsis>y=ode(y0,t0,t,f)
  10             [y,w,iw]=ode([type],y0,t0,t [,rtol [,atol]],f [,jac] [,w,iw])
  11             [y,rd,w,iw]=ode("root",y0,t0,t [,rtol [,atol]],f [,jac],ng,g [,w,iw])
  12             y=ode("discrete",y0,k0,kvect,f)
  13         </synopsis>
  14     </refsynopsisdiv>
  15     <refsection>
  16         <title>Paramètres</title>
  17         <variablelist>
  18             <varlistentry>
  19                 <term>y0</term>
  20                 <listitem>
  21                     <para>vecteur ou matrice réelle (conditions initiales).</para>
  22                 </listitem>
  23             </varlistentry>
  24             <varlistentry>
  25                 <term>t0</term>
  26                 <listitem>
  27                     <para>réel (instant initial).</para>
  28                 </listitem>
  29             </varlistentry>
  30             <varlistentry>
  31                 <term>t</term>
  32                 <listitem>
  33                     <para>vecteur réel (instants où la solution est renvoyée).</para>
  34                 </listitem>
  35             </varlistentry>
  36             <varlistentry>
  37                 <term>f</term>
  38                 <listitem>
  39                     <para>fonction externe (fonction Scilab ou chaîne de caractères ou
  40                         liste).
  41                     </para>
  42                 </listitem>
  43             </varlistentry>
  44             <varlistentry>
  45                 <term>type</term>
  46                 <listitem>
  47                     <para>
  48                         une des chaînes de caractères : <literal>"adams"</literal>,
  49                         <literal>"stiff"</literal>, <literal>"rk"</literal>,
  50                         <literal>"rkf"</literal>, <literal>"fix"</literal>,
  51                         <literal>"discrete"</literal>, <literal>"roots"</literal>
  52                     </para>
  53                 </listitem>
  54             </varlistentry>
  55             <varlistentry>
  56                 <term>rtol, atol</term>
  57                 <listitem>
  58                     <para>constantes ou vecteurs réels de même taille que
  59                         <literal>y</literal>.
  60                     </para>
  61                 </listitem>
  62             </varlistentry>
  63             <varlistentry>
  64                 <term>jac</term>
  65                 <listitem>
  66                     <para>fonction externe (fonction Scilab, chaîne de caractères ou
  67                         liste).
  68                     </para>
  69                 </listitem>
  70             </varlistentry>
  71             <varlistentry>
  72                 <term>w, iw</term>
  73                 <listitem>
  74                     <para>vecteurs réels</para>
  75                 </listitem>
  76             </varlistentry>
  77             <varlistentry>
  78                 <term>ng</term>
  79                 <listitem>
  80                     <para>entier</para>
  81                 </listitem>
  82             </varlistentry>
  83             <varlistentry>
  84                 <term>g</term>
  85                 <listitem>
  86                     <para>fonction externe (fonction Scilab, chaîne de caractères ou
  87                         liste).
  88                     </para>
  89                 </listitem>
  90             </varlistentry>
  91             <varlistentry>
  92                 <term>k0</term>
  93                 <listitem>
  94                     <para>entier (instant initial).</para>
  95                 </listitem>
  96             </varlistentry>
  97             <varlistentry>
  98                 <term>kvect</term>
  99                 <listitem>
 100                     <para>vecteur d'entiers.</para>
 101                 </listitem>
 102             </varlistentry>
 103         </variablelist>
 104     </refsection>
 105     <refsection>
 106         <title>Description</title>
 107         <para>
 108             <literal>ode</literal> est la fonction utilisée pour approcher la
 109             solution d'une équation différentielle ordinaire (EDO) explicite du
 110             premier ordre en temps, définie par : dy/dt=f(t,y) , y(t0)=y0. Il s'agit
 111             d'une interface vers diverses bibliothèques, en particulier ODEPACK. Le type
 112             du problème et la méthode utilisée dépendent de la valeur du premier
 113             argument optionnel <literal>type</literal> qui peut être égal à :
 114         </para>
 115         <variablelist>
 116             <varlistentry>
 117                 <term>&lt;aucun&gt; :</term>
 118                 <listitem>
 119                     <para>
 120                         le solveur <literal>lsoda</literal> du package ODEPACK est
 121                         utilisé par défaut. Il choisit automatiquement entre un schéma
 122                         prédicteur-correcteur d'Adams et un schéma adapté au systèmes raides
 123                         (stiff) de type "Backward Differentiation Formula" (BDF).
 124                         Initialement le schéma adapté aux système non raides est choisi puis
 125                         la méthode adaptée est ensuite choisie dynamiquement.
 126                     </para>
 127                 </listitem>
 128             </varlistentry>
 129             <varlistentry>
 130                 <term>"adams" :</term>
 131                 <listitem>
 132                     <para>
 133                         Problèmes non raides. Le solveur <literal>lsode</literal> du
 134                         package ODEPACK est utilisé (schéma d'Adams).
 135                     </para>
 136                 </listitem>
 137             </varlistentry>
 138             <varlistentry>
 139                 <term>"stiff" :</term>
 140                 <listitem>
 141                     <para>
 142                         Pour les systèmes raides. Le solveur <literal>lsode</literal>
 143                         du package ODEPACK est utilisé avec le schéma BDF.
 144                     </para>
 145                 </listitem>
 146             </varlistentry>
 147             <varlistentry>
 148                 <term>"rk" :</term>
 149                 <listitem>
 150                     <para>Schéma de Runge-Kutta adaptatif d'ordre 4 (RK4).</para>
 151                 </listitem>
 152             </varlistentry>
 153             <varlistentry>
 154                 <term>"rkf" :</term>
 155                 <listitem>
 156                     <para>Formules de Shampine et Watts basées sur les paires de
 157                         Runge-Kutta Fehlberg d'ordre 4 et 5 (RKF45). Bien pour les problèmes
 158                         non raides ou moyennement raides, lorsque le calcul du second membre
 159                         n'est pas trop coûteux. Cette méthode est à éviter si l'on recherche
 160                         une très grande précision.
 161                     </para>
 162                 </listitem>
 163             </varlistentry>
 164             <varlistentry>
 165                 <term>"fix":</term>
 166                 <listitem>
 167                     <para>
 168                         Identique à <literal>"rkf"</literal>, mais l'interface est
 169                         simplifiée, i.e.
 170                         uniquement <literal>rtol</literal> et <literal>atol</literal> sont
 171                         communiqués au solveur.
 172                     </para>
 173                 </listitem>
 174             </varlistentry>
 175             <varlistentry>
 176                 <term>"root":</term>
 177                 <listitem>
 178                     <para>Solveur d'EDO avec recherche de racines. Le solveur
 179                         <literal>lsodar</literal> du package ODEPACK est utilisé. C'est une
 180                         variante de <literal>lsoda</literal> permettant la recherche d'une
 181                         racine d'une fonction vectorielle donnée. Voir <link linkend="ode_root">ode_root</link> pour plus de
 182                         détails.
 183                     </para>
 184                 </listitem>
 185             </varlistentry>
 186             <varlistentry>
 187                 <term>"discrete":</term>
 188                 <listitem>
 189                     <para>
 190                         Simulation en temps discret. Voir <link linkend="ode_discrete">ode_discrete</link> pour plus de
 191                         détails.
 192                     </para>
 193                 </listitem>
 194             </varlistentry>
 195         </variablelist>
 196         <para>
 197             Ici on ne décrit l'usage de <literal>ode</literal> que pour des EDO
 198             explicites.
 199         </para>
 200         <itemizedlist>
 201             <listitem>
 202                 <para>
 203                     L'appel le plus simple de <literal>ode</literal> est du type :
 204                     <literal>y=ode(y0,t0,t,f)</literal> où <literal>y0</literal> est le
 205                     vecteur des conditions initiales, <literal>t0</literal> est le temps
 206                     initial, et <literal>t</literal> est le vecteur des instants où l'on
 207                     veut une approximation de la solution. <literal>y</literal> est
 208                     calculée et <literal>y</literal> est la matrice
 209                     <literal>y=[y(t(1)),y(t(2)),...]</literal>.
 210                 </para>
 211                 <para>
 212                     Le paramètre d'entrée <literal>f</literal> de
 213                     <literal>ode</literal> définit le membre de droite de l'équation
 214                     différentielle du premier ordre dy/dt=f(t,y). C'est un external qui
 215                     peut être :
 216                 </para>
 217                 <itemizedlist>
 218                     <listitem>
 219                         <para>
 220                             Soit une fonction Scilab, sa syntaxe doit être <literal>ydot
 221                                 = f(t,y)
 222                             </literal>
 223                             où <literal>t</literal> est un scalaire (le
 224                             temps), <literal>y</literal> un vecteur (l'état). Cette fonction
 225                             renvoie le second membre de l'équation différentielle
 226                             dy/dt=f(t,y).
 227                         </para>
 228                     </listitem>
 229                     <listitem>
 230                         <para>Soit une chaîne de caractères, elle désigne le nom d'une
 231                             subroutine Fortran ou une procédure C, i.e. si
 232                             <literal>ode(y0,t0,t,"fex")</literal> est la commande, alors la
 233                             procedure <literal>fex</literal> est appelée.
 234                         </para>
 235                         <para>Si c'est une subroutine Fortran, sa liste d'appel doit
 236                             être
 237                         </para>
 238                         <programlisting role=""><![CDATA[
 239 subroutine fex(n,t,y,ydot)
 240 integer n
 241 double precision t,y(*),ydot(*)
 242  ]]></programlisting>
 243                         <para>Si c'est une fonction C son prototype doit être:</para>
 244                         <programlisting role=""><![CDATA[
 245 void fex(int *n,double *t,double *y,double *ydot)
 246  ]]></programlisting>
 247                         <para>
 248                             Cet external peut être compilé par l'utilitaire <link linkend="ilib_for_link">ilib_for_link</link> et chargé
 249                             dynamiquement par la fonction <link linkend="link">link</link>.
 250                         </para>
 251                     </listitem>
 252                     <listitem>
 253                         <para>Soit une liste avec la structure suivante
 254                             <literal>list(vrai_f,u1,u2,...un)</literal> où
 255                             <literal>vrai_f</literal> est une fonction avec la syntaxe
 256                             <literal>ydot = vrai_f(t,y,u1,u2,...,un)</literal>
 257                         </para>
 258                         <para>Cette syntaxe permet de passer des paramètres sous forme
 259                             d'arguments supplémentaires de <literal>vrai_f</literal>.
 260                         </para>
 261                     </listitem>
 262                 </itemizedlist>
 263                 <para>
 264                     La fonction <literal>f</literal> peut renvoyer une matrice
 265                     <literal>p x q</literal> au lieu d'un vecteur. Dans ce cas, on résout
 266                     le système d'EDO <literal>n=p+q</literal>
 267                     <literal>dY/dt=F(t,Y)</literal> où <literal>Y</literal> est une
 268                     matrice <literal>p x q</literal>. La condition initiale
 269                     <literal>Y0</literal> doit aussi être une matrice <literal>p x
 270                         q
 271                     </literal>
 272                     matrix et le résultat renvoyé par <literal>ode</literal>
 273                     est la matrice: <literal>p x q(T+1)</literal> égale à
 274                     <literal>[Y(t_0),Y(t_1),...,Y(t_T)]</literal>.
 275                 </para>
 276             </listitem>
 277             <listitem>
 278                 <para>Des paramètres optionnels contrôlent la tolérance du schéma :
 279                     <literal>rtol</literal> et <literal>atol</literal> sont des valeurs
 280                     seuil sur les erreurs estimées (relative et absolue) L'erreur estimée
 281                     sur <literal>y(i)</literal> est
 282                     <literal>rtol(i)*abs(y(i))+atol(i)</literal>
 283                 </para>
 284                 <para>
 285                     Si <literal>rtol</literal> et/ou <literal>atol</literal> sont
 286                     des constantes <literal>rtol(i)</literal> et/ou
 287                     <literal>atol(i)</literal> prennent ces valeurs. Les valeurs par
 288                     défaut de <literal>rtol</literal> et <literal>atol</literal> sont
 289                     respectivement <literal>rtol=1.d-5</literal> et
 290                     <literal>atol=1.d-7</literal> pour la plupart des solveurs et
 291                     <literal>rtol=1.d-3</literal> et <literal>atol=1.d-4</literal> pour
 292                     <literal>"rfk"</literal> et <literal>"fix"</literal>.
 293                 </para>
 294             </listitem>
 295             <listitem>
 296                 <para>Pour les problèmes raides, il est recommandé de fournir la
 297                     jacobienne du second membre sous forme de l'argument optionnel
 298                     <literal>jac</literal>. Le paramètre <literal>jac</literal> de
 299                     <literal>ode</literal> est par exemple une fonction Scilab, dont la
 300                     syntaxe est imposée, ou le nom d'une subroutine Fortran ou C (chaîne
 301                     de caractères) ou une liste.
 302                 </para>
 303                 <para>
 304                     Si <literal>jac</literal> est une fonction Scilab sa syntaxe
 305                     doit être <literal>J=jac(t,y)</literal>
 306                 </para>
 307                 <para>
 308                     où <literal>t</literal> est un scalaire (le temps) et
 309                     <literal>y</literal> un vecteur (l'état). La matrice
 310                     <literal>J</literal> doit renvoyer df/dx i.e. <literal>J(k,i) = dfk
 311                         /dxi
 312                     </literal>
 313                     avec <literal>fk</literal> = k-ième composante de
 314                     f.
 315                 </para>
 316                 <para>
 317                     Si <literal>f</literal> est une chaîne de caractères, elle
 318                     désigne le nom d'une subroutine Fortran ou C.
 319                 </para>
 320                 <para>En Fortran, Cette routine doit avoir la liste d'appel suivante
 321                     :
 322                 </para>
 323                 <programlisting role=""><![CDATA[
 324 subroutine fex(n,t,y,ml,mu,J,nrpd)
 325 integer n,ml,mu,nrpd
 326 double precision t,y(*),J(*)
 327  ]]></programlisting>
 328                 <para>Si c'est une fonction C son prototype doit être:</para>
 329                 <programlisting role=""><![CDATA[
 330 void fex(int *n,double *t,double *y,int *ml,int *mu,double *J,int *nrpd,)
 331  ]]></programlisting>
 332                 <para>Dans la plupart des cas il n'est pas nécessaire d'utiliser
 333                     <literal>ml</literal>, <literal>mu</literal> et
 334                     <literal>nrpd</literal>, qui sont relatifs à la possibilité de
 335                     stockage "bande" du Jacobien
 336                 </para>
 337                 <para>
 338                     Si <literal>jac</literal> est une liste, les mêmes conventions
 339                     que pour <literal>f</literal> s'appliquent.
 340                 </para>
 341             </listitem>
 342             <listitem>
 343                 <para>
 344                     Les arguments optionnels <literal>w</literal> et
 345                     <literal>iw</literal> sont des vecteurs dans lesquels le solveur stocke des
 346                     informations sur son état (voir <link linkend="ode_optional_output">ode_optional_output</link> pour plus de
 347                     détails) . Lorsque ces paramêtres sont utilisés comme argument
 348                     d'entrée, ils permettent de redémarrer l'intégration au point où elle
 349                     s'était arrêtée à la sortie d'un appel précédent à
 350                     <literal>ode</literal>.
 351                 </para>
 352             </listitem>
 353             <listitem>
 354                 <para>Plus d'options peuvent être passées aux solveurs d'ODEPACK en
 355                     utilisant la variable <literal>%ODEOPTIONS</literal>. Voir <link linkend="odeoptions">odeoptions</link>.
 356                 </para>
 357             </listitem>
 358         </itemizedlist>
 359     </refsection>
 360     <refsection>
 361         <title>Exemples</title>
 362         <programlisting role="example"><![CDATA[
 363 // ---------- EDO Simple  (external : fonction Scilab)
 364 // dy/dt=y^2-y sin(t)+cos(t), y(0)=0
 365 function ydot=f(t,y),ydot=y^2-y*sin(t)+cos(t),endfunction
 366 y0=0;t0=0;t=0:0.1:%pi;
 367 y=ode(y0,t0,t,f)
 368 plot(t,y)
 369
 370 // ---------- EDO Simple  (external : code C)
 371 ccode=['#include <math.h>'
 372        'void myode(int *n,double *t,double *y,double *ydot)'
 373        '{'
 374        '  ydot[0]=y[0]*y[0]-y[0]*sin(*t)+cos(*t);'
 375        '}']
 376 mputl(ccode,TMPDIR+'/myode.c') //create the C file
 377 ilib_for_link('myode','myode.c',[],'c',[],TMPDIR+'/loader.sce');//compile
 378 exec(TMPDIR+'/loader.sce') //incremental linking
 379 y0=0;t0=0;t=0:0.1:%pi;
 380 y=ode(y0,t0,t,'myode');
 381
 382 // ---------- Simulation de dx/dt = A x(t) + B u(t) avec u(t)=sin(omega*t),
 383 // x0=[1;0]
 384 // solution x(t) desired at t=0.1, 0.2, 0.5 ,1.
 385 // A and u function are passed to RHS function in a list.
 386 // B and omega are passed as global variables
 387 function xdot=linear(t,x,A,u),xdot=A*x+B*u(t),endfunction
 388 function ut=u(t),ut=sin(omega*t),endfunction
 389 A=[1 1;0 2];B=[1;1];omega=5;
 390 ode([1;0],0,[0.1,0.2,0.5,1],list(linear,A,u))
 391
 392 // ----------Integration de l'équation différentielle de Riccati (état matriciel)
 393 // Xdot=A'*X + X*A - X'*B*X + C , X(0)=Identity
 394 // Solution at t=[1,2]
 395 function Xdot=ric(t,X),Xdot=A'*X+X*A-X'*B*X+C,endfunction
 396 A=[1,1;0,2]; B=[1,0;0,1]; C=[1,0;0,1];
 397 t0=0;t=0:0.1:%pi;
 398 X=ode(eye(A),0,t,ric)
 399
 400 // ---------- Calcul de exp(A)  (état matriciel)
 401 A=[1,1;0,2];
 402 function xdot=f(t,x),xdot=A*x;,endfunction
 403 ode(eye(A),0,1,f)
 404 ode("adams",eye(A),0,1,f)
 405
 406 // ----------  Calcul de exp(A)  (état matriciel, cas raide, jacobien fourni)
 407 A=[10,0;0,-1];
 408 function xdot=f(t,x),xdot=A*x,endfunction
 409 function J=Jacobian(t,y),J=A,endfunction
 410 ode("stiff",[0;1],0,1,f,Jacobian)
 411  ]]></programlisting>
 412     </refsection>
 413     <refsection role="see also">
 414         <title>Voir aussi</title>
 415         <simplelist type="inline">
 416             <member>
 417                 <link linkend="ode_discrete">ode_discrete</link>
 418             </member>
 419             <member>
 420                 <link linkend="ode_root">ode_root</link>
 421             </member>
 422             <member>
 423                 <link linkend="dassl">dassl</link>
 424             </member>
 425             <member>
 426                 <link linkend="impl">impl</link>
 427             </member>
 428             <member>
 429                 <link linkend="odedc">odedc</link>
 430             </member>
 431             <member>
 432                 <link linkend="odeoptions">odeoptions</link>
 433             </member>
 434             <member>
 435                 <link linkend="csim">csim</link>
 436             </member>
 437             <member>
 438                 <link linkend="ltitr">ltitr</link>
 439             </member>
 440             <member>
 441                 <link linkend="rtitr">rtitr</link>
 442             </member>
 443             <member>
 444                 <link linkend="intg">intg</link>
 445             </member>
 446         </simplelist>
 447     </refsection>
 448     <refsection>
 449         <title>Bibliographie</title>
 450         <para>Alan C. Hindmarsh, lsode and lsodi, two new initial value ordinary
 451             differential equation solvers, acm-signum newsletter, vol. 15, no. 4
 452             (1980), pp. 10-11.
 453         </para>
 454     </refsection>
 455     <refsection>
 456         <title>Fonctions Utilisées</title>
 457         <para>Les sous-programmes associés se trouvent dans le répertoire
 458             SCI/modules/differential_equations/src/fortran/: lsode.f lsoda.f lsodar.f
 459         </para>
 460     </refsection>
 461 </refentry>