Mercurial > octave
changeset 25026:f886561f9696 stable
doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
* NEWS: Announce new functions decic, ode15i, ode15s.
* diffeq.txi: Add decic, ode15i, ode15s to manual. Rewrite table of solvers.
* decic.m, ode15i.m, ode15s.m, ode23.m, ode45.m: Rewrite documentation.
author | Colin Macdonald <cbm@m.fsf.org> |
---|---|
date | Tue, 27 Mar 2018 00:52:29 -0700 |
parents | 7fb82487e828 |
children | c3cc9677db98 |
files | NEWS doc/interpreter/diffeq.txi scripts/ode/decic.m scripts/ode/ode15i.m scripts/ode/ode15s.m scripts/ode/ode23.m scripts/ode/ode45.m |
diffstat | 7 files changed, 141 insertions(+), 118 deletions(-) [+] |
line wrap: on
line diff
--- a/NEWS Tue Mar 27 14:03:52 2018 -0700 +++ b/NEWS Tue Mar 27 00:52:29 2018 -0700 @@ -262,6 +262,7 @@ camzoom corrcoef cosint + decic erase gammaincinv getframe @@ -275,6 +276,8 @@ isgraphics isstring mad + ode15i.m + ode15s.m openvar quad2d repelem
--- a/doc/interpreter/diffeq.txi Tue Mar 27 14:03:52 2018 -0700 +++ b/doc/interpreter/diffeq.txi Tue Mar 27 00:52:29 2018 -0700 @@ -136,25 +136,43 @@ @item @nospell{Runge-Kutta} methods @itemize - @item @code{ode45} Integrates a system of non--stiff ordinary differential equations - (non--stiff ODEs and DAEs) using second order @nospell{Dormand-Prince} - method. This is a fourth--order accurate integrator therefore the local - error normally expected is @math{O(h^5)}. This solver requires six - function evaluations per integration step. + @item @code{ode23} integrates a system of non-stiff ordinary differential + equations (ODEs) or index-1 differential-algebraic equations (DAEs). It + uses the third-order @nospell{Bogacki-Shampine} method and adapts the local + step size in order to satisfy a user-specified tolerance. The solver + requires three function evaluations per integration step. - @item @code{ode23} Integrates a system of non--stiff ordinary differential equations - (non-stiff ODEs and DAEs) using second order @nospell{Bogacki-Shampine} - method. This is a second-order accurate integrator therefore the local - error normally expected is @math{O(h^3)}. This solver requires three - function evaluations per integration step. + @item @code{ode45} integrates a system of non-stiff ODEs (or index-1 DAEs) + using the high-order, variable-step @nospell{Dormand-Prince} method. It + requires six function evaluations per integration step, but may + take larger steps on smooth problems than @code{ode23}: potentially + offering improved efficiency at smaller tolerances. + @end itemize + + @item Linear multistep methods + + @itemize + @item @code{ode15s} integrates a system of stiff ODEs (or index-1 DAEs) + using a variable step, variable order method based on Backward Difference + Formulas (BDF). + + @item @code{ode15i} integrates a system of fully-implicit ODEs (or index-1 + DAEs) using the same variable step, variable order method as @code{ode15s}. + The function @code{decic} can be used to compute consistent initial + conditions. @end itemize @end itemize - @DOCSTRING(ode45) @DOCSTRING(ode23) +@DOCSTRING(ode15s) + +@DOCSTRING(ode15i) + +@DOCSTRING(decic) + @DOCSTRING(odeset) @DOCSTRING(odeget)
--- a/scripts/ode/decic.m Tue Mar 27 14:03:52 2018 -0700 +++ b/scripts/ode/decic.m Tue Mar 27 00:52:29 2018 -0700 @@ -21,22 +21,24 @@ ## @deftypefnx {} {[@var{y0_new}, @var{yp0_new}] =} decic (@var{fun}, @var{t0}, @var{y0}, @var{fixed_y0}, @var{yp0}, @var{fixed_yp0}, @var{options}) ## @deftypefnx {} {[@var{y0_new}, @var{yp0_new}, @var{resnorm}] =} decic (@dots{}) ## -## Compute consistent initial conditions @var{y0_new} and @var{yp0_new}, given -## initial guesses @var{y0} and @var{yp0}. Choose a maximum of -## @code{length(@var{y0})} components between @var{fixed_y0} and -## @var{fixed_yp0} as fixed values. +## Compute consistent implicit ODE initial conditions @var{y0_new} and +## @var{yp0_new} given initial guesses @var{y0} and @var{yp0}. +## +## A maximum of @code{length (@var{y0})} components between @var{fixed_y0} and +## @var{fixed_yp0} may be chosen as fixed values. ## ## @var{fun} is a function handle. The function must accept three inputs where ## the first is time @var{t}, the second is a column vector of unknowns ## @var{y}, and the third is a column vector of unknowns @var{yp}. ## -## @var{t0} is the initial time such that @code{@var{fun}(@var{t0}, -## @var{y0_new}, @var{yp0_new}) = 0}, specified as a scalar. +## @var{t0} is the initial time such that +## @code{@var{fun}(@var{t0}, @var{y0_new}, @var{yp0_new}) = 0}, specified as a +## scalar. ## ## @var{y0} is a vector used as the initial guess for @var{y}. ## ## @var{fixed_y0} is a vector which specifies the components of @var{y0} to -## hold fixed. Choose a maximum of @code{length(@var{y0})} components between +## hold fixed. Choose a maximum of @code{length (@var{y0})} components between ## @var{fixed_y0} and @var{fixed_yp0} as fixed values. ## Set @var{fixed_y0}(i) component to 1 if you want to fix the value of ## @var{y0}(i). @@ -46,41 +48,41 @@ ## @var{yp0} is a vector used as the initial guess for @var{yp}. ## ## @var{fixed_yp0} is a vector which specifies the components of @var{yp0} to -## hold fixed. Choose a maximum of @code{length(@var{yp0})} components +## hold fixed. Choose a maximum of @code{length (@var{yp0})} components ## between @var{fixed_y0} and @var{fixed_yp0} as fixed values. ## Set @var{fixed_yp0}(i) component to 1 if you want to fix the value of ## @var{yp0}(i). ## Set @var{fixed_yp0}(i) component to 0 if you want to allow the value of ## @var{yp0}(i) to change. ## -## The optional seventh argument @var{options} is a structure array. -## Use @code{odeset} to generate this structure. The relevant options are +## The optional seventh argument @var{options} is a structure array. Use +## @code{odeset} to generate this structure. The relevant options are ## @code{RelTol} and @code{AbsTol} which specify the error thresholds used to ## compute the initial conditions. ## ## The function typically returns two outputs. Variable @var{y0_new} is a -## column vector and contains the consistent initial value of y. The -## output @var{yp0_new} is a column vector and contains the consistent initial -## value of yp. +## column vector and contains the consistent initial value of y. The output +## @var{yp0_new} is a column vector and contains the consistent initial value +## of yp. ## ## The optional third output @var{resnorm} is the norm of the vector of ## residuals. If @var{resnorm} is small, @code{decic} has successfully ## computed the initial conditions. If the value of @var{resnorm} is large, ## use @code{RelTol} and @code{AbsTol} to adjust it. ## -## Example: Compute initial conditions of @nospell{Robetson's} equations: +## Example: Compute initial conditions for @nospell{Robertson's} equations: ## -## @example +## @smallexample ## @group -## function r = robertsidae(@var{t}, @var{y}, @var{yp}) -## r = [-(@var{yp}(1) + 0.04*@var{y}(1) - 1e4*@var{y}(2)*@var{y}(3)); -## -(@var{yp}(2) - 0.04*@var{y}(1) + 1e4*@var{y}(2)*@var{y}(3) + 3e7*@var{y}(2)^2); -## @var{y}(1) + @var{y}(2) + @var{y}(3) - 1]; +## function r = robertson_dae (@var{t}, @var{y}, @var{yp}) +## r = [ -(@var{yp}(1) + 0.04*@var{y}(1) - 1e4*@var{y}(2)*@var{y}(3)) +## -(@var{yp}(2) - 0.04*@var{y}(1) + 1e4*@var{y}(2)*@var{y}(3) + 3e7*@var{y}(2)^2) +## @var{y}(1) + @var{y}(2) + @var{y}(3) - 1 ]; ## endfunction ## @end group -## [@var{y0_new},@var{yp0_new}] = decic (@@robertsidae, 0, [1; 0; 0], [1; 1; 0], +## [@var{y0_new},@var{yp0_new}] = decic (@@robertson_dae, 0, [1; 0; 0], [1; 1; 0], ## [-1e-4; 1; 0], [0; 0; 0]); -## @end example +## @end smallexample ## @seealso{ode15i, odeset} ## @end deftypefn
--- a/scripts/ode/ode15i.m Tue Mar 27 14:03:52 2018 -0700 +++ b/scripts/ode/ode15i.m Tue Mar 27 00:52:29 2018 -0700 @@ -23,16 +23,16 @@ ## @deftypefnx {} {@var{solution} =} ode15i (@dots{}) ## @deftypefnx {} {} ode15i (@dots{}) ## -## Solve a set of full-implicit Ordinary Differential Equations and -## Differential Algebraic Equations (DAEs) of index 1, with the variable-step, -## variable order BDF (Backward Differentiation Formula) method, which -## ranges from order 1 to 5. +## Solve a set of fully-implicit Ordinary Differential Equations (ODEs) or +## Differential Algebraic Equations (DAEs) of index 1, with a variable step, +## variable order BDF (Backward Differentiation Formula) method that ranges +## from order 1 to 5. ## ## @var{fun} is a function handle, inline function, or string containing the -## name of the function that defines the ODE: @code{y' = f(t,y,yp)}. The +## name of the function that defines the ODE: @code{0 = f(t,y,yp)}. The ## function must accept three inputs where the first is time @var{t}, the -## second is a column vector of unknowns @var{y}, and the third is a column -## vector of unknowns @var{yp}. +## second is the function value @var{y} (a column vector), and the third +## is the derivative value @var{yp} (a column vector). ## ## @var{trange} specifies the time interval over which the ODE will be ## evaluated. Typically, it is a two-element vector specifying the initial and @@ -46,8 +46,8 @@ ## value in @var{y0} and @var{yp0}. ## ## @var{y0} and @var{yp0} must be consistent initial conditions, meaning that -## @code{f(t,y0,yp0) = 0} is satisfied. You can use function @code{decic} to -## compute consistent initial conditions, given initial guesses. +## @code{f(t,y0,yp0) = 0} is satisfied. The function @code{decic} may be used +## to compute consistent initial conditions given initial guesses. ## ## The optional fifth argument @var{ode_opt} specifies non-default options to ## the ODE solver. It is a structure generated by @code{odeset}. @@ -57,46 +57,46 @@ ## output @var{y} is a matrix in which each column refers to a different ## unknown of the problem and each row corresponds to a time in @var{t}. ## -## The output can also be returned as a structure @var{solution} which -## has a field @var{x} containing a row vector of times where the solution -## was evaluated and a field @var{y} containing the solution matrix such -## that each column corresponds to a time in @var{x}. -## Use @code{fieldnames (@var{solution})} to see the other fields and +## The output can also be returned as a structure @var{solution} which has a +## field @var{x} containing a row vector of times where the solution was +## evaluated and a field @var{y} containing the solution matrix such that each +## column corresponds to a time in @var{x}. Use +## @w{@code{fieldnames (@var{solution})}} to see the other fields and ## additional information returned. ## -## If no output arguments are requested, and no @code{OutputFcn} is -## specified in @var{ode_opt}, then the @code{OutputFcn} is set to -## @code{odeplot} and the results of the solver are plotted immediately. +## If no output arguments are requested, and no @code{OutputFcn} is specified +## in @var{ode_opt}, then the @code{OutputFcn} is set to @code{odeplot} and the +## results of the solver are plotted immediately. ## -## If using the @qcode{"Events"} option then three additional outputs may -## be returned. @var{te} holds the time when an Event function returned a -## zero. @var{ye} holds the value of the solution at time @var{te}. @var{ie} +## If using the @qcode{"Events"} option then three additional outputs may be +## returned. @var{te} holds the time when an Event function returned a zero. +## @var{ye} holds the value of the solution at time @var{te}. @var{ie} ## contains an index indicating which Event function was triggered in the case ## of multiple Event functions. ## -## Example: Solve the @nospell{Robetson's} equations: +## Example: Solve @nospell{Robertson's} equations: ## -## @example +## @smallexample ## @group -## function r = robertsidae (@var{t}, @var{y}, @var{yp}) -## r = [-(@var{yp}(1) + 0.04*@var{y}(1) - 1e4*@var{y}(2)*@var{y}(3)); -## -(@var{yp}(2) - 0.04*@var{y}(1) + 1e4*@var{y}(2)*@var{y}(3) + 3e7*@var{y}(2)^2); -## @var{y}(1) + @var{y}(2) + @var{y}(3) - 1]; +## function r = robertson_dae (@var{t}, @var{y}, @var{yp}) +## r = [ -(@var{yp}(1) + 0.04*@var{y}(1) - 1e4*@var{y}(2)*@var{y}(3)) +## -(@var{yp}(2) - 0.04*@var{y}(1) + 1e4*@var{y}(2)*@var{y}(3) + 3e7*@var{y}(2)^2) +## @var{y}(1) + @var{y}(2) + @var{y}(3) - 1 ]; ## endfunction -## [@var{t},@var{y}] = ode15i (@@robertsidae, [0, 1e3], [1; 0; 0], [-1e-4; 1e-4; 0]); +## [@var{t},@var{y}] = ode15i (@@robertson_dae, [0, 1e3], [1; 0; 0], [-1e-4; 1e-4; 0]); ## @end group -## @end example +## @end smallexample ## @seealso{decic, odeset, odeget} ## @end deftypefn function varargout = ode15i (fun, trange, y0, yp0, varargin) - solver = "ode15i"; - if (nargin < 4) print_usage (); endif + solver = "ode15i"; + n = numel (y0); if (nargin > 4)
--- a/scripts/ode/ode15s.m Tue Mar 27 14:03:52 2018 -0700 +++ b/scripts/ode/ode15s.m Tue Mar 27 00:52:29 2018 -0700 @@ -23,10 +23,10 @@ ## @deftypefnx {} {@var{solution} =} ode15s (@dots{}) ## @deftypefnx {} {} ode15s (@dots{}) ## -## Solve a set of stiff Ordinary Differential Equations and stiff semi-explicit -## Differential Algebraic Equations (DAEs) of index 1, with the variable-step, -## variable order BDF (Backward Differentiation Formula) method, which ranges -## from order 1 to 5. +## Solve a set of stiff Ordinary Differential Equations (ODEs) or stiff +## semi-explicit Differential Algebraic Equations (DAEs) of index 1, with a +## variable step, variable order BDF (Backward Differentiation Formula) method +## that ranges from order 1 to 5. ## ## @var{fun} is a function handle, inline function, or string containing the ## name of the function that defines the ODE: @code{y' = f(t,y)}. The function @@ -51,47 +51,46 @@ ## output @var{y} is a matrix in which each column refers to a different ## unknown of the problem and each row corresponds to a time in @var{t}. ## -## The output can also be returned as a structure @var{solution} which -## has a field @var{x} containing a row vector of times where the solution -## was evaluated and a field @var{y} containing the solution matrix such -## that each column corresponds to a time in @var{x}. -## Use @code{fieldnames (@var{solution})} to see the other fields and +## The output can also be returned as a structure @var{solution} which has a +## field @var{x} containing a row vector of times where the solution was +## evaluated and a field @var{y} containing the solution matrix such that each +## column corresponds to a time in @var{x}. Use +## @w{@code{fieldnames (@var{solution})}} to see the other fields and ## additional information returned. ## -## If no output arguments are requested, and no @code{OutputFcn} is -## specified in @var{ode_opt}, then the @code{OutputFcn} is set to -## @code{odeplot} and the results of the solver are plotted immediately. +## If no output arguments are requested, and no @code{OutputFcn} is specified +## in @var{ode_opt}, then the @code{OutputFcn} is set to @code{odeplot} and the +## results of the solver are plotted immediately. ## -## If using the @qcode{"Events"} option then three additional outputs may -## be returned. @var{te} holds the time when an Event function returned a -## zero. @var{ye} holds the value of the solution at time @var{te}. @var{ie} +## If using the @qcode{"Events"} option then three additional outputs may be +## returned. @var{te} holds the time when an Event function returned a zero. +## @var{ye} holds the value of the solution at time @var{te}. @var{ie} ## contains an index indicating which Event function was triggered in the case ## of multiple Event functions. ## -## Example: Solve the @nospell{Robetson's} equations: +## Example: Solve @nospell{Robertson's} equations: ## -## @example +## @smallexample ## @group -## function r = robertsidae (@var{t}, @var{y}) -## r = [-0.04*@var{y}(1) + 1e4*@var{y}(2)*@var{y}(3); -## 0.04*@var{y}(1) - 1e4*@var{y}(2)*@var{y}(3) - 3e7*@var{y}(2)^2; -## @var{y}(1) + @var{y}(2) + @var{y}(3) - 1]; +## function r = robertson_dae (@var{t}, @var{y}) +## r = [ -0.04*@var{y}(1) + 1e4*@var{y}(2)*@var{y}(3) +## 0.04*@var{y}(1) - 1e4*@var{y}(2)*@var{y}(3) - 3e7*@var{y}(2)^2 +## @var{y}(1) + @var{y}(2) + @var{y}(3) - 1 ]; ## endfunction ## opt = odeset ("Mass", [1 0 0; 0 1 0; 0 0 0], "MStateDependence", "none"); -## [@var{t},@var{y}] = ode15s (@@robertsidae, [0, 1e3], [1; 0; 0], opt); +## [@var{t},@var{y}] = ode15s (@@robertson_dae, [0, 1e3], [1; 0; 0], opt); ## @end group -## @end example -## @seealso{decic, odeset, odeget} +## @end smallexample +## @seealso{decic, odeset, odeget, ode23, ode45} ## @end deftypefn function varargout = ode15s (fun, trange, y0, varargin) - solver = "ode15s"; - if (nargin < 3) print_usage (); endif + solver = "ode15s"; ## Check fun, trange, y0, yp0 fun = check_default_input (fun, trange, solver, y0);
--- a/scripts/ode/ode23.m Tue Mar 27 14:03:52 2018 -0700 +++ b/scripts/ode/ode23.m Tue Mar 27 00:52:29 2018 -0700 @@ -29,8 +29,6 @@ ## ## Solve a set of non-stiff Ordinary Differential Equations (non-stiff ODEs) ## with the well known explicit @nospell{Bogacki-Shampine} method of order 3. -## For the definition of this method see -## @url{http://en.wikipedia.org/wiki/List_of_Runge%E2%80%93Kutta_methods}. ## ## @var{fun} is a function handle, inline function, or string containing the ## name of the function that defines the ODE: @code{y' = f(t,y)}. The function @@ -45,8 +43,8 @@ ## ## By default, @code{ode23} uses an adaptive timestep with the ## @code{integrate_adaptive} algorithm. The tolerance for the timestep -## computation may be changed by using the options @qcode{"RelTol"} -## and @qcode{"AbsTol"}. +## computation may be changed by using the options @qcode{"RelTol"} and +## @qcode{"AbsTol"}. ## ## @var{init} contains the initial value for the unknowns. If it is a row ## vector then the solution @var{y} will be a matrix in which each column is @@ -60,20 +58,20 @@ ## output @var{y} is a matrix in which each column refers to a different ## unknown of the problem and each row corresponds to a time in @var{t}. ## -## The output can also be returned as a structure @var{solution} which -## has a field @var{x} containing a row vector of times where the solution -## was evaluated and a field @var{y} containing the solution matrix such -## that each column corresponds to a time in @var{x}. -## Use @code{fieldnames (@var{solution})} to see the other fields and +## The output can also be returned as a structure @var{solution} which has a +## field @var{x} containing a row vector of times where the solution was +## evaluated and a field @var{y} containing the solution matrix such that each +## column corresponds to a time in @var{x}. Use +## @w{@code{fieldnames (@var{solution})}} to see the other fields and ## additional information returned. ## -## If no output arguments are requested, and no @code{OutputFcn} is -## specified in @var{ode_opt}, then the @code{OutputFcn} is set to -## @code{odeplot} and the results of the solver are plotted immediately. +## If no output arguments are requested, and no @code{OutputFcn} is specified +## in @var{ode_opt}, then the @code{OutputFcn} is set to @code{odeplot} and the +## results of the solver are plotted immediately. ## -## If using the @qcode{"Events"} option then three additional outputs may -## be returned. @var{te} holds the time when an Event function returned a -## zero. @var{ye} holds the value of the solution at time @var{te}. @var{ie} +## If using the @qcode{"Events"} option then three additional outputs may be +## returned. @var{te} holds the time when an Event function returned a zero. +## @var{ye} holds the value of the solution at time @var{te}. @var{ie} ## contains an index indicating which Event function was triggered in the case ## of multiple Event functions. ## @@ -85,7 +83,10 @@ ## [@var{t},@var{y}] = ode23 (fvdp, [0, 20], [2, 0]); ## @end group ## @end example -## @seealso{odeset, odeget, ode45} +## +## Reference: For the definition of this method see +## @url{http://en.wikipedia.org/wiki/List_of_Runge%E2%80%93Kutta_methods}. +## @seealso{odeset, odeget, ode45, ode15s} ## @end deftypefn function varargout = ode23 (fun, trange, init, varargin)
--- a/scripts/ode/ode45.m Tue Mar 27 14:03:52 2018 -0700 +++ b/scripts/ode/ode45.m Tue Mar 27 00:52:29 2018 -0700 @@ -43,8 +43,8 @@ ## ## By default, @code{ode45} uses an adaptive timestep with the ## @code{integrate_adaptive} algorithm. The tolerance for the timestep -## computation may be changed by using the options @qcode{"RelTol"} -## and @qcode{"AbsTol"}. +## computation may be changed by using the options @qcode{"RelTol"} and +## @qcode{"AbsTol"}. ## ## @var{init} contains the initial value for the unknowns. If it is a row ## vector then the solution @var{y} will be a matrix in which each column is @@ -58,20 +58,20 @@ ## output @var{y} is a matrix in which each column refers to a different ## unknown of the problem and each row corresponds to a time in @var{t}. ## -## The output can also be returned as a structure @var{solution} which -## has a field @var{x} containing a row vector of times where the solution -## was evaluated and a field @var{y} containing the solution matrix such -## that each column corresponds to a time in @var{x}. -## Use @code{fieldnames (@var{solution})} to see the other fields and +## The output can also be returned as a structure @var{solution} which has a +## field @var{x} containing a row vector of times where the solution was +## evaluated and a field @var{y} containing the solution matrix such that each +## column corresponds to a time in @var{x}. Use +## @w{@code{fieldnames (@var{solution})}} to see the other fields and ## additional information returned. ## -## If no output arguments are requested, and no @code{OutputFcn} is -## specified in @var{ode_opt}, then the @code{OutputFcn} is set to -## @code{odeplot} and the results of the solver are plotted immediately. +## If no output arguments are requested, and no @code{OutputFcn} is specified +## in @var{ode_opt}, then the @code{OutputFcn} is set to @code{odeplot} and the +## results of the solver are plotted immediately. ## -## If using the @qcode{"Events"} option then three additional outputs may -## be returned. @var{te} holds the time when an Event function returned a -## zero. @var{ye} holds the value of the solution at time @var{te}. @var{ie} +## If using the @qcode{"Events"} option then three additional outputs may be +## returned. @var{te} holds the time when an Event function returned a zero. +## @var{ye} holds the value of the solution at time @var{te}. @var{ie} ## contains an index indicating which Event function was triggered in the case ## of multiple Event functions. ## @@ -83,7 +83,7 @@ ## [@var{t},@var{y}] = ode45 (fvdp, [0, 20], [2, 0]); ## @end group ## @end example -## @seealso{odeset, odeget, ode23} +## @seealso{odeset, odeget, ode23, ode15s} ## @end deftypefn function varargout = ode45 (fun, trange, init, varargin)