octave-nkf: liboctave/DASRT.cc annotate

annotate liboctave/DASRT.cc @ 14101:ff5588774680 stable ss-3-5-91

snapshot 3.5.91 * configure.ac (AC_INIT): Version is now 3.5.91. (OCTAVE_RELEASE_DATE): Now 2011-12-23.

author	John W. Eaton <jwe@octave.org>
date	Fri, 23 Dec 2011 19:24:41 -0500
parents	12df7854fa7c
children	72c96de7a403

rev	line source
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	1 /*
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	2
11523 fd0a3ac60b0e update copyright notices John W. Eaton <jwe@octave.org> parents: 11502 diff changeset	3 Copyright (C) 2002-2011 John W. Eaton
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	4
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	5 This file is part of Octave.
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	6
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	7 Octave is free software; you can redistribute it and/or modify it
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	8 under the terms of the GNU General Public License as published by the
7016 93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5775 diff changeset	9 Free Software Foundation; either version 3 of the License, or (at your
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5775 diff changeset	10 option) any later version.
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	11
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	12 Octave is distributed in the hope that it will be useful, but WITHOUT
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	15 for more details.
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	16
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	17 You should have received a copy of the GNU General Public License
7016 93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5775 diff changeset	18 along with Octave; see the file COPYING. If not, see
93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5775 diff changeset	19 <http://www.gnu.org/licenses/>.
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	20
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	21 */
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	22
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	23 #ifdef HAVE_CONFIG_H
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	24 #include <config.h>
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	25 #endif
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	26
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	27 #include <cfloat>
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	28
5765 7ba9ad1fec11 [project @ 2006-04-17 05:05:15 by jwe] jwe parents: 5307 diff changeset	29 #include <sstream>
7ba9ad1fec11 [project @ 2006-04-17 05:05:15 by jwe] jwe parents: 5307 diff changeset	30
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	31 #include "DASRT.h"
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	32 #include "f77-fcn.h"
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	33 #include "lo-error.h"
7231 2eb392d058bb [project @ 2007-11-30 18:53:29 by jwe] jwe parents: 7017 diff changeset	34 #include "lo-math.h"
4180 84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	35 #include "quit.h"
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	36
11495 8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	37 typedef octave_idx_type (dasrt_fcn_ptr) (const double&, const double,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	38 const double, double,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	39 octave_idx_type&, double*,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	40 octave_idx_type*);
3991 48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	41
11495 8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	42 typedef octave_idx_type (dasrt_jac_ptr) (const double&, const double,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	43 const double, double,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	44 const double&, double*,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	45 octave_idx_type*);
3991 48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	46
11495 8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	47 typedef octave_idx_type (*dasrt_constr_ptr) (const octave_idx_type&,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	48 const double&, const double*,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	49 const octave_idx_type&,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	50 double, double,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	51 octave_idx_type*);
3991 48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	52
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	53 extern "C"
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4548 diff changeset	54 {
6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4548 diff changeset	55 F77_RET_T
11495 8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	56 F77_FUNC (ddasrt, DDASRT) (dasrt_fcn_ptr, const octave_idx_type&,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	57 double&, double, double, const double&,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	58 octave_idx_type, const double,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	59 const double, octave_idx_type&, double,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	60 const octave_idx_type&, octave_idx_type*,
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11574 diff changeset	61 const octave_idx_type&, double*,
11495 8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	62 octave_idx_type*, dasrt_jac_ptr,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	63 dasrt_constr_ptr, const octave_idx_type&,
8a5e980da6aa style fixes John W. Eaton <jwe@octave.org> parents: 10350 diff changeset	64 octave_idx_type*);
4552 6f3382e08a52 [project @ 2003-10-27 20:38:02 by jwe] jwe parents: 4548 diff changeset	65 }
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	66
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	67 static DAEFunc::DAERHSFunc user_fsub;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	68 static DAEFunc::DAEJacFunc user_jsub;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	69 static DAERTFunc::DAERTConstrFunc user_csub;
3993 f23bc69132cc [project @ 2002-07-16 20:18:56 by jwe] jwe parents: 3992 diff changeset	70
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	71 static octave_idx_type nn;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	72
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	73 static octave_idx_type
3991 48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	74 ddasrt_f (const double& t, const double state, const double deriv,
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	75 double delta, octave_idx_type& ires, double , octave_idx_type *)
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	76 {
4180 84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	77 BEGIN_INTERRUPT_WITH_EXCEPTIONS;
84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	78
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	79 ColumnVector tmp_state (nn);
4133 402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	80 ColumnVector tmp_deriv (nn);
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	81
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	82 for (octave_idx_type i = 0; i < nn; i++)
4133 402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	83 {
402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	84 tmp_state(i) = state[i];
402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	85 tmp_deriv(i) = deriv[i];
402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	86 }
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	87
3994 a41827ec5677 [project @ 2002-07-16 23:57:09 by jwe] jwe parents: 3993 diff changeset	88 ColumnVector tmp_fval = (*user_fsub) (tmp_state, tmp_deriv, t, ires);
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	89
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	90 if (tmp_fval.length () == 0)
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	91 ires = -2;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	92 else
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	93 {
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	94 for (octave_idx_type i = 0; i < nn; i++)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	95 delta[i] = tmp_fval(i);
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	96 }
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	97
4180 84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	98 END_INTERRUPT_WITH_EXCEPTIONS;
84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	99
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	100 return 0;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	101 }
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	102
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	103 octave_idx_type
3991 48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	104 ddasrt_j (const double& time, const double state, const double deriv,
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	105 double pd, const double& cj, double , octave_idx_type *)
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	106 {
4180 84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	107 BEGIN_INTERRUPT_WITH_EXCEPTIONS;
84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	108
5775 ace8d8d26933 [project @ 2006-04-24 19:13:06 by jwe] jwe parents: 5765 diff changeset	109 // FIXME -- would be nice to avoid copying the data.
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	110
3991 48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	111 ColumnVector tmp_state (nn);
48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	112 ColumnVector tmp_deriv (nn);
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	113
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	114 for (octave_idx_type i = 0; i < nn; i++)
3991 48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	115 {
48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	116 tmp_deriv.elem (i) = deriv [i];
48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	117 tmp_state.elem (i) = state [i];
48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	118 }
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	119
3994 a41827ec5677 [project @ 2002-07-16 23:57:09 by jwe] jwe parents: 3993 diff changeset	120 Matrix tmp_pd = (*user_jsub) (tmp_state, tmp_deriv, time, cj);
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	121
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	122 for (octave_idx_type j = 0; j < nn; j++)
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	123 for (octave_idx_type i = 0; i < nn; i++)
3991 48d2bc4a3729 [project @ 2002-07-16 17:46:50 by jwe] jwe parents: 3990 diff changeset	124 pd [nn * j + i] = tmp_pd.elem (i, j);
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	125
4180 84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	126 END_INTERRUPT_WITH_EXCEPTIONS;
84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	127
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	128 return 0;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	129 }
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	130
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	131 static octave_idx_type
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	132 ddasrt_g (const octave_idx_type& neq, const double& t, const double *state,
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11574 diff changeset	133 const octave_idx_type& ng, double gout, double , octave_idx_type *)
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	134 {
4180 84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	135 BEGIN_INTERRUPT_WITH_EXCEPTIONS;
84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	136
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	137 octave_idx_type n = neq;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	138
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	139 ColumnVector tmp_state (n);
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	140 for (octave_idx_type i = 0; i < n; i++)
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	141 tmp_state(i) = state[i];
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	142
3994 a41827ec5677 [project @ 2002-07-16 23:57:09 by jwe] jwe parents: 3993 diff changeset	143 ColumnVector tmp_fval = (*user_csub) (tmp_state, t);
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	144
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	145 for (octave_idx_type i = 0; i < ng; i++)
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	146 gout[i] = tmp_fval(i);
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	147
4180 84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	148 END_INTERRUPT_WITH_EXCEPTIONS;
84fe3ca3a246 [project @ 2002-11-15 04:47:01 by jwe] jwe parents: 4133 diff changeset	149
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	150 return 0;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	151 }
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	152
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	153 void
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	154 DASRT::integrate (double tout)
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	155 {
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	156 DASRT_result retval;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	157
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	158 // I suppose this is the safe thing to do. If this is the first
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	159 // call, or if anything about the problem has changed, we should
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	160 // start completely fresh.
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	161
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	162 if (! initialized \|\| restart
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	163 \|\| DAEFunc::reset \|\| DAERTFunc::reset \|\| DASRT_options::reset)
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	164 {
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	165 integration_error = false;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	166
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	167 initialized = true;
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	168
11574 a83bad07f7e3 attempt better backward compatibility for Array resize functions John W. Eaton <jwe@octave.org> parents: 11523 diff changeset	169 info.resize (dim_vector (15, 1));
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	170
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	171 for (octave_idx_type i = 0; i < 15; i++)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	172 info(i) = 0;
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	173
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	174 octave_idx_type n = size ();
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	175
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	176 nn = n;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	177
4133 402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	178 // DAERTFunc
402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	179
402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	180 user_csub = DAERTFunc::constraint_function ();
402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	181
402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	182 if (user_csub)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	183 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	184 ColumnVector tmp = (*user_csub) (x, t);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	185 ng = tmp.length ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	186 }
4133 402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	187 else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	188 ng = 0;
4133 402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	189
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	190 octave_idx_type maxord = maximum_order ();
4133 402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	191 if (maxord >= 0)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	192 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	193 if (maxord > 0 && maxord < 6)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	194 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	195 info(8) = 1;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	196 iwork(2) = maxord;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	197 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	198 else
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	199 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	200 (*current_liboctave_error_handler)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	201 ("dassl: invalid value for maximum order");
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	202 integration_error = true;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	203 return;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	204 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	205 }
4133 402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	206
4428 067160691cc9 [project @ 2003-06-16 19:09:11 by jwe] jwe parents: 4192 diff changeset	207 liw = 21 + n;
4133 402d7b86a0a2 [project @ 2002-10-29 23:57:34 by jwe] jwe parents: 4066 diff changeset	208 lrw = 50 + 9n + nn + 3*ng;
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	209
11574 a83bad07f7e3 attempt better backward compatibility for Array resize functions John W. Eaton <jwe@octave.org> parents: 11523 diff changeset	210 iwork.resize (dim_vector (liw, 1));
a83bad07f7e3 attempt better backward compatibility for Array resize functions John W. Eaton <jwe@octave.org> parents: 11523 diff changeset	211 rwork.resize (dim_vector (lrw, 1));
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	212
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	213 info(0) = 0;
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	214
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	215 if (stop_time_set)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	216 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	217 info(3) = 1;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	218 rwork(0) = stop_time;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	219 }
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	220 else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	221 info(3) = 0;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	222
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	223 restart = false;
3994 a41827ec5677 [project @ 2002-07-16 23:57:09 by jwe] jwe parents: 3993 diff changeset	224
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	225 // DAEFunc
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	226
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	227 user_fsub = DAEFunc::function ();
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	228 user_jsub = DAEFunc::jacobian_function ();
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	229
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	230 if (user_fsub)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	231 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	232 octave_idx_type ires = 0;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	233
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	234 ColumnVector fval = (*user_fsub) (x, xdot, t, ires);
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	235
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	236 if (fval.length () != x.length ())
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	237 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	238 (*current_liboctave_error_handler)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	239 ("dasrt: inconsistent sizes for state and residual vectors");
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	240
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	241 integration_error = true;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	242 return;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	243 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	244 }
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	245 else
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 (*current_liboctave_error_handler)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	248 ("dasrt: no user supplied RHS subroutine!");
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	249
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	250 integration_error = true;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	251 return;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	252 }
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	253
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	254 info(4) = user_jsub ? 1 : 0;
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	255
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	256 DAEFunc::reset = false;
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	257
11574 a83bad07f7e3 attempt better backward compatibility for Array resize functions John W. Eaton <jwe@octave.org> parents: 11523 diff changeset	258 jroot.resize (dim_vector (ng, 1), 1);
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	259
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	260 DAERTFunc::reset = false;
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	261
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	262 // DASRT_options
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	263
4050 6481f41a79f3 [project @ 2002-08-17 02:18:18 by jwe] jwe parents: 4049 diff changeset	264 double mss = maximum_step_size ();
6481f41a79f3 [project @ 2002-08-17 02:18:18 by jwe] jwe parents: 4049 diff changeset	265 if (mss >= 0.0)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	266 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	267 rwork(1) = mss;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	268 info(6) = 1;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	269 }
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	270 else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	271 info(6) = 0;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	272
4050 6481f41a79f3 [project @ 2002-08-17 02:18:18 by jwe] jwe parents: 4049 diff changeset	273 double iss = initial_step_size ();
6481f41a79f3 [project @ 2002-08-17 02:18:18 by jwe] jwe parents: 4049 diff changeset	274 if (iss >= 0.0)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	275 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	276 rwork(2) = iss;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	277 info(7) = 1;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	278 }
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	279 else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	280 info(7) = 0;
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	281
a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	282 if (step_limit () >= 0)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	283 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	284 info(11) = 1;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	285 iwork(20) = step_limit ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	286 }
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	287 else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	288 info(11) = 0;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	289
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	290 abs_tol = absolute_tolerance ();
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	291 rel_tol = relative_tolerance ();
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	292
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	293 octave_idx_type abs_tol_len = abs_tol.length ();
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	294 octave_idx_type rel_tol_len = rel_tol.length ();
3998 f6df65db67f9 [project @ 2002-07-24 18:10:39 by jwe] jwe parents: 3997 diff changeset	295
f6df65db67f9 [project @ 2002-07-24 18:10:39 by jwe] jwe parents: 3997 diff changeset	296 if (abs_tol_len == 1 && rel_tol_len == 1)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	297 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	298 info.elem (1) = 0;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	299 }
3998 f6df65db67f9 [project @ 2002-07-24 18:10:39 by jwe] jwe parents: 3997 diff changeset	300 else if (abs_tol_len == n && rel_tol_len == n)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	301 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	302 info.elem (1) = 1;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	303 }
3998 f6df65db67f9 [project @ 2002-07-24 18:10:39 by jwe] jwe parents: 3997 diff changeset	304 else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	305 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	306 (*current_liboctave_error_handler)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	307 ("dasrt: inconsistent sizes for tolerance arrays");
3998 f6df65db67f9 [project @ 2002-07-24 18:10:39 by jwe] jwe parents: 3997 diff changeset	308
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	309 integration_error = true;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	310 return;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	311 }
3998 f6df65db67f9 [project @ 2002-07-24 18:10:39 by jwe] jwe parents: 3997 diff changeset	312
4049 a35a3c5d4740 [project @ 2002-08-16 08:54:31 by jwe] jwe parents: 4043 diff changeset	313 DASRT_options::reset = false;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	314 }
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	315
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	316 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	317 double *pxdot = xdot.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	318
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	319 octave_idx_type *pinfo = info.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	320
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	321 double *prel_tol = rel_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	322 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	323
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	324 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	325 octave_idx_type *piwork = iwork.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	326
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	327 octave_idx_type *pjroot = jroot.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	328
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	329 double *dummy = 0;
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	330 octave_idx_type *idummy = 0;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	331
4575 c7ba572047df [project @ 2003-10-31 06:05:37 by jwe] jwe parents: 4552 diff changeset	332 F77_XFCN (ddasrt, DDASRT, (ddasrt_f, nn, t, px, pxdot, tout, pinfo,
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	333 prel_tol, pabs_tol, istate, prwork, lrw,
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	334 piwork, liw, dummy, idummy, ddasrt_j,
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	335 ddasrt_g, ng, pjroot));
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	336
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	337 switch (istate)
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	338 {
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	339 case 1: // A step was successfully taken in intermediate-output
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	340 // mode. The code has not yet reached TOUT.
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	341 case 2: // The integration to TOUT was successfully completed
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	342 // (T=TOUT) by stepping exactly to TOUT.
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	343 case 3: // The integration to TOUT was successfully completed
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	344 // (T=TOUT) by stepping past TOUT. Y(*) is obtained by
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	345 // interpolation. YPRIME(*) is obtained by interpolation.
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	346 t = tout;
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	347 break;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	348
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	349 case 4: // The integration was successfully completed
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	350 // by finding one or more roots of G at T.
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	351 break;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	352
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	353 case -1: // A large amount of work has been expended.
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	354 case -2: // The error tolerances are too stringent.
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	355 case -3: // The local error test cannot be satisfied because you
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	356 // specified a zero component in ATOL and the
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	357 // corresponding computed solution component is zero.
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	358 // Thus, a pure relative error test is impossible for
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	359 // this component.
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	360 case -6: // DDASRT had repeated error test failures on the last
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	361 // attempted step.
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	362 case -7: // The corrector could not converge.
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	363 case -8: // The matrix of partial derivatives is singular.
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	364 case -9: // The corrector could not converge. There were repeated
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	365 // error test failures in this step.
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	366 case -10: // The corrector could not converge because IRES was
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	367 // equal to minus one.
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	368 case -11: // IRES equal to -2 was encountered and control is being
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	369 // returned to the calling program.
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	370 case -12: // DASSL failed to compute the initial YPRIME.
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	371 case -33: // The code has encountered trouble from which it cannot
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	372 // recover. A message is printed explaining the trouble
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	373 // and control is returned to the calling program. For
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	374 // example, this occurs when invalid input is detected.
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	375 integration_error = true;
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	376 break;
3996 98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	377
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	378 default:
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	379 integration_error = true;
29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	380 (*current_liboctave_error_handler)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	381 ("unrecognized value of istate (= %d) returned from ddasrt",
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	382 istate);
7482 29980c6b8604 don't check f77_exception_encountered John W. Eaton <jwe@octave.org> parents: 7231 diff changeset	383 break;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	384 }
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	385 }
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	386
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	387 DASRT_result
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	388 DASRT::integrate (const ColumnVector& tout)
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	389 {
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	390 DASRT_result retval;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	391
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	392 Matrix x_out;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	393 Matrix xdot_out;
3994 a41827ec5677 [project @ 2002-07-16 23:57:09 by jwe] jwe parents: 3993 diff changeset	394 ColumnVector t_out = tout;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	395
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	396 octave_idx_type n_out = tout.capacity ();
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	397 octave_idx_type n = size ();
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	398
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	399 if (n_out > 0 && n > 0)
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	400 {
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	401 x_out.resize (n_out, n);
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	402 xdot_out.resize (n_out, n);
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	403
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	404 for (octave_idx_type i = 0; i < n; i++)
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	405 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	406 x_out(0,i) = x(i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	407 xdot_out(0,i) = xdot(i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	408 }
3994 a41827ec5677 [project @ 2002-07-16 23:57:09 by jwe] jwe parents: 3993 diff changeset	409
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	410 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	411 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	412 integrate (tout(j));
3994 a41827ec5677 [project @ 2002-07-16 23:57:09 by jwe] jwe parents: 3993 diff changeset	413
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	414 if (integration_error)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	415 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	416 retval = DASRT_result (x_out, xdot_out, t_out);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	417 return retval;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	418 }
3992 53b4eab68976 [project @ 2002-07-16 19:36:52 by jwe] jwe parents: 3991 diff changeset	419
3997 d4091aff6468 [project @ 2002-07-17 18:00:06 by jwe] jwe parents: 3996 diff changeset	420 if (istate == 4)
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	421 t_out(j) = t;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	422 else
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	423 t_out(j) = tout(j);
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	424
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	425 for (octave_idx_type i = 0; i < n; i++)
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 x_out(j,i) = x(i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	428 xdot_out(j,i) = xdot(i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	429 }
3992 53b4eab68976 [project @ 2002-07-16 19:36:52 by jwe] jwe parents: 3991 diff changeset	430
3997 d4091aff6468 [project @ 2002-07-17 18:00:06 by jwe] jwe parents: 3996 diff changeset	431 if (istate == 4)
10314 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 x_out.resize (j+1, n);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	434 xdot_out.resize (j+1, n);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	435 t_out.resize (j+1);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	436 break;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	437 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	438 }
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	439 }
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	440
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	441 retval = DASRT_result (x_out, xdot_out, t_out);
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	442
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	443 return retval;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	444 }
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	445
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	446 DASRT_result
11586 12df7854fa7c strip trailing whitespace from source files John W. Eaton <jwe@octave.org> parents: 11574 diff changeset	447 DASRT::integrate (const ColumnVector& tout, const ColumnVector& tcrit)
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	448 {
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	449 DASRT_result retval;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	450
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	451 Matrix x_out;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	452 Matrix xdot_out;
3994 a41827ec5677 [project @ 2002-07-16 23:57:09 by jwe] jwe parents: 3993 diff changeset	453 ColumnVector t_outs = tout;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	454
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	455 octave_idx_type n_out = tout.capacity ();
23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	456 octave_idx_type n = size ();
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	457
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	458 if (n_out > 0 && n > 0)
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	459 {
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	460 x_out.resize (n_out, n);
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	461 xdot_out.resize (n_out, n);
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	462
5275 23b37da9fd5b [project @ 2005-04-08 16:07:35 by jwe] jwe parents: 5105 diff changeset	463 octave_idx_type n_crit = tcrit.capacity ();
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	464
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	465 if (n_crit > 0)
10314 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 octave_idx_type i_crit = 0;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	468 octave_idx_type i_out = 1;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	469 double next_crit = tcrit(0);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	470 double next_out;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	471 while (i_out < n_out)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	472 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	473 bool do_restart = false;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	474
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	475 next_out = tout(i_out);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	476 if (i_crit < n_crit)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	477 next_crit = tcrit(i_crit);
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	478
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	479 octave_idx_type save_output;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	480 double t_out;
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	481
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	482 if (next_crit == next_out)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	483 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	484 set_stop_time (next_crit);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	485 t_out = next_out;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	486 save_output = 1;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	487 i_out++;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	488 i_crit++;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	489 do_restart = true;
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 else if (next_crit < next_out)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	492 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	493 if (i_crit < n_crit)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	494 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	495 set_stop_time (next_crit);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	496 t_out = next_crit;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	497 save_output = 0;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	498 i_crit++;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	499 do_restart = true;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	500 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	501 else
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	502 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	503 clear_stop_time ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	504 t_out = next_out;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	505 save_output = 1;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	506 i_out++;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	507 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	508 }
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	509 else
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	510 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	511 set_stop_time (next_crit);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	512 t_out = next_out;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	513 save_output = 1;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	514 i_out++;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	515 }
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	516
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	517 integrate (t_out);
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	518
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	519 if (integration_error)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	520 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	521 retval = DASRT_result (x_out, xdot_out, t_outs);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	522 return retval;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	523 }
3992 53b4eab68976 [project @ 2002-07-16 19:36:52 by jwe] jwe parents: 3991 diff changeset	524
3997 d4091aff6468 [project @ 2002-07-17 18:00:06 by jwe] jwe parents: 3996 diff changeset	525 if (istate == 4)
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	526 t_out = t;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	527
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	528 if (save_output)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	529 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	530 for (octave_idx_type i = 0; i < n; i++)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	531 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	532 x_out(i_out-1,i) = x(i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	533 xdot_out(i_out-1,i) = xdot(i);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	534 }
3994 a41827ec5677 [project @ 2002-07-16 23:57:09 by jwe] jwe parents: 3993 diff changeset	535
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	536 t_outs(i_out-1) = t_out;
3994 a41827ec5677 [project @ 2002-07-16 23:57:09 by jwe] jwe parents: 3993 diff changeset	537
3997 d4091aff6468 [project @ 2002-07-17 18:00:06 by jwe] jwe parents: 3996 diff changeset	538 if (istate == 4)
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	539 {
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	540 x_out.resize (i_out, n);
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	541 xdot_out.resize (i_out, n);
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	542 t_outs.resize (i_out);
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	543 i_out = n_out;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	544 }
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	545 }
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	546
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	547 if (do_restart)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	548 force_restart ();
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	549 }
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	550
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	551 retval = DASRT_result (x_out, xdot_out, t_outs);
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	552 }
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	553 else
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	554 {
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	555 retval = integrate (tout);
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	556
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	557 if (integration_error)
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	558 return retval;
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	559 }
3990 46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	560 }
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	561
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	562 return retval;
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	563 }
46388d6a4e44 [project @ 2002-07-16 06:20:39 by jwe] jwe parents: diff changeset	564
3995 ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	565 std::string
ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	566 DASRT::error_message (void) const
ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	567 {
ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	568 std::string retval;
ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	569
5765 7ba9ad1fec11 [project @ 2006-04-17 05:05:15 by jwe] jwe parents: 5307 diff changeset	570 std::ostringstream buf;
7ba9ad1fec11 [project @ 2006-04-17 05:05:15 by jwe] jwe parents: 5307 diff changeset	571 buf << t;
7ba9ad1fec11 [project @ 2006-04-17 05:05:15 by jwe] jwe parents: 5307 diff changeset	572 std::string t_curr = buf.str ();
4043 6fae69a1796e [project @ 2002-08-15 03:36:27 by jwe] jwe parents: 3998 diff changeset	573
3997 d4091aff6468 [project @ 2002-07-17 18:00:06 by jwe] jwe parents: 3996 diff changeset	574 switch (istate)
3995 ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	575 {
3996 98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	576 case 1:
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	577 retval = "a step was successfully taken in intermediate-output mode.";
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	578 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	579
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	580 case 2:
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	581 retval = "integration completed by stepping exactly to TOUT";
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	582 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	583
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	584 case 3:
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	585 retval = "integration to tout completed by stepping past TOUT";
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	586 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	587
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	588 case 4:
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	589 retval = "integration completed by finding one or more roots of G at T";
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	590 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	591
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	592 case -1:
4043 6fae69a1796e [project @ 2002-08-15 03:36:27 by jwe] jwe parents: 3998 diff changeset	593 retval = std::string ("a large amount of work has been expended (t =")
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	594 + t_curr + ")";
3996 98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	595 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	596
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	597 case -2:
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	598 retval = "the error tolerances are too stringent";
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	599 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	600
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	601 case -3:
4043 6fae69a1796e [project @ 2002-08-15 03:36:27 by jwe] jwe parents: 3998 diff changeset	602 retval = std::string ("error weight became zero during problem. (t = ")
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	603 + t_curr
07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	604 + "; solution component i vanished, and atol or atol(i) == 0)";
3996 98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	605 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	606
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	607 case -6:
4043 6fae69a1796e [project @ 2002-08-15 03:36:27 by jwe] jwe parents: 3998 diff changeset	608 retval = std::string ("repeated error test failures on the last attempted step (t = ")
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	609 + t_curr + ")";
3996 98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	610 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	611
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	612 case -7:
4043 6fae69a1796e [project @ 2002-08-15 03:36:27 by jwe] jwe parents: 3998 diff changeset	613 retval = std::string ("the corrector could not converge (t = ")
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	614 + t_curr + ")";
3996 98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	615 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	616
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	617 case -8:
4043 6fae69a1796e [project @ 2002-08-15 03:36:27 by jwe] jwe parents: 3998 diff changeset	618 retval = std::string ("the matrix of partial derivatives is singular (t = ")
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	619 + t_curr + ")";
3996 98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	620 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	621
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	622 case -9:
4043 6fae69a1796e [project @ 2002-08-15 03:36:27 by jwe] jwe parents: 3998 diff changeset	623 retval = std::string ("the corrector could not converge (t = ")
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	624 + t_curr + "; repeated test failures)";
3996 98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	625 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	626
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	627 case -10:
4043 6fae69a1796e [project @ 2002-08-15 03:36:27 by jwe] jwe parents: 3998 diff changeset	628 retval = std::string ("corrector could not converge because IRES was -1 (t = ")
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	629 + t_curr + ")";
3996 98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	630 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	631
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	632 case -11:
4043 6fae69a1796e [project @ 2002-08-15 03:36:27 by jwe] jwe parents: 3998 diff changeset	633 retval = std::string ("return requested in user-supplied function (t = ")
10314 07ebe522dac2 untabify liboctave C++ sources John W. Eaton <jwe@octave.org> parents: 10158 diff changeset	634 + t_curr + ")";
3996 98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	635 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	636
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	637 case -12:
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	638 retval = "failed to compute consistent initial conditions";
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	639 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	640
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	641 case -33:
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	642 retval = "unrecoverable error (see printed message)";
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	643 break;
98107d72871c [project @ 2002-07-17 16:22:05 by jwe] jwe parents: 3995 diff changeset	644
3995 ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	645 default:
ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	646 retval = "unknown error state";
ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	647 break;
ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	648 }
ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	649
ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	650 return retval;
ee0304212be0 [project @ 2002-07-17 04:32:42 by jwe] jwe parents: 3994 diff changeset	651 }

Mercurial > hg > octave-nkf

annotate liboctave/DASRT.cc @ 14101:ff5588774680 stable ss-3-5-91