Mercurial > octave
annotate liboctave/numeric/LSODE.cc @ 26569:0e77df67b522 stable
Add static compile-time checking of printf functions in liboctave.
* util/lo-error.h: Use macro OCTAVE_FORMAT_PRINTF for printf functions.
* array/Array.h: Use format specifier "z" for type size_t.
* array/Array-util.cc, numeric/oct-norm.cc: Supply a format string of "%s" when
printing a single string.
* array/idx-vector.cc, array/MatrixType.cc, array/Sparse.cc, numeric/DASPK.cc,
numeric/DASRT.cc, numeric/DASSL.cc, numeric/LSODE.cc,
util/lo-array-errwarn.cc, util/sparse-util.cc: Use preprocessor macro
OCTAVE_IDX_TYPE_FORMAT as the format identifier for octave_idx_type.
* util/f77-fcn.c: Cast width argument to int.
* util/lo-utils.cc: Pass missing argument c0.
author | Markus Mützel <markus.muetzel@gmx.de> |
---|---|
date | Thu, 17 Jan 2019 19:25:55 +0100 |
parents | 00f796120a6d |
children | 2e5d8cdc7087 |
rev | line source |
---|---|
3 | 1 /* |
2 | |
26376
00f796120a6d
maint: Update copyright dates in all source files.
John W. Eaton <jwe@octave.org>
parents:
25054
diff
changeset
|
3 Copyright (C) 1993-2019 John W. Eaton |
3 | 4 |
5 This file is part of Octave. | |
6 | |
24534
194eb4bd202b
maint: Update punctuation for GPL v3 license text.
Rik <rik@octave.org>
parents:
23829
diff
changeset
|
7 Octave is free software: you can redistribute it and/or modify it |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22402
diff
changeset
|
8 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:
23829
diff
changeset
|
9 the Free Software Foundation, either version 3 of the License, or |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22402
diff
changeset
|
10 (at your option) any later version. |
3 | 11 |
22755
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22402
diff
changeset
|
12 Octave is distributed in the hope that it will be useful, but |
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22402
diff
changeset
|
13 WITHOUT ANY WARRANTY; without even the implied warranty of |
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22402
diff
changeset
|
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
3a2b891d0b33
maint: Standardize Copyright formatting.
Rik <rik@octave.org>
parents:
22402
diff
changeset
|
15 GNU General Public License for more details. |
3 | 16 |
17 You should have received a copy of the GNU General Public License | |
7016 | 18 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:
23829
diff
changeset
|
19 <https://www.gnu.org/licenses/>. |
3 | 20 |
21 */ | |
22 | |
21724
aba2e6293dd8
use "#if ..." consistently instead of "#ifdef" and "#ifndef"
John W. Eaton <jwe@octave.org>
parents:
21301
diff
changeset
|
23 #if defined (HAVE_CONFIG_H) |
21301
40de9f8f23a6
Use '#include "config.h"' rather than <config.h>.
Rik <rik@octave.org>
parents:
21202
diff
changeset
|
24 # include "config.h" |
3 | 25 #endif |
26 | |
26569
0e77df67b522
Add static compile-time checking of printf functions in liboctave.
Markus Mützel <markus.muetzel@gmx.de>
parents:
26376
diff
changeset
|
27 #include <cinttypes> |
5765 | 28 #include <sstream> |
29 | |
1842 | 30 #include "LSODE.h" |
1847 | 31 #include "f77-fcn.h" |
227 | 32 #include "lo-error.h" |
4180 | 33 #include "quit.h" |
3 | 34 |
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
|
35 typedef F77_INT (*lsode_fcn_ptr) (const F77_INT&, const double&, double*, |
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
|
36 double*, F77_INT&); |
3507 | 37 |
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
|
38 typedef F77_INT (*lsode_jac_ptr) (const F77_INT&, const double&, double*, |
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
|
39 const F77_INT&, const F77_INT&, double*, |
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
|
40 const F77_INT&); |
3507 | 41 |
3 | 42 extern "C" |
4552 | 43 { |
44 F77_RET_T | |
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
|
45 F77_FUNC (dlsode, DLSODE) (lsode_fcn_ptr, F77_INT&, F77_DBLE*, F77_DBLE&, |
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
|
46 F77_DBLE&, F77_INT&, F77_DBLE&, const F77_DBLE*, |
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
|
47 F77_INT&, F77_INT&, F77_INT&, F77_DBLE*, |
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
|
48 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
|
49 F77_INT&); |
4552 | 50 } |
3 | 51 |
532 | 52 static ODEFunc::ODERHSFunc user_fun; |
53 static ODEFunc::ODEJacFunc user_jac; | |
3 | 54 static ColumnVector *tmp_x; |
55 | |
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
|
56 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
|
57 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
|
58 F77_INT& ierr) |
3 | 59 { |
4180 | 60 BEGIN_INTERRUPT_WITH_EXCEPTIONS; |
61 | |
2343 | 62 ColumnVector tmp_deriv; |
3 | 63 |
1360 | 64 // NOTE: this won't work if LSODE passes copies of the state vector. |
65 // In that case we have to create a temporary vector object | |
66 // and copy. | |
67 | |
1251 | 68 tmp_deriv = (*user_fun) (*tmp_x, time); |
3 | 69 |
23577
80c42f4cca13
maint: Deprecate is_empty and replace with isempty.
Rik <rik@octave.org>
parents:
23475
diff
changeset
|
70 if (tmp_deriv.isempty ()) |
1251 | 71 ierr = -1; |
258 | 72 else |
73 { | |
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
|
74 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
|
75 deriv[i] = tmp_deriv.elem (i); |
258 | 76 } |
3 | 77 |
4180 | 78 END_INTERRUPT_WITH_EXCEPTIONS; |
79 | |
3 | 80 return 0; |
81 } | |
82 | |
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
|
83 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
|
84 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
|
85 const F77_INT&, double *pd, const F77_INT& nrowpd) |
3 | 86 { |
4180 | 87 BEGIN_INTERRUPT_WITH_EXCEPTIONS; |
88 | |
1251 | 89 Matrix tmp_jac (neq, neq); |
3 | 90 |
1360 | 91 // NOTE: this won't work if LSODE passes copies of the state vector. |
92 // In that case we have to create a temporary vector object | |
93 // and copy. | |
94 | |
1251 | 95 tmp_jac = (*user_jac) (*tmp_x, time); |
3 | 96 |
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
|
97 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
|
98 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
|
99 pd[nrowpd * j + i] = tmp_jac (i, j); |
3 | 100 |
4180 | 101 END_INTERRUPT_WITH_EXCEPTIONS; |
102 | |
3 | 103 return 0; |
104 } | |
105 | |
106 ColumnVector | |
1842 | 107 LSODE::do_integrate (double tout) |
3 | 108 { |
1945 | 109 ColumnVector retval; |
110 | |
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
|
111 static F77_INT nn = 0; |
4049 | 112 |
113 if (! initialized || restart || ODEFunc::reset || LSODE_options::reset) | |
1945 | 114 { |
4049 | 115 integration_error = false; |
1945 | 116 |
4049 | 117 initialized = true; |
118 | |
119 istate = 1; | |
120 | |
22988
cd33c785e80e
put to_f77_int inside octave namespace
John W. Eaton <jwe@octave.org>
parents:
22960
diff
changeset
|
121 F77_INT n = octave::to_f77_int (size ()); |
4049 | 122 |
123 nn = n; | |
3955 | 124 |
5275 | 125 octave_idx_type max_maxord = 0; |
4231 | 126 |
4049 | 127 if (integration_method () == "stiff") |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
128 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
129 max_maxord = 5; |
4231 | 130 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
131 if (jac) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
132 method_flag = 21; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
133 else |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
134 method_flag = 22; |
3955 | 135 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
136 liw = 20 + n; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
137 lrw = 22 + n * (9 + n); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
138 } |
4049 | 139 else |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
140 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
141 max_maxord = 12; |
4231 | 142 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
143 method_flag = 10; |
3955 | 144 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
145 liw = 20; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
146 lrw = 22 + 16 * n; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
147 } |
4049 | 148 |
11574
a83bad07f7e3
attempt better backward compatibility for Array resize functions
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
149 iwork.resize (dim_vector (liw, 1)); |
5552 | 150 |
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
|
151 for (F77_INT i = 4; i < 9; i++) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
152 iwork(i) = 0; |
5552 | 153 |
11574
a83bad07f7e3
attempt better backward compatibility for Array resize functions
John W. Eaton <jwe@octave.org>
parents:
11523
diff
changeset
|
154 rwork.resize (dim_vector (lrw, 1)); |
5552 | 155 |
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
|
156 for (F77_INT i = 4; i < 9; i++) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
157 rwork(i) = 0; |
5552 | 158 |
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
|
159 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
|
160 |
4231 | 161 if (maxord >= 0) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
162 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
163 if (maxord > 0 && maxord <= max_maxord) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
164 { |
22988
cd33c785e80e
put to_f77_int inside octave namespace
John W. Eaton <jwe@octave.org>
parents:
22960
diff
changeset
|
165 iwork(4) = octave::to_f77_int (maxord); |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
166 iopt = 1; |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11574
diff
changeset
|
167 } |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
168 else |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
169 { |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20974
diff
changeset
|
170 // FIXME: Should this be a warning? |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
171 (*current_liboctave_error_handler) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
172 ("lsode: invalid value for maximum order"); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
173 integration_error = true; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
174 return retval; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
175 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
176 } |
4231 | 177 |
4049 | 178 if (stop_time_set) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
179 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
180 itask = 4; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
181 rwork(0) = stop_time; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
182 iopt = 1; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
183 } |
4049 | 184 else |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
185 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
186 itask = 1; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
187 } |
258 | 188 |
4049 | 189 restart = false; |
190 | |
191 // ODEFunc | |
3 | 192 |
4049 | 193 // NOTE: this won't work if LSODE passes copies of the state vector. |
194 // In that case we have to create a temporary vector object | |
195 // and copy. | |
3 | 196 |
4049 | 197 tmp_x = &x; |
198 | |
199 user_fun = function (); | |
200 user_jac = jacobian_function (); | |
201 | |
2343 | 202 ColumnVector xdot = (*user_fun) (x, t); |
203 | |
20232
a9574e3c6e9e
Deprecate Array::length() and Sparse::length() in favour of ::numel().
Carnë Draug <carandraug@octave.org>
parents:
20229
diff
changeset
|
204 if (x.numel () != xdot.numel ()) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
205 { |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20974
diff
changeset
|
206 // FIXME: Should this be a warning? |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
207 (*current_liboctave_error_handler) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
208 ("lsode: inconsistent sizes for state and derivative vectors"); |
2343 | 209 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
210 integration_error = true; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
211 return retval; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
212 } |
2343 | 213 |
4049 | 214 ODEFunc::reset = false; |
215 | |
216 // LSODE_options | |
217 | |
218 rel_tol = relative_tolerance (); | |
219 abs_tol = absolute_tolerance (); | |
220 | |
22988
cd33c785e80e
put to_f77_int inside octave namespace
John W. Eaton <jwe@octave.org>
parents:
22960
diff
changeset
|
221 F77_INT abs_tol_len = octave::to_f77_int (abs_tol.numel ()); |
4049 | 222 |
223 if (abs_tol_len == 1) | |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
224 itol = 1; |
4049 | 225 else if (abs_tol_len == n) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
226 itol = 2; |
4049 | 227 else |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
228 { |
21136
7cac4e7458f2
maint: clean up code around calls to current_liboctave_error_handler.
Rik <rik@octave.org>
parents:
20974
diff
changeset
|
229 // FIXME: Should this be a warning? |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
230 (*current_liboctave_error_handler) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
231 ("lsode: inconsistent sizes for state and absolute tolerance vectors"); |
4049 | 232 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
233 integration_error = true; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
234 return retval; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
235 } |
2343 | 236 |
4049 | 237 double iss = initial_step_size (); |
238 if (iss >= 0.0) | |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
239 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
240 rwork(4) = iss; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
241 iopt = 1; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
242 } |
4049 | 243 |
244 double maxss = maximum_step_size (); | |
245 if (maxss >= 0.0) | |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
246 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
247 rwork(5) = maxss; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
248 iopt = 1; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
249 } |
4049 | 250 |
251 double minss = minimum_step_size (); | |
252 if (minss >= 0.0) | |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
253 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
254 rwork(6) = minss; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
255 iopt = 1; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
256 } |
4049 | 257 |
22988
cd33c785e80e
put to_f77_int inside octave namespace
John W. Eaton <jwe@octave.org>
parents:
22960
diff
changeset
|
258 F77_INT sl = octave::to_f77_int (step_limit ()); |
4049 | 259 if (sl > 0) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
260 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
261 iwork(5) = sl; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
262 iopt = 1; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
263 } |
4049 | 264 |
265 LSODE_options::reset = false; | |
3 | 266 } |
267 | |
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
|
268 double *px = x.fortran_vec (); |
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
|
269 |
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 double *pabs_tol = abs_tol.fortran_vec (); |
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
|
271 |
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 F77_INT *piwork = iwork.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
|
273 double *prwork = rwork.fortran_vec (); |
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
|
274 |
22988
cd33c785e80e
put to_f77_int inside octave namespace
John W. Eaton <jwe@octave.org>
parents:
22960
diff
changeset
|
275 F77_INT tmp_istate = octave::to_f77_int (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
|
276 |
4583 | 277 F77_XFCN (dlsode, DLSODE, (lsode_f, nn, px, t, tout, itol, rel_tol, |
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 pabs_tol, itask, tmp_istate, iopt, prwork, lrw, |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
279 piwork, liw, lsode_j, method_flag)); |
3 | 280 |
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
|
281 istate = tmp_istate; |
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
|
282 |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
283 switch (istate) |
3 | 284 { |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
285 case 1: // prior to initial integration step. |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
286 case 2: // lsode was successful. |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
287 retval = x; |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
288 t = tout; |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
289 break; |
3996 | 290 |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
291 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
|
292 case -2: // excess accuracy requested (tolerances too small). |
8807
401d54a83690
use 'invalid', not 'illegal'
John W. Eaton <jwe@octave.org>
parents:
7482
diff
changeset
|
293 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
|
294 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
|
295 case -5: // repeated convergence failures (perhaps bad Jacobian |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
296 // 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
|
297 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
|
298 // 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
|
299 case -13: // return requested in user-supplied function. |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
300 integration_error = true; |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
301 break; |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
302 |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
303 default: |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
304 integration_error = true; |
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
305 (*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
|
306 ("unrecognized value of istate (= %" OCTAVE_IDX_TYPE_FORMAT ") " |
0e77df67b522
Add static compile-time checking of printf functions in liboctave.
Markus Mützel <markus.muetzel@gmx.de>
parents:
26376
diff
changeset
|
307 "returned from lsode", istate); |
7482
29980c6b8604
don't check f77_exception_encountered
John W. Eaton <jwe@octave.org>
parents:
7231
diff
changeset
|
308 break; |
3 | 309 } |
310 | |
1945 | 311 return retval; |
3 | 312 } |
313 | |
3959 | 314 std::string |
315 LSODE::error_message (void) const | |
316 { | |
317 std::string retval; | |
318 | |
5765 | 319 std::ostringstream buf; |
320 buf << t; | |
321 std::string t_curr = buf.str (); | |
4042 | 322 |
3959 | 323 switch (istate) |
324 { | |
325 case 1: | |
3996 | 326 retval = "prior to initial integration step"; |
3959 | 327 break; |
328 | |
329 case 2: | |
330 retval = "successful exit"; | |
331 break; | |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11574
diff
changeset
|
332 |
3959 | 333 case 3: |
334 retval = "prior to continuation call with modified parameters"; | |
335 break; | |
11586
12df7854fa7c
strip trailing whitespace from source files
John W. Eaton <jwe@octave.org>
parents:
11574
diff
changeset
|
336 |
3996 | 337 case -1: |
23829
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
338 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
|
339 "; perhaps wrong integration method)"; |
3996 | 340 break; |
341 | |
342 case -2: | |
343 retval = "excess accuracy requested (tolerances too small)"; | |
344 break; | |
345 | |
346 case -3: | |
347 retval = "invalid input detected (see printed message)"; | |
348 break; | |
349 | |
350 case -4: | |
23829
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
351 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
|
352 "; check all inputs)"; |
3996 | 353 break; |
354 | |
355 case -5: | |
23829
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
356 retval = "repeated convergence failures (t = " + t_curr + |
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
357 "; perhaps bad Jacobian supplied or wrong choice of integration method or tolerances)"; |
3996 | 358 break; |
359 | |
360 case -6: | |
23829
01899bdd2a3a
Eliminate unnecessary std::string ("...") constructor calls when "..." suffices.
Rik <rik@octave.org>
parents:
23807
diff
changeset
|
361 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
|
362 "; solution component i vanished, and atol or atol(i) == 0)"; |
3996 | 363 break; |
364 | |
365 case -13: | |
4042 | 366 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
|
367 + t_curr + ')'; |
3996 | 368 break; |
369 | |
3959 | 370 default: |
371 retval = "unknown error state"; | |
372 break; | |
373 } | |
374 | |
375 return retval; | |
376 } | |
377 | |
3 | 378 Matrix |
1842 | 379 LSODE::do_integrate (const ColumnVector& tout) |
3 | 380 { |
381 Matrix retval; | |
4049 | 382 |
20229
5dfaaaae784f
Deprecate Array::capacity() and Sparse::capacity() for numel() and nzmax().
Carnë Draug <carandraug@octave.org>
parents:
19697
diff
changeset
|
383 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
|
384 F77_INT n = octave::to_f77_int (size ()); |
3 | 385 |
386 if (n_out > 0 && n > 0) | |
387 { | |
388 retval.resize (n_out, n); | |
389 | |
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
|
390 for (F77_INT i = 0; i < n; i++) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
391 retval.elem (0, i) = x.elem (i); |
3 | 392 |
5275 | 393 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
|
394 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
395 ColumnVector x_next = do_integrate (tout.elem (j)); |
258 | 396 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
397 if (integration_error) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
398 return retval; |
258 | 399 |
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
|
400 for (F77_INT i = 0; i < n; i++) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
401 retval.elem (j, i) = x_next.elem (i); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
402 } |
3 | 403 } |
404 | |
405 return retval; | |
406 } | |
407 | |
408 Matrix | |
3511 | 409 LSODE::do_integrate (const ColumnVector& tout, const ColumnVector& tcrit) |
3 | 410 { |
411 Matrix retval; | |
4049 | 412 |
20229
5dfaaaae784f
Deprecate Array::capacity() and Sparse::capacity() for numel() and nzmax().
Carnë Draug <carandraug@octave.org>
parents:
19697
diff
changeset
|
413 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
|
414 F77_INT n = octave::to_f77_int (size ()); |
3 | 415 |
416 if (n_out > 0 && n > 0) | |
417 { | |
418 retval.resize (n_out, n); | |
419 | |
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
|
420 for (F77_INT i = 0; i < n; i++) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
421 retval.elem (0, i) = x.elem (i); |
3 | 422 |
20229
5dfaaaae784f
Deprecate Array::capacity() and Sparse::capacity() for numel() and nzmax().
Carnë Draug <carandraug@octave.org>
parents:
19697
diff
changeset
|
423 octave_idx_type n_crit = tcrit.numel (); |
3 | 424 |
425 if (n_crit > 0) | |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
426 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
427 octave_idx_type i_crit = 0; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
428 octave_idx_type i_out = 1; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
429 double next_crit = tcrit.elem (0); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
430 double next_out; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
431 while (i_out < n_out) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
432 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
433 bool do_restart = false; |
3 | 434 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
435 next_out = tout.elem (i_out); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
436 if (i_crit < n_crit) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
437 next_crit = tcrit.elem (i_crit); |
3 | 438 |
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
|
439 bool save_output = false; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
440 double t_out; |
3 | 441 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
442 if (next_crit == next_out) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
443 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
444 set_stop_time (next_crit); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
445 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
|
446 save_output = true; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
447 i_out++; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
448 i_crit++; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
449 do_restart = true; |
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 else if (next_crit < next_out) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
452 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
453 if (i_crit < n_crit) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
454 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
455 set_stop_time (next_crit); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
456 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
|
457 save_output = false; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
458 i_crit++; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
459 do_restart = true; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
460 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
461 else |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
462 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
463 clear_stop_time (); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
464 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
|
465 save_output = true; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
466 i_out++; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
467 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
468 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
469 else |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
470 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
471 set_stop_time (next_crit); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
472 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
|
473 save_output = true; |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
474 i_out++; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
475 } |
3 | 476 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
477 ColumnVector x_next = do_integrate (t_out); |
3 | 478 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
479 if (integration_error) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
480 return retval; |
258 | 481 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
482 if (save_output) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
483 { |
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
|
484 for (F77_INT i = 0; i < n; i++) |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
485 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
|
486 } |
3 | 487 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
488 if (do_restart) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
489 force_restart (); |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
490 } |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
491 } |
3 | 492 else |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
493 { |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
494 retval = do_integrate (tout); |
258 | 495 |
10314
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
496 if (integration_error) |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
497 return retval; |
07ebe522dac2
untabify liboctave C++ sources
John W. Eaton <jwe@octave.org>
parents:
10158
diff
changeset
|
498 } |
3 | 499 } |
500 | |
501 return retval; | |
502 } |