Mercurial > hg > octave-nkf
diff liboctave/DASSL.cc @ 10314:07ebe522dac2
untabify liboctave C++ sources
author | John W. Eaton <jwe@octave.org> |
---|---|
date | Thu, 11 Feb 2010 12:23:32 -0500 |
parents | 4c0cdbe0acca |
children | 12884915a8e4 |
line wrap: on
line diff
--- a/liboctave/DASSL.cc +++ b/liboctave/DASSL.cc @@ -36,20 +36,20 @@ #include "quit.h" typedef octave_idx_type (*dassl_fcn_ptr) (const double&, const double*, const double*, - double*, octave_idx_type&, double*, octave_idx_type*); + double*, octave_idx_type&, double*, octave_idx_type*); typedef octave_idx_type (*dassl_jac_ptr) (const double&, const double*, const double*, - double*, const double&, double*, octave_idx_type*); + double*, const double&, double*, octave_idx_type*); extern "C" { F77_RET_T F77_FUNC (ddassl, DDASSL) (dassl_fcn_ptr, const octave_idx_type&, double&, - double*, double*, double&, const octave_idx_type*, - const double*, const double*, octave_idx_type&, - double*, const octave_idx_type&, octave_idx_type*, const octave_idx_type&, - const double*, const octave_idx_type*, - dassl_jac_ptr); + double*, double*, double&, const octave_idx_type*, + const double*, const double*, octave_idx_type&, + double*, const octave_idx_type&, octave_idx_type*, const octave_idx_type&, + const double*, const octave_idx_type*, + dassl_jac_ptr); } static DAEFunc::DAERHSFunc user_fun; @@ -59,7 +59,7 @@ static octave_idx_type ddassl_f (const double& time, const double *state, const double *deriv, - double *delta, octave_idx_type& ires, double *, octave_idx_type *) + double *delta, octave_idx_type& ires, double *, octave_idx_type *) { BEGIN_INTERRUPT_WITH_EXCEPTIONS; @@ -80,12 +80,12 @@ if (ires >= 0) { if (tmp_delta.length () == 0) - ires = -2; + ires = -2; else - { - for (octave_idx_type i = 0; i < nn; i++) - delta [i] = tmp_delta.elem (i); - } + { + for (octave_idx_type i = 0; i < nn; i++) + delta [i] = tmp_delta.elem (i); + } } END_INTERRUPT_WITH_EXCEPTIONS; @@ -95,7 +95,7 @@ static octave_idx_type ddassl_j (const double& time, const double *state, const double *deriv, - double *pd, const double& cj, double *, octave_idx_type *) + double *pd, const double& cj, double *, octave_idx_type *) { BEGIN_INTERRUPT_WITH_EXCEPTIONS; @@ -135,7 +135,7 @@ info.resize (15); for (octave_idx_type i = 0; i < 15; i++) - info(i) = 0; + info(i) = 0; pinfo = info.fortran_vec (); @@ -152,12 +152,12 @@ info(0) = 0; if (stop_time_set) - { - rwork(0) = stop_time; - info(3) = 1; - } + { + rwork(0) = stop_time; + info(3) = 1; + } else - info(3) = 0; + info(3) = 0; px = x.fortran_vec (); pxdot = xdot.fortran_vec (); @@ -173,28 +173,28 @@ user_jac = DAEFunc::jacobian_function (); if (user_fun) - { - octave_idx_type ires = 0; + { + octave_idx_type ires = 0; - ColumnVector res = (*user_fun) (x, xdot, t, ires); + ColumnVector res = (*user_fun) (x, xdot, t, ires); - if (res.length () != x.length ()) - { - (*current_liboctave_error_handler) - ("dassl: inconsistent sizes for state and residual vectors"); + if (res.length () != x.length ()) + { + (*current_liboctave_error_handler) + ("dassl: inconsistent sizes for state and residual vectors"); - integration_error = true; - return retval; - } - } + integration_error = true; + return retval; + } + } else - { - (*current_liboctave_error_handler) - ("dassl: no user supplied RHS subroutine!"); + { + (*current_liboctave_error_handler) + ("dassl: no user supplied RHS subroutine!"); - integration_error = true; - return retval; - } + integration_error = true; + return retval; + } info(4) = user_jac ? 1 : 0; @@ -204,46 +204,46 @@ double hmax = maximum_step_size (); if (hmax >= 0.0) - { - rwork(1) = hmax; - info(6) = 1; - } + { + rwork(1) = hmax; + info(6) = 1; + } else - info(6) = 0; + info(6) = 0; double h0 = initial_step_size (); if (h0 >= 0.0) - { - rwork(2) = h0; - info(7) = 1; - } + { + rwork(2) = h0; + info(7) = 1; + } else - info(7) = 0; + info(7) = 0; if (step_limit () >= 0) - { - info(11) = 1; - iwork(20) = step_limit (); - } + { + info(11) = 1; + iwork(20) = step_limit (); + } else - info(11) = 0; + info(11) = 0; octave_idx_type maxord = maximum_order (); if (maxord >= 0) - { - if (maxord > 0 && maxord < 6) - { - info(8) = 1; - iwork(2) = maxord; - } - else - { - (*current_liboctave_error_handler) - ("dassl: invalid value for maximum order"); - integration_error = true; - return retval; - } - } + { + if (maxord > 0 && maxord < 6) + { + info(8) = 1; + iwork(2) = maxord; + } + else + { + (*current_liboctave_error_handler) + ("dassl: invalid value for maximum order"); + integration_error = true; + return retval; + } + } octave_idx_type enc = enforce_nonnegativity_constraints (); info(9) = enc ? 1 : 0; @@ -258,21 +258,21 @@ octave_idx_type rel_tol_len = rel_tol.length (); if (abs_tol_len == 1 && rel_tol_len == 1) - { - info(1) = 0; - } + { + info(1) = 0; + } else if (abs_tol_len == n && rel_tol_len == n) - { - info(1) = 1; - } + { + info(1) = 1; + } else - { - (*current_liboctave_error_handler) - ("dassl: inconsistent sizes for tolerance arrays"); + { + (*current_liboctave_error_handler) + ("dassl: inconsistent sizes for tolerance arrays"); - integration_error = true; - return retval; - } + integration_error = true; + return retval; + } pabs_tol = abs_tol.fortran_vec (); prel_tol = rel_tol.fortran_vec (); @@ -284,18 +284,18 @@ static octave_idx_type *idummy = 0; F77_XFCN (ddassl, DDASSL, (ddassl_f, nn, t, px, pxdot, tout, pinfo, - prel_tol, pabs_tol, istate, prwork, lrw, - piwork, liw, dummy, idummy, ddassl_j)); + prel_tol, pabs_tol, istate, prwork, lrw, + piwork, liw, dummy, idummy, ddassl_j)); switch (istate) { case 1: // A step was successfully taken in intermediate-output - // mode. The code has not yet reached TOUT. + // mode. The code has not yet reached TOUT. case 2: // The integration to TSTOP was successfully completed - // (T=TSTOP) by stepping exactly to TSTOP. + // (T=TSTOP) by stepping exactly to TSTOP. case 3: // The integration to TOUT was successfully completed - // (T=TOUT) by stepping past TOUT. Y(*) is obtained by - // interpolation. YPRIME(*) is obtained by interpolation. + // (T=TOUT) by stepping past TOUT. Y(*) is obtained by + // interpolation. YPRIME(*) is obtained by interpolation. retval = x; t = tout; break; @@ -303,33 +303,33 @@ case -1: // A large amount of work has been expended. (~500 steps). case -2: // The error tolerances are too stringent. case -3: // The local error test cannot be satisfied because you - // specified a zero component in ATOL and the - // corresponding computed solution component is zero. - // Thus, a pure relative error test is impossible for - // this component. + // specified a zero component in ATOL and the + // corresponding computed solution component is zero. + // Thus, a pure relative error test is impossible for + // this component. case -6: // DDASSL had repeated error test failures on the last - // attempted step. + // attempted step. case -7: // The corrector could not converge. case -8: // The matrix of partial derivatives is singular. case -9: // The corrector could not converge. There were repeated - // error test failures in this step. + // error test failures in this step. case -10: // The corrector could not converge because IRES was - // equal to minus one. + // equal to minus one. case -11: // IRES equal to -2 was encountered and control is being - // returned to the calling program. + // returned to the calling program. case -12: // DDASSL failed to compute the initial YPRIME. case -33: // The code has encountered trouble from which it cannot - // recover. A message is printed explaining the trouble - // and control is returned to the calling program. For - // example, this occurs when invalid input is detected. + // recover. A message is printed explaining the trouble + // and control is returned to the calling program. For + // example, this occurs when invalid input is detected. integration_error = true; break; default: integration_error = true; (*current_liboctave_error_handler) - ("unrecognized value of istate (= %d) returned from ddassl", - istate); + ("unrecognized value of istate (= %d) returned from ddassl", + istate); break; } @@ -357,24 +357,24 @@ xdot_out.resize (n_out, n); for (octave_idx_type i = 0; i < n; i++) - { - retval.elem (0, i) = x.elem (i); - xdot_out.elem (0, i) = xdot.elem (i); - } + { + retval.elem (0, i) = x.elem (i); + xdot_out.elem (0, i) = xdot.elem (i); + } for (octave_idx_type j = 1; j < n_out; j++) - { - ColumnVector x_next = do_integrate (tout.elem (j)); + { + ColumnVector x_next = do_integrate (tout.elem (j)); - if (integration_error) - return retval; + if (integration_error) + return retval; - for (octave_idx_type i = 0; i < n; i++) - { - retval.elem (j, i) = x_next.elem (i); - xdot_out.elem (j, i) = xdot.elem (i); - } - } + for (octave_idx_type i = 0; i < n; i++) + { + retval.elem (j, i) = x_next.elem (i); + xdot_out.elem (j, i) = xdot.elem (i); + } + } } return retval; @@ -389,7 +389,7 @@ Matrix DASSL::integrate (const ColumnVector& tout, Matrix& xdot_out, - const ColumnVector& tcrit) + const ColumnVector& tcrit) { Matrix retval; @@ -402,90 +402,90 @@ xdot_out.resize (n_out, n); for (octave_idx_type i = 0; i < n; i++) - { - retval.elem (0, i) = x.elem (i); - xdot_out.elem (0, i) = xdot.elem (i); - } + { + retval.elem (0, i) = x.elem (i); + xdot_out.elem (0, i) = xdot.elem (i); + } octave_idx_type n_crit = tcrit.capacity (); if (n_crit > 0) - { - octave_idx_type i_crit = 0; - octave_idx_type i_out = 1; - double next_crit = tcrit.elem (0); - double next_out; - while (i_out < n_out) - { - bool do_restart = false; + { + octave_idx_type i_crit = 0; + octave_idx_type i_out = 1; + double next_crit = tcrit.elem (0); + double next_out; + while (i_out < n_out) + { + bool do_restart = false; - next_out = tout.elem (i_out); - if (i_crit < n_crit) - next_crit = tcrit.elem (i_crit); + next_out = tout.elem (i_out); + if (i_crit < n_crit) + next_crit = tcrit.elem (i_crit); - bool save_output; - double t_out; + bool save_output; + double t_out; - if (next_crit == next_out) - { - set_stop_time (next_crit); - t_out = next_out; - save_output = true; - i_out++; - i_crit++; - do_restart = true; - } - else if (next_crit < next_out) - { - if (i_crit < n_crit) - { - set_stop_time (next_crit); - t_out = next_crit; - save_output = false; - i_crit++; - do_restart = true; - } - else - { - clear_stop_time (); - t_out = next_out; - save_output = true; - i_out++; - } - } - else - { - set_stop_time (next_crit); - t_out = next_out; - save_output = true; - i_out++; - } + if (next_crit == next_out) + { + set_stop_time (next_crit); + t_out = next_out; + save_output = true; + i_out++; + i_crit++; + do_restart = true; + } + else if (next_crit < next_out) + { + if (i_crit < n_crit) + { + set_stop_time (next_crit); + t_out = next_crit; + save_output = false; + i_crit++; + do_restart = true; + } + else + { + clear_stop_time (); + t_out = next_out; + save_output = true; + i_out++; + } + } + else + { + set_stop_time (next_crit); + t_out = next_out; + save_output = true; + i_out++; + } - ColumnVector x_next = do_integrate (t_out); + ColumnVector x_next = do_integrate (t_out); - if (integration_error) - return retval; + if (integration_error) + return retval; - if (save_output) - { - for (octave_idx_type i = 0; i < n; i++) - { - retval.elem (i_out-1, i) = x_next.elem (i); - xdot_out.elem (i_out-1, i) = xdot.elem (i); - } - } + if (save_output) + { + for (octave_idx_type i = 0; i < n; i++) + { + retval.elem (i_out-1, i) = x_next.elem (i); + xdot_out.elem (i_out-1, i) = xdot.elem (i); + } + } - if (do_restart) - force_restart (); - } - } + if (do_restart) + force_restart (); + } + } else - { - retval = integrate (tout, xdot_out); + { + retval = integrate (tout, xdot_out); - if (integration_error) - return retval; - } + if (integration_error) + return retval; + } } return retval; @@ -516,7 +516,7 @@ case -1: retval = std::string ("a large amount of work has been expended (t =") - + t_curr + ")"; + + t_curr + ")"; break; case -2: @@ -525,38 +525,38 @@ case -3: retval = std::string ("error weight became zero during problem. (t = ") - + t_curr - + "; solution component i vanished, and atol or atol(i) == 0)"; + + t_curr + + "; solution component i vanished, and atol or atol(i) == 0)"; break; case -6: retval = std::string ("repeated error test failures on the last attempted step (t = ") - + t_curr + ")"; + + t_curr + ")"; break; case -7: retval = std::string ("the corrector could not converge (t = ") - + t_curr + ")"; + + t_curr + ")"; break; case -8: retval = std::string ("the matrix of partial derivatives is singular (t = ") - + t_curr + ")"; + + t_curr + ")"; break; case -9: retval = std::string ("the corrector could not converge (t = ") - + t_curr + "; repeated test failures)"; + + t_curr + "; repeated test failures)"; break; case -10: retval = std::string ("corrector could not converge because IRES was -1 (t = ") - + t_curr + ")"; + + t_curr + ")"; break; case -11: retval = std::string ("return requested in user-supplied function (t = ") - + t_curr + ")"; + + t_curr + ")"; break; case -12: