23 #if defined (HAVE_CONFIG_H) 36 const double*,
double*,
F77_INT&,
40 const double*,
double*,
const double&,
63 BEGIN_INTERRUPT_WITH_EXCEPTIONS;
73 tmp_deriv.
elem (
i) = deriv[
i];
79 tmp_delta =
user_fun (tmp_state, tmp_deriv, time, tmp_ires);
81 ires = octave::to_f77_int (tmp_ires);
90 delta[
i] = tmp_delta.
elem (
i);
94 END_INTERRUPT_WITH_EXCEPTIONS;
101 double *pd,
const double& cj,
double *,
F77_INT *)
103 BEGIN_INTERRUPT_WITH_EXCEPTIONS;
112 tmp_deriv.
elem (
i) = deriv[
i];
120 pd[
nn * j +
i] = tmp_pd.
elem (
i, j);
122 END_INTERRUPT_WITH_EXCEPTIONS;
146 lrw = 40 + 9*n + n*n;
178 (*current_liboctave_error_handler)
179 (
"dassl: inconsistent sizes for state and residual vectors");
187 (*current_liboctave_error_handler)
188 (
"dassl: no user supplied RHS subroutine!");
200 double hmax = maximum_step_size ();
209 double h0 = initial_step_size ();
218 F77_INT sl = octave::to_f77_int (step_limit ());
228 F77_INT maxord = octave::to_f77_int (maximum_order ());
231 if (maxord > 0 && maxord < 6)
238 (*current_liboctave_error_handler)
239 (
"dassl: invalid value for maximum order: %d", maxord);
245 F77_INT enc = octave::to_f77_int (enforce_nonnegativity_constraints ());
246 info(9) = (enc ? 1 : 0);
248 F77_INT ccic = octave::to_f77_int (compute_consistent_initial_condition ());
249 info(10) = (ccic ? 1 : 0);
251 abs_tol = absolute_tolerance ();
252 rel_tol = relative_tolerance ();
257 if (abs_tol_len == 1 && rel_tol_len == 1)
261 else if (abs_tol_len == n && rel_tol_len == n)
267 (*current_liboctave_error_handler)
268 (
"dassl: inconsistent sizes for tolerance arrays");
274 DASSL_options::reset =
false;
288 double *dummy =
nullptr;
294 prel_tol, pabs_tol, tmp_istate, prwork,
lrw,
339 (*current_liboctave_error_handler)
340 (
"unrecognized value of istate (= %d) returned from ddassl",
363 if (n_out > 0 && n > 0)
366 xdot_out.
resize (n_out, n);
408 if (n_out > 0 && n > 0)
411 xdot_out.
resize (n_out, n);
425 double next_crit = tcrit.
elem (0);
427 while (i_out < n_out)
429 bool do_restart =
false;
431 next_out = tout.
elem (i_out);
433 next_crit = tcrit.
elem (i_crit);
438 if (next_crit == next_out)
447 else if (next_crit < next_out)
508 std::ostringstream buf;
515 retval =
"a step was successfully taken in intermediate-output mode.";
519 retval =
"integration completed by stepping exactly to TOUT";
523 retval =
"integration to tout completed by stepping past TOUT";
527 retval =
"a large amount of work has been expended (t =" + t_curr +
')';
531 retval =
"the error tolerances are too stringent";
535 retval =
"error weight became zero during problem. (t = " + t_curr +
536 "; solution component i vanished, and atol or atol(i) == 0)";
540 retval =
"repeated error test failures on the last attempted step (t = " 545 retval =
"the corrector could not converge (t = " + t_curr +
')';
549 retval =
"the matrix of partial derivatives is singular (t = " + t_curr +
554 retval =
"the corrector could not converge (t = " + t_curr +
555 "; repeated test failures)";
559 retval =
"corrector could not converge because IRES was -1 (t = " 564 retval =
"return requested in user-supplied function (t = " + t_curr +
569 retval =
"failed to compute consistent initial conditions";
573 retval =
"unrecoverable error (see printed message)";
577 retval =
"unknown error state";
ColumnVector do_integrate(double t)
static DAEFunc::DAEJacFunc user_jac
subroutine DDASSL(RES, NEQ, T, Y, YPRIME, TOUT, INFO, RTOL, ATOL, IDID, RWORK, LRW, IWORK, LIW, RPAR, IPAR, JAC)
void resize(octave_idx_type nr, octave_idx_type nc, double rfv=0)
static DAEFunc::DAERHSFunc user_fun
Matrix integrate(const ColumnVector &tout, Matrix &xdot_out)
Array< octave_f77_int_type > info
const T * fortran_vec(void) const
std::string error_message(void) const
F77_INT(* dassl_jac_ptr)(const double &, const double *, const double *, double *, const double &, double *, F77_INT *)
T & elem(octave_idx_type n)
#define F77_XFCN(f, F, args)
F77_RET_T(F77_CONST_CHAR_ARG_DECL, F77_CONST_CHAR_ARG_DECL, F77_CONST_CHAR_ARG_DECL, const F77_INT &, const F77_INT &, const F77_INT &, F77_INT &, F77_INT &, F77_DBLE *, const F77_INT &, F77_DBLE *, const F77_INT &, F77_DBLE *, F77_DBLE *, F77_DBLE *, const F77_INT &, F77_DBLE *, const F77_INT &, F77_DBLE *, const F77_INT &, F77_DBLE *, F77_INT *, F77_INT &F77_CHAR_ARG_LEN_DECL F77_CHAR_ARG_LEN_DECL F77_CHAR_ARG_LEN_DECL)
octave_value resize(const dim_vector &dv, bool fill=false) const
static F77_INT ddassl_f(const double &time, const double *state, const double *deriv, double *delta, F77_INT &ires, double *, F77_INT *)
octave_idx_type size(void) const
Matrix(* DAEJacFunc)(const ColumnVector &x, const ColumnVector &xdot, double t, double cj)
virtual void force_restart(void)
void set_stop_time(double tt)
DAEJacFunc jacobian_function(void) const
void resize(const dim_vector &dv, const T &rfv)
Resizing (with fill).
void clear_stop_time(void)
static uint32_t state[624]
octave_f77_int_type F77_INT
DAERHSFunc function(void) const
ColumnVector(* DAERHSFunc)(const ColumnVector &x, const ColumnVector &xdot, double t, octave_idx_type &ires)
octave_idx_type numel(void) const
Number of elements in the array.
F77_INT(* dassl_fcn_ptr)(const double &, const double *, const double *, double *, F77_INT &, double *, F77_INT *)
Vector representing the dimensions (size) of an Array.
If this string is the system will ring the terminal sometimes it is useful to be able to print the original representation of the string
F77_RET_T F77_FUNC(ddassl, DDASSL)(dassl_fcn_ptr
Array< octave_f77_int_type > iwork
static F77_INT ddassl_j(const double &time, const double *state, const double *deriv, double *pd, const double &cj, double *, F77_INT *)