Mercurial > octave
annotate liboctave/numeric/LSODE.cc @ 30898:51a3d3a69193
maint: Use "fcn" as preferred abbreviation for "function" in liboctave/.
* DAEFunc.h, DASPK.cc, DASSL.cc, LSODE.cc, ODEFunc.h, eigs-base.cc,
eigs-base.h, oct-norm.cc, mx-inlines.cc:
Replace "func", "fun" in documentation and variable names with "fcn".
author | Rik <rik@octave.org> |
---|---|
date | Tue, 05 Apr 2022 13:20:48 -0700 |
parents | 796f54d4ddbf |
children | de6fc38c78c6 |
rev | line source |
---|---|
27923
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
1 //////////////////////////////////////////////////////////////////////// |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
2 // |
30564
796f54d4ddbf
update Octave Project Developers copyright for the new year
John W. Eaton <jwe@octave.org>
parents:
30394
diff
changeset
|
3 // Copyright (C) 1993-2022 The Octave Project Developers |
27923
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
4 // |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
5 // See the file COPYRIGHT.md in the top-level directory of this |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
6 // distribution or <https://octave.org/copyright/>. |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
7 // |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
8 // This file is part of Octave. |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
9 // |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
10 // Octave is free software: you can redistribute it and/or modify it |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
11 // under the terms of the GNU General Public License as published by |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
12 // the Free Software Foundation, either version 3 of the License, or |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
13 // (at your option) any later version. |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
14 // |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
15 // Octave is distributed in the hope that it will be useful, but |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
16 // WITHOUT ANY WARRANTY; without even the implied warranty of |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
18 // GNU General Public License for more details. |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
19 // |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
20 // You should have received a copy of the GNU General Public License |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
21 // along with Octave; see the file COPYING. If not, see |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
22 // <https://www.gnu.org/licenses/>. |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
23 // |
bd51beb6205e
update formatting of copyright notices
John W. Eaton <jwe@octave.org>
parents:
27919
diff
changeset
|
24 //////////////////////////////////////////////////////////////////////// |
3 | 25 |
21724
aba2e6293dd8
use "#if ..." consistently instead of "#ifdef" and "#ifndef"
John W. Eaton <jwe@octave.org>
parents:
21301
diff
changeset
|
26 #if defined (HAVE_CONFIG_H) |
21301
40de9f8f23a6
Use '#include "config.h"' rather than <config.h>.
Rik <rik@octave.org>
parents:
21202
diff
changeset
|
27 # include "config.h" |
3 | 28 #endif |
29 | |
26569
0e77df67b522
Add static compile-time checking of printf functions in liboctave.
Markus Mützel <markus.muetzel@gmx.de>
parents:
26376
diff
changeset
|
30 #include <cinttypes> |
5765 | 31 #include <sstream> |
32 | |
1842 | 33 #include "LSODE.h" |
1847 | 34 #include "f77-fcn.h" |
227 | 35 #include "lo-error.h" |
4180 | 36 #include "quit.h" |
3 | 37 |
30394
f3f3e3793fb5
maint: style check C++ files in liboctave/ ahead of 7.1 release.
Rik <rik@octave.org>
parents:
30067
diff
changeset
|
38 typedef F77_INT (*lsode_fcn_ptr) (const F77_INT&, const double&, double *, |
f3f3e3793fb5
maint: style check C++ files in liboctave/ ahead of 7.1 release.
Rik <rik@octave.org>
parents:
30067
diff
changeset
|
39 double *, F77_INT&); |
3507 | 40 |
30394
f3f3e3793fb5
maint: style check C++ files in liboctave/ ahead of 7.1 release.
Rik <rik@octave.org>
parents:
30067
diff
changeset
|
41 typedef F77_INT (*lsode_jac_ptr) (const F77_INT&, const double&, double *, |
f3f3e3793fb5
maint: style check C++ files in liboctave/ ahead of 7.1 release.
Rik <rik@octave.org>
parents:
30067
diff
changeset
|
42 const F77_INT&, const F77_INT&, double *, |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
43 const F77_INT&); |
3507 | 44 |
3 | 45 extern "C" |
4552 | 46 { |
47 F77_RET_T | |
30394
f3f3e3793fb5
maint: style check C++ files in liboctave/ ahead of 7.1 release.
Rik <rik@octave.org>
parents:
30067
diff
changeset
|
48 F77_FUNC (dlsode, DLSODE) (lsode_fcn_ptr, F77_INT&, F77_DBLE *, F77_DBLE&, |
f3f3e3793fb5
maint: style check C++ files in liboctave/ ahead of 7.1 release.
Rik <rik@octave.org>
parents:
30067
diff
changeset
|
49 F77_DBLE&, F77_INT&, F77_DBLE&, const F77_DBLE *, |
f3f3e3793fb5
maint: style check C++ files in liboctave/ ahead of 7.1 release.
Rik <rik@octave.org>
parents:
30067
diff
changeset
|
50 F77_INT&, F77_INT&, F77_INT&, F77_DBLE *, |
f3f3e3793fb5
maint: style check C++ files in liboctave/ ahead of 7.1 release.
Rik <rik@octave.org>
parents:
30067
diff
changeset
|
51 F77_INT&, F77_INT *, F77_INT&, lsode_jac_ptr, |
22133
59cadee1c74b
new macros for F77 data types
John W. Eaton <jwe@octave.org>
parents:
21724
diff
changeset
|
52 F77_INT&); |
4552 | 53 } |
3 | 54 |
30898
51a3d3a69193
maint: Use "fcn" as preferred abbreviation for "function" in liboctave/.
Rik <rik@octave.org>
parents:
30564
diff
changeset
|
55 static ODEFunc::ODERHSFunc user_fcn; |
532 | 56 static ODEFunc::ODEJacFunc user_jac; |
3 | 57 static ColumnVector *tmp_x; |
58 | |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
59 static F77_INT |
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
60 lsode_f (const F77_INT& neq, const double& time, double *, double *deriv, |
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
61 F77_INT& ierr) |
3 | 62 { |
2343 | 63 ColumnVector tmp_deriv; |
3 | 64 |
1360 | 65 // NOTE: this won't work if LSODE passes copies of the state vector. |
66 // In that case we have to create a temporary vector object | |
67 // and copy. | |
68 | |
30898
51a3d3a69193
maint: Use "fcn" as preferred abbreviation for "function" in liboctave/.
Rik <rik@octave.org>
parents:
30564
diff
changeset
|
69 tmp_deriv = (*user_fcn) (*tmp_x, time); |
3 | 70 |
23577
80c42f4cca13
maint: Deprecate is_empty and replace with isempty.
Rik <rik@octave.org>
parents:
23475
diff
changeset
|
71 if (tmp_deriv.isempty ()) |
1251 | 72 ierr = -1; |
258 | 73 else |
74 { | |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
75 for (F77_INT i = 0; i < neq; i++) |
15018
3d8ace26c5b4
maint: Use Octave coding conventions for cuddled parentheses in liboctave/.
Rik <rik@octave.org>
parents:
14138
diff
changeset
|
76 deriv[i] = tmp_deriv.elem (i); |
258 | 77 } |
3 | 78 |
79 return 0; | |
80 } | |
81 | |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
82 static F77_INT |
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
83 lsode_j (const F77_INT& neq, const double& time, double *, const F77_INT&, |
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
84 const F77_INT&, double *pd, const F77_INT& nrowpd) |
3 | 85 { |
1251 | 86 Matrix tmp_jac (neq, neq); |
3 | 87 |
1360 | 88 // NOTE: this won't work if LSODE passes copies of the state vector. |
89 // In that case we have to create a temporary vector object | |
90 // and copy. | |
91 | |
1251 | 92 tmp_jac = (*user_jac) (*tmp_x, time); |
3 | 93 |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
94 for (F77_INT j = 0; j < neq; j++) |
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
95 for (F77_INT i = 0; i < neq; i++) |
15020
560317fd5977
maint: Cuddle open bracket used for indexing C++ arrays in source code.
Rik <rik@octave.org>
parents:
15018
diff
changeset
|
96 pd[nrowpd * j + i] = tmp_jac (i, j); |
3 | 97 |
98 return 0; | |
99 } | |
100 | |
101 ColumnVector | |
1842 | 102 LSODE::do_integrate (double tout) |
3 | 103 { |
1945 | 104 ColumnVector retval; |
105 | |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
106 static F77_INT nn = 0; |
4049 | 107 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
108 if (! m_initialized || m_restart || ODEFunc::m_reset |
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
109 || LSODE_options::m_reset) |
1945 | 110 { |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
111 m_integration_error = false; |
1945 | 112 |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
113 m_initialized = true; |
4049 | 114 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
115 m_istate = 1; |
4049 | 116 |
22988
cd33c785e80e
put to_f77_int inside octave namespace
John W. Eaton <jwe@octave.org>
parents:
22960
diff
changeset
|
117 F77_INT n = octave::to_f77_int (size ()); |
4049 | 118 |
119 nn = n; | |
3955 | 120 |
5275 | 121 octave_idx_type max_maxord = 0; |
4231 | 122 |
4049 | 123 if (integration_method () == "stiff") |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
124 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
125 max_maxord = 5; |
4231 | 126 |
30047
02e28f7e4d04
maint: use "m_" prefix for member variables DAERTFunc/ODEFunc classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
30046
diff
changeset
|
127 if (m_jac) |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
128 m_method_flag = 21; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
129 else |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
130 m_method_flag = 22; |
3955 | 131 |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
132 m_liw = 20 + n; |
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
133 m_lrw = 22 + n * (9 + n); |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
134 } |
4049 | 135 else |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
136 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
137 max_maxord = 12; |
4231 | 138 |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
139 m_method_flag = 10; |
3955 | 140 |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
141 m_liw = 20; |
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
142 m_lrw = 22 + 16 * n; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
143 } |
4049 | 144 |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
145 m_iwork.resize (dim_vector (m_liw, 1)); |
5552 | 146 |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
147 for (F77_INT i = 4; i < 9; i++) |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
148 m_iwork(i) = 0; |
5552 | 149 |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
150 m_rwork.resize (dim_vector (m_lrw, 1)); |
5552 | 151 |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
152 for (F77_INT i = 4; i < 9; i++) |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
153 m_rwork(i) = 0; |
5552 | 154 |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
155 octave_idx_type maxord = maximum_order (); |
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
156 |
4231 | 157 if (maxord >= 0) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
158 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
159 if (maxord > 0 && maxord <= max_maxord) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
160 { |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
161 m_iwork(4) = octave::to_f77_int (maxord); |
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
162 m_iopt = 1; |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11574
diff
changeset
|
163 } |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
164 else |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
165 { |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20974
diff
changeset
|
166 // FIXME: Should this be a warning? |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
167 (*current_liboctave_error_handler) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
168 ("lsode: invalid value for maximum order"); |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
169 m_integration_error = true; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
170 return retval; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
171 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
172 } |
4231 | 173 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
174 if (m_stop_time_set) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
175 { |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
176 m_itask = 4; |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
177 m_rwork(0) = m_stop_time; |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
178 m_iopt = 1; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
179 } |
4049 | 180 else |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
181 { |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
182 m_itask = 1; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
183 } |
258 | 184 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
185 m_restart = false; |
4049 | 186 |
187 // ODEFunc | |
3 | 188 |
4049 | 189 // NOTE: this won't work if LSODE passes copies of the state vector. |
190 // In that case we have to create a temporary vector object | |
191 // and copy. | |
3 | 192 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
193 tmp_x = &m_x; |
4049 | 194 |
30898
51a3d3a69193
maint: Use "fcn" as preferred abbreviation for "function" in liboctave/.
Rik <rik@octave.org>
parents:
30564
diff
changeset
|
195 user_fcn = function (); |
4049 | 196 user_jac = jacobian_function (); |
197 | |
30898
51a3d3a69193
maint: Use "fcn" as preferred abbreviation for "function" in liboctave/.
Rik <rik@octave.org>
parents:
30564
diff
changeset
|
198 ColumnVector m_xdot = (*user_fcn) (m_x, m_t); |
2343 | 199 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
200 if (m_x.numel () != m_xdot.numel ()) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
201 { |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20974
diff
changeset
|
202 // FIXME: Should this be a warning? |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
203 (*current_liboctave_error_handler) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
204 ("lsode: inconsistent sizes for state and derivative vectors"); |
2343 | 205 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
206 m_integration_error = true; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
207 return retval; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
208 } |
2343 | 209 |
30047
02e28f7e4d04
maint: use "m_" prefix for member variables DAERTFunc/ODEFunc classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
30046
diff
changeset
|
210 ODEFunc::m_reset = false; |
4049 | 211 |
212 // LSODE_options | |
213 | |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
214 m_rel_tol = relative_tolerance (); |
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
215 m_abs_tol = absolute_tolerance (); |
4049 | 216 |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
217 F77_INT abs_tol_len = octave::to_f77_int (m_abs_tol.numel ()); |
4049 | 218 |
219 if (abs_tol_len == 1) | |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
220 m_itol = 1; |
4049 | 221 else if (abs_tol_len == n) |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
222 m_itol = 2; |
4049 | 223 else |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
224 { |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20974
diff
changeset
|
225 // FIXME: Should this be a warning? |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
226 (*current_liboctave_error_handler) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
227 ("lsode: inconsistent sizes for state and absolute tolerance vectors"); |
4049 | 228 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
229 m_integration_error = true; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
230 return retval; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
231 } |
2343 | 232 |
4049 | 233 double iss = initial_step_size (); |
234 if (iss >= 0.0) | |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
235 { |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
236 m_rwork(4) = iss; |
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
237 m_iopt = 1; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
238 } |
4049 | 239 |
240 double maxss = maximum_step_size (); | |
241 if (maxss >= 0.0) | |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
242 { |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
243 m_rwork(5) = maxss; |
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
244 m_iopt = 1; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
245 } |
4049 | 246 |
247 double minss = minimum_step_size (); | |
248 if (minss >= 0.0) | |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
249 { |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
250 m_rwork(6) = minss; |
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
251 m_iopt = 1; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
252 } |
4049 | 253 |
22988
cd33c785e80e
put to_f77_int inside octave namespace
John W. Eaton <jwe@octave.org>
parents:
22960
diff
changeset
|
254 F77_INT sl = octave::to_f77_int (step_limit ()); |
4049 | 255 if (sl > 0) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
256 { |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
257 m_iwork(5) = sl; |
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
258 m_iopt = 1; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
259 } |
4049 | 260 |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
261 LSODE_options::m_reset = false; |
3 | 262 } |
263 | |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
264 double *px = m_x.fortran_vec (); |
11502
4638800cd660
delete data pointer members from liboctave ODE/DAE classes; make destuctors virtual in ODE/DAE base classes
John W. Eaton <jwe@octave.org>
parents:
11495
diff
changeset
|
265 |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
266 double *pabs_tol = m_abs_tol.fortran_vec (); |
11502
4638800cd660
delete data pointer members from liboctave ODE/DAE classes; make destuctors virtual in ODE/DAE base classes
John W. Eaton <jwe@octave.org>
parents:
11495
diff
changeset
|
267 |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
268 F77_INT *piwork = m_iwork.fortran_vec (); |
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
269 double *prwork = m_rwork.fortran_vec (); |
11502
4638800cd660
delete data pointer members from liboctave ODE/DAE classes; make destuctors virtual in ODE/DAE base classes
John W. Eaton <jwe@octave.org>
parents:
11495
diff
changeset
|
270 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
271 F77_INT tmp_istate = octave::to_f77_int (m_istate); |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
272 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
273 F77_XFCN (dlsode, DLSODE, (lsode_f, nn, px, m_t, tout, m_itol, m_rel_tol, |
30046
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
274 pabs_tol, m_itask, tmp_istate, m_iopt, prwork, |
b3717fd85e49
maint: use "m_" prefix for member variables ODE/DAE classes in liboctave/numeric.
Rik <rik@octave.org>
parents:
29358
diff
changeset
|
275 m_lrw, piwork, m_liw, lsode_j, m_method_flag)); |
3 | 276 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
277 m_istate = tmp_istate; |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
278 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
279 switch (m_istate) |
3 | 280 { |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
281 case 1: // prior to initial integration step. |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
282 case 2: // lsode was successful. |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
283 retval = m_x; |
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
284 m_t = tout; |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
285 break; |
3996 | 286 |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
287 case -1: // excess work done on this call (perhaps wrong mf). |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
288 case -2: // excess accuracy requested (tolerances too small). |
8807
401d54a83690
use 'invalid', not 'illegal'
John W. Eaton <jwe@octave.org>
parents:
7482
diff
changeset
|
289 case -3: // invalid input detected (see printed message). |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
290 case -4: // repeated error test failures (check all inputs). |
12600
7c000c70f873
LSODE.cc: Add semicolon to error messages to prevent run-together text.
Rik <octave@nomad.inbox5.com>
parents:
11586
diff
changeset
|
291 case -5: // repeated convergence failures (perhaps bad Jacobian |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
292 // supplied or wrong choice of mf or tolerances). |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
293 case -6: // error weight became zero during problem. (solution |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
294 // component i vanished, and atol or atol(i) = 0.) |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
295 case -13: // return requested in user-supplied function. |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
296 m_integration_error = true; |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
297 break; |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
298 |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
299 default: |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
300 m_integration_error = true; |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
301 (*current_liboctave_error_handler) |
26569
0e77df67b522
Add static compile-time checking of printf functions in liboctave.
Markus Mützel <markus.muetzel@gmx.de>
parents:
26376
diff
changeset
|
302 ("unrecognized value of istate (= %" OCTAVE_IDX_TYPE_FORMAT ") " |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
303 "returned from lsode", m_istate); |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
304 break; |
3 | 305 } |
306 | |
1945 | 307 return retval; |
3 | 308 } |
309 | |
3959 | 310 std::string |
311 LSODE::error_message (void) const | |
312 { | |
313 std::string retval; | |
314 | |
5765 | 315 std::ostringstream buf; |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
316 buf << m_t; |
5765 | 317 std::string t_curr = buf.str (); |
4042 | 318 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
319 switch (m_istate) |
3959 | 320 { |
321 case 1: | |
3996 | 322 retval = "prior to initial integration step"; |
3959 | 323 break; |
324 | |
325 case 2: | |
326 retval = "successful exit"; | |
327 break; | |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11574
diff
changeset
|
328 |
3959 | 329 case 3: |
330 retval = "prior to continuation call with modified parameters"; | |
331 break; | |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11574
diff
changeset
|
332 |
3996 | 333 case -1: |
23829
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
334 retval = "excess work on this call (t = " + t_curr + |
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
335 "; perhaps wrong integration method)"; |
3996 | 336 break; |
337 | |
338 case -2: | |
339 retval = "excess accuracy requested (tolerances too small)"; | |
340 break; | |
341 | |
342 case -3: | |
343 retval = "invalid input detected (see printed message)"; | |
344 break; | |
345 | |
346 case -4: | |
23829
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
347 retval = "repeated error test failures (t = " + t_curr + |
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
348 "; check all inputs)"; |
3996 | 349 break; |
350 | |
351 case -5: | |
23829
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
352 retval = "repeated convergence failures (t = " + t_curr + |
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
353 "; perhaps bad Jacobian supplied or wrong choice of integration method or tolerances)"; |
3996 | 354 break; |
355 | |
356 case -6: | |
23829
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
357 retval = "error weight became zero during problem. (t = " + t_curr + |
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
358 "; solution component i vanished, and atol or atol(i) == 0)"; |
3996 | 359 break; |
360 | |
361 case -13: | |
4042 | 362 retval = "return requested in user-supplied function (t = " |
23807
336f89b6208b
Use character literals 'c' rather than string literals "c" when possible.
Rik <rik@octave.org>
parents:
23577
diff
changeset
|
363 + t_curr + ')'; |
3996 | 364 break; |
365 | |
3959 | 366 default: |
367 retval = "unknown error state"; | |
368 break; | |
369 } | |
370 | |
371 return retval; | |
372 } | |
373 | |
3 | 374 Matrix |
1842 | 375 LSODE::do_integrate (const ColumnVector& tout) |
3 | 376 { |
377 Matrix retval; | |
4049 | 378 |
20229
5dfaaaae784f
Deprecate Array::capacity() and Sparse::capacity() for numel() and nzmax().
Carnë Draug <carandraug@octave.org>
parents:
19697
diff
changeset
|
379 octave_idx_type n_out = tout.numel (); |
22988
cd33c785e80e
put to_f77_int inside octave namespace
John W. Eaton <jwe@octave.org>
parents:
22960
diff
changeset
|
380 F77_INT n = octave::to_f77_int (size ()); |
3 | 381 |
382 if (n_out > 0 && n > 0) | |
383 { | |
384 retval.resize (n_out, n); | |
385 | |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
386 for (F77_INT i = 0; i < n; i++) |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
387 retval.elem (0, i) = m_x.elem (i); |
3 | 388 |
5275 | 389 for (octave_idx_type j = 1; j < n_out; j++) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
390 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
391 ColumnVector x_next = do_integrate (tout.elem (j)); |
258 | 392 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
393 if (m_integration_error) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
394 return retval; |
258 | 395 |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
396 for (F77_INT i = 0; i < n; i++) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
397 retval.elem (j, i) = x_next.elem (i); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
398 } |
3 | 399 } |
400 | |
401 return retval; | |
402 } | |
403 | |
404 Matrix | |
3511 | 405 LSODE::do_integrate (const ColumnVector& tout, const ColumnVector& tcrit) |
3 | 406 { |
407 Matrix retval; | |
4049 | 408 |
20229
5dfaaaae784f
Deprecate Array::capacity() and Sparse::capacity() for numel() and nzmax().
Carnë Draug <carandraug@octave.org>
parents:
19697
diff
changeset
|
409 octave_idx_type n_out = tout.numel (); |
22988
cd33c785e80e
put to_f77_int inside octave namespace
John W. Eaton <jwe@octave.org>
parents:
22960
diff
changeset
|
410 F77_INT n = octave::to_f77_int (size ()); |
3 | 411 |
412 if (n_out > 0 && n > 0) | |
413 { | |
414 retval.resize (n_out, n); | |
415 | |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
416 for (F77_INT i = 0; i < n; i++) |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
417 retval.elem (0, i) = m_x.elem (i); |
3 | 418 |
20229
5dfaaaae784f
Deprecate Array::capacity() and Sparse::capacity() for numel() and nzmax().
Carnë Draug <carandraug@octave.org>
parents:
19697
diff
changeset
|
419 octave_idx_type n_crit = tcrit.numel (); |
3 | 420 |
421 if (n_crit > 0) | |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
422 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
423 octave_idx_type i_crit = 0; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
424 octave_idx_type i_out = 1; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
425 double next_crit = tcrit.elem (0); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
426 double next_out; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
427 while (i_out < n_out) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
428 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
429 bool do_restart = false; |
3 | 430 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
431 next_out = tout.elem (i_out); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
432 if (i_crit < n_crit) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
433 next_crit = tcrit.elem (i_crit); |
3 | 434 |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
435 bool save_output = false; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
436 double t_out; |
3 | 437 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
438 if (next_crit == next_out) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
439 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
440 set_stop_time (next_crit); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
441 t_out = next_out; |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
442 save_output = true; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
443 i_out++; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
444 i_crit++; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
445 do_restart = true; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
446 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
447 else if (next_crit < next_out) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
448 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
449 if (i_crit < n_crit) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
450 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
451 set_stop_time (next_crit); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
452 t_out = next_crit; |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
453 save_output = false; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
454 i_crit++; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
455 do_restart = true; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
456 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
457 else |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
458 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
459 clear_stop_time (); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
460 t_out = next_out; |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
461 save_output = true; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
462 i_out++; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
463 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
464 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
465 else |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
466 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
467 set_stop_time (next_crit); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
468 t_out = next_out; |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
469 save_output = true; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
470 i_out++; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
471 } |
3 | 472 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
473 ColumnVector x_next = do_integrate (t_out); |
3 | 474 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
475 if (m_integration_error) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
476 return retval; |
258 | 477 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
478 if (save_output) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
479 { |
22960
0d1422cb7e93
use F77_INT instead of octave_idx_type for liboctave ODE and DAE classes
John W. Eaton <jwe@octave.org>
parents:
22755
diff
changeset
|
480 for (F77_INT i = 0; i < n; i++) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
481 retval.elem (i_out-1, i) = x_next.elem (i); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
482 } |
3 | 483 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
484 if (do_restart) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
485 force_restart (); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
486 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
487 } |
3 | 488 else |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
489 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
490 retval = do_integrate (tout); |
258 | 491 |
30067
fc509c3c445a
maint: use "m_" prefix for member variables in base-de.h and base-dae.h
Rik <rik@octave.org>
parents:
30047
diff
changeset
|
492 if (m_integration_error) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
493 return retval; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
494 } |
3 | 495 } |
496 | |
497 return retval; | |
498 } |