1 Dies ist odepkg.info, hergestellt von Makeinfo Version 4.8 aus
5 File: odepkg.info, Node: Top, Next: Beginners Guide, Prev: (dir), Up: (dir)
10 Copyright (C) 2006-2012, Thomas Treichl
12 Permission is granted to make and distribute verbatim copies of this
13 manual provided the copyright notice and this permission notice are
14 preserved on all copies.
16 Permission is granted to copy and distribute modified versions of
17 this manual under the conditions for verbatim copying, provided that
18 the entire resulting derived work is distributed under the terms of a
19 permission notice identical to this one.
21 Permission is granted to copy and distribute translations of this
22 manual into another language, under the same conditions as for modified
27 * Beginners Guide:: Manual for users who are completely new to OdePkg
28 * Users Guide:: Manual for users who are already familiar with OdePkg
29 * Programmers Guide:: Manual for users who want to make changes to OdePkg
30 * Function Index:: Reference about all functions from this package
31 * Index:: OdePkg Reference
34 File: odepkg.info, Node: Beginners Guide, Next: Users Guide, Prev: Top, Up: Top
39 The "Beginners Guide" is intended for users who are new to OdePkg and
40 who want to solve differential equations with the Octave language and
41 the package OdePkg. In this section it will be explained what OdePkg is
42 about in *Note About OdePkg:: and how OdePkg grew up from the beginning
43 in *Note OdePkg history and roadmap::. In *Note Installation and
44 deinstallation:: it is explained how OdePkg can be installed in Octave
45 and how it can later be removed from Octave if it is not needed
46 anymore. If you encounter problems while using OdePkg then have a look
47 at *Note Reporting Bugs:: how these bugs can be reported. In the *Note
48 The "foo" example:: a first example is explained.
52 * About OdePkg:: An introduction about OdePkg
53 * OdePkg history and roadmap:: From the first OdePkg release to the future
54 * Installation and deinstallation:: Setting up OdePkg on your system
55 * Reporting Bugs:: Writing comments and bugs to the help list
56 * The "foo" example:: A first example and where to go from here
59 File: odepkg.info, Node: About OdePkg, Next: OdePkg history and roadmap, Prev: Beginners Guide, Up: Beginners Guide
64 OdePkg is part of the Octave Repository (resp. the Octave-Forge
65 project) that was initiated by Matthew W. Roberts and Paul Kienzle in
66 the year 2000 and that is hosted at `http://octave.sourceforge.net'.
67 Since then a lot of contributors joined this project and added a lot
68 more packages and functions to further extend the capabilities of GNU
71 OdePkg includes commands for setting up various options, output
72 functions etc. before solving a set of differential equations with the
73 solver functions that are included. The package formerly was initiated
74 in autumn 2006 to solve ordinary differential equations (ODEs) only but
75 there are already improvements so that differential algebraic equations
76 (DAEs) in explicit form and in implicit form (IDEs) and delay
77 differential equations (DDEs) can also be solved. The goal of OdePkg is
78 to have a package for solving differential equations that is mostly
79 compatible to proprietary solver products.
82 File: odepkg.info, Node: OdePkg history and roadmap, Next: Installation and deinstallation, Prev: About OdePkg, Up: Beginners Guide
84 1.2 OdePkg history and roadmap
85 ==============================
87 OdePkg Version The initial release was a modification of the old
88 0.0.1 "ode" package that is hosted at Octave-Forge and that
89 was written by Marc Compere somewhen between 2000 and
90 2001. The four variable step-size Runge-Kutta
91 algorithms in three solver files and the three fixed
92 step-size solvers have been merged. It was possible
93 to set some options for these solvers. The four
94 output-functions (`odeprint', `odeplot', `odephas2'
95 and `odephas3') have been added along with other
96 examples that initially have not been there.
97 OdePkg Version The major milestone along versions 0.1.x was that
98 0.1.x four stable solvers have been implemented (ie.
99 `ode23', `ode45', `ode54' and `ode78') supporting all
100 options that can be set for these kind of solvers and
101 also all necessary functions for setting their
102 options (eg. `odeset', `odepkg_structure_check,
103 odepkg_event_handle'). Since version 0.1.3 there is
104 also source code available that interfaces the
105 Fortran solver `dopri5.f' (that is written by Ernst
106 Hairer and Gerhard Wanner, cf.
107 `odepkg_mexsolver_dopri5.c' and the helper files
108 `odepkgext.c' and `odepkgmex.c').
109 OdePkg Version The main work along version 0.2.x was to make the
110 0.2.x interface functions for the non-stiff and stiff
111 solvers from Ernst Hairer and Gerhard Wanner enough
112 stable so that they could be compiled and installed
113 by default. Wrapper functions have been added to the
114 package containing a help text and test functions
115 (eg. `ode2r', `ode5r', `oders'). Six testsuite
116 functions have been added to check the performance of
117 the different solvers (eg.
118 `odepkg_testsuite_chemakzo',
119 `odepkg_testsuite_oregonator').
120 OdePkg Version Fixed some minor bugs along version 0.3.x. Thanks to
121 0.3.x Jeff Cash, who released his Fortran `mebdfX' solvers
122 under the GNU GPL V2 after some discussion. The first
123 IDE solver `odebdi' appeared that is an interface
124 function for Cash's `mebdfi' Fortran core solver.
125 With version 0.3.5 of OdePkg a first new interface
126 function was created based on Octave's C++
127 `DEFUN_DLD' interface to achieve the highest
128 performance available. Added more examples and
129 testsuite functions (eg. `odepkg_equations_ilorenz',
130 `odepkg_testsuite_implrober'). Porting all at this
131 time present Mex-file solvers to Octave's C++
132 `DEFUN_DLD' interface. Ongoing work with this manual.
133 OdePkg Version Added a new solver function `odekdi' for the direct
134 0.4.x method (not the Krylov method) of the `daskr.f'
135 solver from the authors Peter N. Brown, Alan C.
136 Hindmarsh, Linda R. Petzold and Clement W. Ulrich
137 that is available under a modified BSD license
138 (without advertising clause). Ongoing work with this
140 OdePkg Version Added new solver functions `ode23d', `ode45d',
141 0.5.x `ode54d' and `ode78d' for solving non-stiff delay
142 differential equations (non-stiff DDEs). These
143 solvers are based on the Runge-Kutta solvers
144 `ode23'..`ode78'. Tests and demos have been included
145 for this type of solvers. Added new functions
146 `odeexamples', `odepkg_examples_ode',
147 `odepkg_examples_dae', `odepkg_examples_ide' and
148 `odepkg_examples_ide'. Ongoing work with this manual.
149 OdePkg Version A lot of compatibility tests, improvements, bugfixes,
151 (current) Version Final releases before version 1.0.0.
153 (future) Version Completed OdePkg release 1.0.0 with M-solvers and
157 File: odepkg.info, Node: Installation and deinstallation, Next: Reporting Bugs, Prev: OdePkg history and roadmap, Up: Beginners Guide
159 1.3 Installation and deinstallation
160 ===================================
162 OdePkg can be installed easily using the `pkg' command in Octave. To
163 install OdePkg download the latest release of OdePkg from the
164 Octave-Forge download site, then get into that directory where the
165 downloaded release of OdePkg has been saved, start Octave and type
166 pkg install odepkg-x.x.x.tar.gz
167 where `x.x.x' in the name of the `*.tar.gz' file is the current
168 release number of OdePkg that is available. If you want to deinstall
169 resp. remove OdePkg then simply type
171 and make sure that OdePkg has been removed completely and does not
172 appear in the list of installed packages anymore with the following
177 File: odepkg.info, Node: Reporting Bugs, Next: The "foo" example, Prev: Installation and deinstallation, Up: Beginners Guide
182 If you encounter problems during the installation process of OdePkg
183 with the `pkg' command or if you have an OdePkg that seems to be broken
184 or if you encounter problems while using OdePkg or if you find bugs in
185 the source codes then please report all of that via email at the
186 Octave-Forge mailing-list using the email address
187 <octave-dev@lists.sourceforge.net>. Not only bugs are welcome but also
188 any kind of comments are welcome (eg. if you think that OdePkg is
189 absolutely useful or even unnecessary).
192 File: odepkg.info, Node: The "foo" example, Prev: Reporting Bugs, Up: Beginners Guide
194 1.5 The "foo" example
195 =====================
197 Have a look at the first ordinary differential equation with the name
198 "`foo'". The `foo' equation of second order may be of the form y"(t) +
199 C_1 y'(t) + C_2 y(t) = C_3. With the substitutions y_1(t) = y(t) and
200 y_2(t) = y'(t) this differential equation of second order can be split
201 into two differential equations of first order, ie. y'_1(t) = y_2(t)
202 and y'_2(t) = - C_1 y_2(t) - C_2 y_1(t) + C_3. Next the numerical
203 values for the constants need to be defined, ie. C_1 = 2.0, C_2 = 5.0,
204 C_3 = 10.0. This set of ordinary differential equations can then be
205 written as an Octave M-file function like
206 function vdy = foo (vt, vy, varargin)
208 vdy(2,1) = - 2.0 * vy(2) - 5.0 * vy(1) + 10.0;
210 It can be seen that this ODEs do not depend on time, nevertheless
211 the first input argument of this function needs to be defined as the
212 time argument VT followed by a solution array argument `vy' as the
213 second input argument and a variable size input argument `varargin'
214 that can be used to set up user defined constants or control variables.
216 As it is known that `foo' is a set of ordinary differential
217 equations we can choose one of the four Runge-Kutta solvers (cf. *Note
218 Solver families::). It is also known that the time period of interest
219 may be between t_0 = 0.0 and t_e = 5.0 as well as that the initial
220 values of the ODEs are y_1(t=0) = 0.0 and y_2(t=0) = 0.0. Solving this
221 set of ODEs can be done by typing the following commands in Octave
222 ode45 (@foo, [0 5], [0 0]);
223 A figure window opens and it can be seen how this ODEs are solved
224 over time. For some of the solvers that come with OdePkg it is possible
225 to define exact time stamps for which an solution is required. Then the
226 example can be called eg.
227 ode45 (@foo, [0:0.1:5], [0 0]);
228 If it is not wanted that a figure window is opened while solving
229 then output arguments have to be used to catch the results of the
230 solving process and to not pass the results to the figure window, eg.
231 [t, y] = ode45 (@foo, [0 5], [0 0]);
232 Results can also be obtained in form of an Octave structure if one
233 output argument is used like in the following example. Then the results
234 are stored in the fields `S.x' and `S.y'.
235 S = ode45 (@foo, [0 5], [0 0]);
236 As noticed before, a function for the ordinary differential
237 equations must not be rewritten all the time if some of the parameters
238 are going to change. That's what the input argument `varargin' can be
239 used for. So rewrite the function `foo' into `newfoo' the following way
240 function vdy = newfoo (vt, vy, varargin)
242 vdy(2,1) = -varargin{1}*vy(2)-varargin{2}*vy(1)+varargin{3};
244 There is nothing said anymore about the constant values but if using
245 the following caller routine in the Octave window then the same results
246 can be obtained with the new function `newfoo' as before with the
247 function `foo' (ie. the parameters are directly feed through from the
248 caller routine `ode45' to the function `newfoo')
249 ode45 (@newfoo, [0 5], [0 0], 2.0, 5.0, 10.0);
250 OdePkg can do much more while solving differential equations, eg.
251 setting up other output functions instead of the function `odeplot' or
252 setting up other tolerances for the solving process etc. As a last
253 example in this beginning chapter it is shown how this can be done, ie.
254 with the command `odeset'
255 A = odeset ('OutputFcn', @odeprint);
256 ode45 (@newfoo, [0 5], [0 0], A, 2.0, 5.0, 10.0);
258 A = odeset ('OutputFcn', @odeprint, 'AbsTol', 1e-5, \
259 'RelTol', 1e-5, 'NormControl', 'on');
260 ode45 (@newfoo, [0 5], [0 0], A, 2.0, 5.0, 10.0);
261 The options structure `A' that can be set up with with the command
262 `odeset' must always be the fourth input argument when using the ODE
263 solvers and the DAE solvers but if you are using an IDE solver then `A'
264 must be the fifth input argument (cf. *Note Solver families::). The
265 various options that can be set with the command `odeset' are described
266 in *Note ODE/DAE/IDE/DDE options::.
268 Further examples have also been implemented. These example files and
269 functions are of the form `odepkg_examples_*'. Different testsuite
270 examples have been added that are stored in files with filenames
271 `odepkg_testsuite_*'. Before reading the next chapter note that nearly
272 every function that comes with OdePkg has its own help text and its own
273 examples. Look for yourself how the different functions, options and
274 combinations can be used. If you want to have a look at the help
275 description of a special function then type
277 in the Octave window where `fcnname' is the name of the function for
278 the help text to be viewed. Type
280 in the Octave window where `fcnname' is the name of the function for
281 the demo to run. Finally write
283 for opening this manual in the texinfo reader of the Octave window.
286 File: odepkg.info, Node: Users Guide, Next: Programmers Guide, Prev: Beginners Guide, Up: Top
291 The "Users Guide" is intended for trained users who already know in
292 principal how to solve differential equations with the Octave language
293 and OdePkg. In this chapter it will be explained which equations can be
294 solved with OdePkg in *Note Differential Equations::. It will be
295 explained which solvers can be used for the different kind of equations
296 in *Note Solver families:: and which options can be set for the
297 optimization of the solving process in *Note ODE/DAE/IDE/DDE options::.
298 The help text of all M-file functions and all Oct-file functions have
299 been extracted and are displayed in the sections *Note M-File Function
300 Reference:: and *Note Oct-File Function Reference::.
304 * Differential Equations:: The different kind of problems that can be solved with OdePkg
305 * Solver families:: The different kind of solvers within OdePkg
306 * ODE/DAE/IDE/DDE options:: The options that can be set for a solving process
307 * M-File Function Reference:: The description about all `*.m'-file functions
308 * Oct-File Function Reference:: The description about all DLD-functions from `*.oct'-files
311 File: odepkg.info, Node: Differential Equations, Next: Solver families, Prev: Users Guide, Up: Users Guide
313 2.1 Differential Equations
314 ==========================
316 In this section the different kind of differential equations that can
317 be solved with OdePkg are explained. The formulation of ordinary
318 differential equations is described in section *Note ODE equations::
319 followed by the description of explicetly formulated differential
320 algebraic equations in section *Note DAE equations::, implicetely
321 formulated differential algebraic equations in section *Note IDE
322 equations:: and delay differential algebraic equations in section *Note
327 * ODE equations:: Ordinary differential equations
328 * DAE equations:: Differential algebraic equations in explicit form
329 * IDE equations:: Differential algebraic equations in implicit form
330 * DDE equations:: Delay differential equations
333 File: odepkg.info, Node: ODE equations, Next: DAE equations, Prev: Differential Equations, Up: Differential Equations
338 ODE equations in general are of the form y'(t) = f(t,y) where y'(t) may
339 be a scalar or vector of derivatives. The variable t always is a scalar
340 describing one point of time and the variable y(t) is a scalar or
341 vector of solutions from the last time step of the set of ordinary
342 differential equations. If the equation is non-stiff then the *Note
343 Runge-Kutta solvers:: can be used to solve such kind of differential
344 equations but if the equation is stiff then it is recommended to use
345 the *Note Hairer-Wanner solvers::. An ODE equation definition in Octave
347 function [dy] = ODEequation (t, y, varargin)
350 File: odepkg.info, Node: DAE equations, Next: IDE equations, Prev: ODE equations, Up: Differential Equations
355 DAE equations in general are of the form M(t,y) \cdot y'(t) = f(t,y)
356 where y'(t) may be a scalar or vector of derivatives. The variable t
357 always is a scalar describing one point of time and the variable y(t)
358 is a scalar or vector of solutions from the set of differential
359 algebraic equations. The variable M(t,y) is the squared singular mass
360 matrix that may depend on y and t. If M(t,y) is not singular then the
361 set of equations from above can normally also be written as an ODE
362 equation. If it does not depend on time then it can be defined as a
363 constant matrix or a function. If it does depend on time then it must
364 be defined as a function. Use the command `odeset' to pass the mass
365 matrix information to the solver function (cf. *Note ODE/DAE/IDE/DDE
366 options::). If the equation is non-stiff then the *Note Runge-Kutta
367 solvers:: can be used to solve such kind of differential equations but
368 if the equation is stiff then it is recommended to use the *Note
369 Hairer-Wanner solvers::. A DAE equation definition in Octave must look
371 function [dy] = DAEequation (t, y, varargin)
372 and the mass matrix definition can either be a constant mass matrix
373 or a valid function handle to a mass matrix calculation function that
374 can be set with the command `odeset' (cf. option `Mass' of section
375 *Note ODE/DAE/IDE/DDE options::).
378 File: odepkg.info, Node: IDE equations, Next: DDE equations, Prev: DAE equations, Up: Differential Equations
383 IDE equations in general are of the form y'(t) + f(t,y) = 0 where y'(t)
384 may be a scalar or vector of derivatives. The variable t always is a
385 scalar describing one point of time and the variable y(t) is a scalar
386 or vector of solutions from the set of implicit differential equations.
387 Only IDE solvers from section *Note Cash modified BDF solvers:: or
388 section *Note DDaskr direct method solver:: can be used to solve such
389 kind of differential equations. A DAE equation definition in Octave
391 function [residual] = IDEequation (t, y, dy, varargin)
394 File: odepkg.info, Node: DDE equations, Prev: IDE equations, Up: Differential Equations
399 DDE equations in general are of the form y'(t) =
400 f(t,y(t),y(t-\tau_1),...,y(t-\tau_n)) where y'(t) may be a scalar or
401 vector of derivatives. The variable t always is a scalar describing one
402 point of time and the variables y(t-\tau_i) are scalars or vectors from
403 the past. Only DDE solvers from section *Note Modified Runge-Kutta
404 solvers:: can be used to solve such kind of differential equations. A
405 DDE equation definition in Octave must look like
406 function [dy] = DDEequation (t, y, z, varargin)
407 NOTE: Only DDEs with constant delays y(t-\tau_i) can be solved with
411 File: odepkg.info, Node: Solver families, Next: ODE/DAE/IDE/DDE options, Prev: Differential Equations, Up: Users Guide
416 In this section the different kind of solvers are introduced that have
417 been implemented in OdePkg. This section starts with the basic
418 Runge-Kutta solvers in section *Note Runge-Kutta solvers:: and is
419 continued with the Mex-file Hairer-Wanner solvers in section *Note
420 Hairer-Wanner solvers::. Performance tests have also been added to the
421 OdePkg. Some of these performance results have been added to section
422 *Note ODE solver performances::.
426 * Runge-Kutta solvers:: ODE solvers written as `*.m' files
427 * Hairer-Wanner solvers:: DAE solvers interfaced by `*.cc' files
428 * Cash modified BDF solvers:: A DAE and an IDE solver interfaced by `*.cc' files
429 * DDaskr direct method solver:: An IDE solver interfaced by a `*.cc' file
430 * Modified Runge-Kutta solvers:: DDE solvers written as `*.m' files
431 * ODE solver performances:: Cross math-engine performance tests
434 File: odepkg.info, Node: Runge-Kutta solvers, Next: Hairer-Wanner solvers, Prev: Solver families, Up: Solver families
436 2.2.1 Runge-Kutta solvers
437 -------------------------
439 The Runge-Kutta solvers are written in the Octave language and that are
440 saved as `m'-files. There have been implemented four different solvers
441 with a very similiar structure, ie. `ode23', `ode45', `ode54' and
442 `ode78'(1). The Runge-Kutta solvers have been added to the OdePkg to
443 solve non-stiff ODEs and DAEs, stiff equations of that form cannot be
444 solved with these solvers.
446 The order of all of the following Runge-Kutta methods is the order
447 of the local truncation error, which is the principle error term in the
448 portion of the Taylor series expansion that gets dropped, or
449 intentionally truncated. This is different from the local error which
450 is the difference between the estimated solution and the actual, or
451 true solution. The local error is used in stepsize selection and may be
452 approximated by the difference between two estimates of different order,
453 l(h) = x(O(h+1)) - x(O(h)). With this definition, the local error will
454 be as large as the error in the lower order method. The local
455 truncation error is within the group of terms that gets multipled by h
456 when solving for a solution from the general Runge-Kutta method.
457 Therefore, the order-p solution created by the Runge-Kunge method will
458 be roughly accurate to O(h^(p+1)) since the local truncation error
459 shows up in the solution as e = h\cdot d which is h-times an
460 O(h^p)-term, or rather O(h^(p+1)).
462 `ode23'Integrates a system of non-stiff ordinary differential equations
463 (non-stiff ODEs and DAEs) using second and third order Runge-Kutta
464 formulas. This particular third order method reduces to Simpson's
465 1/3 rule and uses the third order estimation for the output
466 solutions. Third order accurate Runge-Kutta methods have local and
467 global errors of O(h^4) and O(h^3) respectively and yield exact
468 results when the solution is a cubic (the variable h is the step
469 size from one integration step to another integration step). This
470 solver requires three function evaluations per integration step.
471 `ode45'Integrates a system of non-stiff ordinary differential equations
472 (non-stiff ODEs and DAEs) using fourth and fifth order embedded
473 formulas from Fehlberg. This is a fourth-order accurate integrator
474 therefore the local error normally expected is O(h^5). However,
475 because this particular implementation uses the fifth-order
476 estimate for x_out (ie. local extrapolation) moving forward with
477 the fifth-order estimate should yield local error of O(h^6). This
478 solver requires six function evaluations per integration step.
479 `ode54'Integrates a system of non-stiff ordinary differential equations
480 (non-stiff ODEs and DAEs) using fifth and fourth order Runge-Kutta
481 formulas. The Fehlberg 4(5) of the `ode45' pair is established and
482 works well, however, the Dormand-Prince 5(4) pair minimizes the
483 local truncation error in the fifth-order estimate which is what
484 is used to step forward (local extrapolation). Generally it
485 produces more accurate results and costs roughly the same
486 computationally. This solver requires seven function evaluations
487 per integration step.
488 `ode78'Integrates a system of non-stiff ordinary differential equations
489 (non-stiff ODEs and DAEs) using seventh and eighth order
490 Runge-Kutta formulas. This is a seventh-order accurate integrator
491 therefore the local error normally expected is O(h^8). However,
492 because this particular implementation uses the eighth-order
493 estimate for x_out moving forward with the eighth-order estimate
494 will yield errors on the order of O(h^9). This solver requires
495 thirteen function evaluations per integration step.
497 ---------- Footnotes ----------
499 (1) The descriptions for these Runge-Kutta solvers have been taken
500 from the help texts of the initial Runge-Kutta solvers that were
501 written by Marc Compere, he also pointed out that "a relevant
502 discussion on step size choice can be found on page 90ff in U.M.
503 Ascher, L.R. Petzold, Computer Methods for Ordinary Differential
504 Equations and Differential-Agebraic Equations, Society for Industrial
505 and Applied Mathematics (SIAM), Philadelphia, 1998".
508 File: odepkg.info, Node: Hairer-Wanner solvers, Next: Cash modified BDF solvers, Prev: Runge-Kutta solvers, Up: Solver families
510 2.2.2 Hairer-Wanner solvers
511 ---------------------------
513 The Hairer-Wanner solvers have been written by Ernst Hairer and Gerhard
514 Wanner. They are written in the Fortran language (hosted at
515 `http://www.unige.ch/~hairer') and that have been added to the OdePkg
516 as a compressed file with the name `hairer.tgz'. Papers and other
517 details about these solvers can be found at the adress given before.
518 The licence of these solvers is a modified BSD license (without
519 advertising clause and therefore are GPL compatible) and can be found
520 as `licence.txt' file in the `hairer.tgz' package. The Hairer-Wanner
521 solvers have been added to the OdePkg to solve non-stiff and stiff ODEs
522 and DAEs that cannot be solved with any of the Runge-Kutta solvers.
524 Interface functions for these solvers have been created and have
525 been added to the OdePkg. Their names are `odepkg_octsolver_XXX.cc'
526 where `XXX' is the name of the Fortran file that is interfaced. The
527 file `dldsolver.oct' is created automatically when installing OdePkg
528 with the command `pkg', but developers can also build each solver
529 manually with the instructions given as a preamble of every
530 `odepkg_octsolver_XXX.cc' file.
532 To provide a short name and to circumvent from the syntax of the
533 original solver function wrapper functions have been added, eg. the
534 command `ode2r' calls the solver `radau' from the Fortran file
535 `radau.f'. The other wrapper functions for the Hairer-Wanner solvers
536 are `ode5r' for the `radau5' solver, `oders' for the `rodas' solver and
537 `odesx' for the `seulex' solver. The help text of all these solver
538 functions can be diaplyed by calling `help wrapper' where wrapper is
539 one of `ode2r', `ode5r', `oders' or `odesx'.
542 File: odepkg.info, Node: Cash modified BDF solvers, Next: DDaskr direct method solver, Prev: Hairer-Wanner solvers, Up: Solver families
544 2.2.3 Cash modified BDF solvers
545 -------------------------------
547 The backward differentiation algorithm solvers have been written by
548 Jeff Cash in the Fortran language and that are hosted at
549 `http://pitagora.dm.uniba.it/~testset'. They have been added to the
550 OdePkg as a compressed file with the name `cash.tgz'. The license of
551 these solvers is a General Public License V2 that can be found as a
552 preamble of each Fortran solver source file. Papers and other details
553 about these solvers can be found at the host adress given before and
554 also at Jeff Cash's homepage at `http://www.ma.ic.ac.uk/~jcash'. Jeff
555 Cash's modified BDF solvers have been added to the OdePkg to solve
556 non-stiff and stiff ODEs and DAEs and also IDEs that cannot be solved
557 with any of the Runge-Kutta solvers.
559 Interface functions for these solvers have been created and have
560 been added to the OdePkg. Their names are `odepkg_octsolver_XXX.cc'
561 where `XXX' is the name of the Fortran file that is interfaced. The
562 file `dldsolver.oct' is created automatically when installing OdePkg
563 with the command `pkg', but developers can also build each solver
564 manually with the instructions given as a preamble of every
565 `odepkg_octsolver_XXX.cc' file.
567 To provide a short name and to circumvent from the syntax of the
568 original solver function wrapper functions have been added. The command
569 `odebda' calls the solver `mebdfdae' from the Fortran file `mebdf.f'
570 and the `odebdi' calls the solver `mebdfi' from the Fortran file
574 File: odepkg.info, Node: DDaskr direct method solver, Next: Modified Runge-Kutta solvers, Prev: Cash modified BDF solvers, Up: Solver families
576 2.2.4 DDaskr direct method solver
577 ---------------------------------
579 The direct method from the Krylov solver file `ddaskr.f' has been
580 written by Peter N. Brown, Alan C. Hindmarsh, Linda R. Petzold and
581 Clement W. Ulrich in the Fortran language and that is hosted at
582 `http://www.netlib.org'. The Krylov method has not been implemented
583 within OdePkg, only the direct method has been implemented. The solver
584 and further files for the interface have been added to the OdePkg as a
585 compressed package with the name `ddaskr.tgz'. The license of these
586 solvers is a modfied BSD license (without advertising clause) that can
587 be found inside of the compressed package. Other details about this
588 solver can be found as a preamble in the source file `ddaskr.f'. The
589 direct method solver of the file `ddaskr.f' has been added to the
590 OdePkg to solve non-stiff and stiff IDEs.
592 An interface function for this solver has been created and has been
593 added to the OdePkg. The source file name is
594 `odepkg_octsolver_ddaskr.cc'. The binary function can be found in the
595 file `dldsolver.oct' that is created automatically when installing
596 OdePkg with the command `pkg', but developers can also build the solver
597 wrapper manually with the instructions given as a preamble of the
598 `odepkg_octsolver_ddaskr.cc' file.
600 To provide a short name and to circumvent from the syntax of the
601 original solver function a wrapper function has been added. The command
602 `odekdi' calls the direct method of the solver `ddaskr' from the
603 Fortran file `ddaskr.f'.
606 File: odepkg.info, Node: Modified Runge-Kutta solvers, Next: ODE solver performances, Prev: DDaskr direct method solver, Up: Solver families
608 2.2.5 Modified Runge-Kutta solvers
609 ----------------------------------
611 The modified Runge-Kutta solvers are written in the Octave language and
612 that are saved as m-files. There have been implemented four different
613 solvers that do have a very similiar structure to that solvers found in
614 section *Note Runge-Kutta solvers::. Their names are `ode23d',
615 `ode45d', `ode54d' and `ode78d'. The modified Runge-Kutta solvers have
616 been added to the OdePkg to solve non-stiff DDEs with constant delays
617 only, stiff equations of that form cannot be solved with these solvers.
618 For further information about the error estimation of these solvers cf.
619 section *Note Runge-Kutta solvers::.
621 Note: The four DDE solvers of OdePkg are not syntax compatible to
622 propietary solvers. The reason is that the input arguments of the
623 propietary DDE-solvers are completely mixed up in comparison to ODE,
624 DAE and IDE propietary solvers. The DDE solvers of OdePkg have been
625 implemented in form of a syntax compatible way to the other family
626 solvers, eg. propietary solver calls look like
628 ode23 (@fode, vt, vy) %# for solving an ODE
629 ode15i (@fide, vt, vy, vdy) %# for solving an IDE
630 dde23 (@fdde, vlag, vhist, vt) %# for solving a DDE
631 whereas in OdePkg the same examples would look like
632 ode23 (@fode, vt, vy) %# for solving an ODE
633 odebdi (@fide, vt, vy, vdy) %# for solving an IDE
634 ode23d (@fdde, vt, vy, vlag, vhist) %# for solving a DDE
636 Further, the commands `ddeset' and `ddeget' have not been
637 implemented in OdePkg. Use the functions `odeset' and `odeget' for
638 setting and returning DDE options instead.
641 File: odepkg.info, Node: ODE solver performances, Prev: Modified Runge-Kutta solvers, Up: Solver families
643 2.2.6 ODE solver performances
644 -----------------------------
646 The following tables give an overview about the performance of the
647 OdePkg ODE/DAE solvers in comparison to propietary solvers when running
648 the HIRES function from the OdePkg testsuite (non-stiff ODE test).
650 >> odepkg ('odepkg_performance_mathires');
651 -----------------------------------------------------------------------------------------
652 Solver RelTol AbsTol Init Mescd Scd Steps Accept FEval JEval LUdec Time
653 -----------------------------------------------------------------------------------------
654 ode113 1e-007 1e-007 1e-009 7.57 5.37 24317 21442 45760 11.697
655 ode23 1e-007 1e-007 1e-009 7.23 5.03 13876 13862 41629 2.634
656 ode45 1e-007 1e-007 1e-009 7.91 5.70 11017 10412 66103 2.994
657 ode15s 1e-007 1e-007 1e-009 7.15 4.95 290 273 534 8 59 0.070
658 ode23s 1e-007 1e-007 1e-009 6.24 4.03 702 702 2107 702 702 0.161
659 ode23t 1e-007 1e-007 1e-009 6.00 3.79 892 886 1103 5 72 0.180
660 ode23tb 1e-007 1e-007 1e-009 5.85 3.65 735 731 2011 5 66 0.230
661 -----------------------------------------------------------------------------------------
663 octave:1> odepkg ('odepkg_performance_octavehires');
664 -----------------------------------------------------------------------------------------
665 Solver RelTol AbsTol Init Mescd Scd Steps Accept FEval JEval LUdec Time
666 -----------------------------------------------------------------------------------------
667 ode23 1e-07 1e-07 1e-09 7.86 5.44 17112 13369 51333 138.071
668 ode45 1e-07 1e-07 1e-09 8.05 5.63 9397 9393 56376 92.065
669 ode54 1e-07 1e-07 1e-09 8.25 5.83 9300 7758 65093 84.319
670 ode78 1e-07 1e-07 1e-09 8.54 6.12 7290 6615 94757 97.746
671 ode2r 1e-07 1e-07 1e-09 7.69 5.27 50 50 849 50 59 0.624
672 ode5r 1e-07 1e-07 1e-09 7.55 5.13 71 71 671 71 81 0.447
673 oders 1e-07 1e-07 1e-09 7.08 4.66 138 138 828 138 138 0.661
674 odesx 1e-07 1e-07 1e-09 6.56 4.13 30 26 1808 26 205 1.057
675 odebda 1e-07 1e-07 1e-09 6.53 4.11 401 400 582 42 42 0.378
676 -----------------------------------------------------------------------------------------
678 The following tables give an overview about the performance of the
679 OdePkg ODE/DAE solvers in comparison to propietary solvers when running
680 the chemical AKZO-NOBEL function from the OdePkg testsuite (non-stiff
683 >> odepkg ('odepkg_performance_matchemakzo');
684 -----------------------------------------------------------------------------------------
685 Solver RelTol AbsTol Init Mescd Scd Steps Accept FEval JEval LUdec Time
686 -----------------------------------------------------------------------------------------
687 ode113 1e-007 1e-007 1e-007 NaN Inf - - - - - -
688 ode23 1e-007 1e-007 1e-007 NaN Inf 15 15 47 0.431
689 ode45 1e-007 1e-007 1e-007 NaN Inf 15 15 92 0.170
690 ode15s 1e-007 1e-007 1e-007 7.04 6.20 161 154 4 35 0.521
691 ode23s 1e-007 1e-007 1e-007 7.61 6.77 1676 1676 5029 1676 1677 2.704
692 ode23t 1e-007 1e-007 1e-007 5.95 5.11 406 404 3 39 0.611
693 ode23tb 1e-007 1e-007 1e-007 NaN Inf 607 3036 1 608 6.730
694 -----------------------------------------------------------------------------------------
696 octave:1> odepkg ('odepkg_performance_octavechemakzo');
697 -----------------------------------------------------------------------------------------
698 Solver RelTol AbsTol Init Mescd Scd Steps Accept FEval JEval LUdec Time
699 -----------------------------------------------------------------------------------------
700 ode23 1e-07 1e-07 1e-07 2.95 2.06 424 385 1269 1.270
701 ode45 1e-07 1e-07 1e-07 2.95 2.06 256 218 1530 1.281
702 ode54 1e-07 1e-07 1e-07 2.95 2.06 197 195 1372 1.094
703 ode78 1e-07 1e-07 1e-07 2.95 2.06 184 156 2379 1.933
704 ode2r 1e-07 1e-07 1e-07 8.50 7.57 39 39 372 39 43 0.280
705 ode5r 1e-07 1e-07 1e-07 8.50 7.57 39 39 372 39 43 0.238
706 oders 1e-07 1e-07 1e-07 7.92 7.04 67 66 401 66 67 0.336
707 odesx 1e-07 1e-07 1e-07 7.19 6.26 19 19 457 19 82 0.248
708 odebda 1e-07 1e-07 1e-07 7.47 6.54 203 203 307 25 25 0.182
709 -----------------------------------------------------------------------------------------
711 Other testsuite functions have been added to the OdePkg that can be
712 taken for further performance tests and syntax checks on your own
713 hardware. These functions all have a name `odepkg_testsuite_XXX.m' with
714 `XXX' being the name of the testsuite equation that has been
718 File: odepkg.info, Node: ODE/DAE/IDE/DDE options, Next: M-File Function Reference, Prev: Solver families, Up: Users Guide
720 2.3 ODE/DAE/IDE/DDE options
721 ===========================
723 The default values of an OdePkg options structure can be displayed with
724 the command `odeset'. If `odeset' is called without any input argument
725 and one output argument then a OdePkg options structure with default
726 values is created, eg.
729 There also is an command `odeget' which extracts one or more options
730 from an OdePkg options structure. Other values than default values can
731 also be set with the command `odeset'. The function description of the
732 commands `odeset' and `odeget' can be found in the *Note M-File
733 Function Reference::. The values that can be set with the `odeset'
737 The option `RelTol' is used to set the relative error tolerance
738 for the error estimation of the solver that is used while solving.
739 It can either be a positive scalar or a vector with every element
740 of the vector being a positive scalar (this depends on the solver
741 that is used if both variants are supported). The definite error
742 estimation equation also depends on the solver that is used but
743 generalized (eg. for the solvers `ode23', `ode45', `ode54' and
744 `ode78') it may be a form like e(t) = max (r_tol^T y(t), a_tol).
745 Run the following example to illustrate the effect if this option
747 function yd = fvanderpol (vt, vy, varargin)
748 mu = 1; ## Set mu > 10 for higher stiffness
749 yd = [vy(2); mu * (1 - vy(1)^2) * vy(2) - vy(1)];
752 A = odeset ("RelTol", 1, "OutputFcn", @odeplot);
753 ode78 (@fvanderpol, [0 20], [2 0], A);
754 B = odeset (A, "RelTol", 1e-10);
755 ode78 (@fvanderpol, [0 20], [2 0], B);
758 The option `AbsTol' is used to set the absolute error tolerance
759 for the error estimation of the solver that is used while solving.
760 It can either be a positive scalar or a vector with every element
761 of the vector being a positive scalar (it depends on the solver
762 that is used if both variants are supported). The definite error
763 estimation equation also depends on the solver that is used but
764 generalized (eg. for the solvers `ode23', `ode45', `ode54' and
765 `ode78') it may be a form like e(t) = max (r_tol^T y(t), a_tol).
766 Run the following example to illustrate the effect if this option
768 ## An anonymous implementation of the Van der Pol equation
769 fvdb = @(vt,vy) [vy(2); (1 - vy(1)^2) * vy(2) - vy(1)];
771 A = odeset ("AbsTol", 1e-3, "OutputFcn", @odeplot);
772 ode54 (fvdb, [0 10], [2 0], A);
773 B = odeset (A, "AbsTol", 1e-10);
774 ode54 (fvdb, [0 10], [2 0], B);
777 The option `NormControl' is used to set the type of error
778 tolerance calculation of the solver that is used while solving. It
779 can either be the string `"on"' or `"off"'. At the time the solver
780 starts solving a warning message may be displayed if the solver
781 will ignore the `"on"' setting of this option because of an
782 unhandled resp. missing implementation. If set `"on"' then the
783 definite error estimation equation also depends on the solver that
784 is used but generalized (eg. for the solvers `ode23', `ode45',
785 `ode54' and `ode78') it may be a form like e(t) = max (r_tol^T max
786 (norm (y(t), \infty)), a_tol). Run the following example to
787 illustrate the effect if this option is used
788 function yd = fvanderpol (vt, vy, varargin)
789 mu = 1; ## Set mu > 10 for higher stiffness
790 yd = [vy(2); mu * (1 - vy(1)^2) * vy(2) - vy(1)];
793 A = odeset ("NormControl", "on", "OutputFcn", @odeplot);
794 ode78 (@fvanderpol, [0 20], [2 0], A);
795 B = odeset (A, "NormControl", "off");
796 ode78 (@fvanderpol, [0 20], [2 0], B);
799 The option `MaxStep' is used to set the maximum step size for the
800 solver that is used while solving. It can only be a positive
801 scalar. By default this value is set internally by every solver
802 and also may differ when using different solvers. Run the
803 following example to illustrate the effect if this option is used
804 function yd = fvanderpol (vt, vy, varargin)
805 mu = 1; ## Set mu > 10 for higher stiffness
806 yd = [vy(2); mu * (1 - vy(1)^2) * vy(2) - vy(1)];
809 A = odeset ("MaxStep", 10, "OutputFcn", @odeprint);
810 ode78 (@fvanderpol, [0 20], [2 0], A);
811 B = odeset (A, "MaxStep", 1e-1);
812 ode78 (@fvanderpol, [0 20], [2 0], B);
815 The option `InitialStep' is used to set the initial first step
816 size for the solver. It can only be a positive scalar. By default
817 this value is set internally by every solver and also may be
818 different when using different solvers. Run the following example
819 to illustrate the effect if this option is used
820 function yd = fvanderpol (vt, vy, varargin)
821 mu = 1; ## Set mu > 10 for higher stiffness
822 yd = [vy(2); mu * (1 - vy(1)^2) * vy(2) - vy(1)];
825 A = odeset ("InitialStep", 1, "OutputFcn", @odeprint);
826 ode78 (@fvanderpol, [0 1], [2 0], A);
827 B = odeset (A, "InitialStep", 1e-5);
828 ode78 (@fvanderpol, [0 1], [2 0], B);
831 The option `InitialSlope' is not handled by any of the solvers by
834 The option `OutputFcn' can be used to set up an output function
835 for displaying the results of the solver while solving. It must be
836 a function handle to a valid function. There are four predefined
837 output functions available with OdePkg. `odeprint' prints the
838 actual time values and results in the Octave window while solving,
839 `odeplot' plots the results over time in a new figure window while
840 solving, `odephas2' plots the first result over the second result
841 as a two-dimensional plot while solving and `odephas3' plots the
842 first result over the second result over the third result as a
843 three-dimensional plot while solving. Run the following example to
844 illustrate the effect if this option is used
845 function yd = fvanderpol (vt, vy, varargin)
846 mu = 1; ## Set mu > 10 for higher stiffness
847 yd = [vy(2); mu * (1 - vy(1)^2) * vy(2) - vy(1)];
850 A = odeset ("OutputFcn", @odeprint);
851 ode78 (@fvanderpol, [0 2], [2 0], A);
852 User defined output functions can also be used. A typical
853 framework for a self-made output function may then be of the form
854 function [vret] = odeoutput (vt, vy, vdeci, varargin)
857 ## Do everything needed to intialize output function
859 ## Do everything needed to create output
861 ## Do everything needed to clean up output function
864 The output function `odeplot' is also set automatically if the
865 solver calculation routine is called without any output argument.
866 Run the following example to illustrate the effect if this option
867 is not used and no output argument is given
868 ## An anonymous implementation of the Van der Pol equation
869 fvdb = @(vt,vy) [vy(2); (1 - vy(1)^2) * vy(2) - vy(1)];
871 ode78 (fvdb, [0 20], [2 0]);
874 The option `Refine' is used to set the interpolation factor that
875 is used to increase the quality for the output values if an output
876 function is also set with the option `OutputFcn'. It can only be a
877 integer value 0<=`Refine'<=5. Run the following example to
878 illustrate the effect if this option is used
879 function yd = fvanderpol (vt, vy, varargin)
880 mu = 1; ## Set mu > 10 for higher stiffness
881 yd = [vy(2); mu * (1 - vy(1)^2) * vy(2) - vy(1)];
884 A = odeset ("Refine", 0, "OutputFcn", @odeplot);
885 ode78 (@fvanderpol, [0 20], [2 0], A);
886 B = odeset (A, "Refine", 3);
887 ode78 (@fvanderpol, [0 20], [2 0], B);
890 The option `OutputSel' is used to set the components for which
891 output has to be performed if an output function is also set with
892 the option `OutputFcn'. It can only be a vector of integer values.
893 Run the following example to illustrate the effect if this option
895 function yd = fvanderpol (vt, vy, varargin)
896 mu = 1; ## Set mu > 10 for higher stiffness
897 yd = [vy(2); mu * (1 - vy(1)^2) * vy(2) - vy(1)];
900 A = odeset ("OutputSel", [1, 2], "OutputFcn", @odeplot);
901 ode78 (@fvanderpol, [0 20], [2 0], A);
902 B = odeset (A, "OutputSel", [2]);
903 ode78 (@fvanderpol, [0 20], [2 0], B);
906 The option `Stats' is used to print cost statistics about the
907 solving process after solving has been finished. It can either be
908 the string `"on"' or `"off"'. Run the following example to
909 illustrate the effect if this option is used
910 function yd = fvanderpol (vt, vy, varargin)
911 mu = 1; ## Set mu > 10 for higher stiffness
912 yd = [vy(2); mu * (1 - vy(1)^2) * vy(2) - vy(1)];
915 A = odeset ("Stats", "off");
916 [a, b] = ode78 (@fvanderpol, [0 2], [2 0], A);
917 B = odeset ("Stats", "on");
918 [c, d] = ode78 (@fvanderpol, [0 2], [2 0], B);
919 The cost statistics can also be obtained if the solver calculation
920 routine is called with one output argument. The cost statistics
921 then are in the field `stats' of the output arguemnt structure.
922 Run the following example to illustrate the effect if this option
924 S = ode78 (@fvanderpol, [0 2], [2 0], B);
928 The option `Jacobian' can be used to set up an external Jacobian
929 function or Jacobian matrix for DAE solvers to achieve faster and
930 better results (ODE Runge-Kutta solvers do not need to handle a
931 Jacobian function handle or Jacobian matrix). It must either be a
932 function handle to a valid function or a full constant matrix of
933 size squared the dimension of the set of differential equations.
934 User defined Jacobian functions must have the form `function
935 [vjac] = fjac (vt, vy, varargin)'. Run the following example to
936 illustrate the effect if this option is used
937 function vdy = fpol (vt, vy, varargin)
938 vdy = [vy(2); (1 - vy(1)^2) * vy(2) - vy(1)];
941 function vr = fjac (vt, vy, varargin)
943 -1-2*vy(1)*vy(2), 1-vy(1)^2];
946 A = odeset ("Stats", "on");
947 B = ode5r (@fpol, [0 20], [2 0], A);
948 C = odeset (A, "Jacobian", @fjac);
949 D = ode5r (@fpol, [0 20], [2 0], C);
950 Note: The function definition for Jacobian calculations of IDE
951 equations must have the form `function [vjac, vdjc] = fjac (vt,
952 vy, vyd, varargin)'. Run the following example to illustrate the
953 effect if this option is used
954 function [vres] = fvanderpol (vt, vy, vyd, varargin)
955 vres = [vy(2) - vyd(1);
956 (1 - vy(1)^2) * vy(2) - vy(1) - vyd(2)];
959 function [vjac, vdjc] = fjacobian (vt, vy, vyd, varargin)
960 vjac = [0, 1; -1 - 2 * vy(1) * vy(2), 1 - vy(1)^2];
961 vdjc = [-1, 0; 0, -1];
964 vopt = odeset ("Jacobian", @fjacobian, "Stats", "on");
965 vsol = odebdi (@fvanderpol, [0, 20], [2; 0], [0; -2], vopt, 10);
968 The option `JPattern' is not handled by any of the solvers by now.
970 The option `Vectorized' is not handled by any of the solvers by
973 The option `Mass' can be used to set up an external Mass function
974 or Mass matrix for solving DAE equations. It depends on the solver
975 that is used if `Mass' is supported or not. It must either be a
976 function handle to a valid function or a full constant matrix of
977 size squared the dimension of the set of differential equations.
978 User defined Jacobian functions must have the form `function vmas
979 = fmas (vt, vy, varargin)'. Run the following example to
980 illustrate the effect if this option is used
981 function vdy = frob (t, y, varargin)
982 vdy(1,1) = -0.04*y(1)+1e4*y(2)*y(3);
983 vdy(2,1) = 0.04*y(1)-1e4*y(2)*y(3)-3e7*y(2)^2;
984 vdy(3,1) = y(1)+y(2)+y(3)-1;
987 function vmas = fmas (vt, vy, varargin)
988 vmas = [1, 0, 0; 0, 1, 0; 0, 0, 0];
991 A = odeset ("Mass", @fmas);
992 B = ode5r (@frob, [0 1e8], [1 0 0], A);
993 Note: The function definition for Mass calculations of DDE
994 equations must have the form `function vmas = fmas (vt, vy, vz,
997 The option `MStateDependence' can be used to set up the type of
998 the external Mass function for solving DAE equations if a Mass
999 function handle is set with the option `Mass'. It depends on the
1000 solver that is used if `MStateDependence' is supported or not. It
1001 must be a string of the form `"none"', `"weak"' or `"strong"'. Run
1002 the following example to illustrate the effect if this option is
1004 function vdy = frob (vt, vy, varargin)
1005 vdy(1,1) = -0.04*vy(1)+1e4*vy(2)*vy(3);
1006 vdy(2,1) = 0.04*vy(1)-1e4*vy(2)*vy(3)-3e7*vy(2)^2;
1007 vdy(3,1) = vy(1)+vy(2)+vy(3)-1;
1010 function vmas = fmas (vt, varargin)
1011 vmas = [1, 0, 0; 0, 1, 0; 0, 0, 0];
1014 A = odeset ("Mass", @fmas, "MStateDependence", "none");
1015 B = ode5r (@frob, [0 1e8], [1 0 0], A);
1016 User defined Mass functions must have the form as described before
1017 (ie. `function vmas = fmas (vt, varargin)' if the option
1018 `MStateDependence' was set to `"none"', otherwise the user defined
1019 Mass function must have the form `function vmas = fmas (vt, vy,
1020 varargin)' if the option `MStateDependence' was set to either
1021 `"weak"' or `"strong"'.
1023 The option `MvPattern' is not handled by any of the solvers by now.
1025 The option `MassSingular' is not handled by any of the solvers by
1028 The option `NonNegative' can be used to set solution variables to
1029 zero even if their real solution would be a negative value. It
1030 must be a vector describing the positions in the solution vector
1031 for which the option `NonNegative' should be used. Run the
1032 following example to illustrate the effect if this option is used
1033 vfun = @(vt,vy) -abs(vy);
1034 vopt = odeset ("NonNegative", [1]);
1036 [vt1, vy1] = ode78 (vfun, [0 100], [1]);
1037 [vt2, vy2] = ode78 (vfun, [0 100], [1], vopt);
1039 subplot (2,1,1); plot (vt1, vy1);
1040 subplot (2,1,2); plot (vt2, vy2);
1043 The option `Events' can be used to set up an Event function, ie.
1044 the Event function can be used to find zero crossings in one of
1045 the results. It must either be a function handle to a valid
1046 function. Run the following example to illustrate the effect if
1048 function vdy = fbal (vt, vy, varargin)
1050 vdy(2,1) = -9.81; ## m/s²
1053 function [veve, vterm, vdir] = feve (vt, vy, varargin)
1054 veve = vy(1); ## Which event component should be tread
1055 vterm = 1; ## Terminate if an event is found
1056 vdir = -1; ## In which direction, -1 for falling
1059 A = odeset ("Events", @feve);
1060 B = ode78 (@fbal, [0 1.5], [1 3], A);
1061 plot (B.x, B.y(:,1));
1062 Note: The function definition for Events calculations of DDE
1063 equations must have the form `function [veve, vterm, vdir] = feve
1064 (vt, vy, vz, varargin)' and the function definition for Events
1065 calculations of IDE equations must have the form `function [veve,
1066 vterm, vdir] = feve (vt, vy, vyd, varargin)'.
1068 The option `MaxOrder' can be used to set the maximum order of the
1069 backward differentiation algorithm of the `odebdi' and `odebda'
1070 solvers. It must be a scalar integer value between 1 and 7. Run
1071 the following example to illustrate the effect if this option is
1073 function res = fwei (t, y, yp, varargin)
1074 res = t*y^2*yp^3 - y^3*yp^2 + t*yp*(t^2 + 1) - t^2*y;
1077 function [dy, dyp] = fjac (t, y, yp, varargin)
1078 dy = 2*t*y*yp^3 - 3*y^2*yp^2 - t^2;
1079 dyp = 3*t*y^2*yp^2 - 2*y^3*yp + t*(t^2 + 1);
1082 A = odeset ("AbsTol", 1e-6, "RelTol", 1e-6, "Jacobian", @fjac, ...
1083 "Stats", "on", "MaxOrder", 1, "BDF", "on")
1084 B = odeset (A, "MaxOrder", 5)
1085 C = odebdi (@fwei, [1 10], 1.2257, 0.8165, A);
1086 D = odebdi (@fwei, [1 10], 1.2257, 0.8165, B);
1087 plot (C.x, C.y, "bo-", D.x, D.y, "rx:");
1090 The option `BDF' is only supported by the `odebdi' and `odebda'
1091 solvers. Using these solvers the option `BDF' will automatically
1092 be set `"on"' (even if it was set `"off"' before) because the
1093 `odebdi' and `odebda' solvers all use the backward differentiation
1094 algorithm to solve the different kind of equations.
1099 `MaxNewtonIterations'
1103 File: odepkg.info, Node: M-File Function Reference, Next: Oct-File Function Reference, Prev: ODE/DAE/IDE/DDE options, Up: Users Guide
1105 2.4 M-File Function Reference
1106 =============================
1108 The help texts of this section are autogenerated and refer to commands
1109 that all can be found in the files `*.m'. All commands that are listed
1110 below are loaded automatically everytime you launch Octave.
1111 -- Function File: [] = ode23 (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2,
1113 -- Command: [SOL] = ode23 (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
1114 -- Command: [T, Y, [XE, YE, IE]] = ode23 (@FUN, SLOT, INIT, [OPT],
1116 This function file can be used to solve a set of non-stiff
1117 ordinary differential equations (non-stiff ODEs) or non-stiff
1118 differential algebraic equations (non-stiff DAEs) with the well
1119 known explicit Runge-Kutta method of order (2,3).
1121 If this function is called with no return argument then plot the
1122 solution over time in a figure window while solving the set of
1123 ODEs that are defined in a function and specified by the function
1124 handle @FUN. The second input argument SLOT is a double vector
1125 that defines the time slot, INIT is a double vector that defines
1126 the initial values of the states, OPT can optionally be a
1127 structure array that keeps the options created with the command
1128 `odeset' and PAR1, PAR2, ... can optionally be other input
1129 arguments of any type that have to be passed to the function
1132 If this function is called with one return argument then return
1133 the solution SOL of type structure array after solving the set of
1134 ODEs. The solution SOL has the fields X of type double column
1135 vector for the steps chosen by the solver, Y of type double column
1136 vector for the solutions at each time step of X, SOLVER of type
1137 string for the solver name and optionally the extended time stamp
1138 information XE, the extended solution information YE and the
1139 extended index information IE all of type double column vector
1140 that keep the informations of the event function if an event
1141 function handle is set in the option argument OPT.
1143 If this function is called with more than one return argument then
1144 return the time stamps T, the solution values Y and optionally the
1145 extended time stamp information XE, the extended solution
1146 information YE and the extended index information IE all of type
1147 double column vector.
1149 For example, solve an anonymous implementation of the Van der Pol
1152 fvdb = @(vt,vy) [vy(2); (1 - vy(1)^2) * vy(2) - vy(1)];
1154 vopt = odeset ("RelTol", 1e-3, "AbsTol", 1e-3, \
1155 "NormControl", "on", "OutputFcn", @odeplot);
1156 ode23 (fvdb, [0 20], [2 0], vopt);
1158 -- Function File: [] = ode23d (@FUN, SLOT, INIT, LAGS, HIST, [OPT],
1160 -- Command: [SOL] = ode23d (@FUN, SLOT, INIT, LAGS, HIST, [OPT],
1162 -- Command: [T, Y, [XE, YE, IE]] = ode23d (@FUN, SLOT, INIT, LAGS,
1163 HIST, [OPT], [PAR1, PAR2, ...])
1164 This function file can be used to solve a set of non-stiff delay
1165 differential equations (non-stiff DDEs) with a modified version of
1166 the well known explicit Runge-Kutta method of order (2,3).
1168 If this function is called with no return argument then plot the
1169 solution over time in a figure window while solving the set of
1170 DDEs that are defined in a function and specified by the function
1171 handle @FUN. The second input argument SLOT is a double vector
1172 that defines the time slot, INIT is a double vector that defines
1173 the initial values of the states, LAGS is a double vector that
1174 describes the lags of time, HIST is a double matrix and describes
1175 the history of the DDEs, OPT can optionally be a structure array
1176 that keeps the options created with the command `odeset' and PAR1,
1177 PAR2, ... can optionally be other input arguments of any type that
1178 have to be passed to the function defined by @FUN.
1180 In other words, this function will solve a problem of the form
1181 dy/dt = fun (t, y(t), y(t-lags(1), y(t-lags(2), ...)))
1183 y(slot(1)-lags(1)) = hist(1), y(slot(1)-lags(2)) = hist(2), ...
1185 If this function is called with one return argument then return
1186 the solution SOL of type structure array after solving the set of
1187 DDEs. The solution SOL has the fields X of type double column
1188 vector for the steps chosen by the solver, Y of type double column
1189 vector for the solutions at each time step of X, SOLVER of type
1190 string for the solver name and optionally the extended time stamp
1191 information XE, the extended solution information YE and the
1192 extended index information IE all of type double column vector
1193 that keep the informations of the event function if an event
1194 function handle is set in the option argument OPT.
1196 If this function is called with more than one return argument then
1197 return the time stamps T, the solution values Y and optionally the
1198 extended time stamp information XE, the extended solution
1199 information YE and the extended index information IE all of type
1200 double column vector.
1203 - the following code solves an anonymous implementation of a
1206 fcao = @(vt, vy, vz) [2 * vz / (1 + vz^9.65) - vy];
1208 vopt = odeset ("NormControl", "on", "RelTol", 1e-3);
1209 vsol = ode23d (fcao, [0, 100], 0.5, 2, 0.5, vopt);
1211 vlag = interp1 (vsol.x, vsol.y, vsol.x - 2);
1212 plot (vsol.y, vlag); legend ("fcao (t,y,z)");
1214 - to solve the following problem with two delayed state
1218 d y2(t)/dt = -y2(t) + y1(t-5)
1219 d y3(t)/dt = -y3(t) + y2(t-10)*y1(t-10)
1221 one might do the following
1223 function f = fun (t, y, yd)
1224 f(1) = -y(1); %% y1' = -y1(t)
1225 f(2) = -y(2) + yd(1,1); %% y2' = -y2(t) + y1(t-lags(1))
1226 f(3) = -y(3) + yd(2,2)*yd(1,2); %% y3' = -y3(t) + y2(t-lags(2))*y1(t-lags(2))
1229 res = ode23d (@fun, T, [1;1;1], [5, 10], ones (3,2));
1232 -- Function File: [] = ode45 (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2,
1234 -- Command: [SOL] = ode45 (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
1235 -- Command: [T, Y, [XE, YE, IE]] = ode45 (@FUN, SLOT, INIT, [OPT],
1237 This function file can be used to solve a set of non-stiff
1238 ordinary differential equations (non-stiff ODEs) or non-stiff
1239 differential algebraic equations (non-stiff DAEs) with the well
1240 known explicit Runge-Kutta method of order (4,5).
1242 If this function is called with no return argument then plot the
1243 solution over time in a figure window while solving the set of
1244 ODEs that are defined in a function and specified by the function
1245 handle @FUN. The second input argument SLOT is a double vector
1246 that defines the time slot, INIT is a double vector that defines
1247 the initial values of the states, OPT can optionally be a
1248 structure array that keeps the options created with the command
1249 `odeset' and PAR1, PAR2, ... can optionally be other input
1250 arguments of any type that have to be passed to the function
1253 If this function is called with one return argument then return
1254 the solution SOL of type structure array after solving the set of
1255 ODEs. The solution SOL has the fields X of type double column
1256 vector for the steps chosen by the solver, Y of type double column
1257 vector for the solutions at each time step of X, SOLVER of type
1258 string for the solver name and optionally the extended time stamp
1259 information XE, the extended solution information YE and the
1260 extended index information IE all of type double column vector
1261 that keep the informations of the event function if an event
1262 function handle is set in the option argument OPT.
1264 If this function is called with more than one return argument then
1265 return the time stamps T, the solution values Y and optionally the
1266 extended time stamp information XE, the extended solution
1267 information YE and the extended index information IE all of type
1268 double column vector.
1270 For example, solve an anonymous implementation of the Van der Pol
1273 fvdb = @(vt,vy) [vy(2); (1 - vy(1)^2) * vy(2) - vy(1)];
1275 vopt = odeset ("RelTol", 1e-3, "AbsTol", 1e-3, \
1276 "NormControl", "on", "OutputFcn", @odeplot);
1277 ode45 (fvdb, [0 20], [2 0], vopt);
1279 -- Function File: [] = ode45d (@FUN, SLOT, INIT, LAGS, HIST, [OPT],
1281 -- Command: [SOL] = ode45d (@FUN, SLOT, INIT, LAGS, HIST, [OPT],
1283 -- Command: [T, Y, [XE, YE, IE]] = ode45d (@FUN, SLOT, INIT, LAGS,
1284 HIST, [OPT], [PAR1, PAR2, ...])
1285 This function file can be used to solve a set of non-stiff delay
1286 differential equations (non-stiff DDEs) with a modified version of
1287 the well known explicit Runge-Kutta method of order (4,5).
1289 If this function is called with no return argument then plot the
1290 solution over time in a figure window while solving the set of
1291 DDEs that are defined in a function and specified by the function
1292 handle @FUN. The second input argument SLOT is a double vector
1293 that defines the time slot, INIT is a double vector that defines
1294 the initial values of the states, LAGS is a double vector that
1295 describes the lags of time, HIST is a double matrix and describes
1296 the history of the DDEs, OPT can optionally be a structure array
1297 that keeps the options created with the command `odeset' and PAR1,
1298 PAR2, ... can optionally be other input arguments of any type that
1299 have to be passed to the function defined by @FUN.
1301 In other words, this function will solve a problem of the form
1302 dy/dt = fun (t, y(t), y(t-lags(1), y(t-lags(2), ...)))
1304 y(slot(1)-lags(1)) = hist(1), y(slot(1)-lags(2)) = hist(2), ...
1306 If this function is called with one return argument then return
1307 the solution SOL of type structure array after solving the set of
1308 DDEs. The solution SOL has the fields X of type double column
1309 vector for the steps chosen by the solver, Y of type double column
1310 vector for the solutions at each time step of X, SOLVER of type
1311 string for the solver name and optionally the extended time stamp
1312 information XE, the extended solution information YE and the
1313 extended index information IE all of type double column vector
1314 that keep the informations of the event function if an event
1315 function handle is set in the option argument OPT.
1317 If this function is called with more than one return argument then
1318 return the time stamps T, the solution values Y and optionally the
1319 extended time stamp information XE, the extended solution
1320 information YE and the extended index information IE all of type
1321 double column vector.
1324 - the following code solves an anonymous implementation of a
1327 fcao = @(vt, vy, vz) [2 * vz / (1 + vz^9.65) - vy];
1329 vopt = odeset ("NormControl", "on", "RelTol", 1e-3);
1330 vsol = ode45d (fcao, [0, 100], 0.5, 2, 0.5, vopt);
1332 vlag = interp1 (vsol.x, vsol.y, vsol.x - 2);
1333 plot (vsol.y, vlag); legend ("fcao (t,y,z)");
1335 - to solve the following problem with two delayed state
1339 d y2(t)/dt = -y2(t) + y1(t-5)
1340 d y3(t)/dt = -y3(t) + y2(t-10)*y1(t-10)
1342 one might do the following
1344 function f = fun (t, y, yd)
1345 f(1) = -y(1); %% y1' = -y1(t)
1346 f(2) = -y(2) + yd(1,1); %% y2' = -y2(t) + y1(t-lags(1))
1347 f(3) = -y(3) + yd(2,2)*yd(1,2); %% y3' = -y3(t) + y2(t-lags(2))*y1(t-lags(2))
1350 res = ode45d (@fun, T, [1;1;1], [5, 10], ones (3,2));
1353 -- Function File: [] = ode54 (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2,
1355 -- Command: [SOL] = ode54 (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
1356 -- Command: [T, Y, [XE, YE, IE]] = ode54 (@FUN, SLOT, INIT, [OPT],
1358 This function file can be used to solve a set of non-stiff
1359 ordinary differential equations (non-stiff ODEs) or non-stiff
1360 differential algebraic equations (non-stiff DAEs) with the well
1361 known explicit Runge-Kutta method of order (5,4).
1363 If this function is called with no return argument then plot the
1364 solution over time in a figure window while solving the set of
1365 ODEs that are defined in a function and specified by the function
1366 handle @FUN. The second input argument SLOT is a double vector
1367 that defines the time slot, INIT is a double vector that defines
1368 the initial values of the states, OPT can optionally be a
1369 structure array that keeps the options created with the command
1370 `odeset' and PAR1, PAR2, ... can optionally be other input
1371 arguments of any type that have to be passed to the function
1374 If this function is called with one return argument then return
1375 the solution SOL of type structure array after solving the set of
1376 ODEs. The solution SOL has the fields X of type double column
1377 vector for the steps chosen by the solver, Y of type double column
1378 vector for the solutions at each time step of X, SOLVER of type
1379 string for the solver name and optionally the extended time stamp
1380 information XE, the extended solution information YE and the
1381 extended index information IE all of type double column vector
1382 that keep the informations of the event function if an event
1383 function handle is set in the option argument OPT.
1385 If this function is called with more than one return argument then
1386 return the time stamps T, the solution values Y and optionally the
1387 extended time stamp information XE, the extended solution
1388 information YE and the extended index information IE all of type
1389 double column vector.
1391 For example, solve an anonymous implementation of the Van der Pol
1394 fvdb = @(vt,vy) [vy(2); (1 - vy(1)^2) * vy(2) - vy(1)];
1396 vopt = odeset ("RelTol", 1e-3, "AbsTol", 1e-3, \
1397 "NormControl", "on", "OutputFcn", @odeplot);
1398 ode54 (fvdb, [0 20], [2 0], vopt);
1400 -- Function File: [] = ode54d (@FUN, SLOT, INIT, LAGS, HIST, [OPT],
1402 -- Command: [SOL] = ode54d (@FUN, SLOT, INIT, LAGS, HIST, [OPT],
1404 -- Command: [T, Y, [XE, YE, IE]] = ode54d (@FUN, SLOT, INIT, LAGS,
1405 HIST, [OPT], [PAR1, PAR2, ...])
1406 This function file can be used to solve a set of non-stiff delay
1407 differential equations (non-stiff DDEs) with a modified version of
1408 the well known explicit Runge-Kutta method of order (2,3).
1410 If this function is called with no return argument then plot the
1411 solution over time in a figure window while solving the set of
1412 DDEs that are defined in a function and specified by the function
1413 handle @FUN. The second input argument SLOT is a double vector
1414 that defines the time slot, INIT is a double vector that defines
1415 the initial values of the states, LAGS is a double vector that
1416 describes the lags of time, HIST is a double matrix and describes
1417 the history of the DDEs, OPT can optionally be a structure array
1418 that keeps the options created with the command `odeset' and PAR1,
1419 PAR2, ... can optionally be other input arguments of any type that
1420 have to be passed to the function defined by @FUN.
1422 In other words, this function will solve a problem of the form
1423 dy/dt = fun (t, y(t), y(t-lags(1), y(t-lags(2), ...)))
1425 y(slot(1)-lags(1)) = hist(1), y(slot(1)-lags(2)) = hist(2), ...
1427 If this function is called with one return argument then return
1428 the solution SOL of type structure array after solving the set of
1429 DDEs. The solution SOL has the fields X of type double column
1430 vector for the steps chosen by the solver, Y of type double column
1431 vector for the solutions at each time step of X, SOLVER of type
1432 string for the solver name and optionally the extended time stamp
1433 information XE, the extended solution information YE and the
1434 extended index information IE all of type double column vector
1435 that keep the informations of the event function if an event
1436 function handle is set in the option argument OPT.
1438 If this function is called with more than one return argument then
1439 return the time stamps T, the solution values Y and optionally the
1440 extended time stamp information XE, the extended solution
1441 information YE and the extended index information IE all of type
1442 double column vector.
1445 - the following code solves an anonymous implementation of a
1448 fcao = @(vt, vy, vz) [2 * vz / (1 + vz^9.65) - vy];
1450 vopt = odeset ("NormControl", "on", "RelTol", 1e-3);
1451 vsol = ode54d (fcao, [0, 100], 0.5, 2, 0.5, vopt);
1453 vlag = interp1 (vsol.x, vsol.y, vsol.x - 2);
1454 plot (vsol.y, vlag); legend ("fcao (t,y,z)");
1456 - to solve the following problem with two delayed state
1460 d y2(t)/dt = -y2(t) + y1(t-5)
1461 d y3(t)/dt = -y3(t) + y2(t-10)*y1(t-10)
1463 one might do the following
1465 function f = fun (t, y, yd)
1466 f(1) = -y(1); %% y1' = -y1(t)
1467 f(2) = -y(2) + yd(1,1); %% y2' = -y2(t) + y1(t-lags(1))
1468 f(3) = -y(3) + yd(2,2)*yd(1,2); %% y3' = -y3(t) + y2(t-lags(2))*y1(t-lags(2))
1471 res = ode54d (@fun, T, [1;1;1], [5, 10], ones (3,2));
1474 -- Function File: [] = ode78 (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2,
1476 -- Command: [SOL] = ode78 (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
1477 -- Command: [T, Y, [XE, YE, IE]] = ode78 (@FUN, SLOT, INIT, [OPT],
1479 This function file can be used to solve a set of non-stiff
1480 ordinary differential equations (non-stiff ODEs) or non-stiff
1481 differential algebraic equations (non-stiff DAEs) with the well
1482 known explicit Runge-Kutta method of order (7,8).
1484 If this function is called with no return argument then plot the
1485 solution over time in a figure window while solving the set of
1486 ODEs that are defined in a function and specified by the function
1487 handle @FUN. The second input argument SLOT is a double vector
1488 that defines the time slot, INIT is a double vector that defines
1489 the initial values of the states, OPT can optionally be a
1490 structure array that keeps the options created with the command
1491 `odeset' and PAR1, PAR2, ... can optionally be other input
1492 arguments of any type that have to be passed to the function
1495 If this function is called with one return argument then return
1496 the solution SOL of type structure array after solving the set of
1497 ODEs. The solution SOL has the fields X of type double column
1498 vector for the steps chosen by the solver, Y of type double column
1499 vector for the solutions at each time step of X, SOLVER of type
1500 string for the solver name and optionally the extended time stamp
1501 information XE, the extended solution information YE and the
1502 extended index information IE all of type double column vector
1503 that keep the informations of the event function if an event
1504 function handle is set in the option argument OPT.
1506 If this function is called with more than one return argument then
1507 return the time stamps T, the solution values Y and optionally the
1508 extended time stamp information XE, the extended solution
1509 information YE and the extended index information IE all of type
1510 double column vector.
1512 For example, solve an anonymous implementation of the Van der Pol
1515 fvdb = @(vt,vy) [vy(2); (1 - vy(1)^2) * vy(2) - vy(1)];
1517 vopt = odeset ("RelTol", 1e-3, "AbsTol", 1e-3, \
1518 "NormControl", "on", "OutputFcn", @odeplot);
1519 ode78 (fvdb, [0 20], [2 0], vopt);
1521 -- Function File: [] = ode78d (@FUN, SLOT, INIT, LAGS, HIST, [OPT],
1523 -- Command: [SOL] = ode78d (@FUN, SLOT, INIT, LAGS, HIST, [OPT],
1525 -- Command: [T, Y, [XE, YE, IE]] = ode78d (@FUN, SLOT, INIT, LAGS,
1526 HIST, [OPT], [PAR1, PAR2, ...])
1527 This function file can be used to solve a set of non-stiff delay
1528 differential equations (non-stiff DDEs) with a modified version of
1529 the well known explicit Runge-Kutta method of order (7,8).
1531 If this function is called with no return argument then plot the
1532 solution over time in a figure window while solving the set of
1533 DDEs that are defined in a function and specified by the function
1534 handle @FUN. The second input argument SLOT is a double vector
1535 that defines the time slot, INIT is a double vector that defines
1536 the initial values of the states, LAGS is a double vector that
1537 describes the lags of time, HIST is a double matrix and describes
1538 the history of the DDEs, OPT can optionally be a structure array
1539 that keeps the options created with the command `odeset' and PAR1,
1540 PAR2, ... can optionally be other input arguments of any type that
1541 have to be passed to the function defined by @FUN.
1543 In other words, this function will solve a problem of the form
1544 dy/dt = fun (t, y(t), y(t-lags(1), y(t-lags(2), ...)))
1546 y(slot(1)-lags(1)) = hist(1), y(slot(1)-lags(2)) = hist(2), ...
1548 If this function is called with one return argument then return
1549 the solution SOL of type structure array after solving the set of
1550 DDEs. The solution SOL has the fields X of type double column
1551 vector for the steps chosen by the solver, Y of type double column
1552 vector for the solutions at each time step of X, SOLVER of type
1553 string for the solver name and optionally the extended time stamp
1554 information XE, the extended solution information YE and the
1555 extended index information IE all of type double column vector
1556 that keep the informations of the event function if an event
1557 function handle is set in the option argument OPT.
1559 If this function is called with more than one return argument then
1560 return the time stamps T, the solution values Y and optionally the
1561 extended time stamp information XE, the extended solution
1562 information YE and the extended index information IE all of type
1563 double column vector.
1566 - the following code solves an anonymous implementation of a
1569 fcao = @(vt, vy, vz) [2 * vz / (1 + vz^9.65) - vy];
1571 vopt = odeset ("NormControl", "on", "RelTol", 1e-3);
1572 vsol = ode78d (fcao, [0, 100], 0.5, 2, 0.5, vopt);
1574 vlag = interp1 (vsol.x, vsol.y, vsol.x - 2);
1575 plot (vsol.y, vlag); legend ("fcao (t,y,z)");
1577 - to solve the following problem with two delayed state
1581 d y2(t)/dt = -y2(t) + y1(t-5)
1582 d y3(t)/dt = -y3(t) + y2(t-10)*y1(t-10)
1584 one might do the following
1586 function f = fun (t, y, yd)
1587 f(1) = -y(1); %% y1' = -y1(t)
1588 f(2) = -y(2) + yd(1,1); %% y2' = -y2(t) + y1(t-lags(1))
1589 f(3) = -y(3) + yd(2,2)*yd(1,2); %% y3' = -y3(t) + y2(t-lags(2))*y1(t-lags(2))
1592 res = ode78d (@fun, T, [1;1;1], [5, 10], ones (3,2));
1595 -- Function File: [] = odebwe (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2,
1597 -- Command: [SOL] = odebwe (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
1598 -- Command: [T, Y, [XE, YE, IE]] = odebwe (@FUN, SLOT, INIT, [OPT],
1600 This function file can be used to solve a set of stiff ordinary
1601 differential equations (stiff ODEs) or stiff differential
1602 algebraic equations (stiff DAEs) with the Backward Euler method.
1604 If this function is called with no return argument then plot the
1605 solution over time in a figure window while solving the set of
1606 ODEs that are defined in a function and specified by the function
1607 handle @FUN. The second input argument SLOT is a double vector
1608 that defines the time slot, INIT is a double vector that defines
1609 the initial values of the states, OPT can optionally be a
1610 structure array that keeps the options created with the command
1611 `odeset' and PAR1, PAR2, ... can optionally be other input
1612 arguments of any type that have to be passed to the function
1615 If this function is called with one return argument then return
1616 the solution SOL of type structure array after solving the set of
1617 ODEs. The solution SOL has the fields X of type double column
1618 vector for the steps chosen by the solver, Y of type double column
1619 vector for the solutions at each time step of X, SOLVER of type
1620 string for the solver name and optionally the extended time stamp
1621 information XE, the extended solution information YE and the
1622 extended index information IE all of type double column vector
1623 that keep the informations of the event function if an event
1624 function handle is set in the option argument OPT.
1626 If this function is called with more than one return argument then
1627 return the time stamps T, the solution values Y and optionally the
1628 extended time stamp information XE, the extended solution
1629 information YE and the extended index information IE all of type
1630 double column vector.
1632 For example, solve an anonymous implementation of the Van der Pol
1635 fvdb = @(vt,vy) [vy(2); (1 - vy(1)^2) * vy(2) - vy(1)];
1636 vjac = @(vt,vy) [0, 1; -1 - 2 * vy(1) * vy(2), 1 - vy(1)^2];
1637 vopt = odeset ("RelTol", 1e-3, "AbsTol", 1e-3, \
1638 "NormControl", "on", "OutputFcn", @odeplot, \
1640 odebwe (fvdb, [0 20], [2 0], vopt);
1642 -- Function File: [] = odeexamples ()
1643 Open the differential equations examples menu and allow the user
1644 to select a submenu of ODE, DAE, IDE or DDE examples.
1646 -- Function File: [VALUE] = odeget (ODESTRUCT, OPTION, [DEFAULT])
1647 -- Command: [VALUES] = odeget (ODESTRUCT, {OPT1, OPT2, ...}, [{DEF1,
1649 If this function is called with two input arguments and the first
1650 input argument ODESTRUCT is of type structure array and the second
1651 input argument OPTION is of type string then return the option
1652 value VALUE that is specified by the option name OPTION in the
1653 OdePkg option structure ODESTRUCT. Optionally if this function is
1654 called with a third input argument then return the default value
1655 DEFAULT if OPTION is not set in the structure ODESTRUCT.
1657 If this function is called with two input arguments and the first
1658 input argument ODESTRUCT is of type structure array and the second
1659 input argument OPTION is of type cell array of strings then return
1660 the option values VALUES that are specified by the option names
1661 OPT1, OPT2, ... in the OdePkg option structure ODESTRUCT.
1662 Optionally if this function is called with a third input argument
1663 of type cell array then return the default value DEF1 if OPT1 is
1664 not set in the structure ODESTRUCT, DEF2 if OPT2 is not set in the
1665 structure ODESTRUCT, ...
1667 Run examples with the command
1670 -- Function File: [RET] = odephas2 (T, Y, FLAG)
1671 Open a new figure window and plot the first result from the
1672 variable Y that is of type double column vector over the second
1673 result from the variable Y while solving. The types and the values
1674 of the input parameter T and the output parameter RET depend on
1675 the input value FLAG that is of type string. If FLAG is
1677 then T must be a double column vector of length 2 with the
1678 first and the last time step and nothing is returned from
1682 then T must be a double scalar specifying the actual time
1683 step and the return value is false (resp. value 0) for 'not
1687 then T must be a double scalar specifying the last time step
1688 and nothing is returned from this function.
1690 This function is called by a OdePkg solver function if it was
1691 specified in an OdePkg options structure with the `odeset'. This
1692 function is an OdePkg internal helper function therefore it should
1693 never be necessary that this function is called directly by a
1694 user. There is only little error detection implemented in this
1695 function file to achieve the highest performance.
1697 For example, solve an anonymous implementation of the "Van der
1698 Pol" equation and display the results while solving in a 2D plane
1699 fvdb = @(vt,vy) [vy(2); (1 - vy(1)^2) * vy(2) - vy(1)];
1701 vopt = odeset ('OutputFcn', @odephas2, 'RelTol', 1e-6);
1702 vsol = ode45 (fvdb, [0 20], [2 0], vopt);
1704 -- Function File: [RET] = odephas3 (T, Y, FLAG)
1705 Open a new figure window and plot the first result from the
1706 variable Y that is of type double column vector over the second
1707 and the third result from the variable Y while solving. The types
1708 and the values of the input parameter T and the output parameter
1709 RET depend on the input value FLAG that is of type string. If FLAG
1712 then T must be a double column vector of length 2 with the
1713 first and the last time step and nothing is returned from
1717 then T must be a double scalar specifying the actual time
1718 step and the return value is false (resp. value 0) for 'not
1722 then T must be a double scalar specifying the last time step
1723 and nothing is returned from this function.
1725 This function is called by a OdePkg solver function if it was
1726 specified in an OdePkg options structure with the `odeset'. This
1727 function is an OdePkg internal helper function therefore it should
1728 never be necessary that this function is called directly by a
1729 user. There is only little error detection implemented in this
1730 function file to achieve the highest performance.
1732 For example, solve the "Lorenz attractor" and display the results
1733 while solving in a 3D plane
1734 function vyd = florenz (vt, vx)
1735 vyd = [10 * (vx(2) - vx(1));
1736 vx(1) * (28 - vx(3));
1737 vx(1) * vx(2) - 8/3 * vx(3)];
1740 vopt = odeset ('OutputFcn', @odephas3);
1741 vsol = ode23 (@florenz, [0:0.01:7.5], [3 15 1], vopt);
1743 -- Function File: [] = odepkg ()
1744 OdePkg is part of the GNU Octave Repository (the Octave-Forge
1745 project). The package includes commands for setting up various
1746 options, output functions etc. before solving a set of
1747 differential equations with the solver functions that are also
1748 included. At this time OdePkg is under development with the main
1749 target to make a package that is mostly compatible to proprietary
1752 If this function is called without any input argument then open
1753 the OdePkg tutorial in the Octave window. The tutorial can also be
1754 opened with the following command
1758 -- Function File: [SOL] = odepkg_event_handle (@FUN, TIME, Y, FLAG,
1760 Return the solution of the event function that is specified as the
1761 first input argument @FUN in form of a function handle. The second
1762 input argument TIME is of type double scalar and specifies the
1763 time of the event evaluation, the third input argument Y either is
1764 of type double column vector (for ODEs and DAEs) and specifies the
1765 solutions or is of type cell array (for IDEs and DDEs) and
1766 specifies the derivatives or the history values, the third input
1767 argument FLAG is of type string and can be of the form
1769 then initialize internal persistent variables of the function
1770 `odepkg_event_handle' and return an empty cell array of size
1774 then do the evaluation of the event function and return the
1775 solution SOL as type cell array of size 4,
1778 then cleanup internal variables of the function
1779 `odepkg_event_handle' and return an empty cell array of size
1781 Optionally if further input arguments PAR1, PAR2, ... of any type
1782 are given then pass these parameters through `odepkg_event_handle'
1783 to the event function.
1785 This function is an OdePkg internal helper function therefore it
1786 should never be necessary that this function is called directly by
1787 a user. There is only little error detection implemented in this
1788 function file to achieve the highest performance.
1790 -- Function File: [] = odepkg_examples_dae ()
1791 Open the DAE examples menu and allow the user to select a demo
1792 that will be evaluated.
1794 -- Function File: [] = odepkg_examples_dde ()
1795 Open the DDE examples menu and allow the user to select a demo
1796 that will be evaluated.
1798 -- Function File: [] = odepkg_examples_ide ()
1799 Open the IDE examples menu and allow the user to select a demo
1800 that will be evaluated.
1802 -- Function File: [] = odepkg_examples_ode ()
1803 Open the ODE examples menu and allow the user to select a demo
1804 that will be evaluated.
1806 -- Function File: [NEWSTRUCT] = odepkg_structure_check (OLDSTRUCT,
1808 If this function is called with one input argument of type
1809 structure array then check the field names and the field values of
1810 the OdePkg structure OLDSTRUCT and return the structure as
1811 NEWSTRUCT if no error is found. Optionally if this function is
1812 called with a second input argument "SOLVER" of type string taht
1813 specifies the name of a valid OdePkg solver then a higher level
1814 error detection is performed. The function does not modify any of
1815 the field names or field values but terminates with an error if an
1816 invalid option or value is found.
1818 This function is an OdePkg internal helper function therefore it
1819 should never be necessary that this function is called directly by
1820 a user. There is only little error detection implemented in this
1821 function file to achieve the highest performance.
1823 Run examples with the command
1824 demo odepkg_structure_check
1826 -- Function File: [MESCD] = odepkg_testsuite_calcmescd (SOLUTION,
1827 REFERENCE, ABSTOL, RELTOL)
1828 If this function is called with four input arguments of type
1829 double scalar or column vector then return a normalized value for
1830 the minimum number of correct digits MESCD that is calculated from
1831 the solution at the end of an integration interval SOLUTION and a
1832 set of reference values REFERENCE. The input arguments ABSTOL and
1833 RELTOL are used to calculate a reference solution that depends on
1834 the relative and absolute error tolerances.
1836 Run examples with the command
1837 demo odepkg_testsuite_calcmescd
1839 This function has been ported from the "Test Set for IVP solvers"
1840 which is developed by the INdAM Bari unit project group "Codes and
1841 Test Problems for Differential Equations", coordinator F. Mazzia.
1843 -- Function File: [SCD] = odepkg_testsuite_calcscd (SOLUTION,
1844 REFERENCE, ABSTOL, RELTOL)
1845 If this function is called with four input arguments of type
1846 double scalar or column vector then return a normalized value for
1847 the minimum number of correct digits SCD that is calculated from
1848 the solution at the end of an integration interval SOLUTION and a
1849 set of reference values REFERENCE. The input arguments ABSTOL and
1850 RELTOL are unused but present because of compatibility to the
1851 function `odepkg_testsuite_calcmescd'.
1853 Run examples with the command
1854 demo odepkg_testsuite_calcscd
1856 This function has been ported from the "Test Set for IVP solvers"
1857 which is developed by the INdAM Bari unit project group "Codes and
1858 Test Problems for Differential Equations", coordinator F. Mazzia.
1860 -- Function File: [SOLUTION] = odepkg_testsuite_chemakzo (@SOLVER,
1862 If this function is called with two input arguments and the first
1863 input argument @SOLVER is a function handle describing an OdePkg
1864 solver and the second input argument RELTOL is a double scalar
1865 describing the relative error tolerance then return a cell array
1866 SOLUTION with performance informations about the chemical AKZO
1867 Nobel testsuite of differential algebraic equations after solving
1870 Run examples with the command
1871 demo odepkg_testsuite_chemakzo
1873 This function has been ported from the "Test Set for IVP solvers"
1874 which is developed by the INdAM Bari unit project group "Codes and
1875 Test Problems for Differential Equations", coordinator F. Mazzia.
1877 -- Function File: [SOLUTION] = odepkg_testsuite_hires (@SOLVER, RELTOL)
1878 If this function is called with two input arguments and the first
1879 input argument @SOLVER is a function handle describing an OdePkg
1880 solver and the second input argument RELTOL is a double scalar
1881 describing the relative error tolerance then return a cell array
1882 SOLUTION with performance informations about the HIRES testsuite
1883 of ordinary differential equations after solving (ODE-test).
1885 Run examples with the command
1886 demo odepkg_testsuite_hires
1888 This function has been ported from the "Test Set for IVP solvers"
1889 which is developed by the INdAM Bari unit project group "Codes and
1890 Test Problems for Differential Equations", coordinator F. Mazzia.
1892 -- Function File: [SOLUTION] = odepkg_testsuite_implakzo (@SOLVER,
1894 If this function is called with two input arguments and the first
1895 input argument @SOLVER is a function handle describing an OdePkg
1896 solver and the second input argument RELTOL is a double scalar
1897 describing the relative error tolerance then return a cell array
1898 SOLUTION with performance informations about the chemical AKZO
1899 Nobel testsuite of implicit differential algebraic equations after
1902 Run examples with the command
1903 demo odepkg_testsuite_implakzo
1905 This function has been ported from the "Test Set for IVP solvers"
1906 which is developed by the INdAM Bari unit project group "Codes and
1907 Test Problems for Differential Equations", coordinator F. Mazzia.
1909 -- Function File: [SOLUTION] = odepkg_testsuite_implrober (@SOLVER,
1911 If this function is called with two input arguments and the first
1912 input argument @SOLVER is a function handle describing an OdePkg
1913 solver and the second input argument RELTOL is a double scalar
1914 describing the relative error tolerance then return a cell array
1915 SOLUTION with performance informations about the implicit form of
1916 the modified ROBERTSON testsuite of implicit differential
1917 algebraic equations after solving (IDE-test).
1919 Run examples with the command
1920 demo odepkg_testsuite_implrober
1922 This function has been ported from the "Test Set for IVP solvers"
1923 which is developed by the INdAM Bari unit project group "Codes and
1924 Test Problems for Differential Equations", coordinator F. Mazzia.
1926 -- Function File: [SOLUTION] = odepkg_testsuite_oregonator (@SOLVER,
1928 If this function is called with two input arguments and the first
1929 input argument @SOLVER is a function handle describing an OdePkg
1930 solver and the second input argument RELTOL is a double scalar
1931 describing the relative error tolerance then return a cell array
1932 SOLUTION with performance informations about the OREGONATOR
1933 testsuite of ordinary differential equations after solving
1936 Run examples with the command
1937 demo odepkg_testsuite_oregonator
1939 This function has been ported from the "Test Set for IVP solvers"
1940 which is developed by the INdAM Bari unit project group "Codes and
1941 Test Problems for Differential Equations", coordinator F. Mazzia.
1943 -- Function File: [SOLUTION] = odepkg_testsuite_pollution (@SOLVER,
1945 If this function is called with two input arguments and the first
1946 input argument @SOLVER is a function handle describing an OdePkg
1947 solver and the second input argument RELTOL is a double scalar
1948 describing the relative error tolerance then return the cell array
1949 SOLUTION with performance informations about the POLLUTION
1950 testsuite of ordinary differential equations after solving
1953 Run examples with the command
1954 demo odepkg_testsuite_pollution
1956 This function has been ported from the "Test Set for IVP solvers"
1957 which is developed by the INdAM Bari unit project group "Codes and
1958 Test Problems for Differential Equations", coordinator F. Mazzia.
1960 -- Function File: [SOLUTION] = odepkg_testsuite_robertson (@SOLVER,
1962 If this function is called with two input arguments and the first
1963 input argument @SOLVER is a function handle describing an OdePkg
1964 solver and the second input argument RELTOL is a double scalar
1965 describing the relative error tolerance then return a cell array
1966 SOLUTION with performance informations about the modified
1967 ROBERTSON testsuite of differential algebraic equations after
1970 Run examples with the command
1971 demo odepkg_testsuite_robertson
1973 This function has been ported from the "Test Set for IVP solvers"
1974 which is developed by the INdAM Bari unit project group "Codes and
1975 Test Problems for Differential Equations", coordinator F. Mazzia.
1977 -- Function File: [SOLUTION] = odepkg_testsuite_transistor (@SOLVER,
1979 If this function is called with two input arguments and the first
1980 input argument @SOLVER is a function handle describing an OdePkg
1981 solver and the second input argument RELTOL is a double scalar
1982 describing the relative error tolerance then return the cell array
1983 SOLUTION with performance informations about the TRANSISTOR
1984 testsuite of differential algebraic equations after solving
1987 Run examples with the command
1988 demo odepkg_testsuite_transistor
1990 This function has been ported from the "Test Set for IVP solvers"
1991 which is developed by the INdAM Bari unit project group "Codes and
1992 Test Problems for Differential Equations", coordinator F. Mazzia.
1994 -- Function File: [RET] = odeplot (T, Y, FLAG)
1995 Open a new figure window and plot the results from the variable Y
1996 of type column vector over time while solving. The types and the
1997 values of the input parameter T and the output parameter RET
1998 depend on the input value FLAG that is of type string. If FLAG is
2000 then T must be a double column vector of length 2 with the
2001 first and the last time step and nothing is returned from
2005 then T must be a double scalar specifying the actual time
2006 step and the return value is false (resp. value 0) for 'not
2010 then T must be a double scalar specifying the last time step
2011 and nothing is returned from this function.
2013 This function is called by a OdePkg solver function if it was
2014 specified in an OdePkg options structure with the `odeset'. This
2015 function is an OdePkg internal helper function therefore it should
2016 never be necessary that this function is called directly by a
2017 user. There is only little error detection implemented in this
2018 function file to achieve the highest performance.
2020 For example, solve an anonymous implementation of the "Van der
2021 Pol" equation and display the results while solving
2022 fvdb = @(vt,vy) [vy(2); (1 - vy(1)^2) * vy(2) - vy(1)];
2024 vopt = odeset ('OutputFcn', @odeplot, 'RelTol', 1e-6);
2025 vsol = ode45 (fvdb, [0 20], [2 0], vopt);
2027 -- Function File: [RET] = odeprint (T, Y, FLAG)
2028 Display the results of the set of differential equations in the
2029 Octave window while solving. The first column of the screen output
2030 shows the actual time stamp that is given with the input arguemtn
2031 T, the following columns show the results from the function
2032 evaluation that are given by the column vector Y. The types and
2033 the values of the input parameter T and the output parameter RET
2034 depend on the input value FLAG that is of type string. If FLAG is
2036 then T must be a double column vector of length 2 with the
2037 first and the last time step and nothing is returned from
2041 then T must be a double scalar specifying the actual time
2042 step and the return value is false (resp. value 0) for 'not
2046 then T must be a double scalar specifying the last time step
2047 and nothing is returned from this function.
2049 This function is called by a OdePkg solver function if it was
2050 specified in an OdePkg options structure with the `odeset'. This
2051 function is an OdePkg internal helper function therefore it should
2052 never be necessary that this function is called directly by a
2053 user. There is only little error detection implemented in this
2054 function file to achieve the highest performance.
2056 For example, solve an anonymous implementation of the "Van der
2057 Pol" equation and print the results while solving
2058 fvdb = @(vt,vy) [vy(2); (1 - vy(1)^2) * vy(2) - vy(1)];
2060 vopt = odeset ('OutputFcn', @odeprint, 'RelTol', 1e-6);
2061 vsol = ode45 (fvdb, [0 20], [2 0], vopt);
2063 -- Function File: [ODESTRUCT] = odeset ()
2064 -- Command: [ODESTRUCT] = odeset ("FIELD1", VALUE1, "FIELD2", VALUE2,
2066 -- Command: [ODESTRUCT] = odeset (OLDSTRUCT, "FIELD1", VALUE1,
2067 "FIELD2", VALUE2, ...)
2068 -- Command: [ODESTRUCT] = odeset (OLDSTRUCT, NEWSTRUCT)
2069 If this function is called without an input argument then return a
2070 new OdePkg options structure array that contains all the necessary
2071 fields and sets the values of all fields to default values.
2073 If this function is called with string input arguments "FIELD1",
2074 "FIELD2", ... identifying valid OdePkg options then return a new
2075 OdePkg options structure with all necessary fields and set the
2076 values of the fields "FIELD1", "FIELD2", ... to the values VALUE1,
2079 If this function is called with a first input argument OLDSTRUCT
2080 of type structure array then overwrite all values of the options
2081 "FIELD1", "FIELD2", ... of the structure OLDSTRUCT with new values
2082 VALUE1, VALUE2, ... and return the modified structure array.
2084 If this function is called with two input argumnets OLDSTRUCT and
2085 NEWSTRUCT of type structure array then overwrite all values in the
2086 fields from the structure OLDSTRUCT with new values of the fields
2087 from the structure NEWSTRUCT. Empty values of NEWSTRUCT will not
2088 overwrite values in OLDSTRUCT.
2090 For a detailed explanation about valid fields and field values in
2091 an OdePkg structure aaray have a look at the `odepkg.pdf', Section
2092 'ODE/DAE/IDE/DDE options' or run the command `doc odepkg' to open
2095 Run examples with the command
2099 File: odepkg.info, Node: Oct-File Function Reference, Prev: M-File Function Reference, Up: Users Guide
2101 2.5 Oct-File Function Reference
2102 ===============================
2104 The help texts of this section are autogenerated and refer to commands
2105 that all can be found in the file `dldsolver.oct'. The file
2106 `dldsolver.oct' is generated automatically if you install OdePkg with
2107 the command `pkg'. All commands that are listed below are loaded
2108 automatically everytime you launch Octave.
2109 -- Command: [] = odebda (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
2110 -- Command: [SOL] = odebda (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
2111 -- Command: [T, Y, [XE, YE, IE]] = odebda (@FUN, SLOT, INIT, [OPT],
2113 This function file can be used to solve a set of non-stiff or
2114 stiff ordinary differential equations (ODEs) and non-stiff or
2115 stiff differential algebraic equations (DAEs). This function file
2116 is a wrapper file that uses Jeff Cash's Fortran solver
2119 If this function is called with no return argument then plot the
2120 solution over time in a figure window while solving the set of
2121 ODEs that are defined in a function and specified by the function
2122 handle @FUN. The second input argument SLOT is a double vector
2123 that defines the time slot, INIT is a double vector that defines
2124 the initial values of the states, OPT can optionally be a
2125 structure array that keeps the options created with the command
2126 `odeset' and PAR1, PAR2, ... can optionally be other input
2127 arguments of any type that have to be passed to the function
2130 If this function is called with one return argument then return
2131 the solution SOL of type structure array after solving the set of
2132 ODEs. The solution SOL has the fields X of type double column
2133 vector for the steps chosen by the solver, Y of type double column
2134 vector for the solutions at each time step of X, SOLVER of type
2135 string for the solver name and optionally the extended time stamp
2136 information XE, the extended solution information YE and the
2137 extended index information IE all of type double column vector
2138 that keep the informations of the event function if an event
2139 function handle is set in the option argument OPT.
2141 If this function is called with more than one return argument then
2142 return the time stamps T, the solution values Y and optionally the
2143 extended time stamp information XE, the extended solution
2144 information YE and the extended index information IE all of type
2145 double column vector.
2148 function y = odepkg_equations_lorenz (t, x)
2149 y = [10 * (x(2) - x(1));
2151 x(1) * x(2) - 8/3 * x(3)];
2154 vopt = odeset ("InitialStep", 1e-3, "MaxStep", 1e-1, \\
2155 "OutputFcn", @odephas3, "Refine", 5);
2156 odebda (@odepkg_equations_lorenz, [0, 25], [3 15 1], vopt);
2158 -- Command: [] = odebdi (@FUN, SLOT, Y0, DY0, [OPT], [P1, P2, ...])
2159 -- Command: [SOL] = odebdi (@FUN, SLOT, Y0, DY0, [OPT], [P1, P2, ...])
2160 -- Command: [T, Y, [XE, YE, IE]] = odebdi (@FUN, SLOT, Y0, DY0, [OPT],
2162 This function file can be used to solve a set of non-stiff and
2163 stiff implicit differential equations (IDEs). This function file
2164 is a wrapper file that uses Jeff Cash's Fortran solver `mebdfi.f'.
2166 If this function is called with no return argument then plot the
2167 solution over time in a figure window while solving the set of
2168 IDEs that are defined in a function and specified by the function
2169 handle @FUN. The second input argument SLOT is a double vector
2170 that defines the time slot, Y0 is a double vector that defines the
2171 initial values of the states, DY0 is a double vector that defines
2172 the initial values of the derivatives, OPT can optionally be a
2173 structure array that keeps the options created with the command
2174 `odeset' and PAR1, PAR2, ... can optionally be other input
2175 arguments of any type that have to be passed to the function
2178 If this function is called with one return argument then return
2179 the solution SOL of type structure array after solving the set of
2180 IDEs. The solution SOL has the fields X of type double column
2181 vector for the steps chosen by the solver, Y of type double column
2182 vector for the solutions at each time step of X, SOLVER of type
2183 string for the solver name and optionally the extended time stamp
2184 information XE, the extended solution information YE and the
2185 extended index information IE all of type double column vector
2186 that keep the informations of the event function if an event
2187 function handle is set in the option argument OPT.
2189 If this function is called with more than one return argument then
2190 return the time stamps T, the solution values Y and optionally the
2191 extended time stamp information XE, the extended solution
2192 information YE and the extended index information IE all of type
2193 double column vector.
2196 function res = odepkg_equations_ilorenz (t, y, yd)
2197 res = [10 * (y(2) - y(1)) - yd(1);
2198 y(1) * (28 - y(3)) - yd(2);
2199 y(1) * y(2) - 8/3 * y(3) - yd(3)];
2202 vopt = odeset ("InitialStep", 1e-3, "MaxStep", 1e-1, \\
2203 "OutputFcn", @odephas3, "Refine", 5);
2204 odebdi (@odepkg_equations_ilorenz, [0, 25], [3 15 1], \\
2205 [120 81 42.333333], vopt);
2207 -- Command: [] = odekdi (@FUN, SLOT, Y0, DY0, [OPT], [P1, P2, ...])
2208 -- Command: [SOL] = odekdi (@FUN, SLOT, Y0, DY0, [OPT], [P1, P2, ...])
2209 -- Command: [T, Y, [XE, YE, IE]] = odekdi (@FUN, SLOT, Y0, DY0, [OPT],
2211 This function file can be used to solve a set of non-stiff or
2212 stiff implicit differential equations (IDEs). This function file
2213 is a wrapper file that uses the direct method (not the Krylov
2214 method) of Petzold's, Brown's, Hindmarsh's and Ulrich's Fortran
2217 If this function is called with no return argument then plot the
2218 solution over time in a figure window while solving the set of
2219 IDEs that are defined in a function and specified by the function
2220 handle @FUN. The second input argument SLOT is a double vector
2221 that defines the time slot, Y0 is a double vector that defines the
2222 initial values of the states, DY0 is a double vector that defines
2223 the initial values of the derivatives, OPT can optionally be a
2224 structure array that keeps the options created with the command
2225 `odeset' and PAR1, PAR2, ... can optionally be other input
2226 arguments of any type that have to be passed to the function
2229 If this function is called with one return argument then return
2230 the solution SOL of type structure array after solving the set of
2231 IDEs. The solution SOL has the fields X of type double column
2232 vector for the steps chosen by the solver, Y of type double column
2233 vector for the solutions at each time step of X, SOLVER of type
2234 string for the solver name and optionally the extended time stamp
2235 information XE, the extended solution information YE and the
2236 extended index information IE all of type double column vector
2237 that keep the informations of the event function if an event
2238 function handle is set in the option argument OPT.
2240 If this function is called with more than one return argument then
2241 return the time stamps T, the solution values Y and optionally the
2242 extended time stamp information XE, the extended solution
2243 information YE and the extended index information IE all of type
2244 double column vector.
2247 function res = odepkg_equations_ilorenz (t, y, yd)
2248 res = [10 * (y(2) - y(1)) - yd(1);
2249 y(1) * (28 - y(3)) - yd(2);
2250 y(1) * y(2) - 8/3 * y(3) - yd(3)];
2253 vopt = odeset ("InitialStep", 1e-3, "MaxStep", 1e-1, \\
2254 "OutputFcn", @odephas3, "Refine", 5);
2255 odekdi (@odepkg_equations_ilorenz, [0, 25], [3 15 1], \\
2256 [120 81 42.333333], vopt);
2258 -- Command: [] = ode2r (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
2259 -- Command: [SOL] = ode2r (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
2260 -- Command: [T, Y, [XE, YE, IE]] = ode2r (@FUN, SLOT, INIT, [OPT],
2262 This function file can be used to solve a set of non-stiff or
2263 stiff ordinary differential equations (ODEs) and non-stiff or
2264 stiff differential algebraic equations (DAEs). This function file
2265 is a wrapper to Hairer's and Wanner's Fortran solver `radau.f'.
2267 If this function is called with no return argument then plot the
2268 solution over time in a figure window while solving the set of
2269 ODEs that are defined in a function and specified by the function
2270 handle @FUN. The second input argument SLOT is a double vector
2271 that defines the time slot, INIT is a double vector that defines
2272 the initial values of the states, OPT can optionally be a
2273 structure array that keeps the options created with the command
2274 `odeset' and PAR1, PAR2, ... can optionally be other input
2275 arguments of any type that have to be passed to the function
2278 If this function is called with one return argument then return
2279 the solution SOL of type structure array after solving the set of
2280 ODEs. The solution SOL has the fields X of type double column
2281 vector for the steps chosen by the solver, Y of type double column
2282 vector for the solutions at each time step of X, SOLVER of type
2283 string for the solver name and optionally the extended time stamp
2284 information XE, the extended solution information YE and the
2285 extended index information IE all of type double column vector
2286 that keep the informations of the event function if an event
2287 function handle is set in the option argument OPT.
2289 If this function is called with more than one return argument then
2290 return the time stamps T, the solution values Y and optionally the
2291 extended time stamp information XE, the extended solution
2292 information YE and the extended index information IE all of type
2293 double column vector.
2296 function y = odepkg_equations_lorenz (t, x)
2297 y = [10 * (x(2) - x(1));
2299 x(1) * x(2) - 8/3 * x(3)];
2302 vopt = odeset ("InitialStep", 1e-3, "MaxStep", 1e-1, \\
2303 "OutputFcn", @odephas3, "Refine", 5);
2304 ode2r (@odepkg_equations_lorenz, [0, 25], [3 15 1], vopt);
2306 -- Command: [] = ode5r (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
2307 -- Command: [SOL] = ode5r (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
2308 -- Command: [T, Y, [XE, YE, IE]] = ode5r (@FUN, SLOT, INIT, [OPT],
2310 This function file can be used to solve a set of non-stiff or
2311 stiff ordinary differential equations (ODEs) and non-stiff or
2312 stiff differential algebraic equations (DAEs). This function file
2313 is a wrapper to Hairer's and Wanner's Fortran solver `radau5.f'.
2315 If this function is called with no return argument then plot the
2316 solution over time in a figure window while solving the set of
2317 ODEs that are defined in a function and specified by the function
2318 handle @FUN. The second input argument SLOT is a double vector
2319 that defines the time slot, INIT is a double vector that defines
2320 the initial values of the states, OPT can optionally be a
2321 structure array that keeps the options created with the command
2322 `odeset' and PAR1, PAR2, ... can optionally be other input
2323 arguments of any type that have to be passed to the function
2326 If this function is called with one return argument then return
2327 the solution SOL of type structure array after solving the set of
2328 ODEs. The solution SOL has the fields X of type double column
2329 vector for the steps chosen by the solver, Y of type double column
2330 vector for the solutions at each time step of X, SOLVER of type
2331 string for the solver name and optionally the extended time stamp
2332 information XE, the extended solution information YE and the
2333 extended index information IE all of type double column vector
2334 that keep the informations of the event function if an event
2335 function handle is set in the option argument OPT.
2337 If this function is called with more than one return argument then
2338 return the time stamps T, the solution values Y and optionally the
2339 extended time stamp information XE, the extended solution
2340 information YE and the extended index information IE all of type
2341 double column vector.
2344 function y = odepkg_equations_lorenz (t, x)
2345 y = [10 * (x(2) - x(1));
2347 x(1) * x(2) - 8/3 * x(3)];
2350 vopt = odeset ("InitialStep", 1e-3, "MaxStep", 1e-1, \\
2351 "OutputFcn", @odephas3, "Refine", 5);
2352 ode5r (@odepkg_equations_lorenz, [0, 25], [3 15 1], vopt);
2354 -- Function File: [] = oders (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2,
2356 -- Command: [SOL] = oders (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
2357 -- Command: [T, Y, [XE, YE, IE]] = oders (@FUN, SLOT, INIT, [OPT],
2359 This function file can be used to solve a set of non-stiff or
2360 stiff ordinary differential equations (ODEs) and non-stiff or
2361 stiff differential algebraic equations (DAEs). This function file
2362 is a wrapper to Hairer's and Wanner's Fortran solver `rodas.f'.
2364 If this function is called with no return argument then plot the
2365 solution over time in a figure window while solving the set of
2366 ODEs that are defined in a function and specified by the function
2367 handle @FUN. The second input argument SLOT is a double vector
2368 that defines the time slot, INIT is a double vector that defines
2369 the initial values of the states, OPT can optionally be a
2370 structure array that keeps the options created with the command
2371 `odeset' and PAR1, PAR2, ... can optionally be other input
2372 arguments of any type that have to be passed to the function
2375 If this function is called with one return argument then return
2376 the solution SOL of type structure array after solving the set of
2377 ODEs. The solution SOL has the fields X of type double column
2378 vector for the steps chosen by the solver, Y of type double column
2379 vector for the solutions at each time step of X, SOLVER of type
2380 string for the solver name and optionally the extended time stamp
2381 information XE, the extended solution information YE and the
2382 extended index information IE all of type double column vector
2383 that keep the informations of the event function if an event
2384 function handle is set in the option argument OPT.
2386 If this function is called with more than one return argument then
2387 return the time stamps T, the solution values Y and optionally the
2388 extended time stamp information XE, the extended solution
2389 information YE and the extended index information IE all of type
2390 double column vector.
2393 function y = odepkg_equations_lorenz (t, x)
2394 y = [10 * (x(2) - x(1));
2396 x(1) * x(2) - 8/3 * x(3)];
2399 vopt = odeset ("InitialStep", 1e-3, "MaxStep", 1e-1, \\
2400 "OutputFcn", @odephas3, "Refine", 5);
2401 oders (@odepkg_equations_lorenz, [0, 25], [3 15 1], vopt);
2403 -- Command: [] = odesx (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
2404 -- Command: [SOL] = odesx (@FUN, SLOT, INIT, [OPT], [PAR1, PAR2, ...])
2405 -- Command: [T, Y, [XE, YE, IE]] = odesx (@FUN, SLOT, INIT, [OPT],
2407 This function file can be used to solve a set of stiff or
2408 non-stiff ordinary differential equations (ODEs) and non-stiff or
2409 stiff differential algebraic equations (DAEs). This function file
2410 is a wrapper to Hairer's and Wanner's Fortran solver `seulex.f'.
2412 If this function is called with no return argument then plot the
2413 solution over time in a figure window while solving the set of
2414 ODEs that are defined in a function and specified by the function
2415 handle @FUN. The second input argument SLOT is a double vector
2416 that defines the time slot, INIT is a double vector that defines
2417 the initial values of the states, OPT can optionally be a
2418 structure array that keeps the options created with the command
2419 `odeset' and PAR1, PAR2, ... can optionally be other input
2420 arguments of any type that have to be passed to the function
2423 If this function is called with one return argument then return
2424 the solution SOL of type structure array after solving the set of
2425 ODEs. The solution SOL has the fields X of type double column
2426 vector for the steps chosen by the solver, Y of type double column
2427 vector for the solutions at each time step of X, SOLVER of type
2428 string for the solver name and optionally the extended time stamp
2429 information XE, the extended solution information YE and the
2430 extended index information IE all of type double column vector
2431 that keep the informations of the event function if an event
2432 function handle is set in the option argument OPT.
2434 If this function is called with more than one return argument then
2435 return the time stamps T, the solution values Y and optionally the
2436 extended time stamp information XE, the extended solution
2437 information YE and the extended index information IE all of type
2438 double column vector.
2441 function y = odepkg_equations_lorenz (t, x)
2442 y = [10 * (x(2) - x(1));
2444 x(1) * x(2) - 8/3 * x(3)];
2447 vopt = odeset ("InitialStep", 1e-3, "MaxStep", 1e-1, \\
2448 "OutputFcn", @odephas3, "Refine", 5);
2449 odesx (@odepkg_equations_lorenz, [0, 25], [3 15 1], vopt);
2452 File: odepkg.info, Node: Programmers Guide, Next: Function Index, Prev: Users Guide, Up: Top
2459 * Missing features:: The TODO-list for missing features
2462 File: odepkg.info, Node: Missing features, Prev: Programmers Guide, Up: Programmers Guide
2464 3.1 Missing features
2465 ====================
2467 If somebody want to help improving OdePkg then please contact the
2468 Octave-Forge developer team sending your modifications via the
2469 mailing-list <octave-dev@lists.sourceforge.net>. Here is a TODO-list
2470 about missing features:
2471 * Partial Derivative equations (PDEs) and Boundary Value Equations
2472 (BVPs) cannot be solved with the solvers of the OdePkg. The wish
2473 for solving PDEs and BVPs definitely is there (maybe you'd like to
2474 create another package and call that package PdePkg, which is just
2477 * Some options that can be set with propietary solver products are
2478 not available within OdePkg. Have a look at section *Note
2479 ODE/DAE/IDE/DDE options:: about which options can be set and which
2480 options are not supported and help improving the command `odeset',
2481 `odepkg_structure_check' and the solvers that have to deal with
2484 * OdePkg currently is missing the command `decic' which computes an
2485 initial constraint for IDEs before solving. The command `deval'
2486 also is missing that interpolates the results that can be obtained
2487 from the solvers after solving and then plots the results in a new
2488 figure. However, instead of using `deval' any of the commands
2491 * If you want to include your own solver within OdePkg then either
2492 code in `*.m' or `*.cc'. Solvers in `*.m' are preferred. Choose a
2493 GPL compatible license for your solver and send your solver file
2494 to the mailing list.
2496 * Before interfacing other solvers make sure that the core solver
2497 file is available under a GPL-compatible license (if you'd like to
2498 redistribute your wrapper with OdePkg). There can be found a lot
2499 of solver files at `http://www.netlib.org' but most of them are
2500 not GPL-compatible. Here is a list about authors and their solvers
2501 that do have a GPL compatible license so that their codes can be
2502 redistributed with OdePkg:
2504 * Cecilia Magherini and Luigi Brugnano from the University of
2505 Bari created the BIMD solver that is available at
2506 `http://pitagora.dm.uniba.it/~testset/solvers/bimd.php'. This
2507 solver can be used to solve stiff DAE equations. The
2508 Fortran77 file has been released under the GNU GPL V2.
2510 * Francesca Mazzia and Felice Iavernaro from the University of
2511 Bari created the GAMD solver that is available at
2512 `http://pitagora.dm.uniba.it/~testset/solvers/gamd.php'. This
2513 solver can be used to solve stiff DAE equations. The
2514 Fortran90 file has been released under the GNU GPL V2 but for
2515 OdePkg a Fortran77 implementation would be preferred.
2517 * Ernst Hairer and Gerhard Wanner have been written more
2518 solvers that are released under a modified BSD license than
2519 have been interfaced by OdePkg. Notable solvers that can be
2520 found at `http://www.unige.ch/~hairer/software.html' are
2521 explicit Runge-Kutta methods `dopri5' and `dop853' and
2522 extrapolation methods `odex' and `odex2' for solving ODEs,
2523 `retard' and `radar5' for solving DDEs.
2525 * Jeff Cash has released some more Fortran77 solvers for
2526 different kinds of differential equation problems than are
2527 interfaced by OdePkg, check his website at
2528 `http://www.ma.ic.ac.uk/~jcash'.
2532 File: odepkg.info, Node: Function Index, Next: Index, Prev: Programmers Guide, Up: Top
2540 * ode23: M-File Function Reference.
2542 * ode23d: M-File Function Reference.
2544 * ode2r: Oct-File Function Reference.
2546 * ode45: M-File Function Reference.
2548 * ode45d: M-File Function Reference.
2550 * ode54: M-File Function Reference.
2552 * ode54d: M-File Function Reference.
2554 * ode5r: Oct-File Function Reference.
2556 * ode78: M-File Function Reference.
2558 * ode78d: M-File Function Reference.
2560 * odebda: Oct-File Function Reference.
2562 * odebdi: Oct-File Function Reference.
2564 * odebwe: M-File Function Reference.
2566 * odeexamples: M-File Function Reference.
2568 * odeget: M-File Function Reference.
2570 * odekdi: Oct-File Function Reference.
2572 * odephas2: M-File Function Reference.
2574 * odephas3: M-File Function Reference.
2576 * odepkg: M-File Function Reference.
2578 * odepkg_event_handle: M-File Function Reference.
2580 * odepkg_examples_dae: M-File Function Reference.
2582 * odepkg_examples_dde: M-File Function Reference.
2584 * odepkg_examples_ide: M-File Function Reference.
2586 * odepkg_examples_ode: M-File Function Reference.
2588 * odepkg_structure_check: M-File Function Reference.
2590 * odepkg_testsuite_calcmescd: M-File Function Reference.
2592 * odepkg_testsuite_calcscd: M-File Function Reference.
2594 * odepkg_testsuite_chemakzo: M-File Function Reference.
2596 * odepkg_testsuite_hires: M-File Function Reference.
2598 * odepkg_testsuite_implakzo: M-File Function Reference.
2600 * odepkg_testsuite_implrober: M-File Function Reference.
2602 * odepkg_testsuite_oregonator: M-File Function Reference.
2604 * odepkg_testsuite_pollution: M-File Function Reference.
2606 * odepkg_testsuite_robertson: M-File Function Reference.
2608 * odepkg_testsuite_transistor: M-File Function Reference.
2610 * odeplot: M-File Function Reference.
2612 * odeprint: M-File Function Reference.
2614 * oders: Oct-File Function Reference.
2616 * odeset: M-File Function Reference.
2618 * odesx: Oct-File Function Reference.
2622 File: odepkg.info, Node: Index, Prev: Function Index, Up: Top
2630 * About OdePkg: About OdePkg. (line 6)
2631 * AbsTol option: ODE/DAE/IDE/DDE options.
2633 * BDF option: ODE/DAE/IDE/DDE options.
2635 * BDF solver: Cash modified BDF solvers.
2637 * Beginners guide: Beginners Guide. (line 6)
2638 * bugs: Reporting Bugs. (line 6)
2639 * Cash modified BDF: Cash modified BDF solvers.
2641 * dae equations: DAE equations. (line 6)
2642 * dae options: ODE/DAE/IDE/DDE options.
2644 * ddaskr solver: DDaskr direct method solver.
2646 * dde equations: DDE equations. (line 6)
2647 * dde options: ODE/DAE/IDE/DDE options.
2649 * decic: Missing features. (line 28)
2650 * deinstallation: Installation and deinstallation.
2652 * deval: Missing features. (line 28)
2653 * differential equations: Differential Equations.
2655 * Events option: ODE/DAE/IDE/DDE options.
2657 * foo example: The "foo" example. (line 6)
2658 * Hairer-Wanner: Hairer-Wanner solvers.
2660 * history: OdePkg history and roadmap.
2662 * ide equations: IDE equations. (line 6)
2663 * ide options: ODE/DAE/IDE/DDE options.
2665 * InitialSlope option: ODE/DAE/IDE/DDE options.
2667 * InitialStep option: ODE/DAE/IDE/DDE options.
2669 * installation: Installation and deinstallation.
2671 * Jacobian option: ODE/DAE/IDE/DDE options.
2673 * JPattern option: ODE/DAE/IDE/DDE options.
2675 * m-file reference: M-File Function Reference.
2677 * Mass option: ODE/DAE/IDE/DDE options.
2679 * MassSingular option: ODE/DAE/IDE/DDE options.
2681 * MaxNewtonIterations option: ODE/DAE/IDE/DDE options.
2683 * MaxOrder option: ODE/DAE/IDE/DDE options.
2685 * MaxStep option: ODE/DAE/IDE/DDE options.
2687 * missing features: Missing features. (line 6)
2688 * MStateDependence option: ODE/DAE/IDE/DDE options.
2690 * MvPattern option: ODE/DAE/IDE/DDE options.
2692 * NewtonTol option: ODE/DAE/IDE/DDE options.
2694 * NonNegative option: ODE/DAE/IDE/DDE options.
2696 * NormControl option: ODE/DAE/IDE/DDE options.
2698 * oct-file reference: Oct-File Function Reference.
2700 * ode equations: ODE equations. (line 6)
2701 * ode options: ODE/DAE/IDE/DDE options.
2703 * OutputFcn option: ODE/DAE/IDE/DDE options.
2705 * OutputSel option: ODE/DAE/IDE/DDE options.
2707 * performance: ODE solver performances.
2709 * Programmers guide: Programmers Guide. (line 6)
2710 * Refine option: ODE/DAE/IDE/DDE options.
2712 * RelTol option: ODE/DAE/IDE/DDE options.
2714 * roadmap: OdePkg history and roadmap.
2716 * Runge-Kutta: Runge-Kutta solvers. (line 6)
2717 * Runge-Kutta modified: Modified Runge-Kutta solvers.
2719 * solver families: Solver families. (line 6)
2720 * Stats option: ODE/DAE/IDE/DDE options.
2722 * Users guide: Users Guide. (line 6)
2723 * Vectorized option: ODE/DAE/IDE/DDE options.
2730 Node: Beginners Guide
\7f1175
2731 Node: About OdePkg
\7f2340
2732 Node: OdePkg history and roadmap
\7f3448
2733 Node: Installation and deinstallation
\7f8050
2734 Node: Reporting Bugs
\7f8926
2735 Node: The "foo" example
\7f9602
2736 Node: Users Guide
\7f14647
2737 Node: Differential Equations
\7f15889
2738 Node: ODE equations
\7f16848
2739 Node: DAE equations
\7f17635
2740 Node: IDE equations
\7f19119
2741 Node: DDE equations
\7f19840
2742 Node: Solver families
\7f20556
2743 Node: Runge-Kutta solvers
\7f21606
2744 Ref: Runge-Kutta solvers-Footnote-1
\7f25521
2745 Node: Hairer-Wanner solvers
\7f25973
2746 Node: Cash modified BDF solvers
\7f27800
2747 Node: DDaskr direct method solver
\7f29441
2748 Node: Modified Runge-Kutta solvers
\7f31110
2749 Node: ODE solver performances
\7f32910
2750 Node: ODE/DAE/IDE/DDE options
\7f38580
2751 Node: M-File Function Reference
\7f56436
2752 Node: Oct-File Function Reference
\7f105239
2753 Node: Programmers Guide
\7f123795
2754 Node: Missing features
\7f124013
2755 Node: Function Index
\7f127528
2756 Node: Index
\7f133329
git://git.creatis.insa-lyon.fr / CreaPhase.git / blob