1*3a2a065fSHong Zhang typedef struct _ARKTableau *ARKTableau; 2*3a2a065fSHong Zhang struct _ARKTableau { 3*3a2a065fSHong Zhang char *name; 4*3a2a065fSHong Zhang PetscBool additive; /* If False, it is a DIRK method */ 5*3a2a065fSHong Zhang PetscInt order; /* Classical approximation order of the method */ 6*3a2a065fSHong Zhang PetscInt s; /* Number of stages */ 7*3a2a065fSHong Zhang PetscBool stiffly_accurate; /* The implicit part is stiffly accurate */ 8*3a2a065fSHong Zhang PetscBool FSAL_implicit; /* The implicit part is FSAL */ 9*3a2a065fSHong Zhang PetscBool explicit_first_stage; /* The implicit part has an explicit first stage */ 10*3a2a065fSHong Zhang PetscInt pinterp; /* Interpolation order */ 11*3a2a065fSHong Zhang PetscReal *At, *bt, *ct; /* Stiff tableau */ 12*3a2a065fSHong Zhang PetscReal *A, *b, *c; /* Non-stiff tableau */ 13*3a2a065fSHong Zhang PetscReal *bembedt, *bembed; /* Embedded formula of order one less (order-1) */ 14*3a2a065fSHong Zhang PetscReal *binterpt, *binterp; /* Dense output formula */ 15*3a2a065fSHong Zhang PetscReal ccfl; /* Placeholder for CFL coefficient relative to forward Euler */ 16*3a2a065fSHong Zhang }; 17*3a2a065fSHong Zhang typedef struct _ARKTableauLink *ARKTableauLink; 18*3a2a065fSHong Zhang struct _ARKTableauLink { 19*3a2a065fSHong Zhang struct _ARKTableau tab; 20*3a2a065fSHong Zhang ARKTableauLink next; 21*3a2a065fSHong Zhang }; 22*3a2a065fSHong Zhang 23*3a2a065fSHong Zhang typedef struct { 24*3a2a065fSHong Zhang ARKTableau tableau; 25*3a2a065fSHong Zhang Vec *Y; /* States computed during the step */ 26*3a2a065fSHong Zhang Vec *YdotI; /* Time derivatives for the stiff part */ 27*3a2a065fSHong Zhang Vec *YdotRHS; /* Function evaluations for the non-stiff part */ 28*3a2a065fSHong Zhang Vec *Y_prev; /* States computed during the previous time step */ 29*3a2a065fSHong Zhang Vec *YdotI_prev; /* Time derivatives for the stiff part for the previous time step*/ 30*3a2a065fSHong Zhang Vec *YdotRHS_prev; /* Function evaluations for the non-stiff part for the previous time step*/ 31*3a2a065fSHong Zhang Vec Ydot0; /* Holds the slope from the previous step in FSAL case */ 32*3a2a065fSHong Zhang Vec Ydot; /* Work vector holding Ydot during residual evaluation */ 33*3a2a065fSHong Zhang Vec Z; /* Ydot = shift(Y-Z) */ 34*3a2a065fSHong Zhang IS alg_is; /* Index set for algebraic variables, needed when restarting with DIRK */ 35*3a2a065fSHong Zhang PetscScalar *work; /* Scalar work */ 36*3a2a065fSHong Zhang PetscReal scoeff; /* shift = scoeff/dt */ 37*3a2a065fSHong Zhang PetscReal stage_time; 38*3a2a065fSHong Zhang PetscBool imex; 39*3a2a065fSHong Zhang PetscBool extrapolate; /* Extrapolate initial guess from previous time-step stage values */ 40*3a2a065fSHong Zhang TSStepStatus status; 41*3a2a065fSHong Zhang 42*3a2a065fSHong Zhang /* context for fast-slow split */ 43*3a2a065fSHong Zhang Vec Y_snes; /* Work vector for SNES */ 44*3a2a065fSHong Zhang Vec *YdotI_fast; /* Function evaluations for the fast components in YdotI */ 45*3a2a065fSHong Zhang Vec *YdotRHS_fast; /* Function evaluations for the fast components in YdotRHS */ 46*3a2a065fSHong Zhang Vec *YdotRHS_slow; /* Function evaluations for the slow components in YdotRHS */ 47*3a2a065fSHong Zhang IS is_slow, is_fast; 48*3a2a065fSHong Zhang TS subts_slow, subts_fast; 49*3a2a065fSHong Zhang PetscBool fastslowsplit; 50*3a2a065fSHong Zhang 51*3a2a065fSHong Zhang /* context for sensitivity analysis */ 52*3a2a065fSHong Zhang Vec *VecsDeltaLam; /* Increment of the adjoint sensitivity w.r.t IC at stage */ 53*3a2a065fSHong Zhang Vec *VecsSensiTemp; /* Vectors to be multiplied with Jacobian transpose */ 54*3a2a065fSHong Zhang Vec *VecsSensiPTemp; /* Temporary Vectors to store JacobianP-transpose-vector product */ 55*3a2a065fSHong Zhang } TS_ARKIMEX; 56