scilab/modules/differential_equations/help/en_US/dasrt.xml

   1 <?xml version="1.0" encoding="UTF-8"?>
   2 <!--
   3  * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
   4  * Copyright (C) 2008 - INRIA
   5  * ...
   6  *
   7  * Copyright (C) 2012 - 2016 - Scilab Enterprises
   8  *
   9  * This file is hereby licensed under the terms of the GNU GPL v2.0,
  10  * pursuant to article 5.3.4 of the CeCILL v.2.1.
  11  * This file was originally licensed under the terms of the CeCILL v2.1,
  12  * and continues to be available under such terms.
  13  * For more information, see the COPYING file which you should have received
  14  * along with this program.
  15  *
  16  -->
  17 <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="dasrt" xml:lang="en">
  18     <refnamediv>
  19         <refname>dasrt</refname>
  20         <refpurpose>DAE solver with zero crossing</refpurpose>
  21     </refnamediv>
  22     <refsynopsisdiv>
  23         <title>Syntax</title>
  24         <synopsis>[r,nn,[,hd]] = dasrt(x0,t0,t [,atol,[rtol]],res [,jac],ng, surf [,info] [,hd])</synopsis>
  25     </refsynopsisdiv>
  26     <refsection>
  27         <title>Arguments</title>
  28         <variablelist>
  29             <varlistentry>
  30                 <term>x0</term>
  31                 <listitem>
  32                     <para>
  33                         is either <literal>y0</literal> (<literal>ydot0</literal> is
  34                         estimated by <literal>dassl</literal> with zero as first estimate)
  35                         or the matrix <literal>[y0 ydot0]</literal>.
  36                         <literal>g(t,y0,ydot0)</literal> must be equal to zero. If you only
  37                         know an estimate of <literal>ydot0</literal> set
  38                         <literal>info(7)=1</literal>.
  39                     </para>
  40                     <variablelist>
  41                         <varlistentry>
  42                             <term>y0</term>
  43                             <listitem>
  44                                 <para>a real column vector of initial conditions.</para>
  45                             </listitem>
  46                         </varlistentry>
  47                         <varlistentry>
  48                             <term>ydot0</term>
  49                             <listitem>
  50                                 <para>a real column vector of the time derivative of
  51                                     <literal>y</literal> at <literal>t0</literal> (may be an
  52                                     estimate).
  53                                 </para>
  54                             </listitem>
  55                         </varlistentry>
  56                     </variablelist>
  57                 </listitem>
  58             </varlistentry>
  59             <varlistentry>
  60                 <term>t0</term>
  61                 <listitem>
  62                     <para>a real number, it is the initial instant.</para>
  63                 </listitem>
  64             </varlistentry>
  65             <varlistentry>
  66                 <term>t</term>
  67                 <listitem>
  68                     <para>a real scalar or vector. Gives instants for which you want the
  69                         solution. Note that you can get solution at each dassl's step point
  70                         by setting <literal>info(2)=1</literal>.
  71                     </para>
  72                 </listitem>
  73             </varlistentry>
  74             <varlistentry>
  75                 <term>nn</term>
  76                 <listitem>
  77                     <para>
  78                         a vector with two entries <literal>[times num]</literal>
  79                         <literal>times</literal> is the value of the time at which the
  80                         surface is crossed, <literal>num</literal> is the number of the
  81                         crossed surface.
  82                     </para>
  83                 </listitem>
  84             </varlistentry>
  85             <varlistentry>
  86                 <term>atol, rtol</term>
  87                 <listitem>
  88                     <para>real scalars or column vectors of same size as
  89                         <literal>y</literal> or both of size 1. <literal>atol, rtol</literal> give respectively
  90                         absolute and relative error tolerances of solution. If vectors, the
  91                         tolerances are specified for each component of
  92                         <literal>y</literal>.
  93                     </para>
  94                 </listitem>
  95             </varlistentry>
  96             <varlistentry>
  97                 <term>res</term>
  98                 <listitem>
  99                     <para>
 100                         <link linkend="external">external</link> (function or list or string). Computes the value of
 101                         <literal>g(t,y,ydot)</literal>. It may be :
 102                     </para>
 103                     <itemizedlist>
 104                         <listitem>
 105                             <para>A Scilab function.</para>
 106                             <para>Its syntax must be
 107                                 <literal>[r,ires]=res(t,y,ydot)</literal> and
 108                                 <literal>res</literal> must return the residue
 109                                 <literal>r=g(t,y,ydot)</literal> and error flag
 110                                 <literal>ires</literal>. <literal>ires = 0</literal> if
 111                                 <literal>res</literal> succeeds to compute <literal>r</literal>,
 112                                 <literal>=-1</literal> if residue is locally not defined for
 113                                 <literal>(t,y,ydot)</literal>, <literal>=-2</literal> if
 114                                 parameters are out of admissible range.
 115                             </para>
 116                         </listitem>
 117                         <listitem>
 118                             <para>A list.</para>
 119                             <para>This form allows to pass parameters other than t,y,ydot to
 120                                 the function. It must be as follows:
 121                             </para>
 122                             <screen><![CDATA[
 123 list(res,x1,x2,...)
 124 ]]></screen>
 125                             <para>where the syntax of the function
 126                                 <literal>res</literal> is now
 127                             </para>
 128                             <screen><![CDATA[
 129 r = res(t,y,ydot,x1,x2,...)
 130 ]]></screen>
 131                             <para>
 132                                 <literal>res</literal> still returns
 133                                 <literal>r=g(t,y,ydot)</literal> as a function of
 134                                 <literal>(t,y,ydot,x1,x2,...)</literal>.
 135                             </para>
 136                             <para>
 137                                 <warning>
 138                                     Warning: this form must not be used if there is no extra
 139                                     argument to pass to the function.
 140                                 </warning>
 141                             </para>
 142                         </listitem>
 143                         <listitem>
 144                             <para>A string.</para>
 145                             <para>It must refer to the name of a C or Fortran subroutine
 146                                 linked with Scilab.
 147                             </para>
 148                             <para>In C The syntax must be:</para>
 149                             <screen><![CDATA[
 150 void res(double *t, double y[], double yd[], double r[],
 151          int *ires, double rpar[], int ipar[])
 152 ]]></screen>
 153                             <para>In Fortran it must be:</para>
 154                             <screen><![CDATA[
 155 subroutine res(t,y,yd,r,ires,rpar,ipar)
 156 double precision t, y(*),yd(*),r(*),rpar(*)
 157 integer ires,ipar(*)
 158 ]]></screen>
 159                             <para>
 160                                 The <literal>rpar</literal> and <literal>ipar</literal> arrays must be present but cannot be
 161                                 used.
 162                             </para>
 163                         </listitem>
 164                     </itemizedlist>
 165                 </listitem>
 166             </varlistentry>
 167             <varlistentry>
 168                 <term>jac</term>
 169                 <listitem>
 170                     <para>
 171                         <link linkend="external">external</link> (function or list or string). Computes the value of
 172                         <literal>dg/dy + cj*dg/dydot</literal> for a given value of parameter
 173                         <literal>cj</literal>.
 174                     </para>
 175                     <itemizedlist>
 176                         <listitem>
 177                             <para>A Scilab function.</para>
 178                             <para>Its syntax must be
 179                                 <literal>r=jac(t,y,ydot,cj)</literal> and the
 180                                 <literal>jac</literal> function must return
 181                                 <literal>r=dg(t,y,ydot)/dy+cj*dg(t,y,ydot)/dydot</literal> where
 182                                 <literal>cj</literal> is a real scalar.
 183                             </para>
 184                         </listitem>
 185                         <listitem>
 186                             <para>A list.</para>
 187                             <para>It must be as follows</para>
 188                             <screen><![CDATA[
 189 list(jac,x1,x2,...)
 190 ]]></screen>
 191                             <para>where the syntax of the function
 192                                 <literal>jac</literal> is now
 193                             </para>
 194                             <screen><![CDATA[
 195 r = jac(t,y,ydot,cj,x1,x2,...)
 196 ]]></screen>
 197                             <para>
 198                                 <literal>jac</literal> still returns
 199                                 <literal>dg/dy + cj*dg/dydot</literal> as a function of
 200                                 <literal>(t,y,ydot,cj,x1,x2,...)</literal>.
 201                             </para>
 202                         </listitem>
 203                         <listitem>
 204                             <para>A character string.</para>
 205                             <para>It must refer to the name of a Fortran subroutine linked
 206                                 with scilab
 207                             </para>
 208                             <para>In C The syntax must be:</para>
 209                             <screen><![CDATA[
 210 void jac(double *t, double y[], double yd[], double pd[],
 211          double *cj, double rpar[], int ipar[])
 212 ]]></screen>
 213                             <para>In Fortran it must be:</para>
 214                             <screen><![CDATA[
 215 subroutine jac(t,y,yd,pd,cj,rpar,ipar)
 216 double precision t, y(*),yd(*),pd(*),cj,rpar(*)
 217 integer ipar(*)
 218 ]]></screen>
 219                         </listitem>
 220                     </itemizedlist>
 221                 </listitem>
 222             </varlistentry>
 223             <varlistentry>
 224                 <term>surf</term>
 225                 <listitem>
 226                     <para>
 227                         <link linkend="external">external</link> (function or list or string). Computes the value of
 228                         the column vector <literal>surf(t,y)</literal> with
 229                         <literal>ng</literal> components. Each component defines a surface.
 230                         It may be defined by:
 231                     </para>
 232                     <itemizedlist>
 233                         <listitem>
 234                             <para>A Scilab function.</para>
 235                             <para>Its syntax must be
 236                                 <literal>surf(t,y)</literal>
 237                             </para>
 238                         </listitem>
 239                         <listitem>
 240                             <para>A list.</para>
 241                             <para>It must be as follows</para>
 242                             <screen><![CDATA[
 243 list(surf,x1,x2,...)
 244 ]]></screen>
 245                             <para>where the syntax of the function
 246                                 <literal>surf</literal> is now
 247                             </para>
 248                             <screen><![CDATA[
 249 r = surf(t,y,x1,x2,...)
 250 ]]></screen>
 251                         </listitem>
 252                         <listitem>
 253                             <para>A character string.</para>
 254                             <para>It must refer to the name of a Fortran subroutine linked
 255                                 with scilab.
 256                             </para>
 257                             <para>In C the syntax must be:</para>
 258                             <screen><![CDATA[
 259 void surf(int *ny, double *t, double y[], int *ng, double gout[])
 260 ]]></screen>
 261                             <para>In Fortran it must be:</para>
 262                             <screen><![CDATA[
 263 subroutine surf(ny,t,y,ng,gout)
 264 double precision t, y(*),gout(*)
 265 integer ny,ng
 266 ]]></screen>
 267                         </listitem>
 268                     </itemizedlist>
 269                 </listitem>
 270             </varlistentry>
 271             <varlistentry>
 272                 <term>info</term>
 273                 <listitem>
 274                     <para>
 275                         list which contains <literal>7</literal> elements. Default
 276                         value is <literal>list([],0,[],[],[],0,0)</literal>.
 277                     </para>
 278                     <variablelist>
 279                         <varlistentry>
 280                             <term>info(1)</term>
 281                             <listitem>
 282                                 <para>real scalar which gives the maximum time for which
 283                                     <literal>g</literal> is allowed to be evaluated or an empty
 284                                     matrix <literal>[]</literal> if no limits imposed for
 285                                     time.
 286                                 </para>
 287                             </listitem>
 288                         </varlistentry>
 289                         <varlistentry>
 290                             <term>info(2)</term>
 291                             <listitem>
 292                                 <para>
 293                                     flag which indicates if <literal>dassl</literal> returns
 294                                     its intermediate computed values (<literal>flag=1</literal>)
 295                                     or only the user specified time point values
 296                                     (<literal>flag=0</literal>).
 297                                 </para>
 298                             </listitem>
 299                         </varlistentry>
 300                         <varlistentry>
 301                             <term>info(3)</term>
 302                             <listitem>
 303                                 <para>
 304                                     <literal>2</literal> components vector which give the
 305                                     definition <literal>[ml,mu]</literal> of band matrix computed
 306                                     by <literal>jac</literal>; <literal>r(i - j + ml + mu + 1,j) =
 307                                         "dg(i)/dy(j)+cj*dg(i)/dydot(j)"
 308                                     </literal>
 309                                     .If <literal>jac</literal> returns a full matrix set
 310                                     <literal>info(3)=[]</literal>.
 311                                 </para>
 312                             </listitem>
 313                         </varlistentry>
 314                         <varlistentry>
 315                             <term>info(4)</term>
 316                             <listitem>
 317                                 <para>real scalar which gives the maximum step size. Set
 318                                     <literal>info(4)=[]</literal> if no limitation.
 319                                 </para>
 320                             </listitem>
 321                         </varlistentry>
 322                         <varlistentry>
 323                             <term>info(5)</term>
 324                             <listitem>
 325                                 <para>real scalar which gives the initial step size. Set
 326                                     <literal>info(5)=[]</literal> if not specified.
 327                                 </para>
 328                             </listitem>
 329                         </varlistentry>
 330                         <varlistentry>
 331                             <term>info(6)</term>
 332                             <listitem>
 333                                 <para>
 334                                     set <literal>info(6)=1</literal> if the solution is
 335                                     known to be non negative, else set
 336                                     <literal>info(6)=0</literal>.
 337                                 </para>
 338                             </listitem>
 339                         </varlistentry>
 340                         <varlistentry>
 341                             <term>info(7)</term>
 342                             <listitem>
 343                                 <para>
 344                                     set <literal>info(7)=1</literal> if
 345                                     <literal>ydot0</literal> is just an estimation,
 346                                     <literal>info(7)=0</literal> if
 347                                     <literal>g(t0,y0,ydot0)=0</literal>.
 348                                 </para>
 349                             </listitem>
 350                         </varlistentry>
 351                     </variablelist>
 352                 </listitem>
 353             </varlistentry>
 354             <varlistentry>
 355                 <term>hd</term>
 356                 <listitem>
 357                     <para>
 358                         real vector which allows to store the <literal>dassl</literal>
 359                         context and to resume integration.
 360                     </para>
 361                 </listitem>
 362             </varlistentry>
 363             <varlistentry>
 364                 <term>r</term>
 365                 <listitem>
 366                     <para>
 367                         real matrix . Each column is the vector <literal>[t;x(t);xdot(t)]</literal> where
 368                         <literal>t</literal> is time index for which the solution had been computed.
 369                     </para>
 370                 </listitem>
 371             </varlistentry>
 372         </variablelist>
 373     </refsection>
 374     <refsection>
 375         <title>Description</title>
 376         <para>Solution of the implicit differential equation.</para>
 377         <screen><![CDATA[
 378 g(t,y,ydot) = 0
 379 y(t0) = y0  and   ydot(t0) = ydot0
 380 ]]></screen>
 381         <para>Returns the surface crossing instants and the number of the surface
 382             reached in <literal>nn</literal>.
 383         </para>
 384         <para>Detailed examples can be found in SCIDIR/tests/unit_tests/dassldasrt.tst</para>
 385     </refsection>
 386     <refsection>
 387         <title>Examples</title>
 388         <programlisting role="example"><![CDATA[
 389 // dy/dt = ((2*log(y)+8)/t -5)*y,  y(1) = 1,  1<=t<=6
 390 // g1 = ((2*log(y)+8)/t - 5)*y
 391 // g2 = log(y) - 2.2491
 392 y0 = 1; t = 2:6; t0 = 1; y0d = 3;
 393 atol = 1.d-6; rtol = 0; ng = 2;
 394
 395 deff('[delta,ires] = res1(t,y,ydot)', 'ires=0; delta=ydot-((2*log(y)+8)/t-5)*y')
 396 deff('rts = gr1(t,y)', 'rts=[((2*log(y)+8)/t-5)*y;log(y)-2.2491]')
 397
 398 [yy,nn] = dasrt([y0,y0d],t0,t,atol,rtol,res1,ng,gr1);
 399 //(Should return nn=[2.4698972 2])
 400      ]]></programlisting>
 401     </refsection>
 402     <refsection role="see also">
 403         <title>See also</title>
 404         <simplelist type="inline">
 405             <member>
 406                 <link linkend="ode">ode</link>
 407             </member>
 408             <member>
 409                 <link linkend="dassl">dassl</link>
 410             </member>
 411             <member>
 412                 <link linkend="daskr">daskr</link>
 413             </member>
 414             <member>
 415                 <link linkend="impl">impl</link>
 416             </member>
 417             <member>
 418                 <link linkend="call">call</link>
 419             </member>
 420             <member>
 421                 <link linkend="link">link</link>
 422             </member>
 423             <member>
 424                 <link linkend="external">external</link>
 425             </member>
 426         </simplelist>
 427     </refsection>
 428 </refentry>