annotate scripts/ode/ode15s.m @ 27918:b442ec6dda5c

use centralized file for copyright info for individual contributors * COPYRIGHT.md: New file. * In most other files, use "Copyright (C) YYYY-YYYY The Octave Project Developers" instead of tracking individual names in separate source files. The motivation is to reduce the effort required to update the notices each year. Until now, the Octave source files contained copyright notices that list individual contributors. I adopted these file-scope copyright notices because that is what everyone was doing 30 years ago in the days before distributed version control systems. But now, with many contributors and modern version control systems, having these file-scope copyright notices causes trouble when we update copyright years or refactor code. Over time, the file-scope copyright notices may become outdated as new contributions are made or code is moved from one file to another. Sometimes people contribute significant patches but do not add a line claiming copyright. Other times, people add a copyright notice for their contribution but then a later refactoring moves part or all of their contribution to another file and the notice is not moved with the code. As a practical matter, moving such notices is difficult -- determining what parts are due to a particular contributor requires a time-consuming search through the project history. Even managing the yearly update of copyright years is problematic. We have some contributors who are no longer living. Should we update the copyright dates for their contributions when we release new versions? Probably not, but we do still want to claim copyright for the project as a whole. To minimize the difficulty of maintaining the copyright notices, I would like to change Octave's sources to use what is described here: https://softwarefreedom.org/resources/2012/ManagingCopyrightInformation.html in the section "Maintaining centralized copyright notices": The centralized notice approach consolidates all copyright notices in a single location, usually a top-level file. This file should contain all of the copyright notices provided project contributors, unless the contribution was clearly insignificant. It may also credit -- without a copyright notice -- anyone who helped with the project but did not contribute code or other copyrighted material. This approach captures less information about contributions within individual files, recognizing that the DVCS is better equipped to record those details. As we mentioned before, it does have one disadvantage as compared to the file-scope approach: if a single file is separated from the distribution, the recipient won't see the contributors' copyright notices. But this can be easily remedied by including a single copyright notice in each file's header, pointing to the top-level file: Copyright YYYY-YYYY The Octave Project Developers See the COPYRIGHT file at the top-level directory of this distribution or at https://octave.org/COPYRIGHT.html. followed by the usual GPL copyright statement. For more background, see the discussion here: https://lists.gnu.org/archive/html/octave-maintainers/2020-01/msg00009.html Most files in the following directories have been skipped intentinally in this changeset: doc libgui/qterminal liboctave/external m4
author John W. Eaton <jwe@octave.org>
date Mon, 06 Jan 2020 15:38:17 -0500
parents ee6300e77c92
children 1891570abac8
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
27918
b442ec6dda5c use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents: 26894
diff changeset
1 ## Copyright (C) 2016-2019 The Octave Project Developers
b442ec6dda5c use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents: 26894
diff changeset
2 ##
b442ec6dda5c use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents: 26894
diff changeset
3 ## See the file COPYRIGHT.md in the top-level directory of this distribution
b442ec6dda5c use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents: 26894
diff changeset
4 ## or <https://octave.org/COPYRIGHT.html/>.
b442ec6dda5c use centralized file for copyright info for individual contributors
John W. Eaton <jwe@octave.org>
parents: 26894
diff changeset
5 ##
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
6 ##
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
7 ## This file is part of Octave.
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
8 ##
24534
194eb4bd202b maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents: 24012
diff changeset
9 ## Octave is free software: you can redistribute it and/or modify it
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
10 ## under the terms of the GNU General Public License as published by
24534
194eb4bd202b maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents: 24012
diff changeset
11 ## the Free Software Foundation, either version 3 of the License, or
24007
e8a74d95b4f3 maint: Use same format for Copyright statement throught code base.
Rik <rik@octave.org>
parents: 23023
diff changeset
12 ## (at your option) any later version.
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
13 ##
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
14 ## Octave is distributed in the hope that it will be useful, but
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
15 ## WITHOUT ANY WARRANTY; without even the implied warranty of
24007
e8a74d95b4f3 maint: Use same format for Copyright statement throught code base.
Rik <rik@octave.org>
parents: 23023
diff changeset
16 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
e8a74d95b4f3 maint: Use same format for Copyright statement throught code base.
Rik <rik@octave.org>
parents: 23023
diff changeset
17 ## GNU General Public License for more details.
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
18 ##
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
19 ## You should have received a copy of the GNU General Public License
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
20 ## along with Octave; see the file COPYING. If not, see
24534
194eb4bd202b maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents: 24012
diff changeset
21 ## <https://www.gnu.org/licenses/>.
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
22
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
23 ## -*- texinfo -*-
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
24 ## @deftypefn {} {[@var{t}, @var{y}] =} ode15s (@var{fun}, @var{trange}, @var{y0})
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
25 ## @deftypefnx {} {[@var{t}, @var{y}] =} ode15s (@var{fun}, @var{trange}, @var{y0}, @var{ode_opt})
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
26 ## @deftypefnx {} {[@var{t}, @var{y}, @var{te}, @var{ye}, @var{ie}] =} ode15s (@dots{})
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
27 ## @deftypefnx {} {@var{solution} =} ode15s (@dots{})
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
28 ## @deftypefnx {} {} ode15s (@dots{})
25113
476fc012d09a doc: Shorten very long first sentences of docstrings (bug #53388).
John W. Eaton <jwe@octave.org>
parents: 25054
diff changeset
29 ## Solve a set of stiff Ordinary Differential Equations (ODEs) or stiff
476fc012d09a doc: Shorten very long first sentences of docstrings (bug #53388).
John W. Eaton <jwe@octave.org>
parents: 25054
diff changeset
30 ## semi-explicit index 1 Differential Algebraic Equations (DAEs).
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
31 ##
25113
476fc012d09a doc: Shorten very long first sentences of docstrings (bug #53388).
John W. Eaton <jwe@octave.org>
parents: 25054
diff changeset
32 ## @code{ode15s} uses a variable step, variable order BDF (Backward
476fc012d09a doc: Shorten very long first sentences of docstrings (bug #53388).
John W. Eaton <jwe@octave.org>
parents: 25054
diff changeset
33 ## Differentiation Formula) method that ranges from order 1 to 5.
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
34 ##
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
35 ## @var{fun} is a function handle, inline function, or string containing the
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
36 ## name of the function that defines the ODE: @code{y' = f(t,y)}. The function
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
37 ## must accept two inputs where the first is time @var{t} and the second is a
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
38 ## column vector of unknowns @var{y}.
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
39 ##
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
40 ## @var{trange} specifies the time interval over which the ODE will be
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
41 ## evaluated. Typically, it is a two-element vector specifying the initial and
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
42 ## final times (@code{[tinit, tfinal]}). If there are more than two elements
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
43 ## then the solution will also be evaluated at these intermediate time
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
44 ## instances.
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
45 ##
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
46 ## @var{init} contains the initial value for the unknowns. If it is a row
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
47 ## vector then the solution @var{y} will be a matrix in which each column is
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
48 ## the solution for the corresponding initial value in @var{init}.
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
49 ##
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
50 ## The optional fourth argument @var{ode_opt} specifies non-default options to
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
51 ## the ODE solver. It is a structure generated by @code{odeset}.
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
52 ##
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
53 ## The function typically returns two outputs. Variable @var{t} is a
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
54 ## column vector and contains the times where the solution was found. The
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
55 ## output @var{y} is a matrix in which each column refers to a different
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
56 ## unknown of the problem and each row corresponds to a time in @var{t}.
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
57 ##
25026
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
58 ## The output can also be returned as a structure @var{solution} which has a
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
59 ## field @var{x} containing a row vector of times where the solution was
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
60 ## evaluated and a field @var{y} containing the solution matrix such that each
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
61 ## column corresponds to a time in @var{x}. Use
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
62 ## @w{@code{fieldnames (@var{solution})}} to see the other fields and
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
63 ## additional information returned.
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
64 ##
25026
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
65 ## If no output arguments are requested, and no @code{OutputFcn} is specified
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
66 ## in @var{ode_opt}, then the @code{OutputFcn} is set to @code{odeplot} and the
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
67 ## results of the solver are plotted immediately.
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
68 ##
25026
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
69 ## If using the @qcode{"Events"} option then three additional outputs may be
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
70 ## returned. @var{te} holds the time when an Event function returned a zero.
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
71 ## @var{ye} holds the value of the solution at time @var{te}. @var{ie}
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
72 ## contains an index indicating which Event function was triggered in the case
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
73 ## of multiple Event functions.
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
74 ##
25026
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
75 ## Example: Solve @nospell{Robertson's} equations:
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
76 ##
25026
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
77 ## @smallexample
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
78 ## @group
25026
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
79 ## function r = robertson_dae (@var{t}, @var{y})
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
80 ## r = [ -0.04*@var{y}(1) + 1e4*@var{y}(2)*@var{y}(3)
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
81 ## 0.04*@var{y}(1) - 1e4*@var{y}(2)*@var{y}(3) - 3e7*@var{y}(2)^2
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
82 ## @var{y}(1) + @var{y}(2) + @var{y}(3) - 1 ];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
83 ## endfunction
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
84 ## opt = odeset ("Mass", [1 0 0; 0 1 0; 0 0 0], "MStateDependence", "none");
25026
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
85 ## [@var{t},@var{y}] = ode15s (@@robertson_dae, [0, 1e3], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
86 ## @end group
25026
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
87 ## @end smallexample
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
88 ## @seealso{decic, odeset, odeget, ode23, ode45}
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
89 ## @end deftypefn
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
90
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
91 function varargout = ode15s (fun, trange, y0, varargin)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
92
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
93 if (nargin < 3)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
94 print_usage ();
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
95 endif
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
96
25026
f886561f9696 doc: improve differential eqtn docs and mention ode15i/ode15s (bug #51965).
Colin Macdonald <cbm@m.fsf.org>
parents: 24534
diff changeset
97 solver = "ode15s";
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
98 ## Check fun, trange, y0, yp0
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
99 fun = check_default_input (fun, trange, solver, y0);
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
100
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
101 n = numel (y0);
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
102
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
103 if (nargin > 3)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
104 options = varargin{1};
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
105 else
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
106 options = odeset ();
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
107 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
108
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
109 if (! isempty (options.Mass))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
110 if (ischar (options.Mass))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
111 try
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
112 options.Mass = str2func (options.Mass);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
113 catch
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
114 warning (lasterr);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
115 end_try_catch
25803
23483673ba43 Use is_function_handle instead of isa (x, "function_handle").
Rik <rik@octave.org>
parents: 25113
diff changeset
116 if (! is_function_handle (options.Mass))
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
117 error ("Octave:invalid-input-arg",
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
118 [solver ": invalid value assigned to field 'Mass'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
119 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
120 endif
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
121 endif
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
122
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
123 if (! isempty (options.Jacobian))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
124 if (ischar (options.Jacobian))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
125 try
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
126 options.Jacobian = str2func (options.Jacobian);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
127 catch
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
128 warning (lasterr);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
129 end_try_catch
25803
23483673ba43 Use is_function_handle instead of isa (x, "function_handle").
Rik <rik@octave.org>
parents: 25113
diff changeset
130 if (! is_function_handle (options.Jacobian))
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
131 error ("Octave:invalid-input-arg",
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
132 [solver ": invalid value assigned to field 'Jacobian'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
133 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
134 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
135 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
136
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
137 if (! isempty (options.OutputFcn))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
138 if (ischar (options.OutputFcn))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
139 try
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
140 options.OutputFcn = str2func (options.OutputFcn);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
141 catch
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
142 warning (lasterr);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
143 end_try_catch
25803
23483673ba43 Use is_function_handle instead of isa (x, "function_handle").
Rik <rik@octave.org>
parents: 25113
diff changeset
144 if (! is_function_handle (options.OutputFcn))
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
145 error ("Octave:invalid-input-arg",
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
146 [solver ": invalid value assigned to field '%s'"], "OutputFcn");
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
147 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
148 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
149 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
150
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
151 if (! isempty (options.Events))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
152 if (ischar (options.Events))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
153 try
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
154 options.Events = str2func (options.Events);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
155 catch
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
156 warning (lasterr);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
157 end_try_catch
25803
23483673ba43 Use is_function_handle instead of isa (x, "function_handle").
Rik <rik@octave.org>
parents: 25113
diff changeset
158 if (! is_function_handle (options.Events))
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
159 error ("Octave:invalid-input-arg",
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
160 [solver ": invalid value assigned to field 'Events'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
161 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
162 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
163 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
164
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
165 [defaults, classes, attributes] = odedefaults (n, trange(1), trange(end));
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
166
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
167 classes = odeset (classes, "Vectorized", {});
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
168 attributes = odeset (attributes, "Jacobian", {}, "Vectorized", {});
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
169
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
170 options = odemergeopts ("ode15s", options, defaults,
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
171 classes, attributes, solver);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
172
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
173 ## Mass
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
174 options.havemassfun = false;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
175 options.havestatedep = false;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
176 options.havetimedep = false;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
177 options.havemasssparse = false;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
178
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
179 if (! isempty (options.Mass))
25803
23483673ba43 Use is_function_handle instead of isa (x, "function_handle").
Rik <rik@octave.org>
parents: 25113
diff changeset
180 if (is_function_handle (options.Mass))
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
181 options.havemassfun = true;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
182 if (nargin (options.Mass) == 2)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
183 options.havestatedep = true;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
184 M = options.Mass (trange(1), y0);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
185 options.havemasssparse = issparse (M);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
186 if (any (size (M) != [n n]) || ! isnumeric (M) || ! isreal (M))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
187 error ("Octave:invalid-input-arg",
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
188 [solver ": invalid value assigned to field 'Mass'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
189 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
190 elseif (nargin (options.Mass) == 1)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
191 options.havetimedep = true;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
192 M = options.Mass (trange(1));
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
193 options.havemasssparse = issparse (M);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
194 if (any (size (M) != [n n]) || ! isnumeric (M) || ! isreal (M))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
195 error ("Octave:invalid-input-arg",
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
196 [solver ": invalid value assigned to field 'Mass'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
197 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
198 else
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
199 error ("Octave:invalid-input-arg",
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
200 [solver ": invalid value assigned to field 'Mass'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
201 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
202 elseif (ismatrix (options.Mass))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
203 options.havemasssparse = issparse (options.Mass);
25949
c2bf210ac94f style fixes
John W. Eaton <jwe@octave.org>
parents: 25113
diff changeset
204 if (any (size (options.Mass) != [n n])
c2bf210ac94f style fixes
John W. Eaton <jwe@octave.org>
parents: 25113
diff changeset
205 || ! isnumeric (options.Mass) || ! isreal (options.Mass))
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
206 error ("Octave:invalid-input-arg",
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
207 [solver ": invalid value assigned to field 'Mass'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
208 endif
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
209 else
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
210 error ("Octave:invalid-input-arg",
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
211 [solver ": invalid value assigned to field 'Mass'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
212 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
213 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
214
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
215 ## Jacobian
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
216 options.havejac = false;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
217 options.havejacsparse = false;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
218 options.havejacfun = false;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
219
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
220 if (! isempty (options.Jacobian))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
221 options.havejac = true;
25803
23483673ba43 Use is_function_handle instead of isa (x, "function_handle").
Rik <rik@octave.org>
parents: 25113
diff changeset
222 if (is_function_handle (options.Jacobian))
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
223 options.havejacfun = true;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
224 if (nargin (options.Jacobian) == 2)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
225 [A] = options.Jacobian (trange(1), y0);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
226 if (issparse (A))
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
227 options.havejacsparse = true; # Jac is sparse fun
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
228 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
229 if (any (size (A) != [n n]) || ! isnumeric (A) || ! isreal (A))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
230 error ("Octave:invalid-input-arg",
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
231 [solver ": invalid value assigned to field 'Jacobian'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
232 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
233 else
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
234 error ("Octave:invalid-input-arg",
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
235 [solver ": invalid value assigned to field 'Jacobian'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
236 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
237 elseif (ismatrix (options.Jacobian))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
238 if (issparse (options.Jacobian))
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
239 options.havejacsparse = true; # Jac is sparse matrix
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
240 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
241 if (! issquare (options.Jacobian))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
242 error ("Octave:invalid-input-arg",
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
243 [solver ": invalid value assigned to field 'Jacobian'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
244 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
245 else
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
246 error ("Octave:invalid-input-arg",
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
247 [solver ": invalid value assigned to field 'Jacobian'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
248 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
249 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
250
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
251 ## Derivative of M(t,y) for implicit problem not implemented yet
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
252 if (! isempty (options.Mass) && ! isempty (options.Jacobian))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
253 if (options.MStateDependence != "none" || options.havestatedep == true)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
254 options.havejac = false;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
255 options.Jacobian = [];
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
256 warning ("ode15s:mass_state_dependent_provided",
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
257 ["with MStateDependence != 'none' an internal", ...
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
258 " approximation of Jacobian Matrix will be used.", ...
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
259 " Set MStateDependence equal to 'none' if you want ", ...
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
260 " to provide a constant or time-dependent Jacobian"]);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
261 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
262 endif
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
263
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
264 ## Use sparse methods only if all matrices are sparse
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
265 if (! options.havemasssparse)
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
266 options.havejacsparse = false;
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
267 endif
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
268
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
269 ## If Mass or Jacobian is fun, then new Jacobian is fun
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
270 if (options.havejac)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
271 if (options.havejacfun || options.havetimedep)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
272 options.Jacobian = @ (t, y, yp) wrapjacfun (t, y, yp,
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
273 options.Jacobian,
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
274 options.Mass,
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
275 options.havetimedep,
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
276 options.havejacfun);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
277 options.havejacfun = true;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
278 else ## All matrices are constant
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
279 options.Jacobian = {[- options.Jacobian], [options.Mass]};
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
280
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
281 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
282 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
283
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
284 ## Abstol and Reltol
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
285 options.haveabstolvec = false;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
286
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
287 if (numel (options.AbsTol) != 1 && numel (options.AbsTol) != n)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
288 error ("Octave:invalid-input-arg",
22941
d92ec2901770 Make ode15i,ode15s doc consistent with other ode functions.
Rik <rik@octave.org>
parents: 22938
diff changeset
289 [solver ": invalid value assigned to field 'AbsTol'"]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
290 elseif (numel (options.AbsTol) == n)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
291 options.haveabstolvec = true;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
292 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
293
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
294 ## Stats
22935
c9344df03da5 Allow case-insensitive option argument 'on' to ode solvers (bug #49918).
Rik <rik@octave.org>
parents: 22933
diff changeset
295 options.havestats = strcmpi (options.Stats, "on");
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
296
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
297 ## Don't use Refine when the output is a structure
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
298 if (nargout == 1)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
299 options.Refine = 1;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
300 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
301
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
302 ## OutputFcn and OutputSel
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
303 if (isempty (options.OutputFcn) && nargout == 0)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
304 options.OutputFcn = @odeplot;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
305 options.haveoutputfunction = true;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
306 else
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
307 options.haveoutputfunction = ! isempty (options.OutputFcn);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
308 endif
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
309
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
310 options.haveoutputselection = ! isempty (options.OutputSel);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
311 if (options.haveoutputselection)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
312 options.OutputSel = options.OutputSel - 1;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
313 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
314
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
315 ## Events
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
316 options.haveeventfunction = ! isempty (options.Events);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
317
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
318 yp0 = options.InitialSlope;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
319
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
320 [t, y, te, ye, ie] = __ode15__ (@ (t, y, yp) wrap (t, y, yp, options.Mass,
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
321 options.havetimedep,
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
322 options.havestatedep,
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
323 fun),
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
324 trange, y0, yp0, options);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
325
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
326 if (nargout == 2)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
327 varargout{1} = t;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
328 varargout{2} = y;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
329 elseif (nargout == 1)
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
330 varargout{1}.x = t; # Time stamps are saved in field x
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
331 varargout{1}.y = y; # Results are saved in field y
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
332 varargout{1}.solver = solver;
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
333 if (options.haveeventfunction)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
334 varargout{1}.xe = te; # Time info when an event occurred
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
335 varargout{1}.ye = ye; # Results when an event occurred
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
336 varargout{1}.ie = ie; # Index info which event occurred
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
337 endif
22938
a88ceac2aa53 Fix undefined return argument for more than 2 outputs from ode15i,ode15s.
Rik <rik@octave.org>
parents: 22936
diff changeset
338 elseif (nargout > 2)
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
339 varargout = cell (1,5);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
340 varargout{1} = t;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
341 varargout{2} = y;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
342 if (options.haveeventfunction)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
343 varargout{3} = te; # Time info when an event occurred
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
344 varargout{4} = ye; # Results when an event occurred
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
345 varargout{5} = ie; # Index info which event occurred
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
346 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
347 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
348
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
349 endfunction
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
350
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
351 function res = wrap (t, y, yp, Mass, havetimedep, havestatedep, fun)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
352
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
353 if (! isempty (Mass) && havestatedep)
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
354 res = Mass (t, y) * yp - fun (t, y);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
355 elseif (! isempty (Mass) && havetimedep)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
356 res = Mass (t) * yp - fun (t, y);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
357 elseif (! isempty (Mass))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
358 res = Mass * yp - fun (t, y);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
359 else
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
360 res = yp - fun (t, y);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
361 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
362
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
363 endfunction
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
364
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
365 function [jac, jact] = wrapjacfun (t, y, yp, Jac, Mass,
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
366 havetimedep, havejacfun)
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
367 if (havejacfun)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
368 jac = - Jac (t, y);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
369 else
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
370 jac = - Jac;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
371 endif
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
372
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
373 if (! isempty (Mass) && havetimedep)
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
374 jact = Mass (t);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
375 elseif (! isempty (Mass))
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
376 jact = Mass;
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
377 else
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
378 jact = speye (numel (y));
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
379 endif
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
380
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
381 endfunction
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
382
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
383
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
384 %!demo
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
385 %! ## Solve Robertson's equations with ode15s
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
386 %! fun = @ (t, y) [-0.04*y(1) + 1e4*y(2).*y(3);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
387 %! 0.04*y(1) - 1e4*y(2).*y(3) - 3e7*y(2).^2;
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
388 %! y(1) + y(2) + y(3) - 1];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
389 %!
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
390 %! y0 = [1; 0; 0];
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
391 %! tspan = [0, 4*logspace(-6, 6)];
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
392 %! M = [1, 0, 0; 0, 1, 0; 0, 0, 0];
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
393 %!
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
394 %! options = odeset ("RelTol", 1e-4, "AbsTol", [1e-6, 1e-10, 1e-6],
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
395 %! "MStateDependence", "none", "Mass", M);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
396 %!
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
397 %! [t, y] = ode15s (fun, tspan, y0, options);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
398 %!
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
399 %! y(:,2) = 1e4 * y(:,2);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
400 %! figure (2);
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
401 %! semilogx (t, y, "o");
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
402 %! xlabel ("time");
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
403 %! ylabel ("species concentration");
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
404 %! title ("Robertson DAE problem with a Conservation Law");
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
405 %! legend ("y1", "y2", "y3");
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
406
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
407 %!function ydot = fpol (t, y) # Van der Pol equation
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
408 %! ydot = [y(2); (1 - y(1)^2) * y(2) - y(1)];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
409 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
410 %!
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
411 %!function ref = fref () # The computed reference sol
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
412 %! ref = [0.32331666704577, -1.83297456798624];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
413 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
414 %!
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
415 %!function jac = fjac (t, y) # its Jacobian
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
416 %! jac = [0, 1; -1 - 2 * y(1) * y(2), 1 - y(1)^2];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
417 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
418 %!
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
419 %!function jac = fjcc (t, y) # sparse type
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
420 %! jac = sparse ([0, 1; -1 - 2 * y(1) * y(2), 1 - y(1)^2]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
421 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
422 %!
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
423 %!function mas = fmas (t, y)
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
424 %! mas = [1, 0; 0, 1]; # Dummy mass matrix for tests
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
425 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
426 %!
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
427 %!function mas = fmsa (t, y)
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
428 %! mas = sparse ([1, 0; 0, 1]); # A sparse dummy matrix
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
429 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
430 %!
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
431 %!function res = rob (t, y)
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
432 %! res = [-0.04*y(1) + 1e4*y(2).*y(3);
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
433 %! 0.04*y(1) - 1e4*y(2).*y(3) - 3e7*y(2).^2;
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
434 %! y(1) + y(2) + y(3) - 1];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
435 %!endfunction
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
436 %!
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
437 %!function refrob = frefrob ()
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
438 %! refrob = [100, 0.617234887614937, 0.000006153591397, 0.382758958793666];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
439 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
440 %!
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
441 %!function [val, isterminal, direction] = feve (t, y)
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
442 %! isterminal = [0, 1];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
443 %! if (t < 1e1)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
444 %! val = [-1, -2];
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
445 %! else
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
446 %! val = [1, 3];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
447 %! endif
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
448 %!
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
449 %! direction = [1, 0];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
450 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
451 %!
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
452 %!function masrob = massdensefunstate (t, y)
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
453 %! masrob = [1, 0, 0; 0, 1, 0; 0, 0, 0];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
454 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
455 %!
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
456 %!function masrob = masssparsefunstate (t, y)
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
457 %! masrob = sparse ([1, 0, 0; 0, 1, 0; 0, 0, 0]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
458 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
459 %!
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
460 %!function masrob = massdensefuntime (t)
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
461 %! masrob = [1, 0, 0; 0, 1, 0; 0, 0, 0];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
462 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
463 %!
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
464 %!function masrob = masssparsefuntime (t)
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
465 %! masrob = sparse ([1, 0, 0; 0, 1, 0; 0, 0, 0]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
466 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
467 %!
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
468 %!function jac = jacfundense (t, y)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
469 %! jac = [-0.04, 1e4*y(3), 1e4*y(2);
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
470 %! 0.04, -1e4*y(3)-6e7*y(2), -1e4*y(2);
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
471 %! 1, 1, 1];
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
472 %!endfunction
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
473 %!
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
474 %!function jac = jacfunsparse (t, y)
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
475 %! jac = sparse([-0.04, 1e4*y(3), 1e4*y(2);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
476 %! 0.04, -1e4*y(3)-6e7*y(2), -1e4*y(2);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
477 %! 1, 1, 1]);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
478 %!endfunction
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
479
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
480 %!testif HAVE_SUNDIALS
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
481 %! opt = odeset ("MStateDependence", "none",
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
482 %! "Mass", [1, 0, 0; 0, 1, 0; 0, 0, 0]);
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
483 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
484 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
485
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
486 %!testif HAVE_SUNDIALS
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
487 %! opt = odeset ("MStateDependence", "none",
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
488 %! "Mass", sparse ([1, 0, 0; 0, 1, 0; 0, 0, 0]));
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
489 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
490 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
491
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
492 %!testif HAVE_SUNDIALS
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
493 %! opt = odeset ("MStateDependence", "none",
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
494 %! "Mass", @massdensefunstate);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
495 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
496 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
497
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
498 %!testif HAVE_SUNDIALS
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
499 %! opt = odeset ("MStateDependence", "none",
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
500 %! "Mass", @masssparsefunstate);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
501 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
502 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
503
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
504 %!testif HAVE_SUNDIALS
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
505 %! opt = odeset ("MStateDependence", "none",
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
506 %! "Mass", "massdensefuntime");
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
507 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
508 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
509
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
510 %!testif HAVE_SUNDIALS
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
511 %! opt = odeset ("MStateDependence", "none",
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
512 %! "Mass", [1, 0, 0; 0, 1, 0; 0, 0, 0],
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
513 %! "Jacobian", "jacfundense");
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
514 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
515 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
516
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
517 %!testif HAVE_SUNDIALS
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
518 %! opt = odeset ("MStateDependence", "none",
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
519 %! "Mass", sparse ([1, 0, 0; 0, 1, 0; 0, 0, 0]),
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
520 %! "Jacobian", @jacfundense);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
521 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
522 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
523
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
524 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
525 %! warning ("off", "ode15s:mass_state_dependent_provided", "local");
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
526 %! opt = odeset ("MStateDependence", "none",
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
527 %! "Mass", @massdensefunstate,
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
528 %! "Jacobian", @jacfundense);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
529 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
530 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
531
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
532 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
533 %! warning ("off", "ode15s:mass_state_dependent_provided", "local");
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
534 %! opt = odeset ("MStateDependence", "none",
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
535 %! "Mass", @masssparsefunstate,
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
536 %! "Jacobian", @jacfundense);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
537 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
538 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
539
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
540 %!testif HAVE_SUNDIALS
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
541 %! opt = odeset ("MStateDependence", "none",
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
542 %! "Mass", @massdensefuntime,
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
543 %! "Jacobian", @jacfundense);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
544 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
545 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
546 %! opt = odeset ("MStateDependence", "none",
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
547 %! "Mass", "masssparsefuntime",
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
548 %! "Jacobian", "jacfundense");
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
549 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
550 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
551
26894
ee6300e77c92 Update detection of sundials in the build system (bug #52475).
Carlo de Falco <carlo.defalco@polimi.it>
parents: 26376
diff changeset
552 %!testif HAVE_SUNDIALS_SUNLINSOL_KLU
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
553 %! opt = odeset ("MStateDependence", "none",
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
554 %! "Mass", [1, 0, 0; 0, 1, 0; 0, 0, 0],
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
555 %! "Jacobian", @jacfunsparse);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
556 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
557 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
558
26894
ee6300e77c92 Update detection of sundials in the build system (bug #52475).
Carlo de Falco <carlo.defalco@polimi.it>
parents: 26376
diff changeset
559 %!testif HAVE_SUNDIALS_SUNLINSOL_KLU
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
560 %! opt = odeset ("MStateDependence", "none",
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
561 %! "Mass", sparse ([1, 0, 0; 0, 1, 0; 0, 0, 0]),
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
562 %! "Jacobian", @jacfunsparse);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
563 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
564 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
565
26894
ee6300e77c92 Update detection of sundials in the build system (bug #52475).
Carlo de Falco <carlo.defalco@polimi.it>
parents: 26376
diff changeset
566 %!testif HAVE_SUNDIALS_SUNLINSOL_KLU
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
567 %! warning ("off", "ode15s:mass_state_dependent_provided", "local");
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
568 %! opt = odeset ("MStateDependence", "none",
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
569 %! "Mass", @massdensefunstate,
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
570 %! "Jacobian", @jacfunsparse);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
571 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
572 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
573
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
574 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
575 %! warning ("off", "ode15s:mass_state_dependent_provided", "local");
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
576 %! opt = odeset ("MStateDependence", "none",
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
577 %! "Mass", @masssparsefunstate,
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
578 %! "Jacobian", @jacfunsparse);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
579 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
580 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
581
26894
ee6300e77c92 Update detection of sundials in the build system (bug #52475).
Carlo de Falco <carlo.defalco@polimi.it>
parents: 26376
diff changeset
582 %!testif HAVE_SUNDIALS_SUNLINSOL_KLU
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
583 %! opt = odeset ("MStateDependence", "none",
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
584 %! "Mass", @massdensefuntime,
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
585 %! "Jacobian", @jacfunsparse);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
586 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
587 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
588
26894
ee6300e77c92 Update detection of sundials in the build system (bug #52475).
Carlo de Falco <carlo.defalco@polimi.it>
parents: 26376
diff changeset
589 %!testif HAVE_SUNDIALS_SUNLINSOL_KLU
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
590 %! opt = odeset ("MStateDependence", "none",
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
591 %! "Mass", @masssparsefuntime,
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
592 %! "Jacobian", @jacfunsparse);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
593 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
594 %! assert ([t(end), y(end,:)], frefrob, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
595
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
596 ## two output arguments
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
597 %!testif HAVE_SUNDIALS
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
598 %! [t, y] = ode15s (@fpol, [0, 2], [2, 0]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
599 %! assert ([t(end), y(end,:)], [2, fref], 1e-2);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
600
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
601 ## anonymous function instead of real function
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
602 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
603 %! fvdb = @(t,y) [y(2); (1 - y(1)^2) * y(2) - y(1)];
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
604 %! [t, y] = ode15s (fvdb, [0, 2], [2, 0]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
605 %! assert ([t(end), y(end,:)], [2, fref], 1e-2);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
606
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
607 ## Solve another anonymous function below zero
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
608 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
609 %! ref = [0, 14.77810590694212];
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
610 %! [t, y] = ode15s (@(t,y) y, [-2, 0], 2);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
611 %! assert ([t(end), y(end,:)], ref, 5e-2);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
612
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
613 ## InitialStep option
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
614 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
615 %! opt = odeset ("InitialStep", 1e-8);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
616 %! [t, y] = ode15s (@fpol, [0, 0.2], [2, 0], opt);
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
617 %! assert (t(2)-t(1), 1e-8, 1e-9);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
618
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
619 ## MaxStep option
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
620 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
621 %! opt = odeset ("MaxStep", 1e-3);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
622 %! sol = ode15s (@fpol, [0, 0.2], [2, 0], opt);
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
623 %! assert (sol.x(5)-sol.x(4), 1e-3, 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
624
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
625 ## Solve in backward direction starting at t=0
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
626 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
627 %! ref = [-1.205364552835178, 0.951542399860817];
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
628 %! sol = ode15s (@fpol, [0 -2], [2, 0]);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
629 %! assert ([sol.x(end), sol.y(end,:)], [-2, ref], 5e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
630
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
631 ## Solve in backward direction starting at t=2
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
632 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
633 %! ref = [-1.205364552835178, 0.951542399860817];
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
634 %! sol = ode15s (@fpol, [2, 0 -2], fref);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
635 %! assert ([sol.x(end), sol.y(end,:)], [-2, ref], 3e-2);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
636
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
637 ## Solve another anonymous function in backward direction
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
638 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
639 %! ref = [-1, 0.367879437558975];
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
640 %! sol = ode15s (@(t,y) y, [0 -1], 1);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
641 %! assert ([sol.x(end), sol.y(end,:)], ref, 1e-2);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
642
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
643 ## Solve another anonymous function below zero
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
644 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
645 %! ref = [0, 14.77810590694212];
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
646 %! sol = ode15s (@(t,y) y, [-2, 0], 2);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
647 %! assert ([sol.x(end), sol.y(end,:)], ref, 5e-2);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
648
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
649 ## Solve in backward direction starting at t=0 with MaxStep option
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
650 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
651 %! ref = [-1.205364552835178, 0.951542399860817];
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
652 %! opt = odeset ("MaxStep", 1e-3);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
653 %! sol = ode15s (@fpol, [0 -2], [2, 0], opt);
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
654 %! assert (abs (sol.x(8)-sol.x(7)), 1e-3, 1e-3);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
655 %! assert ([sol.x(end), sol.y(end,:)], [-2, ref], 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
656
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
657 ## AbsTol option
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
658 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
659 %! opt = odeset ("AbsTol", 1e-5);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
660 %! sol = ode15s (@fpol, [0, 2], [2, 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
661 %! assert ([sol.x(end), sol.y(end,:)], [2, fref], 4e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
662
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
663 ## AbsTol and RelTol option
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
664 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
665 %! opt = odeset ("AbsTol", 1e-8, "RelTol", 1e-8);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
666 %! sol = ode15s (@fpol, [0, 2], [2, 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
667 %! assert ([sol.x(end), sol.y(end,:)], [2, fref], 1e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
668
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
669 ## RelTol option -- higher accuracy
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
670 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
671 %! opt = odeset ("RelTol", 1e-8);
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
672 %! sol = ode15s (@fpol, [0, 2], [2, 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
673 %! assert ([sol.x(end), sol.y(end,:)], [2, fref], 1e-4);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
674
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
675 ## Mass option as function
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
676 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
677 %! opt = odeset ("Mass", @fmas, "MStateDependence", "none");
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
678 %! sol = ode15s (@fpol, [0, 2], [2, 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
679 %! assert ([sol.x(end), sol.y(end,:)], [2, fref], 3e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
680
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
681 ## Mass option as matrix
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
682 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
683 %! opt = odeset ("Mass", eye (2,2), "MStateDependence", "none");
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
684 %! sol = ode15s (@fpol, [0, 2], [2, 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
685 %! assert ([sol.x(end), sol.y(end,:)], [2, fref], 3e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
686
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
687 ## Mass option as sparse matrix
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
688 %!testif HAVE_SUNDIALS
23023
63a12df71848 avoid sparse jacobian tests for ode15i and ode15s if IDAKLU is missing
John W. Eaton <jwe@octave.org>
parents: 22941
diff changeset
689 %! opt = odeset ("Mass", speye (2), "MStateDependence", "none");
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
690 %! sol = ode15s (@fpol, [0, 2], [2, 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
691 %! assert ([sol.x(end), sol.y(end,:)], [2, fref], 3e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
692
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
693 ## Mass option as function and sparse matrix
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
694 %!testif HAVE_SUNDIALS
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
695 %! opt = odeset ("Mass", "fmsa", "MStateDependence", "none");
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
696 %! sol = ode15s (@fpol, [0, 2], [2, 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
697 %! assert ([sol.x(end), sol.y(end,:)], [2, fref], 3e-3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
698
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
699 ## Refine
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
700 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
701 %! opt2 = odeset ("Refine", 3, "Mass", @massdensefunstate,
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
702 %! "MStateDependence", "none");
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
703 %! opt1 = odeset ("Mass", @massdensefunstate, "MStateDependence", "none");
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
704 %! [t, y] = ode15s (@rob, [0, 100], [1; 0; 0], opt1);
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
705 %! [t2, y2] = ode15s (@rob, [0, 100], [1; 0; 0], opt2);
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
706 %! assert (numel (t2), numel (t) * 3, 3);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
707
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
708 ## Refine ignored if numel (trange) > 2
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
709 %!testif HAVE_SUNDIALS
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
710 %! opt2 = odeset ("Refine", 3, "Mass", "massdensefunstate",
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
711 %! "MStateDependence", "none");
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
712 %! opt1 = odeset ("Mass", @massdensefunstate, "MStateDependence", "none");
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
713 %! [t, y] = ode15s ("rob", [0, 10, 100], [1; 0; 0], opt1);
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
714 %! [t2, y2] = ode15s ("rob", [0, 10, 100], [1; 0; 0], opt2);
22933
c3428bb9aca9 ode15i.m, ode15s.m: Follow Octave coding conventions.
Rik <rik@octave.org>
parents: 22928
diff changeset
715 %! assert (numel (t2), numel (t));
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
716
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
717 ## Events option add further elements in sol
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
718 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
719 %! opt = odeset ("Events", @feve, "Mass", @massdensefunstate,
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
720 %! "MStateDependence", "none");
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
721 %! sol = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
722 %! assert (isfield (sol, "ie"));
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
723 %! assert (sol.ie, [0;1]);
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
724 %! assert (isfield (sol, "xe"));
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
725 %! assert (isfield (sol, "ye"));
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
726 %! assert (sol.x(end), 10, 1);
22928
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
727
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
728 ## Events option, five output arguments
dec22bceafa2 use testif HAVE_SUNDIALS for ode15 tests
John W. Eaton <jwe@octave.org>
parents: 22910
diff changeset
729 %!testif HAVE_SUNDIALS
22901
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
730 %! opt = odeset ("Events", @feve, "Mass", @massdensefunstate,
4c56f3ffec04 Add functions ode15i and ode15s
Francesco Faccio <francesco.faccio@mail.polimi.it>
parents:
diff changeset
731 %! "MStateDependence", "none");
22936
f2a1fc90a903 * ode15s.m: More minor style fixes in tests.
John W. Eaton <jwe@octave.org>
parents: 22935
diff changeset
732 %! [t, y, te, ye, ie] = ode15s (@rob, [0, 100], [1; 0; 0], opt);
22910
23847979b91e maint: Strip trailing whitespace from source files.
John W. Eaton <jwe@octave.org>
parents: 22902
diff changeset
733 %! assert ([t(end), te', ie'], [10, 10, 10, 0, 1], [1, 0.5, 0.5, 0, 0]);
24008
6e7bb85e32b8 Allow ode15s and ode15i to work with row initial vectors (bug #50192).
Marco Caliari <marco.caliari@univr.it>
parents: 24007
diff changeset
734
6e7bb85e32b8 Allow ode15s and ode15i to work with row initial vectors (bug #50192).
Marco Caliari <marco.caliari@univr.it>
parents: 24007
diff changeset
735 ## Initial solution as row vector
6e7bb85e32b8 Allow ode15s and ode15i to work with row initial vectors (bug #50192).
Marco Caliari <marco.caliari@univr.it>
parents: 24007
diff changeset
736 %!testif HAVE_SUNDIALS
6e7bb85e32b8 Allow ode15s and ode15i to work with row initial vectors (bug #50192).
Marco Caliari <marco.caliari@univr.it>
parents: 24007
diff changeset
737 %! A = zeros (2);
6e7bb85e32b8 Allow ode15s and ode15i to work with row initial vectors (bug #50192).
Marco Caliari <marco.caliari@univr.it>
parents: 24007
diff changeset
738 %! [tout, yout] = ode15s (@(t, y) A * y, [0, 1], [1, 1]);
6e7bb85e32b8 Allow ode15s and ode15i to work with row initial vectors (bug #50192).
Marco Caliari <marco.caliari@univr.it>
parents: 24007
diff changeset
739 %! assert (yout, ones (18, 2))
6e7bb85e32b8 Allow ode15s and ode15i to work with row initial vectors (bug #50192).
Marco Caliari <marco.caliari@univr.it>
parents: 24007
diff changeset
740
6e7bb85e32b8 Allow ode15s and ode15i to work with row initial vectors (bug #50192).
Marco Caliari <marco.caliari@univr.it>
parents: 24007
diff changeset
741 %!testif HAVE_SUNDIALS
6e7bb85e32b8 Allow ode15s and ode15i to work with row initial vectors (bug #50192).
Marco Caliari <marco.caliari@univr.it>
parents: 24007
diff changeset
742 %! A = zeros (2);
6e7bb85e32b8 Allow ode15s and ode15i to work with row initial vectors (bug #50192).
Marco Caliari <marco.caliari@univr.it>
parents: 24007
diff changeset
743 %! fail ("ode15s (@(t, y) A * y, [0, 1], eye (2))",
24012
3d65720cd68a ode15s.m: Fix message in BIST test (bug #50192).
Rik <rik@octave.org>
parents: 24008
diff changeset
744 %! "ode15s: Y0 must be a numeric vector");