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