1 /* ----------------------------------------------------------------- 2 * Programmer(s): Radu Serban @ LLNL 3 * ----------------------------------------------------------------- 4 * SUNDIALS Copyright Start 5 * Copyright (c) 2002-2021, 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 * This is the header file (private version) for the main IDAS solver. 15 * ----------------------------------------------------------------- 16 */ 17 18 #ifndef _IDAS_IMPL_H 19 #define _IDAS_IMPL_H 20 21 #include <stdarg.h> 22 23 #include "idas/idas.h" 24 25 #ifdef __cplusplus /* wrapper to enable C++ usage */ 26 extern "C" { 27 #endif 28 29 /* 30 * ================================================================= 31 * 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 32 * ================================================================= 33 */ 34 35 36 /* Basic IDA constants */ 37 38 #define HMAX_INV_DEFAULT RCONST(0.0) /* hmax_inv default value */ 39 #define MAXORD_DEFAULT 5 /* maxord default value */ 40 #define MXORDP1 6 /* max. number of N_Vectors in phi */ 41 #define MXSTEP_DEFAULT 500 /* mxstep default value */ 42 43 /* Return values for lower level routines used by IDASolve and functions 44 provided to the nonlinear solver */ 45 46 #define IDA_RES_RECVR +1 47 #define IDA_LSETUP_RECVR +2 48 #define IDA_LSOLVE_RECVR +3 49 #define IDA_CONSTR_RECVR +5 50 #define IDA_NLS_SETUP_RECVR +6 51 52 #define IDA_QRHS_RECVR +10 53 #define IDA_SRES_RECVR +11 54 #define IDA_QSRHS_RECVR +12 55 56 /* itol */ 57 #define IDA_NN 0 58 #define IDA_SS 1 59 #define IDA_SV 2 60 #define IDA_WF 3 61 #define IDA_EE 4 62 63 /* 64 * ---------------------------------------------------------------- 65 * Types: struct IDAMemRec, IDAMem 66 * ---------------------------------------------------------------- 67 * The type IDAMem is type pointer to struct IDAMemRec. 68 * This structure contains fields to keep track of problem state. 69 * ---------------------------------------------------------------- 70 */ 71 72 typedef struct IDAMemRec { 73 74 realtype ida_uround; /* machine unit roundoff */ 75 76 /*-------------------------- 77 Problem Specification Data 78 --------------------------*/ 79 80 IDAResFn ida_res; /* F(t,y(t),y'(t))=0; the function F */ 81 void *ida_user_data; /* user pointer passed to res */ 82 83 int ida_itol; /* itol = IDA_SS, IDA_SV, IDA_WF, IDA_NN */ 84 realtype ida_rtol; /* relative tolerance */ 85 realtype ida_Satol; /* scalar absolute tolerance */ 86 N_Vector ida_Vatol; /* vector absolute tolerance */ 87 booleantype ida_atolmin0; /* flag indicating that min(atol) = 0 */ 88 booleantype ida_user_efun; /* SUNTRUE if user provides efun */ 89 IDAEwtFn ida_efun; /* function to set ewt */ 90 void *ida_edata; /* user pointer passed to efun */ 91 92 booleantype ida_constraintsSet; /* constraints vector present: 93 do constraints calc */ 94 booleantype ida_suppressalg; /* SUNTRUE means suppress algebraic vars 95 in local error tests */ 96 97 /*----------------------- 98 Quadrature Related Data 99 -----------------------*/ 100 101 booleantype ida_quadr; 102 103 IDAQuadRhsFn ida_rhsQ; 104 void *ida_user_dataQ; 105 106 booleantype ida_errconQ; 107 108 int ida_itolQ; 109 realtype ida_rtolQ; 110 realtype ida_SatolQ; /* scalar absolute tolerance for quadratures */ 111 N_Vector ida_VatolQ; /* vector absolute tolerance for quadratures */ 112 booleantype ida_atolQmin0; /* flag indicating that min(atolQ) = 0 */ 113 114 /*------------------------ 115 Sensitivity Related Data 116 ------------------------*/ 117 118 booleantype ida_sensi; 119 int ida_Ns; 120 int ida_ism; 121 122 IDASensResFn ida_resS; 123 void *ida_user_dataS; 124 booleantype ida_resSDQ; 125 126 realtype *ida_p; 127 realtype *ida_pbar; 128 int *ida_plist; 129 int ida_DQtype; 130 realtype ida_DQrhomax; 131 132 booleantype ida_errconS; /* SUNTRUE if sensitivities in err. control */ 133 134 int ida_itolS; 135 realtype ida_rtolS; /* relative tolerance for sensitivities */ 136 realtype *ida_SatolS; /* scalar absolute tolerances for sensi. */ 137 N_Vector *ida_VatolS; /* vector absolute tolerances for sensi. */ 138 booleantype *ida_atolSmin0; /* flag indicating that min(atolS[is]) = 0 */ 139 140 /*----------------------------------- 141 Quadrature Sensitivity Related Data 142 -----------------------------------*/ 143 144 booleantype ida_quadr_sensi; /* SUNTRUE if computing sensitivities of quadrs. */ 145 146 IDAQuadSensRhsFn ida_rhsQS; /* fQS = (dfQ/dy)*yS + (dfQ/dp) */ 147 void *ida_user_dataQS; /* data pointer passed to fQS */ 148 booleantype ida_rhsQSDQ; /* SUNTRUE if using internal DQ functions */ 149 150 booleantype ida_errconQS; /* SUNTRUE if yQS are considered in err. con. */ 151 152 int ida_itolQS; 153 realtype ida_rtolQS; /* relative tolerance for yQS */ 154 realtype *ida_SatolQS; /* scalar absolute tolerances for yQS */ 155 N_Vector *ida_VatolQS; /* vector absolute tolerances for yQS */ 156 booleantype *ida_atolQSmin0; /* flag indicating that min(atolQS[is]) = 0 */ 157 158 /*----------------------------------------------- 159 Divided differences array and associated arrays 160 -----------------------------------------------*/ 161 162 N_Vector ida_phi[MXORDP1]; /* phi = (maxord+1) arrays of divided differences */ 163 164 realtype ida_psi[MXORDP1]; /* differences in t (sums of recent step sizes) */ 165 realtype ida_alpha[MXORDP1]; /* ratios of current stepsize to psi values */ 166 realtype ida_beta[MXORDP1]; /* ratios of current to previous product of psi's */ 167 realtype ida_sigma[MXORDP1]; /* product successive alpha values and factorial */ 168 realtype ida_gamma[MXORDP1]; /* sum of reciprocals of psi values */ 169 170 /*------------------------- 171 N_Vectors for integration 172 -------------------------*/ 173 174 N_Vector ida_ewt; /* error weight vector */ 175 N_Vector ida_yy; /* work space for y vector (= user's yret) */ 176 N_Vector ida_yp; /* work space for y' vector (= user's ypret) */ 177 N_Vector ida_yypredict; /* predicted y vector */ 178 N_Vector ida_yppredict; /* predicted y' vector */ 179 N_Vector ida_delta; /* residual vector */ 180 N_Vector ida_id; /* bit vector for diff./algebraic components */ 181 N_Vector ida_constraints; /* vector of inequality constraint options */ 182 N_Vector ida_savres; /* saved residual vector */ 183 N_Vector ida_ee; /* accumulated corrections to y vector, but 184 set equal to estimated local errors upon 185 successful return */ 186 N_Vector ida_tempv1; /* work space vector */ 187 N_Vector ida_tempv2; /* work space vector */ 188 N_Vector ida_tempv3; /* work space vector */ 189 N_Vector ida_ynew; /* work vector for y in IDACalcIC (= tempv2) */ 190 N_Vector ida_ypnew; /* work vector for yp in IDACalcIC (= ee) */ 191 N_Vector ida_delnew; /* work vector for delta in IDACalcIC (= phi[2]) */ 192 N_Vector ida_dtemp; /* work vector in IDACalcIC (= phi[3]) */ 193 194 195 /*---------------------------- 196 Quadrature Related N_Vectors 197 ----------------------------*/ 198 199 N_Vector ida_phiQ[MXORDP1]; 200 N_Vector ida_yyQ; 201 N_Vector ida_ypQ; 202 N_Vector ida_ewtQ; 203 N_Vector ida_eeQ; 204 205 /*--------------------------- 206 Sensitivity Related Vectors 207 ---------------------------*/ 208 209 N_Vector *ida_phiS[MXORDP1]; 210 N_Vector *ida_ewtS; 211 212 N_Vector *ida_eeS; /* cumulative sensitivity corrections */ 213 214 N_Vector *ida_yyS; /* allocated and used for: */ 215 N_Vector *ida_ypS; /* ism = SIMULTANEOUS */ 216 N_Vector *ida_yySpredict; /* ism = STAGGERED */ 217 N_Vector *ida_ypSpredict; 218 N_Vector *ida_deltaS; 219 220 N_Vector ida_tmpS1; /* work space vectors | tmpS1 = tempv1 */ 221 N_Vector ida_tmpS2; /* for resS | tmpS2 = tempv2 */ 222 N_Vector ida_tmpS3; /* | tmpS3 = allocated */ 223 224 N_Vector *ida_savresS; /* work vector in IDACalcIC for stg (= phiS[2]) */ 225 N_Vector *ida_delnewS; /* work vector in IDACalcIC for stg (= phiS[3]) */ 226 227 N_Vector *ida_yyS0; /* initial yS, ypS vectors allocated and */ 228 N_Vector *ida_ypS0; /* deallocated in IDACalcIC function */ 229 230 N_Vector *ida_yyS0new; /* work vector in IDASensLineSrch (= phiS[4]) */ 231 N_Vector *ida_ypS0new; /* work vector in IDASensLineSrch (= eeS) */ 232 233 /*-------------------------------------- 234 Quadrature Sensitivity Related Vectors 235 --------------------------------------*/ 236 237 N_Vector *ida_phiQS[MXORDP1];/* Mod. div. diffs. for quadr. sensitivities */ 238 N_Vector *ida_ewtQS; /* error weight vectors for sensitivities */ 239 240 N_Vector *ida_eeQS; /* cumulative quadr.sensi.corrections */ 241 242 N_Vector *ida_yyQS; /* Unlike yS, yQS is not allocated by the user */ 243 N_Vector *ida_tempvQS; /* temporary storage vector (~ tempv) */ 244 N_Vector ida_savrhsQ; /* saved quadr. rhs (needed for rhsQS calls) */ 245 246 /*------------------------------ 247 Variables for use by IDACalcIC 248 ------------------------------*/ 249 250 realtype ida_t0; /* initial t */ 251 N_Vector ida_yy0; /* initial y vector (user-supplied). */ 252 N_Vector ida_yp0; /* initial y' vector (user-supplied). */ 253 254 int ida_icopt; /* IC calculation user option */ 255 booleantype ida_lsoff; /* IC calculation linesearch turnoff option */ 256 int ida_maxnh; /* max. number of h tries in IC calculation */ 257 int ida_maxnj; /* max. number of J tries in IC calculation */ 258 int ida_maxnit; /* max. number of Netwon iterations in IC calc. */ 259 int ida_nbacktr; /* number of IC linesearch backtrack operations */ 260 int ida_sysindex; /* computed system index (0 or 1) */ 261 int ida_maxbacks; /* max backtracks per Newton step */ 262 realtype ida_epiccon; /* IC nonlinear convergence test constant */ 263 realtype ida_steptol; /* minimum Newton step size in IC calculation */ 264 realtype ida_tscale; /* time scale factor = abs(tout1 - t0) */ 265 266 /* Tstop information */ 267 268 booleantype ida_tstopset; 269 realtype ida_tstop; 270 271 /* Step Data */ 272 273 int ida_kk; /* current BDF method order */ 274 int ida_kused; /* method order used on last successful step */ 275 int ida_knew; /* order for next step from order decrease decision */ 276 int ida_phase; /* flag to trigger step doubling in first few steps */ 277 int ida_ns; /* counts steps at fixed stepsize and order */ 278 279 realtype ida_hin; /* initial step */ 280 realtype ida_h0u; /* actual initial stepsize */ 281 realtype ida_hh; /* current step size h */ 282 realtype ida_hused; /* step size used on last successful step */ 283 realtype ida_rr; /* rr = hnext / hused */ 284 realtype ida_tn; /* current internal value of t */ 285 realtype ida_tretlast; /* value of tret previously returned by IDASolve */ 286 realtype ida_cj; /* current value of scalar (-alphas/hh) in Jacobian */ 287 realtype ida_cjlast; /* cj value saved from last successful step */ 288 realtype ida_cjold; /* cj value saved from last call to lsetup */ 289 realtype ida_cjratio; /* ratio of cj values: cj/cjold */ 290 realtype ida_ss; /* scalar used in Newton iteration convergence test */ 291 realtype ida_oldnrm; /* norm of previous nonlinear solver update */ 292 realtype ida_epsNewt; /* test constant in Newton convergence test */ 293 realtype ida_epcon; /* coeficient of the Newton covergence test */ 294 realtype ida_toldel; /* tolerance in direct test on Newton corrections */ 295 296 realtype ida_ssS; /* scalar ss for staggered sensitivities */ 297 298 /*------ 299 Limits 300 ------*/ 301 302 int ida_maxncf; /* max numer of convergence failures */ 303 int ida_maxnef; /* max number of error test failures */ 304 305 int ida_maxord; /* max value of method order k: */ 306 int ida_maxord_alloc; /* value of maxord used when allocating memory */ 307 long int ida_mxstep; /* max number of internal steps for one user call */ 308 realtype ida_hmax_inv; /* inverse of max. step size hmax (default = 0.0) */ 309 310 /*-------- 311 Counters 312 --------*/ 313 314 long int ida_nst; /* number of internal steps taken */ 315 316 long int ida_nre; /* number of function (res) calls */ 317 long int ida_nrQe; 318 long int ida_nrSe; 319 long int ida_nrQSe; /* number of fQS calls */ 320 long int ida_nreS; 321 long int ida_nrQeS; /* number of fQ calls from sensi DQ */ 322 323 324 long int ida_ncfn; /* number of corrector convergence failures */ 325 long int ida_ncfnQ; 326 long int ida_ncfnS; 327 328 long int ida_netf; /* number of error test failures */ 329 long int ida_netfQ; 330 long int ida_netfS; 331 long int ida_netfQS; /* number of quadr. sensi. error test failures */ 332 333 long int ida_nni; /* number of Newton iterations performed */ 334 long int ida_nniS; 335 336 long int ida_nsetups; /* number of lsetup calls */ 337 long int ida_nsetupsS; 338 339 /*------------------ 340 Space requirements 341 ------------------*/ 342 343 sunindextype ida_lrw1; /* no. of realtype words in 1 N_Vector */ 344 sunindextype ida_liw1; /* no. of integer words in 1 N_Vector */ 345 sunindextype ida_lrw1Q; 346 sunindextype ida_liw1Q; 347 long int ida_lrw; /* number of realtype words in IDA work vectors */ 348 long int ida_liw; /* no. of integer words in IDA work vectors */ 349 350 realtype ida_tolsf; /* tolerance scale factor (saved value) */ 351 352 /*------------------------------------------- 353 Error handler function and error ouput file 354 -------------------------------------------*/ 355 356 IDAErrHandlerFn ida_ehfun; /* Error messages are handled by ehfun */ 357 void *ida_eh_data; /* dats pointer passed to ehfun */ 358 FILE *ida_errfp; /* IDA error messages are sent to errfp */ 359 360 /* Flags to verify correct calling sequence */ 361 362 booleantype ida_SetupDone; /* set to SUNFALSE by IDAMalloc and IDAReInit 363 set to SUNTRUE by IDACalcIC or IDASolve */ 364 365 booleantype ida_VatolMallocDone; 366 booleantype ida_constraintsMallocDone; 367 booleantype ida_idMallocDone; 368 369 booleantype ida_MallocDone; /* set to SUNFALSE by IDACreate 370 set to SUNTRUE by IDAMAlloc 371 tested by IDAReInit and IDASolve */ 372 373 booleantype ida_VatolQMallocDone; 374 booleantype ida_quadMallocDone; 375 376 booleantype ida_VatolSMallocDone; 377 booleantype ida_SatolSMallocDone; 378 booleantype ida_sensMallocDone; 379 380 booleantype ida_VatolQSMallocDone; 381 booleantype ida_SatolQSMallocDone; 382 booleantype ida_quadSensMallocDone; 383 384 /*--------------------- 385 Nonlinear Solver Data 386 ---------------------*/ 387 388 SUNNonlinearSolver NLS; /* nonlinear solver object */ 389 booleantype ownNLS; /* flag indicating NLS ownership */ 390 391 SUNNonlinearSolver NLSsim; /* nonlinear solver object for DAE+Sens solves 392 with the simultaneous corrector option */ 393 booleantype ownNLSsim; /* flag indicating NLS ownership */ 394 395 SUNNonlinearSolver NLSstg; /* nonlinear solver object for DAE+Sens solves 396 with the staggered corrector option */ 397 booleantype ownNLSstg; /* flag indicating NLS ownership */ 398 399 /* The following vectors are NVector wrappers for use with the simultaneous 400 and staggered corrector methods: 401 402 Simult: ypredictSim = [ida_delta, ida_deltaS] 403 ycorSim = [ida_ee, ida_eeS] 404 ewtSim = [ida_ewt, ida_ewtS] 405 406 Stagger: ypredictStg = ida_deltaS 407 ycorStg = ida_eeS 408 ewtStg = ida_ewtS 409 */ 410 N_Vector ypredictSim, ycorSim, ewtSim; 411 N_Vector ypredictStg, ycorStg, ewtStg; 412 413 /* flags indicating if vector wrappers for the simultaneous and staggered 414 correctors have been allocated */ 415 booleantype simMallocDone; 416 booleantype stgMallocDone; 417 418 /*------------------ 419 Linear Solver Data 420 ------------------*/ 421 422 /* Linear Solver functions to be called */ 423 424 int (*ida_linit)(struct IDAMemRec *idamem); 425 426 int (*ida_lsetup)(struct IDAMemRec *idamem, N_Vector yyp, 427 N_Vector ypp, N_Vector resp, 428 N_Vector tempv1, N_Vector tempv2, N_Vector tempv3); 429 430 int (*ida_lsolve)(struct IDAMemRec *idamem, N_Vector b, N_Vector weight, 431 N_Vector ycur, N_Vector ypcur, N_Vector rescur); 432 433 int (*ida_lperf)(struct IDAMemRec *idamem, int perftask); 434 435 int (*ida_lfree)(struct IDAMemRec *idamem); 436 437 /* Linear Solver specific memory */ 438 439 void *ida_lmem; 440 441 /* Flag to request a call to the setup routine */ 442 443 booleantype ida_forceSetup; 444 445 /* Flag to indicate successful ida_linit call */ 446 447 booleantype ida_linitOK; 448 449 /*---------------- 450 Rootfinding Data 451 ----------------*/ 452 453 IDARootFn ida_gfun; /* Function g for roots sought */ 454 int ida_nrtfn; /* number of components of g */ 455 int *ida_iroots; /* array for root information */ 456 int *ida_rootdir; /* array specifying direction of zero-crossing */ 457 realtype ida_tlo; /* nearest endpoint of interval in root search */ 458 realtype ida_thi; /* farthest endpoint of interval in root search */ 459 realtype ida_trout; /* t return value from rootfinder routine */ 460 realtype *ida_glo; /* saved array of g values at t = tlo */ 461 realtype *ida_ghi; /* saved array of g values at t = thi */ 462 realtype *ida_grout; /* array of g values at t = trout */ 463 realtype ida_toutc; /* copy of tout (if NORMAL mode) */ 464 realtype ida_ttol; /* tolerance on root location */ 465 int ida_taskc; /* copy of parameter itask */ 466 int ida_irfnd; /* flag showing whether last step had a root */ 467 long int ida_nge; /* counter for g evaluations */ 468 booleantype *ida_gactive; /* array with active/inactive event functions */ 469 int ida_mxgnull; /* number of warning messages about possible g==0 */ 470 471 /* Arrays for Fused Vector Operations */ 472 473 /* scalar arrays */ 474 realtype* ida_cvals; 475 realtype ida_dvals[MAXORD_DEFAULT]; 476 477 /* vector arrays */ 478 N_Vector* ida_Xvecs; 479 N_Vector* ida_Zvecs; 480 481 /*------------------------ 482 Adjoint sensitivity data 483 ------------------------*/ 484 485 booleantype ida_adj; /* SUNTRUE if performing ASA */ 486 487 struct IDAadjMemRec *ida_adj_mem; /* Pointer to adjoint memory structure */ 488 489 booleantype ida_adjMallocDone; 490 491 } *IDAMem; 492 493 /* 494 * ================================================================= 495 * A D J O I N T M O D U L E M E M O R Y B L O C K 496 * ================================================================= 497 */ 498 499 /* 500 * ----------------------------------------------------------------- 501 * Forward references for pointers to various structures 502 * ----------------------------------------------------------------- 503 */ 504 505 typedef struct IDAadjMemRec *IDAadjMem; 506 typedef struct CkpntMemRec *CkpntMem; 507 typedef struct DtpntMemRec *DtpntMem; 508 typedef struct IDABMemRec *IDABMem; 509 510 /* 511 * ----------------------------------------------------------------- 512 * Types for functions provided by an interpolation module 513 * ----------------------------------------------------------------- 514 * IDAAMMallocFn: Type for a function that initializes the content 515 * field of the structures in the dt array 516 * IDAAMFreeFn: Type for a function that deallocates the content 517 * field of the structures in the dt array 518 * IDAAGetYFn: Function type for a function that returns the 519 * interpolated forward solution. 520 * IDAAStorePnt: Function type for a function that stores a new 521 * point in the structure d 522 * ----------------------------------------------------------------- 523 */ 524 525 typedef booleantype (*IDAAMMallocFn)(IDAMem IDA_mem); 526 typedef void (*IDAAMFreeFn)(IDAMem IDA_mem); 527 typedef int (*IDAAGetYFn)(IDAMem IDA_mem, realtype t, 528 N_Vector yy, N_Vector yp, 529 N_Vector *yyS, N_Vector *ypS); 530 typedef int (*IDAAStorePntFn)(IDAMem IDA_mem, DtpntMem d); 531 532 /* 533 * ----------------------------------------------------------------- 534 * Types : struct CkpntMemRec, CkpntMem 535 * ----------------------------------------------------------------- 536 * The type CkpntMem is type pointer to struct CkpntMemRec. 537 * This structure contains fields to store all information at a 538 * check point that is needed to 'hot' start IDAS. 539 * ----------------------------------------------------------------- 540 */ 541 542 struct CkpntMemRec { 543 544 /* Integration limits */ 545 realtype ck_t0; 546 realtype ck_t1; 547 548 /* Modified divided difference array */ 549 N_Vector ck_phi[MXORDP1]; 550 551 /* Do we need to carry quadratures? */ 552 booleantype ck_quadr; 553 554 /* Modified divided difference array for quadratures */ 555 N_Vector ck_phiQ[MXORDP1]; 556 557 /* Do we need to carry sensitivities? */ 558 booleantype ck_sensi; 559 560 /* number of sensitivities */ 561 int ck_Ns; 562 563 /* Modified divided difference array for sensitivities */ 564 N_Vector *ck_phiS[MXORDP1]; 565 566 /* Do we need to carry quadrature sensitivities? */ 567 booleantype ck_quadr_sensi; 568 569 /* Modified divided difference array for quadrature sensitivities */ 570 N_Vector *ck_phiQS[MXORDP1]; 571 572 573 /* Step data */ 574 long int ck_nst; 575 realtype ck_tretlast; 576 int ck_ns; 577 int ck_kk; 578 int ck_kused; 579 int ck_knew; 580 int ck_phase; 581 582 realtype ck_hh; 583 realtype ck_hused; 584 realtype ck_rr; 585 realtype ck_cj; 586 realtype ck_cjlast; 587 realtype ck_cjold; 588 realtype ck_cjratio; 589 realtype ck_ss; 590 realtype ck_ssS; 591 592 realtype ck_psi[MXORDP1]; 593 realtype ck_alpha[MXORDP1]; 594 realtype ck_beta[MXORDP1]; 595 realtype ck_sigma[MXORDP1]; 596 realtype ck_gamma[MXORDP1]; 597 598 /* How many phi, phiS, phiQ and phiQS were allocated? */ 599 int ck_phi_alloc; 600 601 /* Pointer to next structure in list */ 602 struct CkpntMemRec *ck_next; 603 }; 604 605 /* 606 * ----------------------------------------------------------------- 607 * Type : struct DtpntMemRec 608 * ----------------------------------------------------------------- 609 * This structure contains fields to store all information at a 610 * data point that is needed to interpolate solution of forward 611 * simulations. Its content field is interpType-dependent. 612 * ----------------------------------------------------------------- 613 */ 614 615 struct DtpntMemRec { 616 realtype t; /* time */ 617 void *content; /* interpType-dependent content */ 618 }; 619 620 /* Data for cubic Hermite interpolation */ 621 typedef struct HermiteDataMemRec { 622 N_Vector y; 623 N_Vector yd; 624 N_Vector *yS; 625 N_Vector *ySd; 626 } *HermiteDataMem; 627 628 /* Data for polynomial interpolation */ 629 typedef struct PolynomialDataMemRec { 630 N_Vector y; 631 N_Vector *yS; 632 633 /* yd and ySd store the derivative(s) only for the first dt 634 point. NULL otherwise. */ 635 N_Vector yd; 636 N_Vector *ySd; 637 int order; 638 } *PolynomialDataMem; 639 640 /* 641 * ----------------------------------------------------------------- 642 * Type : struct IDABMemRec 643 * ----------------------------------------------------------------- 644 * The type IDABMemRec is a pointer to a structure which stores all 645 * information for ONE backward problem. 646 * The IDAadjMem struct contains a linked list of IDABMem pointers 647 * ----------------------------------------------------------------- 648 */ 649 struct IDABMemRec { 650 651 /* Index of this backward problem */ 652 int ida_index; 653 654 /* Time at which the backward problem is initialized. */ 655 realtype ida_t0; 656 657 /* Memory for this backward problem */ 658 IDAMem IDA_mem; 659 660 /* Flags to indicate that this backward problem's RHS or quad RHS 661 * require forward sensitivities */ 662 booleantype ida_res_withSensi; 663 booleantype ida_rhsQ_withSensi; 664 665 /* Residual function for backward run */ 666 IDAResFnB ida_res; 667 IDAResFnBS ida_resS; 668 669 /* Right hand side quadrature function (fQB) for backward run */ 670 IDAQuadRhsFnB ida_rhsQ; 671 IDAQuadRhsFnBS ida_rhsQS; 672 673 /* User user_data */ 674 void *ida_user_data; 675 676 /* Linear solver's data and functions */ 677 678 /* Memory block for a linear solver's interface to IDAA */ 679 void *ida_lmem; 680 681 /* Function to free any memory allocated by the linear solver */ 682 int (*ida_lfree)(IDABMem IDAB_mem); 683 684 /* Memory block for a preconditioner's module interface to IDAA */ 685 void *ida_pmem; 686 687 /* Function to free any memory allocated by the preconditioner module */ 688 int (*ida_pfree)(IDABMem IDAB_mem); 689 690 /* Time at which to extract solution / quadratures */ 691 realtype ida_tout; 692 693 /* Workspace Nvectors */ 694 N_Vector ida_yy; 695 N_Vector ida_yp; 696 697 /* Link to next structure in list. */ 698 struct IDABMemRec *ida_next; 699 }; 700 701 702 /* 703 * ----------------------------------------------------------------- 704 * Type : struct IDAadjMemRec 705 * ----------------------------------------------------------------- 706 * The type IDAadjMem is type pointer to struct IDAadjMemRec. 707 * This structure contins fields to store all information 708 * necessary for adjoint sensitivity analysis. 709 * ----------------------------------------------------------------- 710 */ 711 712 struct IDAadjMemRec { 713 714 /* -------------------- 715 * Forward problem data 716 * -------------------- */ 717 718 /* Integration interval */ 719 realtype ia_tinitial, ia_tfinal; 720 721 /* Flag for first call to IDASolveF */ 722 booleantype ia_firstIDAFcall; 723 724 /* Flag if IDASolveF was called with TSTOP */ 725 booleantype ia_tstopIDAFcall; 726 realtype ia_tstopIDAF; 727 728 /* Flag if IDASolveF was called in IDA_NORMAL_MODE and encountered 729 a root after tout */ 730 booleantype ia_rootret; 731 realtype ia_troot; 732 733 /* ---------------------- 734 * Backward problems data 735 * ---------------------- */ 736 737 /* Storage for backward problems */ 738 struct IDABMemRec *IDAB_mem; 739 740 /* Number of backward problems. */ 741 int ia_nbckpbs; 742 743 /* Address of current backward problem (iterator). */ 744 struct IDABMemRec *ia_bckpbCrt; 745 746 /* Flag for first call to IDASolveB */ 747 booleantype ia_firstIDABcall; 748 749 /* ---------------- 750 * Check point data 751 * ---------------- */ 752 753 /* Storage for check point information */ 754 struct CkpntMemRec *ck_mem; 755 756 /* address of the check point structure for which data is available */ 757 struct CkpntMemRec *ia_ckpntData; 758 759 /* Number of checkpoints. */ 760 int ia_nckpnts; 761 762 /* ------------------ 763 * Interpolation data 764 * ------------------ */ 765 766 /* Number of steps between 2 check points */ 767 long int ia_nsteps; 768 769 /* Last index used in IDAAfindIndex */ 770 long int ia_ilast; 771 772 /* Storage for data from forward runs */ 773 struct DtpntMemRec **dt_mem; 774 775 /* Actual number of data points saved in current dt_mem */ 776 /* Commonly, np = nsteps+1 */ 777 long int ia_np; 778 779 /* Interpolation type */ 780 int ia_interpType; 781 782 783 /* Functions set by the interpolation module */ 784 IDAAStorePntFn ia_storePnt; /* store a new interpolation point */ 785 IDAAGetYFn ia_getY; /* interpolate forward solution */ 786 IDAAMMallocFn ia_malloc; /* allocate new data point */ 787 IDAAMFreeFn ia_free; /* destroys data point */ 788 789 /* Flags controlling the interpolation module */ 790 booleantype ia_mallocDone; /* IM initialized? */ 791 booleantype ia_newData; /* new data available in dt_mem? */ 792 booleantype ia_storeSensi; /* store sensitivities? */ 793 booleantype ia_interpSensi; /* interpolate sensitivities? */ 794 795 booleantype ia_noInterp; /* interpolations are temporarly */ 796 /* disabled ( IDACalcICB ) */ 797 798 /* Workspace for polynomial interpolation */ 799 N_Vector ia_Y[MXORDP1]; /* pointers phi[i] */ 800 N_Vector *ia_YS[MXORDP1]; /* pointers phiS[i] */ 801 realtype ia_T[MXORDP1]; 802 803 /* Workspace for wrapper functions */ 804 N_Vector ia_yyTmp, ia_ypTmp; 805 N_Vector *ia_yySTmp, *ia_ypSTmp; 806 807 }; 808 809 810 /* 811 * ================================================================= 812 * I N T E R F A C E T O L I N E A R S O L V E R S 813 * ================================================================= 814 */ 815 816 /* 817 * ----------------------------------------------------------------- 818 * int (*ida_linit)(IDAMem IDA_mem); 819 * ----------------------------------------------------------------- 820 * The purpose of ida_linit is to allocate memory for the 821 * solver-specific fields in the structure *(idamem->ida_lmem) and 822 * perform any needed initializations of solver-specific memory, 823 * such as counters/statistics. An (*ida_linit) should return 824 * 0 if it has successfully initialized the IDA linear solver and 825 * a non-zero value otherwise. If an error does occur, an 826 * appropriate message should be sent to the error handler function. 827 * ---------------------------------------------------------------- 828 */ 829 830 /* 831 * ----------------------------------------------------------------- 832 * int (*ida_lsetup)(IDAMem IDA_mem, N_Vector yyp, N_Vector ypp, 833 * N_Vector resp, N_Vector tempv1, 834 * N_Vector tempv2, N_Vector tempv3); 835 * ----------------------------------------------------------------- 836 * The job of ida_lsetup is to prepare the linear solver for 837 * subsequent calls to ida_lsolve. Its parameters are as follows: 838 * 839 * idamem - problem memory pointer of type IDAMem. See the big 840 * typedef earlier in this file. 841 * 842 * yyp - the predicted y vector for the current IDA internal 843 * step. 844 * 845 * ypp - the predicted y' vector for the current IDA internal 846 * step. 847 * 848 * resp - F(tn, yyp, ypp). 849 * 850 * tempv1, tempv2, tempv3 - temporary N_Vectors provided for use 851 * by ida_lsetup. 852 * 853 * The ida_lsetup routine should return 0 if successful, 854 * a positive value for a recoverable error, and a negative value 855 * for an unrecoverable error. 856 * ----------------------------------------------------------------- 857 */ 858 859 /* 860 * ----------------------------------------------------------------- 861 * int (*ida_lsolve)(IDAMem IDA_mem, N_Vector b, N_Vector weight, 862 * N_Vector ycur, N_Vector ypcur, N_Vector rescur); 863 * ----------------------------------------------------------------- 864 * ida_lsolve must solve the linear equation P x = b, where 865 * P is some approximation to the system Jacobian 866 * J = (dF/dy) + cj (dF/dy') 867 * evaluated at (tn,ycur,ypcur) and the RHS vector b is input. 868 * The N-vector ycur contains the solver's current approximation 869 * to y(tn), ypcur contains that for y'(tn), and the vector rescur 870 * contains the N-vector residual F(tn,ycur,ypcur). 871 * The solution is to be returned in the vector b. 872 * 873 * The ida_lsolve routine should return 0 if successful, 874 * a positive value for a recoverable error, and a negative value 875 * for an unrecoverable error. 876 * ----------------------------------------------------------------- 877 */ 878 879 /* 880 * ----------------------------------------------------------------- 881 * int (*ida_lperf)(IDAMem IDA_mem, int perftask); 882 * ----------------------------------------------------------------- 883 * ida_lperf is called two places in IDAS where linear solver 884 * performance data is required by IDAS. For perftask = 0, an 885 * initialization of performance variables is performed, while for 886 * perftask = 1, the performance is evaluated. 887 * ----------------------------------------------------------------- 888 */ 889 890 /* 891 * ----------------------------------------------------------------- 892 * int (*ida_lfree)(IDAMem IDA_mem); 893 * ----------------------------------------------------------------- 894 * ida_lfree should free up any memory allocated by the linear 895 * solver. This routine is called once a problem has been 896 * completed and the linear solver is no longer needed. It should 897 * return 0 upon success, nonzero on failure. 898 * ----------------------------------------------------------------- 899 */ 900 901 /* 902 * ================================================================= 903 * I N T E R N A L F U N C T I O N S 904 * ================================================================= 905 */ 906 907 /* Prototype of internal ewtSet function */ 908 909 int IDAEwtSet(N_Vector ycur, N_Vector weight, void *data); 910 911 /* High level error handler */ 912 913 void IDAProcessError(IDAMem IDA_mem, 914 int error_code, const char *module, const char *fname, 915 const char *msgfmt, ...); 916 917 /* Prototype of internal errHandler function */ 918 919 void IDAErrHandler(int error_code, const char *module, const char *function, 920 char *msg, void *data); 921 922 /* Norm functions. Also used for IC, so they are global.*/ 923 924 realtype IDAWrmsNorm(IDAMem IDA_mem, N_Vector x, N_Vector w, 925 booleantype mask); 926 927 realtype IDASensWrmsNorm(IDAMem IDA_mem, N_Vector *xS, N_Vector *wS, 928 booleantype mask); 929 930 realtype IDASensWrmsNormUpdate(IDAMem IDA_mem, realtype old_nrm, 931 N_Vector *xS, N_Vector *wS, 932 booleantype mask); 933 934 /* Nonlinear solver initialization */ 935 936 int idaNlsInit(IDAMem IDA_mem); 937 int idaNlsInitSensSim(IDAMem IDA_mem); 938 int idaNlsInitSensStg(IDAMem IDA_mem); 939 940 /* Prototype for internal sensitivity residual DQ function */ 941 942 int IDASensResDQ(int Ns, realtype t, 943 N_Vector yy, N_Vector yp, N_Vector resval, 944 N_Vector *yyS, N_Vector *ypS, N_Vector *resvalS, 945 void *user_dataS, 946 N_Vector ytemp, N_Vector yptemp, N_Vector restemp); 947 948 /* 949 * ================================================================= 950 * E R R O R M E S S A G E S 951 * ================================================================= 952 */ 953 954 #if defined(SUNDIALS_EXTENDED_PRECISION) 955 956 #define MSG_TIME "t = %Lg, " 957 #define MSG_TIME_H "t = %Lg and h = %Lg, " 958 #define MSG_TIME_INT "t = %Lg is not between tcur - hu = %Lg and tcur = %Lg." 959 #define MSG_TIME_TOUT "tout = %Lg" 960 #define MSG_TIME_TSTOP "tstop = %Lg" 961 962 #elif defined(SUNDIALS_DOUBLE_PRECISION) 963 964 #define MSG_TIME "t = %lg, " 965 #define MSG_TIME_H "t = %lg and h = %lg, " 966 #define MSG_TIME_INT "t = %lg is not between tcur - hu = %lg and tcur = %lg." 967 #define MSG_TIME_TOUT "tout = %lg" 968 #define MSG_TIME_TSTOP "tstop = %lg" 969 970 #else 971 972 #define MSG_TIME "t = %g, " 973 #define MSG_TIME_H "t = %g and h = %g, " 974 #define MSG_TIME_INT "t = %g is not between tcur - hu = %g and tcur = %g." 975 #define MSG_TIME_TOUT "tout = %g" 976 #define MSG_TIME_TSTOP "tstop = %g" 977 978 #endif 979 980 /* General errors */ 981 982 #define MSG_MEM_FAIL "A memory request failed." 983 #define MSG_NO_MEM "ida_mem = NULL illegal." 984 #define MSG_NO_MALLOC "Attempt to call before IDAMalloc." 985 #define MSG_BAD_NVECTOR "A required vector operation is not implemented." 986 987 /* Initialization errors */ 988 989 #define MSG_Y0_NULL "y0 = NULL illegal." 990 #define MSG_YP0_NULL "yp0 = NULL illegal." 991 #define MSG_BAD_ITOL "Illegal value for itol. The legal values are IDA_SS, IDA_SV, and IDA_WF." 992 #define MSG_RES_NULL "res = NULL illegal." 993 #define MSG_BAD_RTOL "rtol < 0 illegal." 994 #define MSG_ATOL_NULL "atol = NULL illegal." 995 #define MSG_BAD_ATOL "Some atol component < 0.0 illegal." 996 #define MSG_ROOT_FUNC_NULL "g = NULL illegal." 997 998 #define MSG_MISSING_ID "id = NULL but suppressalg option on." 999 #define MSG_NO_TOLS "No integration tolerances have been specified." 1000 #define MSG_FAIL_EWT "The user-provide EwtSet function failed." 1001 #define MSG_BAD_EWT "Some initial ewt component = 0.0 illegal." 1002 #define MSG_Y0_FAIL_CONSTR "y0 fails to satisfy constraints." 1003 #define MSG_BAD_ISM_CONSTR "Constraints can not be enforced while forward sensitivity is used with simultaneous method." 1004 #define MSG_LSOLVE_NULL "The linear solver's solve routine is NULL." 1005 #define MSG_LINIT_FAIL "The linear solver's init routine failed." 1006 #define MSG_NLS_INIT_FAIL "The nonlinear solver's init routine failed." 1007 1008 #define MSG_NO_QUAD "Illegal attempt to call before calling IDAQuadInit." 1009 #define MSG_BAD_EWTQ "Initial ewtQ has component(s) equal to zero (illegal)." 1010 #define MSG_BAD_ITOLQ "Illegal value for itolQ. The legal values are IDA_SS and IDA_SV." 1011 #define MSG_NO_TOLQ "No integration tolerances for quadrature variables have been specified." 1012 #define MSG_NULL_ATOLQ "atolQ = NULL illegal." 1013 #define MSG_BAD_RTOLQ "rtolQ < 0 illegal." 1014 #define MSG_BAD_ATOLQ "atolQ has negative component(s) (illegal)." 1015 1016 #define MSG_NO_SENSI "Illegal attempt to call before calling IDASensInit." 1017 #define MSG_BAD_EWTS "Initial ewtS has component(s) equal to zero (illegal)." 1018 #define MSG_BAD_ITOLS "Illegal value for itolS. The legal values are IDA_SS, IDA_SV, and IDA_EE." 1019 #define MSG_NULL_ATOLS "atolS = NULL illegal." 1020 #define MSG_BAD_RTOLS "rtolS < 0 illegal." 1021 #define MSG_BAD_ATOLS "atolS has negative component(s) (illegal)." 1022 #define MSG_BAD_PBAR "pbar has zero component(s) (illegal)." 1023 #define MSG_BAD_PLIST "plist has negative component(s) (illegal)." 1024 #define MSG_BAD_NS "NS <= 0 illegal." 1025 #define MSG_NULL_YYS0 "yyS0 = NULL illegal." 1026 #define MSG_NULL_YPS0 "ypS0 = NULL illegal." 1027 #define MSG_BAD_ISM "Illegal value for ism. Legal values are: IDA_SIMULTANEOUS and IDA_STAGGERED." 1028 #define MSG_BAD_IS "Illegal value for is." 1029 #define MSG_NULL_DKYA "dkyA = NULL illegal." 1030 #define MSG_BAD_DQTYPE "Illegal value for DQtype. Legal values are: IDA_CENTERED and IDA_FORWARD." 1031 #define MSG_BAD_DQRHO "DQrhomax < 0 illegal." 1032 1033 #define MSG_NULL_ABSTOLQS "abstolQS = NULL illegal parameter." 1034 #define MSG_BAD_RELTOLQS "reltolQS < 0 illegal parameter." 1035 #define MSG_BAD_ABSTOLQS "abstolQS has negative component(s) (illegal)." 1036 #define MSG_NO_QUADSENSI "Forward sensitivity analysis for quadrature variables was not activated." 1037 #define MSG_NULL_YQS0 "yQS0 = NULL illegal parameter." 1038 1039 1040 /* IDACalcIC error messages */ 1041 1042 #define MSG_IC_BAD_ICOPT "icopt has an illegal value." 1043 #define MSG_IC_BAD_MAXBACKS "maxbacks <= 0 illegal." 1044 #define MSG_IC_MISSING_ID "id = NULL conflicts with icopt." 1045 #define MSG_IC_TOO_CLOSE "tout1 too close to t0 to attempt initial condition calculation." 1046 #define MSG_IC_BAD_ID "id has illegal values." 1047 #define MSG_IC_BAD_EWT "Some initial ewt component = 0.0 illegal." 1048 #define MSG_IC_RES_NONREC "The residual function failed unrecoverably. " 1049 #define MSG_IC_RES_FAIL "The residual function failed at the first call. " 1050 #define MSG_IC_SETUP_FAIL "The linear solver setup failed unrecoverably." 1051 #define MSG_IC_SOLVE_FAIL "The linear solver solve failed unrecoverably." 1052 #define MSG_IC_NO_RECOVERY "The residual routine or the linear setup or solve routine had a recoverable error, but IDACalcIC was unable to recover." 1053 #define MSG_IC_FAIL_CONSTR "Unable to satisfy the inequality constraints." 1054 #define MSG_IC_FAILED_LINS "The linesearch algorithm failed: step too small or too many backtracks." 1055 #define MSG_IC_CONV_FAILED "Newton/Linesearch algorithm failed to converge." 1056 1057 /* IDASolve error messages */ 1058 1059 #define MSG_YRET_NULL "yret = NULL illegal." 1060 #define MSG_YPRET_NULL "ypret = NULL illegal." 1061 #define MSG_TRET_NULL "tret = NULL illegal." 1062 #define MSG_BAD_ITASK "itask has an illegal value." 1063 #define MSG_TOO_CLOSE "tout too close to t0 to start integration." 1064 #define MSG_BAD_HINIT "Initial step is not towards tout." 1065 #define MSG_BAD_TSTOP "The value " MSG_TIME_TSTOP " is behind current " MSG_TIME "in the direction of integration." 1066 #define MSG_CLOSE_ROOTS "Root found at and very near " MSG_TIME "." 1067 #define MSG_MAX_STEPS "At " MSG_TIME ", mxstep steps taken before reaching tout." 1068 #define MSG_EWT_NOW_FAIL "At " MSG_TIME "the user-provide EwtSet function failed." 1069 #define MSG_EWT_NOW_BAD "At " MSG_TIME "some ewt component has become <= 0.0." 1070 #define MSG_TOO_MUCH_ACC "At " MSG_TIME "too much accuracy requested." 1071 1072 #define MSG_BAD_T "Illegal value for t. " MSG_TIME_INT 1073 #define MSG_BAD_TOUT "Trouble interpolating at " MSG_TIME_TOUT ". tout too far back in direction of integration." 1074 1075 #define MSG_BAD_K "Illegal value for k." 1076 #define MSG_NULL_DKY "dky = NULL illegal." 1077 #define MSG_NULL_DKYP "dkyp = NULL illegal." 1078 1079 #define MSG_ERR_FAILS "At " MSG_TIME_H "the error test failed repeatedly or with |h| = hmin." 1080 #define MSG_CONV_FAILS "At " MSG_TIME_H "the corrector convergence failed repeatedly or with |h| = hmin." 1081 #define MSG_SETUP_FAILED "At " MSG_TIME "the linear solver setup failed unrecoverably." 1082 #define MSG_SOLVE_FAILED "At " MSG_TIME "the linear solver solve failed unrecoverably." 1083 #define MSG_REP_RES_ERR "At " MSG_TIME "repeated recoverable residual errors." 1084 #define MSG_RES_NONRECOV "At " MSG_TIME "the residual function failed unrecoverably." 1085 #define MSG_FAILED_CONSTR "At " MSG_TIME "unable to satisfy inequality constraints." 1086 #define MSG_RTFUNC_FAILED "At " MSG_TIME ", the rootfinding routine failed in an unrecoverable manner." 1087 #define MSG_NO_ROOT "Rootfinding was not initialized." 1088 #define MSG_INACTIVE_ROOTS "At the end of the first step, there are still some root functions identically 0. This warning will not be issued again." 1089 #define MSG_NLS_INPUT_NULL "At " MSG_TIME ", the nonlinear solver was passed a NULL input." 1090 #define MSG_NLS_SETUP_FAILED "At " MSG_TIME ", the nonlinear solver setup failed unrecoverably." 1091 #define MSG_NLS_FAIL "At " MSG_TIME ", the nonlinear solver failed in an unrecoverable manner." 1092 1093 1094 #define MSG_EWTQ_NOW_BAD "At " MSG_TIME ", a component of ewtQ has become <= 0." 1095 #define MSG_QRHSFUNC_FAILED "At " MSG_TIME ", the quadrature right-hand side routine failed in an unrecoverable manner." 1096 #define MSG_QRHSFUNC_UNREC "At " MSG_TIME ", the quadrature right-hand side failed in a recoverable manner, but no recovery is possible." 1097 #define MSG_QRHSFUNC_REPTD "At " MSG_TIME "repeated recoverable quadrature right-hand side function errors." 1098 #define MSG_QRHSFUNC_FIRST "The quadrature right-hand side routine failed at the first call." 1099 1100 #define MSG_NULL_P "p = NULL when using internal DQ for sensitivity residual is illegal." 1101 #define MSG_EWTS_NOW_BAD "At " MSG_TIME ", a component of ewtS has become <= 0." 1102 #define MSG_SRES_FAILED "At " MSG_TIME ", the sensitivity residual routine failed in an unrecoverable manner." 1103 #define MSG_SRES_UNREC "At " MSG_TIME ", the sensitivity residual failed in a recoverable manner, but no recovery is possible." 1104 #define MSG_SRES_REPTD "At " MSG_TIME "repeated recoverable sensitivity residual function errors." 1105 1106 #define MSG_NO_TOLQS "No integration tolerances for quadrature sensitivity variables have been specified." 1107 #define MSG_NULL_RHSQ "IDAS is expected to use DQ to evaluate the RHS of quad. sensi., but quadratures were not initialized." 1108 #define MSG_BAD_EWTQS "Initial ewtQS has component(s) equal to zero (illegal)." 1109 #define MSG_EWTQS_NOW_BAD "At " MSG_TIME ", a component of ewtQS has become <= 0." 1110 #define MSG_QSRHSFUNC_FAILED "At " MSG_TIME ", the sensitivity quadrature right-hand side routine failed in an unrecoverable manner." 1111 #define MSG_QSRHSFUNC_REPTD "At " MSG_TIME "repeated recoverable sensitivity quadrature right-hand side function errors." 1112 #define MSG_QSRHSFUNC_FIRST "The quadrature right-hand side routine failed at the first call." 1113 1114 /* IDASet* / IDAGet* error messages */ 1115 #define MSG_NEG_MAXORD "maxord<=0 illegal." 1116 #define MSG_BAD_MAXORD "Illegal attempt to increase maximum order." 1117 #define MSG_NEG_HMAX "hmax < 0 illegal." 1118 #define MSG_NEG_EPCON "epcon <= 0.0 illegal." 1119 #define MSG_BAD_CONSTR "Illegal values in constraints vector." 1120 #define MSG_BAD_EPICCON "epiccon <= 0.0 illegal." 1121 #define MSG_BAD_MAXNH "maxnh <= 0 illegal." 1122 #define MSG_BAD_MAXNJ "maxnj <= 0 illegal." 1123 #define MSG_BAD_MAXNIT "maxnit <= 0 illegal." 1124 #define MSG_BAD_STEPTOL "steptol <= 0.0 illegal." 1125 1126 #define MSG_TOO_LATE "IDAGetConsistentIC can only be called before IDASolve." 1127 1128 /* 1129 * ================================================================= 1130 * I D A A E R R O R M E S S A G E S 1131 * ================================================================= 1132 */ 1133 1134 #define MSGAM_NULL_IDAMEM "ida_mem = NULL illegal." 1135 #define MSGAM_NO_ADJ "Illegal attempt to call before calling IDAadjInit." 1136 #define MSGAM_BAD_INTERP "Illegal value for interp." 1137 #define MSGAM_BAD_STEPS "Steps nonpositive illegal." 1138 #define MSGAM_BAD_WHICH "Illegal value for which." 1139 #define MSGAM_NO_BCK "No backward problems have been defined yet." 1140 #define MSGAM_NO_FWD "Illegal attempt to call before calling IDASolveF." 1141 #define MSGAM_BAD_TB0 "The initial time tB0 is outside the interval over which the forward problem was solved." 1142 #define MSGAM_BAD_SENSI "At least one backward problem requires sensitivities, but they were not stored for interpolation." 1143 #define MSGAM_BAD_ITASKB "Illegal value for itaskB. Legal values are IDA_NORMAL and IDA_ONE_STEP." 1144 #define MSGAM_BAD_TBOUT "The final time tBout is outside the interval over which the forward problem was solved." 1145 #define MSGAM_BACK_ERROR "Error occured while integrating backward problem # %d" 1146 #define MSGAM_BAD_TINTERP "Bad t = %g for interpolation." 1147 #define MSGAM_BAD_T "Bad t for interpolation." 1148 #define MSGAM_WRONG_INTERP "This function cannot be called for the specified interp type." 1149 #define MSGAM_MEM_FAIL "A memory request failed." 1150 #define MSGAM_NO_INITBS "Illegal attempt to call before calling IDAInitBS." 1151 1152 #ifdef __cplusplus 1153 } 1154 #endif 1155 1156 #endif 1157