1 | /* |
---|
2 | * ----------------------------------------------------------------- |
---|
3 | * $Revision$ |
---|
4 | * $Date$ |
---|
5 | * ----------------------------------------------------------------- |
---|
6 | * Programmer(s): Radu Serban @ LLNL |
---|
7 | * ----------------------------------------------------------------- |
---|
8 | * LLNS Copyright Start |
---|
9 | * Copyright (c) 2014, Lawrence Livermore National Security |
---|
10 | * This work was performed under the auspices of the U.S. Department |
---|
11 | * of Energy by Lawrence Livermore National Laboratory in part under |
---|
12 | * Contract W-7405-Eng-48 and in part under Contract DE-AC52-07NA27344. |
---|
13 | * Produced at the Lawrence Livermore National Laboratory. |
---|
14 | * All rights reserved. |
---|
15 | * For details, see the LICENSE file. |
---|
16 | * LLNS Copyright End |
---|
17 | * ----------------------------------------------------------------- |
---|
18 | * Implementation header file for the main CVODES integrator. |
---|
19 | * ----------------------------------------------------------------- |
---|
20 | */ |
---|
21 | |
---|
22 | #ifndef _CVODES_IMPL_H |
---|
23 | #define _CVODES_IMPL_H |
---|
24 | |
---|
25 | #include <stdarg.h> |
---|
26 | |
---|
27 | #include <cvodes/cvodes.h> |
---|
28 | #include <sundials/sundials_nvector.h> |
---|
29 | #include <sundials/sundials_types.h> |
---|
30 | |
---|
31 | #ifdef __cplusplus /* wrapper to enable C++ usage */ |
---|
32 | extern "C" { |
---|
33 | #endif |
---|
34 | |
---|
35 | /* |
---|
36 | * ================================================================= |
---|
37 | * I N T E R N A L C V O D E S C O N S T A N T S |
---|
38 | * ================================================================= |
---|
39 | */ |
---|
40 | |
---|
41 | /* Basic CVODES constants */ |
---|
42 | |
---|
43 | #define ADAMS_Q_MAX 12 /* max value of q for lmm == ADAMS */ |
---|
44 | #define BDF_Q_MAX 5 /* max value of q for lmm == BDF */ |
---|
45 | #define Q_MAX ADAMS_Q_MAX /* max value of q for either lmm */ |
---|
46 | #define L_MAX (Q_MAX+1) /* max value of L for either lmm */ |
---|
47 | #define NUM_TESTS 5 /* number of error test quantities */ |
---|
48 | |
---|
49 | #define HMIN_DEFAULT RCONST(0.0) /* hmin default value */ |
---|
50 | #define HMAX_INV_DEFAULT RCONST(0.0) /* hmax_inv default value */ |
---|
51 | #define MXHNIL_DEFAULT 10 /* mxhnil default value */ |
---|
52 | #define MXSTEP_DEFAULT 500 /* mxstep default value */ |
---|
53 | |
---|
54 | /* |
---|
55 | * ================================================================= |
---|
56 | * F O R W A R D P O I N T E R R E F E R E N C E S |
---|
57 | * ================================================================= |
---|
58 | */ |
---|
59 | |
---|
60 | typedef struct CVadjMemRec *CVadjMem; |
---|
61 | typedef struct CkpntMemRec *CkpntMem; |
---|
62 | typedef struct DtpntMemRec *DtpntMem; |
---|
63 | typedef struct CVodeBMemRec *CVodeBMem; |
---|
64 | |
---|
65 | /* |
---|
66 | * ================================================================= |
---|
67 | * M A I N I N T E G R A T O R M E M O R Y B L O C K |
---|
68 | * ================================================================= |
---|
69 | */ |
---|
70 | |
---|
71 | |
---|
72 | /* |
---|
73 | * ----------------------------------------------------------------- |
---|
74 | * Types: struct CVodeMemRec, CVodeMem |
---|
75 | * ----------------------------------------------------------------- |
---|
76 | * The type CVodeMem is type pointer to struct CVodeMemRec. |
---|
77 | * This structure contains fields to keep track of problem state. |
---|
78 | * ----------------------------------------------------------------- |
---|
79 | */ |
---|
80 | |
---|
81 | typedef struct CVodeMemRec { |
---|
82 | |
---|
83 | realtype cv_uround; /* machine unit roundoff */ |
---|
84 | |
---|
85 | /*-------------------------- |
---|
86 | Problem Specification Data |
---|
87 | --------------------------*/ |
---|
88 | |
---|
89 | CVRhsFn cv_f; /* y' = f(t,y(t)) */ |
---|
90 | void *cv_user_data; /* user pointer passed to f */ |
---|
91 | |
---|
92 | int cv_lmm; /* lmm = ADAMS or BDF */ |
---|
93 | int cv_iter; /* iter = FUNCTIONAL or NEWTON */ |
---|
94 | |
---|
95 | int cv_itol; /* itol = CV_SS, CV_SV, or CV_WF, or CV_NN */ |
---|
96 | realtype cv_reltol; /* relative tolerance */ |
---|
97 | realtype cv_Sabstol; /* scalar absolute tolerance */ |
---|
98 | N_Vector cv_Vabstol; /* vector absolute tolerance */ |
---|
99 | booleantype cv_user_efun; /* SUNTRUE if user sets efun */ |
---|
100 | CVEwtFn cv_efun; /* function to set ewt */ |
---|
101 | void *cv_e_data; /* user pointer passed to efun */ |
---|
102 | |
---|
103 | booleantype cv_constraintsSet; /* constraints vector present: |
---|
104 | do constraints calc */ |
---|
105 | /*----------------------- |
---|
106 | Quadrature Related Data |
---|
107 | -----------------------*/ |
---|
108 | |
---|
109 | booleantype cv_quadr; /* SUNTRUE if integrating quadratures */ |
---|
110 | |
---|
111 | CVQuadRhsFn cv_fQ; /* q' = fQ(t, y(t)) */ |
---|
112 | |
---|
113 | booleantype cv_errconQ; /* SUNTRUE if quadrs. are included in error test */ |
---|
114 | |
---|
115 | int cv_itolQ; /* itolQ = CV_SS or CV_SV */ |
---|
116 | realtype cv_reltolQ; /* relative tolerance for quadratures */ |
---|
117 | realtype cv_SabstolQ; /* scalar absolute tolerance for quadratures */ |
---|
118 | N_Vector cv_VabstolQ; /* vector absolute tolerance for quadratures */ |
---|
119 | |
---|
120 | /*------------------------ |
---|
121 | Sensitivity Related Data |
---|
122 | ------------------------*/ |
---|
123 | |
---|
124 | booleantype cv_sensi; /* SUNTRUE if computing sensitivities */ |
---|
125 | |
---|
126 | int cv_Ns; /* Number of sensitivities */ |
---|
127 | |
---|
128 | int cv_ism; /* ism = SIMULTANEOUS or STAGGERED */ |
---|
129 | |
---|
130 | CVSensRhsFn cv_fS; /* fS = (df/dy)*yS + (df/dp) */ |
---|
131 | CVSensRhs1Fn cv_fS1; /* fS1 = (df/dy)*yS_i + (df/dp) */ |
---|
132 | void *cv_fS_data; /* data pointer passed to fS */ |
---|
133 | booleantype cv_fSDQ; /* SUNTRUE if using internal DQ functions */ |
---|
134 | int cv_ifS; /* ifS = ALLSENS or ONESENS */ |
---|
135 | |
---|
136 | realtype *cv_p; /* parameters in f(t,y,p) */ |
---|
137 | realtype *cv_pbar; /* scale factors for parameters */ |
---|
138 | int *cv_plist; /* list of sensitivities */ |
---|
139 | int cv_DQtype; /* central/forward finite differences */ |
---|
140 | realtype cv_DQrhomax; /* cut-off value for separate/simultaneous FD */ |
---|
141 | |
---|
142 | booleantype cv_errconS; /* SUNTRUE if yS are considered in err. control */ |
---|
143 | |
---|
144 | int cv_itolS; |
---|
145 | realtype cv_reltolS; /* relative tolerance for sensitivities */ |
---|
146 | realtype *cv_SabstolS; /* scalar absolute tolerances for sensi. */ |
---|
147 | N_Vector *cv_VabstolS; /* vector absolute tolerances for sensi. */ |
---|
148 | |
---|
149 | /*----------------------------------- |
---|
150 | Quadrature Sensitivity Related Data |
---|
151 | -----------------------------------*/ |
---|
152 | |
---|
153 | booleantype cv_quadr_sensi; /* SUNTRUE if computing sensitivties of quadrs. */ |
---|
154 | |
---|
155 | CVQuadSensRhsFn cv_fQS; /* fQS = (dfQ/dy)*yS + (dfQ/dp) */ |
---|
156 | void *cv_fQS_data; /* data pointer passed to fQS */ |
---|
157 | booleantype cv_fQSDQ; /* SUNTRUE if using internal DQ functions */ |
---|
158 | |
---|
159 | booleantype cv_errconQS; /* SUNTRUE if yQS are considered in err. con. */ |
---|
160 | |
---|
161 | int cv_itolQS; |
---|
162 | realtype cv_reltolQS; /* relative tolerance for yQS */ |
---|
163 | realtype *cv_SabstolQS; /* scalar absolute tolerances for yQS */ |
---|
164 | N_Vector *cv_VabstolQS; /* vector absolute tolerances for yQS */ |
---|
165 | |
---|
166 | /*----------------------- |
---|
167 | Nordsieck History Array |
---|
168 | -----------------------*/ |
---|
169 | |
---|
170 | N_Vector cv_zn[L_MAX]; /* Nordsieck array, of size N x (q+1). |
---|
171 | zn[j] is a vector of length N (j=0,...,q) |
---|
172 | zn[j] = [1/factorial(j)] * h^j * |
---|
173 | (jth derivative of the interpolating poly.) */ |
---|
174 | |
---|
175 | /*------------------- |
---|
176 | Vectors of length N |
---|
177 | -------------------*/ |
---|
178 | |
---|
179 | N_Vector cv_ewt; /* error weight vector */ |
---|
180 | N_Vector cv_y; /* y is used as temporary storage by the solver. |
---|
181 | The memory is provided by the user to CVode |
---|
182 | where the vector is named yout. */ |
---|
183 | N_Vector cv_acor; /* In the context of the solution of the |
---|
184 | nonlinear equation, acor = y_n(m) - y_n(0). |
---|
185 | On return, this vector is scaled to give |
---|
186 | the estimated local error in y. */ |
---|
187 | N_Vector cv_tempv; /* temporary storage vector */ |
---|
188 | N_Vector cv_ftemp; /* temporary storage vector */ |
---|
189 | |
---|
190 | N_Vector cv_mm; /* mask vector in constraints tests */ |
---|
191 | N_Vector cv_constraints; /* vector of inequality constraint options */ |
---|
192 | |
---|
193 | /*-------------------------- |
---|
194 | Quadrature Related Vectors |
---|
195 | --------------------------*/ |
---|
196 | |
---|
197 | N_Vector cv_znQ[L_MAX]; /* Nordsieck arrays for quadratures */ |
---|
198 | N_Vector cv_ewtQ; /* error weight vector for quadratures */ |
---|
199 | N_Vector cv_yQ; /* Unlike y, yQ is not allocated by the user */ |
---|
200 | N_Vector cv_acorQ; /* acorQ = yQ_n(m) - yQ_n(0) */ |
---|
201 | N_Vector cv_tempvQ; /* temporary storage vector (~ tempv) */ |
---|
202 | |
---|
203 | /*--------------------------- |
---|
204 | Sensitivity Related Vectors |
---|
205 | ---------------------------*/ |
---|
206 | |
---|
207 | N_Vector *cv_znS[L_MAX]; /* Nordsieck arrays for sensitivities */ |
---|
208 | N_Vector *cv_ewtS; /* error weight vectors for sensitivities */ |
---|
209 | N_Vector *cv_yS; /* yS=yS0 (allocated by the user) */ |
---|
210 | N_Vector *cv_acorS; /* acorS = yS_n(m) - yS_n(0) */ |
---|
211 | N_Vector *cv_tempvS; /* temporary storage vector (~ tempv) */ |
---|
212 | N_Vector *cv_ftempS; /* temporary storage vector (~ ftemp) */ |
---|
213 | |
---|
214 | booleantype cv_stgr1alloc; /* Did we allocate ncfS1, ncfnS1, and nniS1? */ |
---|
215 | |
---|
216 | /*-------------------------------------- |
---|
217 | Quadrature Sensitivity Related Vectors |
---|
218 | --------------------------------------*/ |
---|
219 | |
---|
220 | N_Vector *cv_znQS[L_MAX]; /* Nordsieck arrays for quadr. sensitivities */ |
---|
221 | N_Vector *cv_ewtQS; /* error weight vectors for sensitivities */ |
---|
222 | N_Vector *cv_yQS; /* Unlike yS, yQS is not allocated by the user */ |
---|
223 | N_Vector *cv_acorQS; /* acorQS = yQS_n(m) - yQS_n(0) */ |
---|
224 | N_Vector *cv_tempvQS; /* temporary storage vector (~ tempv) */ |
---|
225 | N_Vector cv_ftempQ; /* temporary storage vector (~ ftemp) */ |
---|
226 | |
---|
227 | /*----------------- |
---|
228 | Tstop information |
---|
229 | -----------------*/ |
---|
230 | |
---|
231 | booleantype cv_tstopset; |
---|
232 | realtype cv_tstop; |
---|
233 | |
---|
234 | /*--------- |
---|
235 | Step Data |
---|
236 | ---------*/ |
---|
237 | |
---|
238 | int cv_q; /* current order */ |
---|
239 | int cv_qprime; /* order to be used on the next step |
---|
240 | * qprime = q-1, q, or q+1 */ |
---|
241 | int cv_next_q; /* order to be used on the next step */ |
---|
242 | int cv_qwait; /* number of internal steps to wait before |
---|
243 | * considering a change in q */ |
---|
244 | int cv_L; /* L = q + 1 */ |
---|
245 | |
---|
246 | realtype cv_hin; |
---|
247 | realtype cv_h; /* current step size */ |
---|
248 | realtype cv_hprime; /* step size to be used on the next step */ |
---|
249 | realtype cv_next_h; /* step size to be used on the next step */ |
---|
250 | realtype cv_eta; /* eta = hprime / h */ |
---|
251 | realtype cv_hscale; /* value of h used in zn */ |
---|
252 | realtype cv_tn; /* current internal value of t */ |
---|
253 | realtype cv_tretlast; /* last value of t returned */ |
---|
254 | |
---|
255 | realtype cv_tau[L_MAX+1]; /* array of previous q+1 successful step |
---|
256 | * sizes indexed from 1 to q+1 */ |
---|
257 | realtype cv_tq[NUM_TESTS+1]; /* array of test quantities indexed from |
---|
258 | * 1 to NUM_TESTS(=5) */ |
---|
259 | realtype cv_l[L_MAX]; /* coefficients of l(x) (degree q poly) */ |
---|
260 | |
---|
261 | realtype cv_rl1; /* the scalar 1/l[1] */ |
---|
262 | realtype cv_gamma; /* gamma = h * rl1 */ |
---|
263 | realtype cv_gammap; /* gamma at the last setup call */ |
---|
264 | realtype cv_gamrat; /* gamma / gammap */ |
---|
265 | |
---|
266 | realtype cv_crate; /* est. corrector conv. rate in Nls */ |
---|
267 | realtype cv_crateS; /* est. corrector conv. rate in NlsStgr */ |
---|
268 | realtype cv_acnrm; /* | acor | */ |
---|
269 | realtype cv_acnrmQ; /* | acorQ | */ |
---|
270 | realtype cv_acnrmS; /* | acorS | */ |
---|
271 | realtype cv_acnrmQS; /* | acorQS | */ |
---|
272 | realtype cv_nlscoef; /* coeficient in nonlinear convergence test */ |
---|
273 | int cv_mnewt; /* Newton iteration counter */ |
---|
274 | int *cv_ncfS1; /* Array of Ns local counters for conv. |
---|
275 | * failures (used in CVStep for STAGGERED1) */ |
---|
276 | |
---|
277 | /*------ |
---|
278 | Limits |
---|
279 | ------*/ |
---|
280 | |
---|
281 | int cv_qmax; /* q <= qmax */ |
---|
282 | long int cv_mxstep; /* maximum number of internal steps for one |
---|
283 | user call */ |
---|
284 | int cv_maxcor; /* maximum number of corrector iterations for |
---|
285 | the solution of the nonlinear equation */ |
---|
286 | int cv_maxcorS; |
---|
287 | int cv_mxhnil; /* max. number of warning messages issued to the |
---|
288 | user that t + h == t for the next internal step */ |
---|
289 | int cv_maxnef; /* maximum number of error test failures */ |
---|
290 | int cv_maxncf; /* maximum number of nonlinear conv. failures */ |
---|
291 | |
---|
292 | realtype cv_hmin; /* |h| >= hmin */ |
---|
293 | realtype cv_hmax_inv; /* |h| <= 1/hmax_inv */ |
---|
294 | realtype cv_etamax; /* eta <= etamax */ |
---|
295 | |
---|
296 | /*---------- |
---|
297 | Counters |
---|
298 | ----------*/ |
---|
299 | |
---|
300 | long int cv_nst; /* number of internal steps taken */ |
---|
301 | |
---|
302 | long int cv_nfe; /* number of f calls */ |
---|
303 | long int cv_nfQe; /* number of fQ calls */ |
---|
304 | long int cv_nfSe; /* number of fS calls */ |
---|
305 | long int cv_nfeS; /* number of f calls from sensi DQ */ |
---|
306 | long int cv_nfQSe; /* number of fQS calls */ |
---|
307 | long int cv_nfQeS; /* number of fQ calls from sensi DQ */ |
---|
308 | |
---|
309 | |
---|
310 | long int cv_ncfn; /* number of corrector convergence failures */ |
---|
311 | long int cv_ncfnS; /* number of total sensi. corr. conv. failures */ |
---|
312 | long int *cv_ncfnS1; /* number of sensi. corrector conv. failures */ |
---|
313 | |
---|
314 | long int cv_nni; /* number of nonlinear iterations performed */ |
---|
315 | long int cv_nniS; /* number of total sensi. nonlinear iterations */ |
---|
316 | long int *cv_nniS1; /* number of sensi. nonlinear iterations */ |
---|
317 | |
---|
318 | long int cv_netf; /* number of error test failures */ |
---|
319 | long int cv_netfQ; /* number of quadr. error test failures */ |
---|
320 | long int cv_netfS; /* number of sensi. error test failures */ |
---|
321 | long int cv_netfQS; /* number of quadr. sensi. error test failures */ |
---|
322 | |
---|
323 | long int cv_nsetups; /* number of setup calls */ |
---|
324 | long int cv_nsetupsS; /* number of setup calls due to sensitivities */ |
---|
325 | |
---|
326 | int cv_nhnil; /* number of messages issued to the user that |
---|
327 | t + h == t for the next iternal step */ |
---|
328 | |
---|
329 | /*----------------------------- |
---|
330 | Space requirements for CVODES |
---|
331 | -----------------------------*/ |
---|
332 | |
---|
333 | sunindextype cv_lrw1; /* no. of realtype words in 1 N_Vector y */ |
---|
334 | sunindextype cv_liw1; /* no. of integer words in 1 N_Vector y */ |
---|
335 | sunindextype cv_lrw1Q; /* no. of realtype words in 1 N_Vector yQ */ |
---|
336 | sunindextype cv_liw1Q; /* no. of integer words in 1 N_Vector yQ */ |
---|
337 | long int cv_lrw; /* no. of realtype words in CVODES work vectors */ |
---|
338 | long int cv_liw; /* no. of integer words in CVODES work vectors */ |
---|
339 | |
---|
340 | /*---------------- |
---|
341 | Step size ratios |
---|
342 | ----------------*/ |
---|
343 | |
---|
344 | realtype cv_etaqm1; /* ratio of new to old h for order q-1 */ |
---|
345 | realtype cv_etaq; /* ratio of new to old h for order q */ |
---|
346 | realtype cv_etaqp1; /* ratio of new to old h for order q+1 */ |
---|
347 | |
---|
348 | /*------------------ |
---|
349 | Linear Solver Data |
---|
350 | ------------------*/ |
---|
351 | |
---|
352 | /* Linear Solver functions to be called */ |
---|
353 | |
---|
354 | int (*cv_linit)(struct CVodeMemRec *cv_mem); |
---|
355 | |
---|
356 | int (*cv_lsetup)(struct CVodeMemRec *cv_mem, int convfail, |
---|
357 | N_Vector ypred, N_Vector fpred, booleantype *jcurPtr, |
---|
358 | N_Vector vtemp1, N_Vector vtemp2, N_Vector vtemp3); |
---|
359 | |
---|
360 | int (*cv_lsolve)(struct CVodeMemRec *cv_mem, N_Vector b, N_Vector weight, |
---|
361 | N_Vector ycur, N_Vector fcur); |
---|
362 | |
---|
363 | int (*cv_lfree)(struct CVodeMemRec *cv_mem); |
---|
364 | |
---|
365 | /* Linear Solver specific memory */ |
---|
366 | |
---|
367 | void *cv_lmem; |
---|
368 | |
---|
369 | /* Flag to request a call to the setup routine */ |
---|
370 | |
---|
371 | booleantype cv_forceSetup; |
---|
372 | |
---|
373 | /*------------ |
---|
374 | Saved Values |
---|
375 | ------------*/ |
---|
376 | |
---|
377 | int cv_qu; /* last successful q value used */ |
---|
378 | long int cv_nstlp; /* step number of last setup call */ |
---|
379 | realtype cv_h0u; /* actual initial stepsize */ |
---|
380 | realtype cv_hu; /* last successful h value used */ |
---|
381 | realtype cv_saved_tq5; /* saved value of tq[5] */ |
---|
382 | booleantype cv_jcur; /* is Jacobian info for linear solver current? */ |
---|
383 | int cv_convfail; /* flag storing previous solver failure mode */ |
---|
384 | realtype cv_tolsf; /* tolerance scale factor */ |
---|
385 | int cv_qmax_alloc; /* qmax used when allocating mem */ |
---|
386 | int cv_qmax_allocQ; /* qmax used when allocating quad. mem */ |
---|
387 | int cv_qmax_allocS; /* qmax used when allocating sensi. mem */ |
---|
388 | int cv_qmax_allocQS; /* qmax used when allocating quad. sensi. mem */ |
---|
389 | int cv_indx_acor; /* index of zn vector in which acor is saved */ |
---|
390 | |
---|
391 | /*-------------------------------------------------------------------- |
---|
392 | Flags turned ON by CVodeInit, CVodeSensMalloc, and CVodeQuadMalloc |
---|
393 | and read by CVodeReInit, CVodeSensReInit, and CVodeQuadReInit |
---|
394 | --------------------------------------------------------------------*/ |
---|
395 | |
---|
396 | booleantype cv_VabstolMallocDone; |
---|
397 | booleantype cv_MallocDone; |
---|
398 | booleantype cv_constraintsMallocDone; |
---|
399 | |
---|
400 | booleantype cv_VabstolQMallocDone; |
---|
401 | booleantype cv_QuadMallocDone; |
---|
402 | |
---|
403 | booleantype cv_VabstolSMallocDone; |
---|
404 | booleantype cv_SabstolSMallocDone; |
---|
405 | booleantype cv_SensMallocDone; |
---|
406 | |
---|
407 | booleantype cv_VabstolQSMallocDone; |
---|
408 | booleantype cv_SabstolQSMallocDone; |
---|
409 | booleantype cv_QuadSensMallocDone; |
---|
410 | |
---|
411 | /*------------------------------------------- |
---|
412 | Error handler function and error ouput file |
---|
413 | -------------------------------------------*/ |
---|
414 | |
---|
415 | CVErrHandlerFn cv_ehfun; /* Error messages are handled by ehfun */ |
---|
416 | void *cv_eh_data; /* dats pointer passed to ehfun */ |
---|
417 | FILE *cv_errfp; /* CVODES error messages are sent to errfp */ |
---|
418 | |
---|
419 | /*------------------------- |
---|
420 | Stability Limit Detection |
---|
421 | -------------------------*/ |
---|
422 | |
---|
423 | booleantype cv_sldeton; /* Is Stability Limit Detection on? */ |
---|
424 | realtype cv_ssdat[6][4]; /* scaled data array for STALD */ |
---|
425 | int cv_nscon; /* counter for STALD method */ |
---|
426 | long int cv_nor; /* counter for number of order reductions */ |
---|
427 | |
---|
428 | /*---------------- |
---|
429 | Rootfinding Data |
---|
430 | ----------------*/ |
---|
431 | |
---|
432 | CVRootFn cv_gfun; /* Function g for roots sought */ |
---|
433 | int cv_nrtfn; /* number of components of g */ |
---|
434 | int *cv_iroots; /* array for root information */ |
---|
435 | int *cv_rootdir; /* array specifying direction of zero-crossing */ |
---|
436 | realtype cv_tlo; /* nearest endpoint of interval in root search */ |
---|
437 | realtype cv_thi; /* farthest endpoint of interval in root search */ |
---|
438 | realtype cv_trout; /* t value returned by rootfinding routine */ |
---|
439 | realtype *cv_glo; /* saved array of g values at t = tlo */ |
---|
440 | realtype *cv_ghi; /* saved array of g values at t = thi */ |
---|
441 | realtype *cv_grout; /* array of g values at t = trout */ |
---|
442 | realtype cv_toutc; /* copy of tout (if NORMAL mode) */ |
---|
443 | realtype cv_ttol; /* tolerance on root location trout */ |
---|
444 | int cv_taskc; /* copy of parameter itask */ |
---|
445 | int cv_irfnd; /* flag showing whether last step had a root */ |
---|
446 | long int cv_nge; /* counter for g evaluations */ |
---|
447 | booleantype *cv_gactive; /* array with active/inactive event functions */ |
---|
448 | int cv_mxgnull; /* number of warning messages about possible g==0 */ |
---|
449 | |
---|
450 | /*------------------------ |
---|
451 | Adjoint sensitivity data |
---|
452 | ------------------------*/ |
---|
453 | |
---|
454 | booleantype cv_adj; /* SUNTRUE if performing ASA */ |
---|
455 | |
---|
456 | struct CVadjMemRec *cv_adj_mem; /* Pointer to adjoint memory structure */ |
---|
457 | |
---|
458 | booleantype cv_adjMallocDone; |
---|
459 | |
---|
460 | } *CVodeMem; |
---|
461 | |
---|
462 | |
---|
463 | /* |
---|
464 | * ================================================================= |
---|
465 | * A D J O I N T M O D U L E M E M O R Y B L O C K |
---|
466 | * ================================================================= |
---|
467 | */ |
---|
468 | |
---|
469 | /* |
---|
470 | * ----------------------------------------------------------------- |
---|
471 | * Types : struct CkpntMemRec, CkpntMem |
---|
472 | * ----------------------------------------------------------------- |
---|
473 | * The type CkpntMem is type pointer to struct CkpntMemRec. |
---|
474 | * This structure contains fields to store all information at a |
---|
475 | * check point that is needed to 'hot' start cvodes. |
---|
476 | * ----------------------------------------------------------------- |
---|
477 | */ |
---|
478 | |
---|
479 | struct CkpntMemRec { |
---|
480 | |
---|
481 | /* Integration limits */ |
---|
482 | realtype ck_t0; |
---|
483 | realtype ck_t1; |
---|
484 | |
---|
485 | /* Nordsieck History Array */ |
---|
486 | N_Vector ck_zn[L_MAX]; |
---|
487 | |
---|
488 | /* Do we need to carry quadratures? */ |
---|
489 | booleantype ck_quadr; |
---|
490 | |
---|
491 | /* Nordsieck History Array for quadratures */ |
---|
492 | N_Vector ck_znQ[L_MAX]; |
---|
493 | |
---|
494 | /* Do we need to carry sensitivities? */ |
---|
495 | booleantype ck_sensi; |
---|
496 | |
---|
497 | /* number of sensitivities */ |
---|
498 | int ck_Ns; |
---|
499 | |
---|
500 | /* Nordsieck History Array for sensitivities */ |
---|
501 | N_Vector *ck_znS[L_MAX]; |
---|
502 | |
---|
503 | /* Do we need to carry quadrature sensitivities? */ |
---|
504 | booleantype ck_quadr_sensi; |
---|
505 | |
---|
506 | /* Nordsieck History Array for quadrature sensitivities */ |
---|
507 | N_Vector *ck_znQS[L_MAX]; |
---|
508 | |
---|
509 | /* Was ck_zn[qmax] allocated? |
---|
510 | ck_zqm = 0 - no |
---|
511 | ck_zqm = qmax - yes */ |
---|
512 | int ck_zqm; |
---|
513 | |
---|
514 | /* Step data */ |
---|
515 | long int ck_nst; |
---|
516 | realtype ck_tretlast; |
---|
517 | int ck_q; |
---|
518 | int ck_qprime; |
---|
519 | int ck_qwait; |
---|
520 | int ck_L; |
---|
521 | realtype ck_gammap; |
---|
522 | realtype ck_h; |
---|
523 | realtype ck_hprime; |
---|
524 | realtype ck_hscale; |
---|
525 | realtype ck_eta; |
---|
526 | realtype ck_etamax; |
---|
527 | realtype ck_tau[L_MAX+1]; |
---|
528 | realtype ck_tq[NUM_TESTS+1]; |
---|
529 | realtype ck_l[L_MAX]; |
---|
530 | |
---|
531 | /* Saved values */ |
---|
532 | realtype ck_saved_tq5; |
---|
533 | |
---|
534 | /* Pointer to next structure in list */ |
---|
535 | struct CkpntMemRec *ck_next; |
---|
536 | |
---|
537 | }; |
---|
538 | |
---|
539 | /* |
---|
540 | * ----------------------------------------------------------------- |
---|
541 | * Types for functions provided by an interpolation module |
---|
542 | * ----------------------------------------------------------------- |
---|
543 | * cvaIMMallocFn: Type for a function that initializes the content |
---|
544 | * field of the structures in the dt array |
---|
545 | * cvaIMFreeFn: Type for a function that deallocates the content |
---|
546 | * field of the structures in the dt array |
---|
547 | * cvaIMGetYFn: Type for a function that returns the |
---|
548 | * interpolated forward solution. |
---|
549 | * cvaIMStorePnt: Type for a function that stores a new |
---|
550 | * point in the structure d |
---|
551 | * ----------------------------------------------------------------- |
---|
552 | */ |
---|
553 | |
---|
554 | typedef booleantype (*cvaIMMallocFn)(CVodeMem cv_mem); |
---|
555 | typedef void (*cvaIMFreeFn)(CVodeMem cv_mem); |
---|
556 | typedef int (*cvaIMGetYFn)(CVodeMem cv_mem, realtype t, N_Vector y, N_Vector *yS); |
---|
557 | typedef int (*cvaIMStorePntFn)(CVodeMem cv_mem, DtpntMem d); |
---|
558 | |
---|
559 | /* |
---|
560 | * ----------------------------------------------------------------- |
---|
561 | * Type : struct DtpntMemRec |
---|
562 | * ----------------------------------------------------------------- |
---|
563 | * This structure contains fields to store all information at a |
---|
564 | * data point that is needed to interpolate solution of forward |
---|
565 | * simulations. Its content field depends on IMtype. |
---|
566 | * ----------------------------------------------------------------- |
---|
567 | */ |
---|
568 | |
---|
569 | struct DtpntMemRec { |
---|
570 | realtype t; /* time */ |
---|
571 | void *content; /* IMtype-dependent content */ |
---|
572 | }; |
---|
573 | |
---|
574 | /* Data for cubic Hermite interpolation */ |
---|
575 | typedef struct HermiteDataMemRec { |
---|
576 | N_Vector y; |
---|
577 | N_Vector yd; |
---|
578 | N_Vector *yS; |
---|
579 | N_Vector *ySd; |
---|
580 | } *HermiteDataMem; |
---|
581 | |
---|
582 | /* Data for polynomial interpolation */ |
---|
583 | typedef struct PolynomialDataMemRec { |
---|
584 | N_Vector y; |
---|
585 | N_Vector *yS; |
---|
586 | int order; |
---|
587 | } *PolynomialDataMem; |
---|
588 | |
---|
589 | |
---|
590 | /* |
---|
591 | * ----------------------------------------------------------------- |
---|
592 | * Type : struct CVodeBMemRec |
---|
593 | * ----------------------------------------------------------------- |
---|
594 | * The type CVodeBMem is a pointer to a structure which stores all |
---|
595 | * information for ONE backward problem. |
---|
596 | * The CVadjMem structure contains a linked list of CVodeBMem pointers |
---|
597 | * ----------------------------------------------------------------- |
---|
598 | */ |
---|
599 | |
---|
600 | struct CVodeBMemRec { |
---|
601 | |
---|
602 | /* Index of this backward problem */ |
---|
603 | int cv_index; |
---|
604 | |
---|
605 | /* Time at which the backward problem is initialized */ |
---|
606 | realtype cv_t0; |
---|
607 | |
---|
608 | /* CVODES memory for this backward problem */ |
---|
609 | CVodeMem cv_mem; |
---|
610 | |
---|
611 | /* Flags to indicate that this backward problem's RHS or quad RHS |
---|
612 | * require forward sensitivities */ |
---|
613 | booleantype cv_f_withSensi; |
---|
614 | booleantype cv_fQ_withSensi; |
---|
615 | |
---|
616 | /* Right hand side function for backward run */ |
---|
617 | CVRhsFnB cv_f; |
---|
618 | CVRhsFnBS cv_fs; |
---|
619 | |
---|
620 | /* Right hand side quadrature function for backward run */ |
---|
621 | CVQuadRhsFnB cv_fQ; |
---|
622 | CVQuadRhsFnBS cv_fQs; |
---|
623 | |
---|
624 | /* User user_data */ |
---|
625 | void *cv_user_data; |
---|
626 | |
---|
627 | /* Memory block for a linear solver's interface to CVODEA */ |
---|
628 | void *cv_lmem; |
---|
629 | |
---|
630 | /* Function to free any memory allocated by the linear solver */ |
---|
631 | int (*cv_lfree)(CVodeBMem cvB_mem); |
---|
632 | |
---|
633 | /* Memory block for a preconditioner's module interface to CVODEA */ |
---|
634 | void *cv_pmem; |
---|
635 | |
---|
636 | /* Function to free any memory allocated by the preconditioner module */ |
---|
637 | int (*cv_pfree)(CVodeBMem cvB_mem); |
---|
638 | |
---|
639 | /* Time at which to extract solution / quadratures */ |
---|
640 | realtype cv_tout; |
---|
641 | |
---|
642 | /* Workspace Nvector */ |
---|
643 | N_Vector cv_y; |
---|
644 | |
---|
645 | /* Pointer to next structure in list */ |
---|
646 | struct CVodeBMemRec *cv_next; |
---|
647 | |
---|
648 | }; |
---|
649 | |
---|
650 | /* |
---|
651 | * ----------------------------------------------------------------- |
---|
652 | * Type : struct CVadjMemRec |
---|
653 | * ----------------------------------------------------------------- |
---|
654 | * The type CVadjMem is type pointer to struct CVadjMemRec. |
---|
655 | * This structure contins fields to store all information |
---|
656 | * necessary for adjoint sensitivity analysis. |
---|
657 | * ----------------------------------------------------------------- |
---|
658 | */ |
---|
659 | |
---|
660 | struct CVadjMemRec { |
---|
661 | |
---|
662 | /* -------------------- |
---|
663 | * Forward problem data |
---|
664 | * -------------------- */ |
---|
665 | |
---|
666 | /* Integration interval */ |
---|
667 | realtype ca_tinitial, ca_tfinal; |
---|
668 | |
---|
669 | /* Flag for first call to CVodeF */ |
---|
670 | booleantype ca_firstCVodeFcall; |
---|
671 | |
---|
672 | /* Flag if CVodeF was called with TSTOP */ |
---|
673 | booleantype ca_tstopCVodeFcall; |
---|
674 | realtype ca_tstopCVodeF; |
---|
675 | |
---|
676 | /* ---------------------- |
---|
677 | * Backward problems data |
---|
678 | * ---------------------- */ |
---|
679 | |
---|
680 | /* Storage for backward problems */ |
---|
681 | struct CVodeBMemRec *cvB_mem; |
---|
682 | |
---|
683 | /* Number of backward problems */ |
---|
684 | int ca_nbckpbs; |
---|
685 | |
---|
686 | /* Address of current backward problem */ |
---|
687 | struct CVodeBMemRec *ca_bckpbCrt; |
---|
688 | |
---|
689 | /* Flag for first call to CVodeB */ |
---|
690 | booleantype ca_firstCVodeBcall; |
---|
691 | |
---|
692 | /* ---------------- |
---|
693 | * Check point data |
---|
694 | * ---------------- */ |
---|
695 | |
---|
696 | /* Storage for check point information */ |
---|
697 | struct CkpntMemRec *ck_mem; |
---|
698 | |
---|
699 | /* Number of check points */ |
---|
700 | int ca_nckpnts; |
---|
701 | |
---|
702 | /* address of the check point structure for which data is available */ |
---|
703 | struct CkpntMemRec *ca_ckpntData; |
---|
704 | |
---|
705 | /* ------------------ |
---|
706 | * Interpolation data |
---|
707 | * ------------------ */ |
---|
708 | |
---|
709 | /* Number of steps between 2 check points */ |
---|
710 | long int ca_nsteps; |
---|
711 | |
---|
712 | /* Last index used in CVAfindIndex */ |
---|
713 | long int ca_ilast; |
---|
714 | |
---|
715 | /* Storage for data from forward runs */ |
---|
716 | struct DtpntMemRec **dt_mem; |
---|
717 | |
---|
718 | /* Actual number of data points in dt_mem (typically np=nsteps+1) */ |
---|
719 | long int ca_np; |
---|
720 | |
---|
721 | /* Interpolation type */ |
---|
722 | int ca_IMtype; |
---|
723 | |
---|
724 | /* Functions set by the interpolation module */ |
---|
725 | cvaIMMallocFn ca_IMmalloc; |
---|
726 | cvaIMFreeFn ca_IMfree; |
---|
727 | cvaIMStorePntFn ca_IMstore; /* store a new interpolation point */ |
---|
728 | cvaIMGetYFn ca_IMget; /* interpolate forward solution */ |
---|
729 | |
---|
730 | /* Flags controlling the interpolation module */ |
---|
731 | booleantype ca_IMmallocDone; /* IM initialized? */ |
---|
732 | booleantype ca_IMnewData; /* new data available in dt_mem?*/ |
---|
733 | booleantype ca_IMstoreSensi; /* store sensitivities? */ |
---|
734 | booleantype ca_IMinterpSensi; /* interpolate sensitivities? */ |
---|
735 | |
---|
736 | /* Workspace for the interpolation module */ |
---|
737 | N_Vector ca_Y[L_MAX]; /* pointers to zn[i] */ |
---|
738 | N_Vector *ca_YS[L_MAX]; /* pointers to znS[i] */ |
---|
739 | realtype ca_T[L_MAX]; |
---|
740 | |
---|
741 | /* ------------------------------- |
---|
742 | * Workspace for wrapper functions |
---|
743 | * ------------------------------- */ |
---|
744 | |
---|
745 | N_Vector ca_ytmp; |
---|
746 | |
---|
747 | N_Vector *ca_yStmp; |
---|
748 | |
---|
749 | }; |
---|
750 | |
---|
751 | |
---|
752 | /* |
---|
753 | * ================================================================= |
---|
754 | * I N T E R F A C E T O L I N E A R S O L V E R S |
---|
755 | * ================================================================= |
---|
756 | */ |
---|
757 | |
---|
758 | /* |
---|
759 | * ----------------------------------------------------------------- |
---|
760 | * Communication between CVODE and a CVODE Linear Solver |
---|
761 | * ----------------------------------------------------------------- |
---|
762 | * convfail (input to cv_lsetup) |
---|
763 | * |
---|
764 | * CV_NO_FAILURES : Either this is the first cv_setup call for this |
---|
765 | * step, or the local error test failed on the |
---|
766 | * previous attempt at this step (but the Newton |
---|
767 | * iteration converged). |
---|
768 | * |
---|
769 | * CV_FAIL_BAD_J : This value is passed to cv_lsetup if |
---|
770 | * |
---|
771 | * (a) The previous Newton corrector iteration |
---|
772 | * did not converge and the linear solver's |
---|
773 | * setup routine indicated that its Jacobian- |
---|
774 | * related data is not current |
---|
775 | * or |
---|
776 | * (b) During the previous Newton corrector |
---|
777 | * iteration, the linear solver's solve routine |
---|
778 | * failed in a recoverable manner and the |
---|
779 | * linear solver's setup routine indicated that |
---|
780 | * its Jacobian-related data is not current. |
---|
781 | * |
---|
782 | * CV_FAIL_OTHER : During the current internal step try, the |
---|
783 | * previous Newton iteration failed to converge |
---|
784 | * even though the linear solver was using current |
---|
785 | * Jacobian-related data. |
---|
786 | * ----------------------------------------------------------------- |
---|
787 | */ |
---|
788 | |
---|
789 | /* Constants for convfail (input to cv_lsetup) */ |
---|
790 | |
---|
791 | #define CV_NO_FAILURES 0 |
---|
792 | #define CV_FAIL_BAD_J 1 |
---|
793 | #define CV_FAIL_OTHER 2 |
---|
794 | |
---|
795 | /* |
---|
796 | * ----------------------------------------------------------------- |
---|
797 | * int (*cv_linit)(CVodeMem cv_mem); |
---|
798 | * ----------------------------------------------------------------- |
---|
799 | * The purpose of cv_linit is to complete initializations for a |
---|
800 | * specific linear solver, such as counters and statistics. |
---|
801 | * An LInitFn should return 0 if it has successfully initialized the |
---|
802 | * CVODE linear solver and a negative value otherwise. |
---|
803 | * If an error does occur, an appropriate message should be sent to |
---|
804 | * the error handler function. |
---|
805 | * ----------------------------------------------------------------- |
---|
806 | */ |
---|
807 | |
---|
808 | /* |
---|
809 | * ----------------------------------------------------------------- |
---|
810 | * int (*cv_lsetup)(CVodeMem cv_mem, int convfail, N_Vector ypred, |
---|
811 | * N_Vector fpred, booleantype *jcurPtr, |
---|
812 | * N_Vector vtemp1, N_Vector vtemp2, |
---|
813 | * N_Vector vtemp3); |
---|
814 | * ----------------------------------------------------------------- |
---|
815 | * The job of cv_lsetup is to prepare the linear solver for |
---|
816 | * subsequent calls to cv_lsolve. It may recompute Jacobian- |
---|
817 | * related data is it deems necessary. Its parameters are as |
---|
818 | * follows: |
---|
819 | * |
---|
820 | * cv_mem - problem memory pointer of type CVodeMem. See the |
---|
821 | * typedef earlier in this file. |
---|
822 | * |
---|
823 | * convfail - a flag to indicate any problem that occurred during |
---|
824 | * the solution of the nonlinear equation on the |
---|
825 | * current time step for which the linear solver is |
---|
826 | * being used. This flag can be used to help decide |
---|
827 | * whether the Jacobian data kept by a CVODE linear |
---|
828 | * solver needs to be updated or not. |
---|
829 | * Its possible values have been documented above. |
---|
830 | * |
---|
831 | * ypred - the predicted y vector for the current CVODE internal |
---|
832 | * step. |
---|
833 | * |
---|
834 | * fpred - f(tn, ypred). |
---|
835 | * |
---|
836 | * jcurPtr - a pointer to a boolean to be filled in by cv_lsetup. |
---|
837 | * The function should set *jcurPtr=SUNTRUE if its Jacobian |
---|
838 | * data is current after the call and should set |
---|
839 | * *jcurPtr=SUNFALSE if its Jacobian data is not current. |
---|
840 | * Note: If cv_lsetup calls for re-evaluation of |
---|
841 | * Jacobian data (based on convfail and CVODE state |
---|
842 | * data), it should return *jcurPtr=SUNTRUE always; |
---|
843 | * otherwise an infinite loop can result. |
---|
844 | * |
---|
845 | * vtemp1 - temporary N_Vector provided for use by cv_lsetup. |
---|
846 | * |
---|
847 | * vtemp3 - temporary N_Vector provided for use by cv_lsetup. |
---|
848 | * |
---|
849 | * vtemp3 - temporary N_Vector provided for use by cv_lsetup. |
---|
850 | * |
---|
851 | * The cv_lsetup routine should return 0 if successful, a positive |
---|
852 | * value for a recoverable error, and a negative value for an |
---|
853 | * unrecoverable error. |
---|
854 | * ----------------------------------------------------------------- |
---|
855 | */ |
---|
856 | |
---|
857 | /* |
---|
858 | * ----------------------------------------------------------------- |
---|
859 | * int (*cv_lsolve)(CVodeMem cv_mem, N_Vector b, N_Vector weight, |
---|
860 | * N_Vector ycur, N_Vector fcur); |
---|
861 | * ----------------------------------------------------------------- |
---|
862 | * cv_lsolve must solve the linear equation P x = b, where |
---|
863 | * P is some approximation to (I - gamma J), J = (df/dy)(tn,ycur) |
---|
864 | * and the RHS vector b is input. The N-vector ycur contains |
---|
865 | * the solver's current approximation to y(tn) and the vector |
---|
866 | * fcur contains the N_Vector f(tn,ycur). The solution is to be |
---|
867 | * returned in the vector b. cv_lsolve returns a positive value |
---|
868 | * for a recoverable error and a negative value for an |
---|
869 | * unrecoverable error. Success is indicated by a 0 return value. |
---|
870 | * ----------------------------------------------------------------- |
---|
871 | */ |
---|
872 | |
---|
873 | /* |
---|
874 | * ----------------------------------------------------------------- |
---|
875 | * int (*cv_lfree)(CVodeMem cv_mem); |
---|
876 | * ----------------------------------------------------------------- |
---|
877 | * cv_lfree should free up any memory allocated by the linear |
---|
878 | * solver. This routine is called once a problem has been |
---|
879 | * completed and the linear solver is no longer needed. It should |
---|
880 | * return 0 upon success, nonzero on failure. |
---|
881 | * ----------------------------------------------------------------- |
---|
882 | */ |
---|
883 | |
---|
884 | /* |
---|
885 | * ================================================================= |
---|
886 | * C V O D E S I N T E R N A L F U N C T I O N S |
---|
887 | * ================================================================= |
---|
888 | */ |
---|
889 | |
---|
890 | /* Prototype of internal ewtSet function */ |
---|
891 | |
---|
892 | int cvEwtSet(N_Vector ycur, N_Vector weight, void *data); |
---|
893 | |
---|
894 | /* High level error handler */ |
---|
895 | |
---|
896 | void cvProcessError(CVodeMem cv_mem, |
---|
897 | int error_code, const char *module, const char *fname, |
---|
898 | const char *msgfmt, ...); |
---|
899 | |
---|
900 | /* Prototype of internal errHandler function */ |
---|
901 | |
---|
902 | void cvErrHandler(int error_code, const char *module, const char *function, |
---|
903 | char *msg, void *data); |
---|
904 | |
---|
905 | /* Prototypes for internal sensitivity rhs wrappers */ |
---|
906 | |
---|
907 | int cvSensRhsWrapper(CVodeMem cv_mem, realtype time, |
---|
908 | N_Vector ycur, N_Vector fcur, |
---|
909 | N_Vector *yScur, N_Vector *fScur, |
---|
910 | N_Vector temp1, N_Vector temp2); |
---|
911 | |
---|
912 | int cvSensRhs1Wrapper(CVodeMem cv_mem, realtype time, |
---|
913 | N_Vector ycur, N_Vector fcur, |
---|
914 | int is, N_Vector yScur, N_Vector fScur, |
---|
915 | N_Vector temp1, N_Vector temp2); |
---|
916 | |
---|
917 | /* Prototypes for internal sensitivity rhs DQ functions */ |
---|
918 | |
---|
919 | int cvSensRhsInternalDQ(int Ns, realtype t, |
---|
920 | N_Vector y, N_Vector ydot, |
---|
921 | N_Vector *yS, N_Vector *ySdot, |
---|
922 | void *fS_data, |
---|
923 | N_Vector tempv, N_Vector ftemp); |
---|
924 | |
---|
925 | int cvSensRhs1InternalDQ(int Ns, realtype t, |
---|
926 | N_Vector y, N_Vector ydot, |
---|
927 | int is, N_Vector yS, N_Vector ySdot, |
---|
928 | void *fS_data, |
---|
929 | N_Vector tempv, N_Vector ftemp); |
---|
930 | |
---|
931 | /* |
---|
932 | * ================================================================= |
---|
933 | * C V O D E S E R R O R M E S S A G E S |
---|
934 | * ================================================================= |
---|
935 | */ |
---|
936 | |
---|
937 | #if defined(SUNDIALS_EXTENDED_PRECISION) |
---|
938 | |
---|
939 | #define MSG_TIME "t = %Lg" |
---|
940 | #define MSG_TIME_H "t = %Lg and h = %Lg" |
---|
941 | #define MSG_TIME_INT "t = %Lg is not between tcur - hu = %Lg and tcur = %Lg." |
---|
942 | #define MSG_TIME_TOUT "tout = %Lg" |
---|
943 | #define MSG_TIME_TSTOP "tstop = %Lg" |
---|
944 | |
---|
945 | #elif defined(SUNDIALS_DOUBLE_PRECISION) |
---|
946 | |
---|
947 | #define MSG_TIME "t = %lg" |
---|
948 | #define MSG_TIME_H "t = %lg and h = %lg" |
---|
949 | #define MSG_TIME_INT "t = %lg is not between tcur - hu = %lg and tcur = %lg." |
---|
950 | #define MSG_TIME_TOUT "tout = %lg" |
---|
951 | #define MSG_TIME_TSTOP "tstop = %lg" |
---|
952 | |
---|
953 | #else |
---|
954 | |
---|
955 | #define MSG_TIME "t = %g" |
---|
956 | #define MSG_TIME_H "t = %g and h = %g" |
---|
957 | #define MSG_TIME_INT "t = %g is not between tcur - hu = %g and tcur = %g." |
---|
958 | #define MSG_TIME_TOUT "tout = %g" |
---|
959 | #define MSG_TIME_TSTOP "tstop = %g" |
---|
960 | |
---|
961 | #endif |
---|
962 | |
---|
963 | |
---|
964 | /* Initialization and I/O error messages */ |
---|
965 | |
---|
966 | #define MSGCV_NO_MEM "cvode_mem = NULL illegal." |
---|
967 | #define MSGCV_CVMEM_FAIL "Allocation of cvode_mem failed." |
---|
968 | #define MSGCV_MEM_FAIL "A memory request failed." |
---|
969 | #define MSGCV_BAD_LMM "Illegal value for lmm. The legal values are CV_ADAMS and CV_BDF." |
---|
970 | #define MSGCV_BAD_ITER "Illegal value for iter. The legal values are CV_FUNCTIONAL and CV_NEWTON." |
---|
971 | #define MSGCV_NO_MALLOC "Attempt to call before CVodeInit." |
---|
972 | #define MSGCV_NEG_MAXORD "maxord <= 0 illegal." |
---|
973 | #define MSGCV_BAD_MAXORD "Illegal attempt to increase maximum method order." |
---|
974 | #define MSGCV_SET_SLDET "Attempt to use stability limit detection with the CV_ADAMS method illegal." |
---|
975 | #define MSGCV_NEG_HMIN "hmin < 0 illegal." |
---|
976 | #define MSGCV_NEG_HMAX "hmax < 0 illegal." |
---|
977 | #define MSGCV_BAD_HMIN_HMAX "Inconsistent step size limits: hmin > hmax." |
---|
978 | #define MSGCV_BAD_RELTOL "reltol < 0 illegal." |
---|
979 | #define MSGCV_BAD_ABSTOL "abstol has negative component(s) (illegal)." |
---|
980 | #define MSGCV_NULL_ABSTOL "abstol = NULL illegal." |
---|
981 | #define MSGCV_NULL_Y0 "y0 = NULL illegal." |
---|
982 | #define MSGCV_Y0_FAIL_CONSTR "y0 fails to satisfy constraints." |
---|
983 | #define MSGCV_BAD_ISM_CONSTR "Constraints can not be enforced while forward sensitivity is used with simultaneous method" |
---|
984 | #define MSGCV_NULL_F "f = NULL illegal." |
---|
985 | #define MSGCV_NULL_G "g = NULL illegal." |
---|
986 | #define MSGCV_BAD_NVECTOR "A required vector operation is not implemented." |
---|
987 | #define MSGCV_BAD_CONSTR "Illegal values in constraints vector." |
---|
988 | #define MSGCV_BAD_K "Illegal value for k." |
---|
989 | #define MSGCV_NULL_DKY "dky = NULL illegal." |
---|
990 | #define MSGCV_BAD_T "Illegal value for t." MSG_TIME_INT |
---|
991 | #define MSGCV_NO_ROOT "Rootfinding was not initialized." |
---|
992 | |
---|
993 | #define MSGCV_NO_QUAD "Quadrature integration not activated." |
---|
994 | #define MSGCV_BAD_ITOLQ "Illegal value for itolQ. The legal values are CV_SS and CV_SV." |
---|
995 | #define MSGCV_NULL_ABSTOLQ "abstolQ = NULL illegal." |
---|
996 | #define MSGCV_BAD_RELTOLQ "reltolQ < 0 illegal." |
---|
997 | #define MSGCV_BAD_ABSTOLQ "abstolQ has negative component(s) (illegal)." |
---|
998 | |
---|
999 | #define MSGCV_SENSINIT_2 "Sensitivity analysis already initialized." |
---|
1000 | #define MSGCV_NO_SENSI "Forward sensitivity analysis not activated." |
---|
1001 | #define MSGCV_BAD_ITOLS "Illegal value for itolS. The legal values are CV_SS, CV_SV, and CV_EE." |
---|
1002 | #define MSGCV_NULL_ABSTOLS "abstolS = NULL illegal." |
---|
1003 | #define MSGCV_BAD_RELTOLS "reltolS < 0 illegal." |
---|
1004 | #define MSGCV_BAD_ABSTOLS "abstolS has negative component(s) (illegal)." |
---|
1005 | #define MSGCV_BAD_PBAR "pbar has zero component(s) (illegal)." |
---|
1006 | #define MSGCV_BAD_PLIST "plist has negative component(s) (illegal)." |
---|
1007 | #define MSGCV_BAD_NS "NS <= 0 illegal." |
---|
1008 | #define MSGCV_NULL_YS0 "yS0 = NULL illegal." |
---|
1009 | #define MSGCV_BAD_ISM "Illegal value for ism. Legal values are: CV_SIMULTANEOUS, CV_STAGGERED and CV_STAGGERED1." |
---|
1010 | #define MSGCV_BAD_IFS "Illegal value for ifS. Legal values are: CV_ALLSENS and CV_ONESENS." |
---|
1011 | #define MSGCV_BAD_ISM_IFS "Illegal ism = CV_STAGGERED1 for CVodeSensInit." |
---|
1012 | #define MSGCV_BAD_IS "Illegal value for is." |
---|
1013 | #define MSGCV_NULL_DKYA "dkyA = NULL illegal." |
---|
1014 | #define MSGCV_BAD_DQTYPE "Illegal value for DQtype. Legal values are: CV_CENTERED and CV_FORWARD." |
---|
1015 | #define MSGCV_BAD_DQRHO "DQrhomax < 0 illegal." |
---|
1016 | |
---|
1017 | #define MSGCV_BAD_ITOLQS "Illegal value for itolQS. The legal values are CV_SS, CV_SV, and CV_EE." |
---|
1018 | #define MSGCV_NULL_ABSTOLQS "abstolQS = NULL illegal." |
---|
1019 | #define MSGCV_BAD_RELTOLQS "reltolQS < 0 illegal." |
---|
1020 | #define MSGCV_BAD_ABSTOLQS "abstolQS has negative component(s) (illegal)." |
---|
1021 | #define MSGCV_NO_QUADSENSI "Forward sensitivity analysis for quadrature variables not activated." |
---|
1022 | #define MSGCV_NULL_YQS0 "yQS0 = NULL illegal." |
---|
1023 | |
---|
1024 | /* CVode Error Messages */ |
---|
1025 | |
---|
1026 | #define MSGCV_NO_TOL "No integration tolerances have been specified." |
---|
1027 | #define MSGCV_LSOLVE_NULL "The linear solver's solve routine is NULL." |
---|
1028 | #define MSGCV_YOUT_NULL "yout = NULL illegal." |
---|
1029 | #define MSGCV_TRET_NULL "tret = NULL illegal." |
---|
1030 | #define MSGCV_BAD_EWT "Initial ewt has component(s) equal to zero (illegal)." |
---|
1031 | #define MSGCV_EWT_NOW_BAD "At " MSG_TIME ", a component of ewt has become <= 0." |
---|
1032 | #define MSGCV_BAD_ITASK "Illegal value for itask." |
---|
1033 | #define MSGCV_BAD_H0 "h0 and tout - t0 inconsistent." |
---|
1034 | #define MSGCV_BAD_TOUT "Trouble interpolating at " MSG_TIME_TOUT ". tout too far back in direction of integration" |
---|
1035 | #define MSGCV_EWT_FAIL "The user-provide EwtSet function failed." |
---|
1036 | #define MSGCV_EWT_NOW_FAIL "At " MSG_TIME ", the user-provide EwtSet function failed." |
---|
1037 | #define MSGCV_LINIT_FAIL "The linear solver's init routine failed." |
---|
1038 | #define MSGCV_HNIL_DONE "The above warning has been issued mxhnil times and will not be issued again for this problem." |
---|
1039 | #define MSGCV_TOO_CLOSE "tout too close to t0 to start integration." |
---|
1040 | #define MSGCV_MAX_STEPS "At " MSG_TIME ", mxstep steps taken before reaching tout." |
---|
1041 | #define MSGCV_TOO_MUCH_ACC "At " MSG_TIME ", too much accuracy requested." |
---|
1042 | #define MSGCV_HNIL "Internal " MSG_TIME_H " are such that t + h = t on the next step. The solver will continue anyway." |
---|
1043 | #define MSGCV_ERR_FAILS "At " MSG_TIME_H ", the error test failed repeatedly or with |h| = hmin." |
---|
1044 | #define MSGCV_CONV_FAILS "At " MSG_TIME_H ", the corrector convergence test failed repeatedly or with |h| = hmin." |
---|
1045 | #define MSGCV_SETUP_FAILED "At " MSG_TIME ", the setup routine failed in an unrecoverable manner." |
---|
1046 | #define MSGCV_SOLVE_FAILED "At " MSG_TIME ", the solve routine failed in an unrecoverable manner." |
---|
1047 | #define MSGCV_FAILED_CONSTR "At " MSG_TIME ", unable to satisfy inequality constraints." |
---|
1048 | #define MSGCV_RHSFUNC_FAILED "At " MSG_TIME ", the right-hand side routine failed in an unrecoverable manner." |
---|
1049 | #define MSGCV_RHSFUNC_UNREC "At " MSG_TIME ", the right-hand side failed in a recoverable manner, but no recovery is possible." |
---|
1050 | #define MSGCV_RHSFUNC_REPTD "At " MSG_TIME " repeated recoverable right-hand side function errors." |
---|
1051 | #define MSGCV_RHSFUNC_FIRST "The right-hand side routine failed at the first call." |
---|
1052 | #define MSGCV_RTFUNC_FAILED "At " MSG_TIME ", the rootfinding routine failed in an unrecoverable manner." |
---|
1053 | #define MSGCV_CLOSE_ROOTS "Root found at and very near " MSG_TIME "." |
---|
1054 | #define MSGCV_BAD_TSTOP "The value " MSG_TIME_TSTOP " is behind current " MSG_TIME " in the direction of integration." |
---|
1055 | #define MSGCV_INACTIVE_ROOTS "At the end of the first step, there are still some root functions identically 0. This warning will not be issued again." |
---|
1056 | |
---|
1057 | #define MSGCV_NO_TOLQ "No integration tolerances for quadrature variables have been specified." |
---|
1058 | #define MSGCV_BAD_EWTQ "Initial ewtQ has component(s) equal to zero (illegal)." |
---|
1059 | #define MSGCV_EWTQ_NOW_BAD "At " MSG_TIME ", a component of ewtQ has become <= 0." |
---|
1060 | #define MSGCV_QRHSFUNC_FAILED "At " MSG_TIME ", the quadrature right-hand side routine failed in an unrecoverable manner." |
---|
1061 | #define MSGCV_QRHSFUNC_UNREC "At " MSG_TIME ", the quadrature right-hand side failed in a recoverable manner, but no recovery is possible." |
---|
1062 | #define MSGCV_QRHSFUNC_REPTD "At " MSG_TIME " repeated recoverable quadrature right-hand side function errors." |
---|
1063 | #define MSGCV_QRHSFUNC_FIRST "The quadrature right-hand side routine failed at the first call." |
---|
1064 | |
---|
1065 | #define MSGCV_NO_TOLS "No integration tolerances for sensitivity variables have been specified." |
---|
1066 | #define MSGCV_NULL_P "p = NULL when using internal DQ for sensitivity RHS illegal." |
---|
1067 | #define MSGCV_BAD_EWTS "Initial ewtS has component(s) equal to zero (illegal)." |
---|
1068 | #define MSGCV_EWTS_NOW_BAD "At " MSG_TIME ", a component of ewtS has become <= 0." |
---|
1069 | #define MSGCV_SRHSFUNC_FAILED "At " MSG_TIME ", the sensitivity right-hand side routine failed in an unrecoverable manner." |
---|
1070 | #define MSGCV_SRHSFUNC_UNREC "At " MSG_TIME ", the sensitivity right-hand side failed in a recoverable manner, but no recovery is possible." |
---|
1071 | #define MSGCV_SRHSFUNC_REPTD "At " MSG_TIME " repeated recoverable sensitivity right-hand side function errors." |
---|
1072 | #define MSGCV_SRHSFUNC_FIRST "The sensitivity right-hand side routine failed at the first call." |
---|
1073 | |
---|
1074 | #define MSGCV_NULL_FQ "CVODES is expected to use DQ to evaluate the RHS of quad. sensi., but quadratures were not initialized." |
---|
1075 | #define MSGCV_NO_TOLQS "No integration tolerances for quadrature sensitivity variables have been specified." |
---|
1076 | #define MSGCV_BAD_EWTQS "Initial ewtQS has component(s) equal to zero (illegal)." |
---|
1077 | #define MSGCV_EWTQS_NOW_BAD "At " MSG_TIME ", a component of ewtQS has become <= 0." |
---|
1078 | #define MSGCV_QSRHSFUNC_FAILED "At " MSG_TIME ", the quadrature sensitivity right-hand side routine failed in an unrecoverable manner." |
---|
1079 | #define MSGCV_QSRHSFUNC_UNREC "At " MSG_TIME ", the quadrature sensitivity right-hand side failed in a recoverable manner, but no recovery is possible." |
---|
1080 | #define MSGCV_QSRHSFUNC_REPTD "At " MSG_TIME " repeated recoverable quadrature sensitivity right-hand side function errors." |
---|
1081 | #define MSGCV_QSRHSFUNC_FIRST "The quadrature sensitivity right-hand side routine failed at the first call." |
---|
1082 | |
---|
1083 | /* |
---|
1084 | * ================================================================= |
---|
1085 | * C V O D E A E R R O R M E S S A G E S |
---|
1086 | * ================================================================= |
---|
1087 | */ |
---|
1088 | |
---|
1089 | #define MSGCV_NO_ADJ "Illegal attempt to call before calling CVodeAdjMalloc." |
---|
1090 | #define MSGCV_BAD_STEPS "Steps nonpositive illegal." |
---|
1091 | #define MSGCV_BAD_INTERP "Illegal value for interp." |
---|
1092 | #define MSGCV_BAD_WHICH "Illegal value for which." |
---|
1093 | #define MSGCV_NO_BCK "No backward problems have been defined yet." |
---|
1094 | #define MSGCV_NO_FWD "Illegal attempt to call before calling CVodeF." |
---|
1095 | #define MSGCV_BAD_TB0 "The initial time tB0 for problem %d is outside the interval over which the forward problem was solved." |
---|
1096 | #define MSGCV_BAD_SENSI "At least one backward problem requires sensitivities, but they were not stored for interpolation." |
---|
1097 | #define MSGCV_BAD_ITASKB "Illegal value for itaskB. Legal values are CV_NORMAL and CV_ONE_STEP." |
---|
1098 | #define MSGCV_BAD_TBOUT "The final time tBout is outside the interval over which the forward problem was solved." |
---|
1099 | #define MSGCV_BACK_ERROR "Error occured while integrating backward problem # %d" |
---|
1100 | #define MSGCV_BAD_TINTERP "Bad t = %g for interpolation." |
---|
1101 | #define MSGCV_WRONG_INTERP "This function cannot be called for the specified interp type." |
---|
1102 | |
---|
1103 | #ifdef __cplusplus |
---|
1104 | } |
---|
1105 | #endif |
---|
1106 | |
---|
1107 | #endif |
---|