Xcos Solvers : Scilab Documentation 42/9542/33
Paul BIGNIER [Fri, 26 Oct 2012 14:40:56 +0000 (16:40 +0200)]
Added help files to describe the functioning of the solvers (CVode, IDA, Runge-Kutta)

Change-Id: Ib5d1510b3084187ca40524f13c4897a231c2bc75

scilab/modules/xcos/examples/solvers/IDA_Example.xcos [new file with mode: 0644]
scilab/modules/xcos/examples/solvers/ODE_Example.xcos [new file with mode: 0644]
scilab/modules/xcos/examples/solvers/integ.sce [new file with mode: 0644]
scilab/modules/xcos/examples/solvers/integRK.sce [new file with mode: 0644]
scilab/modules/xcos/help/en_US/solvers/CHAPTER [new file with mode: 0644]
scilab/modules/xcos/help/en_US/solvers/CVode.xml [new file with mode: 0644]
scilab/modules/xcos/help/en_US/solvers/IDA.xml [new file with mode: 0644]
scilab/modules/xcos/help/en_US/solvers/Runge-Kutta.xml [new file with mode: 0644]

diff --git a/scilab/modules/xcos/examples/solvers/IDA_Example.xcos b/scilab/modules/xcos/examples/solvers/IDA_Example.xcos
new file mode 100644 (file)
index 0000000..cd547a9
Binary files /dev/null and b/scilab/modules/xcos/examples/solvers/IDA_Example.xcos differ
diff --git a/scilab/modules/xcos/examples/solvers/ODE_Example.xcos b/scilab/modules/xcos/examples/solvers/ODE_Example.xcos
new file mode 100644 (file)
index 0000000..3c925b1
Binary files /dev/null and b/scilab/modules/xcos/examples/solvers/ODE_Example.xcos differ
diff --git a/scilab/modules/xcos/examples/solvers/integ.sce b/scilab/modules/xcos/examples/solvers/integ.sce
new file mode 100644 (file)
index 0000000..a8a82ad
--- /dev/null
@@ -0,0 +1,50 @@
+// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+// Copyright (C) 2009 - Paul Bignier
+//
+// This file must be used under the terms of the CeCILL.
+// This source file is licensed as described in the file COPYING, which
+// you should have received as part of this distribution.  The terms
+// are also available at
+// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+
+// Run with exec("SCI/modules/xcos/help/en_US/solvers/integ.sce");
+
+// Import the diagram and augment the ending time
+importXcosDiagram("SCI/modules/xcos/examples/solvers/ODE_Example.xcos");
+scs_m.props.tf = 3500;
+
+// BDF / Newton
+// Select the solver
+scs_m.props.tol(6) = 0;
+// Start the timer, launch the simulation and display time
+timer();
+xcos_simulate(scs_m, 4);
+t = timer();
+disp(t, "Time for BDF / Newton :");
+
+// BDF / Functional
+// Select the solver
+scs_m.props.tol(6) = 1;
+// Start the timer, launch the simulation and display time
+timer();
+xcos_simulate(scs_m, 4);
+t = timer();
+disp(t, "Time for BDF / Functional :");
+
+// Adams / Functional
+// Select the solver
+scs_m.props.tol(6) = 3;
+// Start the timer, launch the simulation and display time
+timer();
+xcos_simulate(scs_m, 4);
+t = timer();
+disp(t, "Time for Adams / Functional :");
+
+// Adams / Newton
+// Select the solver
+scs_m.props.tol(6) = 2;
+// Start the timer, launch the simulation and display time
+timer();
+xcos_simulate(scs_m, 4);
+t = timer();
+disp(t, "Time for Adams / Newton :");
diff --git a/scilab/modules/xcos/examples/solvers/integRK.sce b/scilab/modules/xcos/examples/solvers/integRK.sce
new file mode 100644 (file)
index 0000000..608eda5
--- /dev/null
@@ -0,0 +1,60 @@
+// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+// Copyright (C) 2009 - Paul Bignier
+//
+// This file must be used under the terms of the CeCILL.
+// This source file is licensed as described in the file COPYING, which
+// you should have received as part of this distribution.  The terms
+// are also available at
+// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+
+// Run with exec("SCI/modules/xcos/help/en_US/solvers/integRK.sce");
+
+// Import the diagram and augment the ending time
+importXcosDiagram("SCI/modules/xcos/examples/solvers/ODE_Example.xcos");
+scs_m.props.tf = 3500;
+
+// BDF / Newton
+// Select the solver
+scs_m.props.tol(6) = 0;
+// Start the timer, launch the simulation and display time
+timer();
+xcos_simulate(scs_m, 4);
+t = timer();
+disp(t, "Time for BDF / Newton :");
+
+// BDF / Functional
+// Select the solver
+scs_m.props.tol(6) = 1;
+// Start the timer, launch the simulation and display time
+timer();
+xcos_simulate(scs_m, 4);
+t = timer();
+disp(t, "Time for BDF / Functional :");
+
+// Adams / Functional
+// Select the solver
+scs_m.props.tol(6) = 3;
+// Start the timer, launch the simulation and display time
+timer();
+xcos_simulate(scs_m, 4);
+t = timer();
+disp(t, "Time for Adams / Functional :");
+
+// Adams / Newton
+// Select the solver
+scs_m.props.tol(6) = 2;
+// Start the timer, launch the simulation and display time
+timer();
+xcos_simulate(scs_m, 4);
+t = timer();
+disp(t, "Time for Adams / Newton :");
+
+// Runge-Kutta
+// Select the solver and set abstol to 10^-2
+scs_m.props.tol(6) = 4;
+scs_m.props.tol(1) = 0.01;
+// Start the timer, launch the simulation and display time
+timer();
+xcos_simulate(scs_m, 4);
+t = timer();
+disp(t, "Time for Runge-Kutta :");
diff --git a/scilab/modules/xcos/help/en_US/solvers/CHAPTER b/scilab/modules/xcos/help/en_US/solvers/CHAPTER
new file mode 100644 (file)
index 0000000..0b8fad9
--- /dev/null
@@ -0,0 +1 @@
+title = Solvers
diff --git a/scilab/modules/xcos/help/en_US/solvers/CVode.xml b/scilab/modules/xcos/help/en_US/solvers/CVode.xml
new file mode 100644 (file)
index 0000000..2a782c7
--- /dev/null
@@ -0,0 +1,276 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<!--
+ * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ * Copyright (C) Scilab Enterprises - 2012 - Paul Bignier
+ *
+ * This file must be used under the terms of the CeCILL.
+ * This source file is licensed as described in the file COPYING, which
+ * you should have received as part of this distribution.
+ * The terms are also available at
+ * http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+ -->
+<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg"  xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:lang="en_US" xml:id="CVode">
+    <refnamediv>
+        <refname>CVode</refname>
+        <refpurpose>
+            <emphasis>CVode</emphasis> is a numerical solver providing an efficient and stable method to solve ODE Initial Value Problems. Called by scicos(), it uses either <emphasis>BDF</emphasis> or <emphasis>Adams</emphasis> as implicit integration method, and <emphasis>Newton</emphasis> or <emphasis>Functional</emphasis> iterations
+        </refpurpose>
+    </refnamediv>
+    <refsection>
+        <title>Description</title>
+        <para>
+            <emphasis>CVode</emphasis> is a numerical solver providing an efficient and stable method to solve Initial Value Problems of the form :
+        </para>
+        <para>
+            <latex>
+                \begin{eqnarray}
+                \dot{y} = f(t,y), \hspace{3 mm} y(t_0) = y_0, \hspace{3 mm} y \in R^N
+                \end{eqnarray}
+            </latex>
+        </para>
+        <para>
+            Starting with 
+            <emphasis>
+                y<subscript>0</subscript>
+            </emphasis>
+            ,<emphasis>CVode</emphasis> approximates 
+            <emphasis>
+                y<subscript>n+1</subscript>
+            </emphasis>
+            with the formula :
+        </para>
+        <para>
+            <latex>
+                \begin{eqnarray}
+                \sum_{i=0}^{K_1} \alpha_{n,i} y_{n-i} + h_n\sum_{i=0}^{K_2} \beta_{n,i} \dot{y}_{n-i} = 0,\hspace{10 mm} (1)
+                \end{eqnarray}
+            </latex>
+            <para>
+                &#160; with 
+                <emphasis>
+                    y<subscript>n</subscript>
+                </emphasis>
+                the approximation of 
+                <emphasis>
+                    y(t<subscript>n</subscript>)
+                </emphasis>
+                ,and 
+                <emphasis>
+                    h<subscript>n</subscript>
+                </emphasis>
+                =
+                <emphasis>
+                    t<subscript>n</subscript> - t<subscript>n-1</subscript>
+                </emphasis>
+                the step size.
+            </para>
+        </para>
+        <para>
+            These implicit methods are characterized by their respective order <emphasis>q</emphasis>, which indicates the number of intermediate points required to compute 
+            <emphasis>
+                y<subscript>n+1</subscript>
+            </emphasis>
+            .
+        </para>
+        <para>
+            This is where the difference between <emphasis>BDF</emphasis> and <emphasis>Adams</emphasis> intervenes (<emphasis>Backward Differenciation Formula</emphasis> and <emphasis>Adams-Moulton formula</emphasis>) :
+        </para>
+        <para>
+            <caution>
+                If the problem is stiff, the user should select <emphasis>BDF</emphasis> :
+            </caution>
+        </para>
+        <itemizedlist>
+            <listitem>
+                <emphasis>q</emphasis>, the order of the method, is set between 1 and 5 (automated),
+            </listitem>
+            <listitem>
+                <emphasis>K1 = q</emphasis> and <emphasis>K2 = 0</emphasis>.
+            </listitem>
+        </itemizedlist>
+        <para>
+            In the case of nonstiffness, <emphasis>Adams</emphasis> is preferred :
+        </para>
+        <itemizedlist>
+            <listitem>
+                <emphasis>q</emphasis> is set between 1 and 12 (automated),
+            </listitem>
+            <listitem>
+                <emphasis>K1 = 1</emphasis> and <emphasis>K2 = q</emphasis>.
+            </listitem>
+        </itemizedlist>
+        <para>
+            The coefficients are fixed, uniquely determined by the method type, its order, the history of the step sizes, and the normalization 
+            <emphasis>
+                &#x3B1;<subscript>n, 0</subscript> = -1
+            </emphasis>
+            .
+        </para>
+        <para>
+            For either choice and at each step, injecting this integration in <emphasis>(1)</emphasis> yields the nonlinear system :
+        </para>
+        <para>
+            <latex>
+                G(y_n)\equiv y_n-h_n\beta_{n,0}f(t_n,y_n)-a_n=0, \hspace{2 mm} where \hspace{2 mm} a_n\equiv \sum_{i>0} (\alpha_{n,i} y_{n-i} + h_n\beta_{n,i}\dot{y}_{n-i})
+            </latex>
+        </para>
+        <para>
+            This system can be solved by either <emphasis>Functional</emphasis> or <emphasis>Newton</emphasis> iterations, described hereafter.
+        </para>
+        <para>
+            In both following cases, the initial "predicted" 
+            <emphasis>
+                y<subscript>n(0)</subscript>
+            </emphasis>
+            is explicitly computed from the history data, by adding derivatives.
+        </para>
+        <itemizedlist>
+            <para>  </para>
+            <listitem>
+                <emphasis role="bold">Functional :</emphasis> this method only involves evaluations of <emphasis>f</emphasis>, it simply computes 
+                <emphasis>
+                    y<subscript>n(0)</subscript>
+                </emphasis>
+                by iterating the formula :
+                <para>
+                    <latex>
+                        y_{n(m+1)} = h_n &#x3B2;_{n,0} f(t_n,y_{n(m+1)}) + a_n
+                    </latex>
+                    <latex>
+                        where \hspace{2 mm} a_n\equiv \sum_{i>0} (\alpha_{n,i} y_{n-i} + h_n\beta_{n,i}\dot{y}_{n-i})
+                    </latex>
+                </para>
+            </listitem>
+            <para>  </para>
+            <listitem>
+                <emphasis role="bold">Newton :</emphasis> here, we use an implemented direct dense solver on the linear system :
+                <para>
+                    <latex>
+                        M[y_{n(m+1)}-y_{n(m)}]=-G(y_{n(m)}), \hspace{4 mm} M \approx I-\gamma J, \hspace{2 mm} J=\frac{\partial f}{\partial y}, \hspace{2 mm} and \hspace{2 mm} \gamma = h_n\beta_{n,0}
+                    </latex>
+                </para>
+            </listitem>
+        </itemizedlist>
+        <para>
+            In both situations, <emphasis>CVode</emphasis> uses the history array to control the local error 
+            <emphasis>
+                y<subscript>n(m)</subscript> - y<subscript>n(0)</subscript>
+            </emphasis>
+            and recomputes 
+            <emphasis>
+                h<subscript>n</subscript>
+            </emphasis>
+            if that error is not satisfying.
+        </para>
+        <para>
+            <emphasis role="bold">
+                The recommended choices are <emphasis>BDF / Newton</emphasis> for stiff problems and <emphasis>Adams / Functional</emphasis> for the nonstiff ones.
+            </emphasis>
+        </para>
+        <para>
+            The function is called in between activations, because a discrete activation may change the system.
+        </para>
+        <para>
+            Following the criticality of the event (its effect on the continuous problem), we either relaunch the solver with different start and final times as if nothing happened, or, if the system has been modified, we need to "cold-restart" the problem by reinitializing it anew and relaunching the solver.
+        </para>
+    </refsection>
+    <refsection>
+        <title>Examples</title>
+        <para>
+            <link type="scilab" linkend="scilab.xcos/xcos/examples/solvers/ODE_Example.xcos">
+                <inlinemediaobject>
+                    <imageobject>
+                        <imagedata align="center" fileref="../../../examples/solvers/ODE_Example.xcos" valign="middle"/>
+                    </imageobject>
+                </inlinemediaobject>
+            </link>
+            <scilab:image><![CDATA[
+importXcosDiagram(SCI + "/modules/xcos/examples/solvers/ODE_Example.xcos");
+xcos_simulate(scs_m, 4);
+]]></scilab:image>
+        </para>
+        <para>
+            The integral block returns its continuous state, we can evaluate it with BDF / Newton by running the example :
+        </para>
+        <para>
+            <programlisting language="example"><![CDATA[
+      // Import the diagram and set the ending time
+      importXcosDiagram("SCI/modules/xcos/examples/solvers/ODE_Example.xcos");
+      scs_m.props.tf = 3500;
+
+      // Select the solver BDF / Newton
+      scs_m.props.tol(6) = 0;
+
+      // Start the timer, launch the simulation and display time
+      timer();
+      xcos_simulate(scs_m, 4);
+      t = timer();
+      disp(t, "Time for BDF / Newton :");
+      ]]></programlisting>
+        </para>
+        <para>
+            The Scilab console displays :
+            <screen><![CDATA[
+Time for BDF / Newton :
+ 27.549
+            ]]></screen>
+        </para>
+        <para>
+            Now, in the following script, we compare the time difference between the methods by running the example with the four solvers in turn :
+            <link type="scilab" linkend ="scilab.scinotes/xcos/examples/solvers/integ.sce">
+                Open the script
+            </link>
+        </para>
+        <para>
+            Results :
+            <screen language="example"><![CDATA[
+Time for BDF / Newton :
+ 25.457
+
+Time for BDF / Functional :
+ 24.893
+
+Time for Adams / Functional :
+ 20.049
+
+Time for Adams / Newton :
+ 20.469
+            ]]></screen>
+        </para>
+        <para>
+            The results show that for a simple nonstiff continuous problem, Adams / Functional is fastest.
+        </para>
+    </refsection>
+    <refsection>
+        <title>See Also</title>
+        <simplelist type="inline">
+            <member>
+                <link linkend="IDA">IDA</link>
+            </member>
+            <member>
+                <link linkend="RK">Runge-Kutta 4(5)</link>
+            </member>
+            <member>
+                <link linkend="ode">ode</link>
+            </member>
+            <member>
+                <link linkend="ode_discrete">ode_discrete</link>
+            </member>
+            <member>
+                <link linkend="ode_root">ode_root</link>
+            </member>
+            <member>
+                <link linkend="odedc">odedc</link>
+            </member>
+            <member>
+                <link linkend="impl">impl</link>
+            </member>
+        </simplelist>
+    </refsection>
+    <refsection>
+        <title>Bibliography</title>
+        <para>
+            <ulink url="https://computation.llnl.gov/casc/sundials/documentation/documentation.html">Sundials Documentation</ulink>
+        </para>
+    </refsection>
+</refentry>
diff --git a/scilab/modules/xcos/help/en_US/solvers/IDA.xml b/scilab/modules/xcos/help/en_US/solvers/IDA.xml
new file mode 100644 (file)
index 0000000..7aa39bd
--- /dev/null
@@ -0,0 +1,190 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<!--
+ * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ * Copyright (C) Scilab Enterprises - 2012 - Paul Bignier
+ *
+ * This file must be used under the terms of the CeCILL.
+ * This source file is licensed as described in the file COPYING, which
+ * you should have received as part of this distribution.
+ * The terms are also available at
+ * http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+ -->
+<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg"  xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:lang="en_US" xml:id="IDA">
+    <refnamediv>
+        <refname>IDA</refname>
+        <refpurpose>
+            <emphasis>IDA</emphasis> (Implicit Differential Algebraic) is a numerical solver providing an efficient and stable method to solve DAE Initial Value Problems. Called by <link linkend="scicos">scicos</link>.
+        </refpurpose>
+    </refnamediv>
+    <refsection>
+        <title>Description</title>
+        <para>
+            <emphasis>IDA</emphasis> is a numerical solver providing an efficient and stable method to solve Initial Value Problems of the form :
+        </para>
+        <para>
+            <latex>
+                \begin{eqnarray}
+                F(t,y,\dot{y}) = 0, \hspace{2 mm} y(t_0)=y_0, \hspace{2 mm} \dot{y}(t_0)=\dot{y}_0, \hspace{3 mm} y, \hspace{1.5 mm} \dot{y}  \hspace{1.5 mm} and \hspace{1.5 mm} F \in R^N \hspace{10 mm} (1)
+                \end{eqnarray}
+            </latex>
+        </para>
+        <para>
+        </para>
+        Before solving the problem, <emphasis>IDA</emphasis> runs an implemented routine to find consistent values for 
+        <emphasis>
+            y<subscript>0</subscript>
+        </emphasis>
+        and 
+        <emphasis>
+            yPrime<subscript>0</subscript>
+        </emphasis>
+        .
+        <para>
+            Starting then with those 
+            <emphasis>
+                y<subscript>0</subscript>
+            </emphasis>
+            and 
+            <emphasis>
+                yPrime<subscript>0</subscript>
+            </emphasis>
+            ,<emphasis>IDA</emphasis> approximates 
+            <emphasis>
+                y<subscript>n+1</subscript>
+            </emphasis>
+            with the <emphasis>BDF</emphasis> formula :
+        </para>
+        <para>
+            <latex>
+                \begin{eqnarray}
+                \sum_{i=0}^{q} \alpha_{n,i} y_{n-i} = h_n\dot{y}_{n}
+                \end{eqnarray}
+            </latex>
+            <para>
+                &#160; with, like in <link linkend="CVode">CVode</link>, 
+                <emphasis>
+                    y<subscript>n</subscript>
+                </emphasis>
+                the approximation of 
+                <emphasis>
+                    y(t<subscript>n</subscript>)
+                </emphasis>
+                ,
+                <emphasis>
+                    h<subscript>n</subscript>
+                </emphasis>
+                =
+                <emphasis>
+                    t<subscript>n</subscript> - t<subscript>n-1</subscript>
+                </emphasis>
+                the step size, and the coefficients are fixed, uniquely determined by the method type, its order <emphasis>q</emphasis> ranging from 1 to 5 and the history of the step sizes.
+            </para>
+        </para>
+        <para>
+            Injecting this formula in <emphasis>(1)</emphasis> yields the system :
+        </para>
+        <para>
+            <latex>
+                G(y_n) \equiv F \left( t_n, \hspace{1.5mm} y_n, \hspace{1.5mm} h_n^{-1}\sum_{i=0}^{q} \alpha_{n,i} y_{n-i} \right) = 0
+            </latex>
+        </para>
+        <para>
+            To apply <emphasis>Newton</emphasis> iterations to it, we rewrite it into :
+        </para>
+        <para>
+            <latex>
+                J \left[y_{n(m+1)}-y_{n(m)} \right] = -G(y_{n(m)})
+            </latex>
+        </para>
+        <para>
+            &#160; with <emphasis>J</emphasis> an approximation of the Jacobian :
+        </para>
+        <para>
+            &#160; <latex>
+                J = \frac{\partial{G}}{\partial{y}} = \frac{\partial{F}}{\partial{y}}+\alpha\frac{\partial{F}}{\partial{\dot{y}}}, \hspace{4 mm} \alpha = \frac{\alpha_{n,0}}{h_n},
+            </latex>
+        </para>
+        <para>
+            &#160; <emphasis>&#x3B1;</emphasis> changes whenever the step size or the method order varies.
+        </para>
+        <para>
+            An implemented direct dense solver is used and we go on to the next step.
+        </para>
+        <para>
+            <emphasis>IDA</emphasis> uses the history array to control the local error 
+            <emphasis>
+                y<subscript>n(m)</subscript> - y<subscript>n(0)</subscript>
+            </emphasis>
+            and recomputes 
+            <emphasis>
+                h<subscript>n</subscript>
+            </emphasis>
+            if that error is not satisfying.
+        </para>
+        <para>
+            The function is called in between activations, because a discrete activation may change the system.
+        </para>
+        <para>
+            Following the criticality of the event (its effect on the continuous problem), we either relaunch the solver with different start and final times as if nothing happened, or, if the system has been modified, we need to "cold-restart" the problem by reinitializing it anew and relaunching the solver.
+        </para>
+    </refsection>
+    <refsection id="Example_IDA">
+        <title>Example</title>
+        <para>
+            The 'Modelica Generic' block returns its continuous states, we can evaluate them with IDA by running the example :
+        </para>
+        <para>
+            <link type="scilab" linkend="scilab.xcos/xcos/examples/solvers/IDA_Example.xcos">
+                <inlinemediaobject>
+                    <imageobject>
+                        <imagedata align="center" fileref="../../../examples/solvers/IDA_Example.xcos" valign="middle"/>
+                    </imageobject>
+                </inlinemediaobject>
+            </link>
+            <scilab:image><![CDATA[
+importXcosDiagram(SCI + "/modules/xcos/examples/solvers/IDA_Example.xcos");
+
+// Redefining messagebox() to avoid popup
+function messagebox(msg, title)
+disp(title);
+disp(msg);
+endfunction
+
+// Simulation
+xcos_simulate(scs_m, 4);
+]]></scilab:image>
+        </para>
+    </refsection>
+    <refsection>
+        <title>See Also</title>
+        <simplelist type="inline">
+            <member>
+                <link linkend="CVode">CVode</link>
+            </member>
+            <member>
+                <link linkend="RK">Runge-Kutta 4(5)</link>
+            </member>
+            <member>
+                <link linkend="ode">ode</link>
+            </member>
+            <member>
+                <link linkend="ode_discrete">ode_discrete</link>
+            </member>
+            <member>
+                <link linkend="ode_root">ode_root</link>
+            </member>
+            <member>
+                <link linkend="odedc">odedc</link>
+            </member>
+            <member>
+                <link linkend="impl">impl</link>
+            </member>
+        </simplelist>
+    </refsection>
+    <refsection>
+        <title>Bibliography</title>
+        <para>
+            <ulink url="https://computation.llnl.gov/casc/sundials/documentation/documentation.html">Sundials Documentation</ulink>
+        </para>
+    </refsection>
+</refentry>
diff --git a/scilab/modules/xcos/help/en_US/solvers/Runge-Kutta.xml b/scilab/modules/xcos/help/en_US/solvers/Runge-Kutta.xml
new file mode 100644 (file)
index 0000000..0de65e2
--- /dev/null
@@ -0,0 +1,238 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<!--
+ * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ * Copyright (C) Scilab Enterprises - 2012 - Paul Bignier
+ *
+ * This file must be used under the terms of the CeCILL.
+ * This source file is licensed as described in the file COPYING, which
+ * you should have received as part of this distribution.
+ * The terms are also available at
+ * http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+ -->
+<refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg"  xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook" xmlns:scilab="http://www.scilab.org" xml:lang="en_US" xml:id="RK">
+    <refnamediv>
+        <refname>Runge-Kutta 4(5)</refname>
+        <refpurpose>
+            Runge-Kutta is a numerical solver providing an efficient explicit method to solve ODE Initial Value Problems. Called by <link linkend="scicos">scicos</link>.
+        </refpurpose>
+    </refnamediv>
+    <refsection>
+        <title>Description</title>
+        <para>
+            Runge-Kutta is a numerical solver providing an efficient fixed-size step method to solve Initial Value Problems of the form :
+        </para>
+        <para>
+            <latex>
+                \begin{eqnarray}
+                \dot{y} = f(t,y), \hspace{3 mm} y(t_0) = y_0, \hspace{3 mm} y \in R^N
+                \end{eqnarray}
+            </latex>
+        </para>
+        <para>
+            <emphasis>CVode</emphasis> and <emphasis>IDA</emphasis> use variable-size steps for the integration.
+        </para>
+        <para>
+            A drawback of that is the unpredictable computation time. With Runge-Kutta, we don't adapt to the complexity of the problem, but we guarantee a stable computation time.
+        </para>
+        <para>
+            As of now, this method is explicit, so it is not concerned with Newton or Functional iterations, and not advised for stiff problems.
+        </para>
+        <para>
+            It is an enhancement of the Euler method, which approximates 
+            <emphasis>
+                y<subscript>n+1</subscript>
+            </emphasis>
+            by truncating the Taylor expansion.
+        </para>
+        <para>
+            By convention, to use fixed-size steps, the program first computes a fitting <emphasis>h</emphasis> that approaches the simulation parameter <emphasis>abstol</emphasis>.
+        </para>
+        <para>
+            An important difference of Runge-Kutta with the previous methods is that it computes up to the fourth derivative of <emphasis>y</emphasis>, while the others only use linear combinations of <emphasis>y</emphasis> and <emphasis>y'</emphasis>.
+        </para>
+        <para>
+            Here, the next value is determined by the present value 
+            <emphasis>
+                y<subscript>n</subscript>
+            </emphasis>
+            plus the weighted average of four increments, where each increment is the product of the size of the interval, <emphasis>h</emphasis>, and an estimated slope specified by the function <emphasis>f(t,y)</emphasis> :
+            <itemizedlist>
+                <listitem>
+                    <emphasis>k1</emphasis> is the increment based on the slope at the beginning of the interval, using 
+                    <emphasis>
+                        y<subscript>n</subscript>
+                    </emphasis>
+                    (Euler's method),
+                </listitem>
+                <listitem>
+                    <emphasis>k2</emphasis>is the increment based on the slope at the midpoint of the interval, using 
+                    <emphasis>
+                        y<subscript>n</subscript> + h*k1/2
+                    </emphasis>
+                    ,
+                </listitem>
+                <listitem>
+                    <emphasis>k3</emphasis>is again the increment based on the slope at the midpoint, but now using 
+                    <emphasis>
+                        y<subscript>n</subscript> + h*k2/2
+                    </emphasis>
+                </listitem>
+                <listitem>
+                    <emphasis>k4</emphasis>is the increment based on the slope at the end of the interval, using 
+                    <emphasis>
+                        y<subscript>n</subscript> + h*k3
+                    </emphasis>
+                </listitem>
+            </itemizedlist>
+        </para>
+        <para>
+            We can see that with the <emphasis>ki</emphasis>, we progress in the derivatives of 
+            <emphasis>
+                y<subscript>n</subscript>
+            </emphasis>
+            .So in <emphasis>k4</emphasis>, we are approximating 
+            <emphasis>
+                y<superscript>(4)</superscript><subscript>n</subscript>
+            </emphasis>
+            ,thus making an error in 
+            <emphasis>
+                O(h<superscript>5</superscript>)
+            </emphasis>
+            .
+        </para>
+        <para>
+            So the total error is 
+            <emphasis>
+                number of steps * O(h<superscript>5</superscript>)
+            </emphasis>
+            .And since <emphasis>number of steps = interval size / h</emphasis> by definition, the total error is in 
+            <emphasis>
+                O(h<superscript>4</superscript>)
+            </emphasis>
+            .
+        </para>
+        <para>
+            That error analysis baptized the method <emphasis>Runge-Kutta 4(5)</emphasis>, 
+            <emphasis>
+                O(h<superscript>5</superscript>)
+            </emphasis>
+            per step, 
+            <emphasis>
+                O(h<superscript>4</superscript>)
+            </emphasis>
+            in total.
+        </para>
+        <para>
+            Althought the solver works fine for abstol up to <emphasis>abstol = 10^{-3}</emphasis>, rounding errors sometimes come into play as we approach 
+            <emphasis>
+                4*10<superscript>-4</superscript>
+            </emphasis>
+            .Indeed, the interval splitting can't be done properly and we get capricious results.
+        </para>
+    </refsection>
+    <refsection>
+        <title>Examples</title>
+        <para>
+            <link type="scilab" linkend="scilab.xcos/xcos/examples/solvers/ODE_Example.xcos">
+                <inlinemediaobject>
+                    <imageobject>
+                        <imagedata align="center" fileref="../../../examples/solvers/ODE_Example.xcos" valign="middle"/>
+                    </imageobject>
+                </inlinemediaobject>
+            </link>
+            <scilab:image><![CDATA[
+importXcosDiagram(SCI + "/modules/xcos/examples/solvers/ODE_Example.xcos");
+xcos_simulate(scs_m, 4);
+]]></scilab:image>
+        </para>
+        <para>
+            The integral block returns its continuous state, we can evaluate it with Runge-Kutta by running the example :
+        </para>
+        <para>
+            <programlisting language="example"><![CDATA[
+      // Import the diagram and set the ending time
+      importXcosDiagram("SCI/modules/xcos/examples/solvers/ODE_Example.xcos");
+      scs_m.props.tf = 3500;
+
+      // Select the solver Runge-Kutta and set the precision
+      scs_m.props.tol(6) = 4;
+      scs_m.props.tol(1) = 10^-2;
+
+      // Start the timer, launch the simulation and display time
+      timer();
+      xcos_simulate(scs_m, 4);
+      t = timer();
+      disp(t, "Time for Runge-Kutta :");
+      ]]></programlisting>
+        </para>
+        <para>
+            The Scilab console displays :
+            <screen><![CDATA[
+Time for Runge-Kutta :
+ 18.254
+            ]]></screen>
+        </para>
+        <para>
+            Now, in the following script, we compare the time difference between the methods by running the example with the five solvers in turn :
+            <link type="scilab" linkend ="scilab.scinotes/xcos/examples/solvers/integRK.sce">
+                Open the script
+            </link>
+        </para>
+        <para>
+            <screen><![CDATA[
+Time for BDF / Newton :
+ 25.457
+
+Time for BDF / Functional :
+ 24.893
+
+Time for Adams / Functional :
+ 20.049
+
+Time for Adams / Newton :
+ 20.469
+
+Time for Runge-Kutta :
+ 18.254
+            ]]></screen>
+        </para>
+        <para>
+            These results show that on a nonstiff problem, for relatively same precision required and forcing the same step size, Runge-Kutta is faster.
+        </para>
+        <para>
+            Variable step-size ODE solvers are not appropriate for deterministic real-time applications because the computational overhead of taking a time step varies over the course of an application.
+        </para>
+    </refsection>
+    <refsection>
+        <title>See Also</title>
+        <simplelist type="inline">
+            <member>
+                <link linkend="CVode">Cvode</link>
+            </member>
+            <member>
+                <link linkend="IDA">IDA</link>
+            </member>
+            <member>
+                <link linkend="ode">ode</link>
+            </member>
+            <member>
+                <link linkend="ode_discrete">ode_discrete</link>
+            </member>
+            <member>
+                <link linkend="ode_root">ode_root</link>
+            </member>
+            <member>
+                <link linkend="odedc">odedc</link>
+            </member>
+            <member>
+                <link linkend="impl">impl</link>
+            </member>
+        </simplelist>
+    </refsection>
+    <refsection>
+        <title>Bibliography</title>
+        <para>
+            <ulink url="https://computation.llnl.gov/casc/sundials/documentation/documentation.html">Sundials Documentation</ulink>
+        </para>
+    </refsection>
+</refentry>