Mercurial > octave
diff liboctave/LSODE.cc @ 258:1c9a678906fb
[project @ 1993-12-14 06:24:04 by jwe]
author | jwe |
---|---|
date | Tue, 14 Dec 1993 06:24:21 +0000 |
parents | 780cbbc57b7c |
children | c23f50e61c58 |
line wrap: on
line diff
--- a/liboctave/LSODE.cc Wed Dec 08 23:36:43 1993 +0000 +++ b/liboctave/LSODE.cc Tue Dec 14 06:24:21 1993 +0000 @@ -31,9 +31,12 @@ extern "C" { - int F77_FCN (lsode) (int (*)(), int *, double *, double *, double *, + int F77_FCN (lsode) (int (*)(int*, double*, double*, double*, int*), + int *, double *, double *, double *, int *, double *, double *, int *, int *, int *, - double *, int *, int *, int *, int (*)(), int *); + double *, int *, int *, int *, + int (*)(int*, double*, double*, int*, int*, + double*, int*), int *); } static ColumnVector (*user_fun) (ColumnVector&, double); @@ -51,6 +54,7 @@ stop_time_set = 0; stop_time = 0.0; + integration_error = 0; restart = 1; istate = 1; @@ -84,6 +88,7 @@ stop_time_set = 0; stop_time = 0.0; + integration_error = 0; restart = 1; istate = 1; @@ -118,6 +123,7 @@ stop_time_set = 0; stop_time = 0.0; + integration_error = 0; restart = 1; istate = 1; @@ -147,7 +153,7 @@ } int -lsode_f (int *neq, double *time, double *state, double *deriv) +lsode_f (int *neq, double *time, double *state, double *deriv, int *ierr) { int nn = *neq; ColumnVector tmp_deriv (nn); @@ -159,8 +165,13 @@ */ tmp_deriv = (*user_fun) (*tmp_x, *time); - for (int i = 0; i < nn; i++) - deriv [i] = tmp_deriv.elem (i); + if (tmp_deriv.length () == 0) + *ierr = -1; + else + { + for (int i = 0; i < nn; i++) + deriv [i] = tmp_deriv.elem (i); + } return 0; } @@ -194,6 +205,8 @@ else method_flag = 21; + integration_error = 0; + double *xp = x.fortran_vec (); // NOTE: this won't work if LSODE passes copies of the state vector. @@ -236,17 +249,16 @@ switch (istate) { + case -13: // Return requested in user-supplied function. case -6: // error weight became zero during problem. (solution // component i vanished, and atol or atol(i) = 0.) - break; case -5: // repeated convergence failures (perhaps bad jacobian // supplied or wrong choice of mf or tolerances). - break; case -4: // repeated error test failures (check all inputs). - break; case -3: // illegal input detected (see printed message). - break; case -2: // excess accuracy requested (tolerances too small). + integration_error = 1; + return ColumnVector (); break; case -1: // excess work done on this call (perhaps wrong mf). working_too_hard++; @@ -255,6 +267,7 @@ (*current_liboctave_error_handler) ("giving up after more than %d steps attempted in lsode", iwork[5] * 20); + integration_error = 1; return ColumnVector (); } else @@ -317,6 +330,10 @@ for (int j = 1; j < n_out; j++) { ColumnVector x_next = integrate (tout.elem (j)); + + if (integration_error) + return retval; + for (i = 0; i < n; i++) retval.elem (j, i) = x_next.elem (i); } @@ -394,6 +411,9 @@ ColumnVector x_next = integrate (t_out); + if (integration_error) + return retval; + if (save_output) { for (i = 0; i < n; i++) @@ -405,7 +425,12 @@ } } else - retval = integrate (tout); + { + retval = integrate (tout); + + if (integration_error) + return retval; + } } return retval;