1 #if !defined(PETSCGL_H) 2 #define PETSCGL_H 3 4 #include <petsc/private/tsimpl.h> 5 6 typedef enum {TSGLLEERROR_FORWARD,TSGLLEERROR_BACKWARD} TSGLLEErrorDirection; 7 8 typedef struct _TSGLLEScheme *TSGLLEScheme; 9 struct _TSGLLEScheme { 10 PetscInt p; /* order of the method */ 11 PetscInt q; /* stage-order of the method */ 12 PetscInt r; /* number of items carried between stages */ 13 PetscInt s; /* number of stages */ 14 PetscScalar *c; /* location of the stages */ 15 PetscScalar *a,*b,*u,*v; /* tableau for the method */ 16 17 /* For use in rescale & modify */ 18 PetscScalar *alpha; /* X_n(t_n) - X_{n-1}(t_n) = - alpha^T h^{p+1} x^{(p+1)}(t_n) */ 19 PetscScalar *beta; /* - beta^T h^{p+2} x^{(p+2)}(t_n) */ 20 PetscScalar *gamma; /* - gamma^T h^{p+2} f' x^{(p+1)}(t_n) + O(h^{p+3}) */ 21 22 /* Error estimates */ 23 /* h^{p+1}x^{(p+1)}(t_n) ~= phi[0]*h*Ydot + psi[0]*X[1:] */ 24 /* h^{p+2}x^{(p+2)}(t_n) ~= phi[1]*h*Ydot + psi[1]*X[1:] */ 25 /* h^{p+2}f' x^{(p+1)}(t_n) ~= phi[2]*h*Ydot + psi[2]*X[1:] */ 26 PetscScalar *phi; /* dim=[3][s] for estimating higher moments, see B,J,W 2007 */ 27 PetscScalar *psi; /* dim=[3][r-1], [0 psi^T] of B,J,W 2007 */ 28 PetscScalar *stage_error; 29 30 /* Desirable properties which enable extra optimizations */ 31 PetscBool stiffly_accurate; /* Last row of [A U] is equal t first row of [B V]? */ 32 PetscBool fsal; /* First Same As Last: X[1] = h*Ydot[s-1] (and stiffly accurate) */ 33 }; 34 35 typedef struct TS_GLLE { 36 TSGLLEAcceptFunction Accept; /* Decides whether to accept a given time step, given estimates of local truncation error */ 37 TSGLLEAdapt adapt; 38 39 /* These names are only stored so that they can be printed in TSView_GLLE() without making these schemes full-blown 40 objects (the implementations I'm thinking of do not have state and I'm lazy). */ 41 char accept_name[256]; 42 43 /* specific to the family of GL method */ 44 PetscErrorCode (*EstimateHigherMoments)(TSGLLEScheme,PetscReal,Vec*,Vec*,Vec*); /* Provide local error estimates */ 45 PetscErrorCode (*CompleteStep)(TSGLLEScheme,PetscReal,TSGLLEScheme,PetscReal,Vec*,Vec*,Vec*); 46 PetscErrorCode (*Destroy)(struct TS_GLLE*); 47 PetscErrorCode (*View)(struct TS_GLLE*,PetscViewer); 48 char type_name[256]; 49 PetscInt nschemes; 50 TSGLLEScheme *schemes; 51 52 Vec *X; /* Items to carry between steps */ 53 Vec *Xold; /* Values of these items at the last step */ 54 Vec W; /* = 1/(atol+rtol*|X0|), used for WRMS norm */ 55 Vec *himom; /* len=3, Estimates of h^{p+1}x^{(p+1)}, h^{p+2}x^{(p+2)}, h^{p+2}(df/dx) x^{(p+1)} */ 56 PetscReal wrms_atol,wrms_rtol; 57 58 /* Stages (Y,Ydot) are computed sequentially */ 59 Vec *Ydot; /* Derivatives of stage vectors, must be stored */ 60 Vec Y; /* Stage vector, only used while solving the stage so we don't need to store it */ 61 Vec Z; /* Affine vector */ 62 PetscReal scoeff; /* Ydot = Z + shift*Y; shift = scoeff/ts->time_step */ 63 PetscReal stage_time; /* time at current stage */ 64 PetscInt stage; /* index of the stage we are currently solving for */ 65 66 /* Runtime options */ 67 PetscInt current_scheme; 68 PetscInt max_order,min_order,start_order; 69 PetscBool extrapolate; /* use extrapolation to produce initial Newton iterate? */ 70 TSGLLEErrorDirection error_direction; /* TSGLLEERROR_FORWARD or TSGLLEERROR_BACKWARD */ 71 72 PetscInt max_step_rejections; 73 74 PetscBool setupcalled; 75 void *data; 76 } TS_GLLE; 77 78 #endif 79