1 /* 2 * ----------------------------------------------------------------- 3 * Programmer(s): Daniel Reynolds @ SMU 4 * Scott D. Cohen, Alan C. Hindmarsh, Radu Serban 5 * and Dan Shumaker @ LLNL 6 * ----------------------------------------------------------------- 7 * SUNDIALS Copyright Start 8 * Copyright (c) 2002-2020, Lawrence Livermore National Security 9 * and Southern Methodist University. 10 * All rights reserved. 11 * 12 * See the top-level LICENSE and NOTICE files for details. 13 * 14 * SPDX-License-Identifier: BSD-3-Clause 15 * SUNDIALS Copyright End 16 * ----------------------------------------------------------------- 17 * Implementation header file for the main CVODE integrator. 18 * ----------------------------------------------------------------- 19 */ 20 21 #ifndef _CVODE_IMPL_H 22 #define _CVODE_IMPL_H 23 24 #include <stdarg.h> 25 #include <cvode/cvode.h> 26 27 #include "cvode_proj_impl.h" 28 29 #ifdef __cplusplus /* wrapper to enable C++ usage */ 30 extern "C" { 31 #endif 32 33 /* 34 * ================================================================= 35 * 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 36 * ================================================================= 37 */ 38 39 /* Basic CVODE constants */ 40 41 #define ADAMS_Q_MAX 12 /* max value of q for lmm == ADAMS */ 42 #define BDF_Q_MAX 5 /* max value of q for lmm == BDF */ 43 #define Q_MAX ADAMS_Q_MAX /* max value of q for either lmm */ 44 #define L_MAX (Q_MAX+1) /* max value of L for either lmm */ 45 #define NUM_TESTS 5 /* number of error test quantities */ 46 47 #define HMIN_DEFAULT RCONST(0.0) /* hmin default value */ 48 #define HMAX_INV_DEFAULT RCONST(0.0) /* hmax_inv default value */ 49 #define MXHNIL_DEFAULT 10 /* mxhnil default value */ 50 #define MXSTEP_DEFAULT 500 /* mxstep default value */ 51 52 /* Control constants for lower-level functions used by cvStep 53 * ---------------------------------------------------------- 54 * 55 * cvHin return values: 56 * CV_SUCCESS 57 * CV_RHSFUNC_FAIL 58 * CV_TOO_CLOSE 59 * 60 * cvStep control constants: 61 * DO_ERROR_TEST 62 * PREDICT_AGAIN 63 * 64 * cvStep return values: 65 * CV_SUCCESS, 66 * CV_LSETUP_FAIL, CV_LSOLVE_FAIL, 67 * CV_RHSFUNC_FAIL, CV_RTFUNC_FAIL 68 * CV_CONV_FAILURE, CV_ERR_FAILURE, 69 * CV_FIRST_RHSFUNC_ERR 70 * 71 * cvNls input nflag values: 72 * FIRST_CALL 73 * PREV_CONV_FAIL 74 * PREV_PROJ_FAIL 75 * PREV_ERR_FAIL 76 * 77 * cvNls return values: 78 * CV_SUCCESS, 79 * CV_LSETUP_FAIL, CV_LSOLVE_FAIL, CV_RHSFUNC_FAIL, 80 * CONV_FAIL, RHSFUNC_RECVR 81 * 82 * cvNewtonIteration return values: 83 * CV_SUCCESS, 84 * CV_LSOLVE_FAIL, CV_RHSFUNC_FAIL 85 * CONV_FAIL, RHSFUNC_RECVR, 86 * TRY_AGAIN 87 * 88 */ 89 90 #define DO_ERROR_TEST +2 91 #define PREDICT_AGAIN +3 92 93 #define TRY_AGAIN +5 94 #define FIRST_CALL +6 95 #define PREV_CONV_FAIL +7 96 #define PREV_PROJ_FAIL +8 97 #define PREV_ERR_FAIL +9 98 99 #define RHSFUNC_RECVR +10 100 #define CONSTR_RECVR +11 101 #define CONSTRFUNC_RECVR +12 102 #define PROJFUNC_RECVR +13 103 104 /* 105 * ----------------------------------------------------------------- 106 * Types : struct CVodeMemRec, CVodeMem 107 * ----------------------------------------------------------------- 108 * The type CVodeMem is type pointer to struct CVodeMemRec. 109 * This structure contains fields to keep track of problem state. 110 * ----------------------------------------------------------------- 111 */ 112 113 typedef struct CVodeMemRec { 114 115 realtype cv_uround; /* machine unit roundoff */ 116 117 /*-------------------------- 118 Problem Specification Data 119 --------------------------*/ 120 121 CVRhsFn cv_f; /* y' = f(t,y(t)) */ 122 void *cv_user_data; /* user pointer passed to f */ 123 int cv_lmm; /* lmm = CV_ADAMS or CV_BDF */ 124 int cv_itol; /* itol = CV_SS, CV_SV, CV_WF, CV_NN */ 125 126 realtype cv_reltol; /* relative tolerance */ 127 realtype cv_Sabstol; /* scalar absolute tolerance */ 128 N_Vector cv_Vabstol; /* vector absolute tolerance */ 129 booleantype cv_atolmin0; /* flag indicating that min(abstol) = 0 */ 130 booleantype cv_user_efun; /* SUNTRUE if user sets efun */ 131 CVEwtFn cv_efun; /* function to set ewt */ 132 void *cv_e_data; /* user pointer passed to efun */ 133 134 booleantype cv_constraintsSet; /* constraints vector present: 135 do constraints calc */ 136 137 /*----------------------- 138 Nordsieck History Array 139 -----------------------*/ 140 141 N_Vector cv_zn[L_MAX]; /* Nordsieck array, of size N x (q+1). 142 zn[j] is a vector of length N (j=0,...,q) 143 zn[j] = [1/factorial(j)] * h^j * (jth 144 derivative of the interpolating polynomial */ 145 146 /*-------------------------- 147 other vectors of length N 148 -------------------------*/ 149 150 N_Vector cv_ewt; /* error weight vector */ 151 N_Vector cv_y; /* y is used as temporary storage by the solver 152 The memory is provided by the user to CVode 153 where the vector is named yout. */ 154 N_Vector cv_acor; /* In the context of the solution of the nonlinear 155 equation, acor = y_n(m) - y_n(0). On return, 156 this vector is scaled to give the est. local err. */ 157 N_Vector cv_tempv; /* temporary storage vector */ 158 N_Vector cv_ftemp; /* temporary storage vector */ 159 N_Vector cv_vtemp1; /* temporary storage vector */ 160 N_Vector cv_vtemp2; /* temporary storage vector */ 161 N_Vector cv_vtemp3; /* temporary storage vector */ 162 163 N_Vector cv_constraints; /* vector of inequality constraint options */ 164 165 /*----------------- 166 Tstop information 167 -----------------*/ 168 169 booleantype cv_tstopset; 170 realtype cv_tstop; 171 172 /*--------- 173 Step Data 174 ---------*/ 175 176 int cv_q; /* current order */ 177 int cv_qprime; /* order to be used on the next step 178 = q-1, q, or q+1 */ 179 int cv_next_q; /* order to be used on the next step */ 180 int cv_qwait; /* number of internal steps to wait before 181 considering a change in q */ 182 int cv_L; /* L = q + 1 */ 183 184 realtype cv_hin; /* initial step size */ 185 realtype cv_h; /* current step size */ 186 realtype cv_hprime; /* step size to be used on the next step */ 187 realtype cv_next_h; /* step size to be used on the next step */ 188 realtype cv_eta; /* eta = hprime / h */ 189 realtype cv_hscale; /* value of h used in zn */ 190 realtype cv_tn; /* current internal value of t */ 191 realtype cv_tretlast; /* value of tret last returned by CVode */ 192 193 realtype cv_tau[L_MAX+1]; /* array of previous q+1 successful step 194 sizes indexed from 1 to q+1 */ 195 realtype cv_tq[NUM_TESTS+1]; /* array of test quantities indexed from 196 1 to NUM_TESTS(=5) */ 197 realtype cv_l[L_MAX]; /* coefficients of l(x) (degree q poly) */ 198 199 realtype cv_rl1; /* the scalar 1/l[1] */ 200 realtype cv_gamma; /* gamma = h * rl1 */ 201 realtype cv_gammap; /* gamma at the last setup call */ 202 realtype cv_gamrat; /* gamma / gammap */ 203 204 realtype cv_crate; /* estimated corrector convergence rate */ 205 realtype cv_delp; /* norm of previous nonlinear solver update */ 206 realtype cv_acnrm; /* | acor | wrms */ 207 booleantype cv_acnrmcur; /* is | acor | wrms current? */ 208 realtype cv_nlscoef; /* coeficient in nonlinear convergence test */ 209 210 /*------ 211 Limits 212 ------*/ 213 214 int cv_qmax; /* q <= qmax */ 215 long int cv_mxstep; /* maximum number of internal steps for one user call */ 216 int cv_maxcor; /* maximum number of corrector iterations for the 217 solution of the nonlinear equation */ 218 int cv_mxhnil; /* maximum number of warning messages issued to the 219 user that t + h == t for the next internal step */ 220 int cv_maxnef; /* maximum number of error test failures */ 221 int cv_maxncf; /* maximum number of nonlinear convergence failures */ 222 223 realtype cv_hmin; /* |h| >= hmin */ 224 realtype cv_hmax_inv; /* |h| <= 1/hmax_inv */ 225 realtype cv_etamax; /* eta <= etamax */ 226 227 /*-------- 228 Counters 229 --------*/ 230 231 long int cv_nst; /* number of internal steps taken */ 232 long int cv_nfe; /* number of f calls */ 233 long int cv_ncfn; /* number of corrector convergence failures */ 234 long int cv_netf; /* number of error test failures */ 235 long int cv_nni; /* number of Newton iterations performed */ 236 long int cv_nsetups; /* number of setup calls */ 237 int cv_nhnil; /* number of messages issued to the user that 238 t + h == t for the next iternal step */ 239 240 realtype cv_etaqm1; /* ratio of new to old h for order q-1 */ 241 realtype cv_etaq; /* ratio of new to old h for order q */ 242 realtype cv_etaqp1; /* ratio of new to old h for order q+1 */ 243 244 /*---------------------------- 245 Space requirements for CVODE 246 ----------------------------*/ 247 248 sunindextype cv_lrw1; /* no. of realtype words in 1 N_Vector */ 249 sunindextype cv_liw1; /* no. of integer words in 1 N_Vector */ 250 long int cv_lrw; /* no. of realtype words in CVODE work vectors */ 251 long int cv_liw; /* no. of integer words in CVODE work vectors */ 252 253 /*--------------------- 254 Nonlinear Solver Data 255 ---------------------*/ 256 257 SUNNonlinearSolver NLS; /* Sundials generic nonlinear solver object */ 258 booleantype ownNLS; /* flag indicating if CVODE created the nonlinear 259 solver object */ 260 int convfail; /* flag to indicate when a Jacbian update may 261 be needed */ 262 263 /*------------------ 264 Linear Solver Data 265 ------------------*/ 266 267 /* Linear Solver functions to be called */ 268 269 int (*cv_linit)(struct CVodeMemRec *cv_mem); 270 271 int (*cv_lsetup)(struct CVodeMemRec *cv_mem, int convfail, N_Vector ypred, 272 N_Vector fpred, booleantype *jcurPtr, N_Vector vtemp1, 273 N_Vector vtemp2, N_Vector vtemp3); 274 275 int (*cv_lsolve)(struct CVodeMemRec *cv_mem, N_Vector b, N_Vector weight, 276 N_Vector ycur, N_Vector fcur); 277 278 int (*cv_lfree)(struct CVodeMemRec *cv_mem); 279 280 /* Linear Solver specific memory */ 281 282 void *cv_lmem; 283 284 /*------------ 285 Saved Values 286 ------------*/ 287 288 int cv_qu; /* last successful q value used */ 289 long int cv_nstlp; /* step number of last setup call */ 290 realtype cv_h0u; /* actual initial stepsize */ 291 realtype cv_hu; /* last successful h value used */ 292 realtype cv_saved_tq5; /* saved value of tq[5] */ 293 booleantype cv_jcur; /* is Jacobian info. for lin. solver current? */ 294 realtype cv_tolsf; /* tolerance scale factor */ 295 int cv_qmax_alloc; /* value of qmax used when allocating memory */ 296 int cv_indx_acor; /* index of the zn vector with saved acor */ 297 298 booleantype cv_VabstolMallocDone; 299 booleantype cv_MallocDone; 300 booleantype cv_constraintsMallocDone; 301 302 /*------------------------------------------- 303 Error handler function and error ouput file 304 -------------------------------------------*/ 305 306 CVErrHandlerFn cv_ehfun; /* error messages are handled by ehfun */ 307 void *cv_eh_data; /* data pointer passed to ehfun */ 308 FILE *cv_errfp; /* CVODE error messages are sent to errfp */ 309 310 /*------------------------------------------- 311 User access function 312 -------------------------------------------*/ 313 CVMonitorFn cv_monitorfun; /* func called with CVODE mem and user data */ 314 long int cv_monitor_interval; /* step interval to call cv_monitorfun */ 315 316 /*------------------------- 317 Stability Limit Detection 318 -------------------------*/ 319 320 booleantype cv_sldeton; /* is Stability Limit Detection on? */ 321 realtype cv_ssdat[6][4]; /* scaled data array for STALD */ 322 int cv_nscon; /* counter for STALD method */ 323 long int cv_nor; /* counter for number of order reductions */ 324 325 /*---------------- 326 Rootfinding Data 327 ----------------*/ 328 329 CVRootFn cv_gfun; /* function g for roots sought */ 330 int cv_nrtfn; /* number of components of g */ 331 int *cv_iroots; /* array for root information */ 332 int *cv_rootdir; /* array specifying direction of zero-crossing */ 333 realtype cv_tlo; /* nearest endpoint of interval in root search */ 334 realtype cv_thi; /* farthest endpoint of interval in root search */ 335 realtype cv_trout; /* t value returned by rootfinding routine */ 336 realtype *cv_glo; /* saved array of g values at t = tlo */ 337 realtype *cv_ghi; /* saved array of g values at t = thi */ 338 realtype *cv_grout; /* array of g values at t = trout */ 339 realtype cv_toutc; /* copy of tout (if NORMAL mode) */ 340 realtype cv_ttol; /* tolerance on root location */ 341 int cv_taskc; /* copy of parameter itask */ 342 int cv_irfnd; /* flag showing whether last step had a root */ 343 long int cv_nge; /* counter for g evaluations */ 344 booleantype *cv_gactive; /* array with active/inactive event functions */ 345 int cv_mxgnull; /* number of warning messages about possible g==0 */ 346 347 /*---------------- 348 Projection Data 349 ----------------*/ 350 351 CVodeProjMem proj_mem; /* projection memory structure */ 352 booleantype proj_enabled; /* flag indicating if projection is enabled */ 353 booleantype proj_applied; /* flag indicating if projection was applied */ 354 realtype cv_p[L_MAX]; /* coefficients of p(x) (degree q poly) */ 355 356 /*----------------------- 357 Fused Vector Operations 358 -----------------------*/ 359 360 realtype cv_cvals[L_MAX]; /* array of scalars */ 361 N_Vector cv_Xvecs[L_MAX]; /* array of vectors */ 362 363 booleantype cv_usefused; /* flag indicating if CVODE specific fused kernels should be used */ 364 365 } *CVodeMem; 366 367 /* 368 * ================================================================= 369 * I N T E R F A C E T O L I N E A R S O L V E R S 370 * ================================================================= 371 */ 372 373 /* 374 * ----------------------------------------------------------------- 375 * Communication between CVODE and a CVODE Linear Solver 376 * ----------------------------------------------------------------- 377 * convfail (input to cv_lsetup) 378 * 379 * CV_NO_FAILURES : Either this is the first cv_setup call for this 380 * step, or the local error test failed on the 381 * previous attempt at this step (but the Newton 382 * iteration converged). 383 * 384 * CV_FAIL_BAD_J : This value is passed to cv_lsetup if 385 * 386 * (a) The previous Newton corrector iteration 387 * did not converge and the linear solver's 388 * setup routine indicated that its Jacobian- 389 * related data is not current 390 * or 391 * (b) During the previous Newton corrector 392 * iteration, the linear solver's solve routine 393 * failed in a recoverable manner and the 394 * linear solver's setup routine indicated that 395 * its Jacobian-related data is not current. 396 * 397 * CV_FAIL_OTHER : During the current internal step try, the 398 * previous Newton iteration failed to converge 399 * even though the linear solver was using current 400 * Jacobian-related data. 401 * ----------------------------------------------------------------- 402 */ 403 404 /* Constants for convfail (input to cv_lsetup) */ 405 406 #define CV_NO_FAILURES 0 407 #define CV_FAIL_BAD_J 1 408 #define CV_FAIL_OTHER 2 409 410 /* 411 * ----------------------------------------------------------------- 412 * int (*cv_linit)(CVodeMem cv_mem); 413 * ----------------------------------------------------------------- 414 * The purpose of cv_linit is to complete initializations for a 415 * specific linear solver, such as counters and statistics. 416 * An LInitFn should return 0 if it has successfully initialized the 417 * CVODE linear solver and a negative value otherwise. 418 * If an error does occur, an appropriate message should be sent to 419 * the error handler function. 420 * ----------------------------------------------------------------- 421 */ 422 423 /* 424 * ----------------------------------------------------------------- 425 * int (*cv_lsetup)(CVodeMem cv_mem, int convfail, N_Vector ypred, 426 * N_Vector fpred, booleantype *jcurPtr, 427 * N_Vector vtemp1, N_Vector vtemp2, 428 * N_Vector vtemp3); 429 * ----------------------------------------------------------------- 430 * The job of cv_lsetup is to prepare the linear solver for 431 * subsequent calls to cv_lsolve. It may recompute Jacobian- 432 * related data is it deems necessary. Its parameters are as 433 * follows: 434 * 435 * cv_mem - problem memory pointer of type CVodeMem. See the 436 * typedef earlier in this file. 437 * 438 * convfail - a flag to indicate any problem that occurred during 439 * the solution of the nonlinear equation on the 440 * current time step for which the linear solver is 441 * being used. This flag can be used to help decide 442 * whether the Jacobian data kept by a CVODE linear 443 * solver needs to be updated or not. 444 * Its possible values have been documented above. 445 * 446 * ypred - the predicted y vector for the current CVODE internal 447 * step. 448 * 449 * fpred - f(tn, ypred). 450 * 451 * jcurPtr - a pointer to a boolean to be filled in by cv_lsetup. 452 * The function should set *jcurPtr=SUNTRUE if its Jacobian 453 * data is current after the call and should set 454 * *jcurPtr=SUNFALSE if its Jacobian data is not current. 455 * Note: If cv_lsetup calls for re-evaluation of 456 * Jacobian data (based on convfail and CVODE state 457 * data), it should return *jcurPtr=SUNTRUE always; 458 * otherwise an infinite loop can result. 459 * 460 * vtemp1 - temporary N_Vector provided for use by cv_lsetup. 461 * 462 * vtemp3 - temporary N_Vector provided for use by cv_lsetup. 463 * 464 * vtemp3 - temporary N_Vector provided for use by cv_lsetup. 465 * 466 * The cv_lsetup routine should return 0 if successful, a positive 467 * value for a recoverable error, and a negative value for an 468 * unrecoverable error. 469 * ----------------------------------------------------------------- 470 */ 471 472 /* 473 * ----------------------------------------------------------------- 474 * int (*cv_lsolve)(CVodeMem cv_mem, N_Vector b, N_Vector weight, 475 * N_Vector ycur, N_Vector fcur); 476 * ----------------------------------------------------------------- 477 * cv_lsolve must solve the linear equation P x = b, where 478 * P is some approximation to (I - gamma J), J = (df/dy)(tn,ycur) 479 * and the RHS vector b is input. The N-vector ycur contains 480 * the solver's current approximation to y(tn) and the vector 481 * fcur contains the N_Vector f(tn,ycur). The solution is to be 482 * returned in the vector b. cv_lsolve returns a positive value 483 * for a recoverable error and a negative value for an 484 * unrecoverable error. Success is indicated by a 0 return value. 485 * ----------------------------------------------------------------- 486 */ 487 488 /* 489 * ----------------------------------------------------------------- 490 * int (*cv_lfree)(CVodeMem cv_mem); 491 * ----------------------------------------------------------------- 492 * cv_lfree should free up any memory allocated by the linear 493 * solver. This routine is called once a problem has been 494 * completed and the linear solver is no longer needed. It should 495 * return 0 upon success, nonzero on failure. 496 * ----------------------------------------------------------------- 497 */ 498 499 /* 500 * ================================================================= 501 * C V O D E I N T E R N A L F U N C T I O N S 502 * ================================================================= 503 */ 504 505 /* Prototype of internal ewtSet function */ 506 507 int cvEwtSet(N_Vector ycur, N_Vector weight, void *data); 508 509 /* High level error handler */ 510 511 void cvProcessError(CVodeMem cv_mem, 512 int error_code, const char *module, const char *fname, 513 const char *msgfmt, ...); 514 515 /* Prototype of internal ErrHandler function */ 516 517 void cvErrHandler(int error_code, const char *module, const char *function, 518 char *msg, void *data); 519 520 /* Nonlinear solver initializtion function */ 521 522 int cvNlsInit(CVodeMem cv_mem); 523 524 /* Projection functions */ 525 526 int cvDoProjection(CVodeMem cv_mem, int *nflagPtr, realtype saved_t, 527 int *npfPtr); 528 int cvProjInit(CVodeProjMem proj_mem); 529 int cvProjFree(CVodeProjMem *proj_mem); 530 531 /* Restore tn and undo prediction to reattempt a step */ 532 533 void cvRestore(CVodeMem cv_mem, realtype saved_t); 534 535 /* Reset h and rescale history array to prepare for a step */ 536 537 void cvRescale(CVodeMem cv_mem); 538 539 /* 540 * ================================================================= 541 * C V O D E E R R O R M E S S A G E S 542 * ================================================================= 543 */ 544 545 #if defined(SUNDIALS_EXTENDED_PRECISION) 546 547 #define MSG_TIME "t = %Lg" 548 #define MSG_TIME_H "t = %Lg and h = %Lg" 549 #define MSG_TIME_INT "t = %Lg is not between tcur - hu = %Lg and tcur = %Lg." 550 #define MSG_TIME_TOUT "tout = %Lg" 551 #define MSG_TIME_TSTOP "tstop = %Lg" 552 553 #elif defined(SUNDIALS_DOUBLE_PRECISION) 554 555 #define MSG_TIME "t = %lg" 556 #define MSG_TIME_H "t = %lg and h = %lg" 557 #define MSG_TIME_INT "t = %lg is not between tcur - hu = %lg and tcur = %lg." 558 #define MSG_TIME_TOUT "tout = %lg" 559 #define MSG_TIME_TSTOP "tstop = %lg" 560 561 #else 562 563 #define MSG_TIME "t = %g" 564 #define MSG_TIME_H "t = %g and h = %g" 565 #define MSG_TIME_INT "t = %g is not between tcur - hu = %g and tcur = %g." 566 #define MSG_TIME_TOUT "tout = %g" 567 #define MSG_TIME_TSTOP "tstop = %g" 568 569 #endif 570 571 /* Initialization and I/O error messages */ 572 573 #define MSGCV_NO_MEM "cvode_mem = NULL illegal." 574 #define MSGCV_CVMEM_FAIL "Allocation of cvode_mem failed." 575 #define MSGCV_MEM_FAIL "A memory request failed." 576 #define MSGCV_BAD_LMM "Illegal value for lmm. The legal values are CV_ADAMS and CV_BDF." 577 #define MSGCV_NO_MALLOC "Attempt to call before CVodeInit." 578 #define MSGCV_NEG_MAXORD "maxord <= 0 illegal." 579 #define MSGCV_BAD_MAXORD "Illegal attempt to increase maximum method order." 580 #define MSGCV_SET_SLDET "Attempt to use stability limit detection with the CV_ADAMS method illegal." 581 #define MSGCV_NEG_HMIN "hmin < 0 illegal." 582 #define MSGCV_NEG_HMAX "hmax < 0 illegal." 583 #define MSGCV_BAD_HMIN_HMAX "Inconsistent step size limits: hmin > hmax." 584 #define MSGCV_BAD_RELTOL "reltol < 0 illegal." 585 #define MSGCV_BAD_ABSTOL "abstol has negative component(s) (illegal)." 586 #define MSGCV_NULL_ABSTOL "abstol = NULL illegal." 587 #define MSGCV_NULL_Y0 "y0 = NULL illegal." 588 #define MSGCV_Y0_FAIL_CONSTR "y0 fails to satisfy constraints." 589 #define MSGCV_NULL_F "f = NULL illegal." 590 #define MSGCV_NULL_G "g = NULL illegal." 591 #define MSGCV_BAD_NVECTOR "A required vector operation is not implemented." 592 #define MSGCV_BAD_CONSTR "Illegal values in constraints vector." 593 #define MSGCV_BAD_K "Illegal value for k." 594 #define MSGCV_NULL_DKY "dky = NULL illegal." 595 #define MSGCV_BAD_T "Illegal value for t." MSG_TIME_INT 596 #define MSGCV_NO_ROOT "Rootfinding was not initialized." 597 #define MSGCV_NLS_INIT_FAIL "The nonlinear solver's init routine failed." 598 599 /* CVode Error Messages */ 600 601 #define MSGCV_NO_TOL "No integration tolerances have been specified." 602 #define MSGCV_LSOLVE_NULL "The linear solver's solve routine is NULL." 603 #define MSGCV_YOUT_NULL "yout = NULL illegal." 604 #define MSGCV_TRET_NULL "tret = NULL illegal." 605 #define MSGCV_BAD_EWT "Initial ewt has component(s) equal to zero (illegal)." 606 #define MSGCV_EWT_NOW_BAD "At " MSG_TIME ", a component of ewt has become <= 0." 607 #define MSGCV_BAD_ITASK "Illegal value for itask." 608 #define MSGCV_BAD_H0 "h0 and tout - t0 inconsistent." 609 #define MSGCV_BAD_TOUT "Trouble interpolating at " MSG_TIME_TOUT ". tout too far back in direction of integration" 610 #define MSGCV_EWT_FAIL "The user-provide EwtSet function failed." 611 #define MSGCV_EWT_NOW_FAIL "At " MSG_TIME ", the user-provide EwtSet function failed." 612 #define MSGCV_LINIT_FAIL "The linear solver's init routine failed." 613 #define MSGCV_HNIL_DONE "The above warning has been issued mxhnil times and will not be issued again for this problem." 614 #define MSGCV_TOO_CLOSE "tout too close to t0 to start integration." 615 #define MSGCV_MAX_STEPS "At " MSG_TIME ", mxstep steps taken before reaching tout." 616 #define MSGCV_TOO_MUCH_ACC "At " MSG_TIME ", too much accuracy requested." 617 #define MSGCV_HNIL "Internal " MSG_TIME_H " are such that t + h = t on the next step. The solver will continue anyway." 618 #define MSGCV_ERR_FAILS "At " MSG_TIME_H ", the error test failed repeatedly or with |h| = hmin." 619 #define MSGCV_CONV_FAILS "At " MSG_TIME_H ", the corrector convergence test failed repeatedly or with |h| = hmin." 620 #define MSGCV_SETUP_FAILED "At " MSG_TIME ", the setup routine failed in an unrecoverable manner." 621 #define MSGCV_SOLVE_FAILED "At " MSG_TIME ", the solve routine failed in an unrecoverable manner." 622 #define MSGCV_FAILED_CONSTR "At " MSG_TIME ", unable to satisfy inequality constraints." 623 #define MSGCV_RHSFUNC_FAILED "At " MSG_TIME ", the right-hand side routine failed in an unrecoverable manner." 624 #define MSGCV_RHSFUNC_UNREC "At " MSG_TIME ", the right-hand side failed in a recoverable manner, but no recovery is possible." 625 #define MSGCV_RHSFUNC_REPTD "At " MSG_TIME " repeated recoverable right-hand side function errors." 626 #define MSGCV_RHSFUNC_FIRST "The right-hand side routine failed at the first call." 627 #define MSGCV_RTFUNC_FAILED "At " MSG_TIME ", the rootfinding routine failed in an unrecoverable manner." 628 #define MSGCV_CLOSE_ROOTS "Root found at and very near " MSG_TIME "." 629 #define MSGCV_BAD_TSTOP "The value " MSG_TIME_TSTOP " is behind current " MSG_TIME " in the direction of integration." 630 #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." 631 #define MSGCV_NLS_SETUP_FAILED "At " MSG_TIME ", the nonlinear solver setup failed unrecoverably." 632 #define MSGCV_NLS_INPUT_NULL "At " MSG_TIME ", the nonlinear solver was passed a NULL input." 633 #define MSGCV_NLS_FAIL "At " MSG_TIME ", the nonlinear solver failed in an unrecoverable manner." 634 635 /* CVode Projection Error Messages */ 636 637 #define MSG_CV_MEM_NULL "cvode_mem = NULL illegal." 638 #define MSG_CV_MEM_FAIL "A memory request failed." 639 640 #define MSG_CV_PROJ_MEM_NULL "proj_mem = NULL illegal." 641 #define MSG_CV_PROJFUNC_FAIL "At " MSG_TIME " the projection function failed with an unrecoverable error." 642 #define MSG_CV_REPTD_PROJFUNC_ERR "At " MSG_TIME " the projection function had repeated recoverable errors." 643 644 #ifdef __cplusplus 645 } 646 #endif 647 648 #endif 649