Actual source code: tsimpl.h

  2: #ifndef __TSIMPL_H

 5:  #include petscts.h

  7: /*
  8:     Timesteping context.
  9:       General DAE: F(t,U,U_t) = 0, required Jacobian is G'(U) where G(U) = F(t,U,U0+a*U)
 10:       General ODE: U_t = F(t,U) <-- the right-hand-side function
 11:       Linear  ODE: U_t = A(t) U <-- the right-hand-side matrix
 12:       Linear (no time) ODE: U_t = A U <-- the right-hand-side matrix
 13: */

 15: /*
 16:      Maximum number of monitors you can run with a single TS
 17: */
 18: #define MAXTSMONITORS 5 

 20: struct _TSOps {
 21:   PetscErrorCode (*rhsmatrix)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*);
 22:   PetscErrorCode (*lhsmatrix)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*);
 23:   PetscErrorCode (*rhsfunction)(TS,PetscReal,Vec,Vec,void*);
 24:   PetscErrorCode (*rhsjacobian)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*);
 25:   PetscErrorCode (*ifunction)(TS,PetscReal,Vec,Vec,Vec,void*);
 26:   PetscErrorCode (*ijacobian)(TS,PetscReal,Vec,Vec,PetscReal,Mat*,Mat*,MatStructure*,void*);
 27:   PetscErrorCode (*prestep)(TS);
 28:   PetscErrorCode (*poststep)(TS);
 29:   PetscErrorCode (*setup)(TS);
 30:   PetscErrorCode (*step)(TS,PetscInt*,PetscReal*);
 31:   PetscErrorCode (*setfromoptions)(TS);
 32:   PetscErrorCode (*destroy)(TS);
 33:   PetscErrorCode (*view)(TS,PetscViewer);
 34: };

 36: struct _p_TS {
 37:   PETSCHEADER(struct _TSOps);
 38:   TSProblemType problem_type;
 39:   Vec           vec_sol,vec_sol_always;

 41:   /* ---------------- User (or PETSc) Provided stuff ---------------------*/
 42:   PetscErrorCode (*monitor[MAXTSMONITORS])(TS,PetscInt,PetscReal,Vec,void*); /* returns control to user after */
 43:   PetscErrorCode (*mdestroy[MAXTSMONITORS])(void*);
 44:   void *monitorcontext[MAXTSMONITORS];                 /* residual calculation, allows user */
 45:   PetscInt  numbermonitors;                                 /* to, for instance, print residual norm, etc. */

 47:   /* ---------------------Linear Iteration---------------------------------*/
 48:   KSP ksp;
 49:   Mat A,B;           /* internel matrix and preconditioner used for KSPSolve() */
 50:   Mat Arhs,Alhs;     /* user provided right/left hand side matrix and preconditioner */
 51:   MatStructure matflg; /* flag indicating the matrix structure of Arhs and Alhs */

 53:   /* ---------------------Nonlinear Iteration------------------------------*/
 54:   SNES  snes;
 55:   void *funP;
 56:   void *jacP,*jacPlhs;
 57:   void *bcP;


 60:   /* --- Data that is unique to each particular solver --- */
 61:   PetscInt setupcalled;            /* true if setup has been called */
 62:   void     *data;                   /* implementationspecific data */
 63:   void     *user;                   /* user context */

 65:   /* ------------------  Parameters -------------------------------------- */
 66:   PetscInt  max_steps;              /* max number of steps */
 67:   PetscReal max_time;               /* max time allowed */
 68:   PetscReal time_step;              /* current time increment */
 69:   PetscReal time_step_old;          /* previous time increment */
 70:   PetscReal initial_time_step;      /* initial time increment */
 71:   PetscInt  steps;                  /* steps taken so far */
 72:   PetscReal ptime;                  /* time taken so far */
 73:   PetscInt  linear_its;             /* total number of linear solver iterations */
 74:   PetscInt  nonlinear_its;          /* total number of nonlinear solver iterations */

 76:   /* ------------------- Default work-area management ------------------ */
 77:   PetscInt nwork;
 78:   Vec      *work;
 79: };

 81: EXTERN PetscErrorCode TSMonitor(TS,PetscInt,PetscReal,Vec);
 82: EXTERN PetscErrorCode TSScaleShiftMatrices(TS,Mat,Mat,MatStructure);


 86: #endif