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