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