Actual source code: gmresimpl.h
1: /*
2: Private data structure used by the GMRES method. This data structure
3: must be identical to the beginning of the KSP_FGMRES data structure
4: so if you CHANGE anything here you must also change it there.
5: */
9: #include <private/kspimpl.h> /*I "petscksp.h" I*/
11: typedef struct {
12: /* Hessenberg matrix and orthogonalization information. Hes holds
13: the original (unmodified) hessenberg matrix which may be used
14: to estimate the Singular Values of the matrix */
15: PetscScalar *hh_origin,*hes_origin,*cc_origin,*ss_origin,*rs_origin;
17: PetscScalar *orthogwork; /* holds dot products computed in orthogonalization */
19: /* Work space for computing eigenvalues/singular values */
20: PetscReal *Dsvd;
21: PetscScalar *Rsvd;
22:
24: PetscReal haptol; /* tolerance for happy ending */
25: PetscInt max_k; /* number of vectors in Krylov space, restart size */
27: PetscErrorCode (*orthog)(KSP,PetscInt);
28: KSPGMRESCGSRefinementType cgstype;
29:
30: Vec *vecs; /* the work vectors */
31: PetscInt q_preallocate,delta_allocate;
32: PetscInt vv_allocated; /* number of allocated gmres direction vectors */
33: PetscInt vecs_allocated; /* total number of vecs available */
34: /* Since we may call the user "obtain_work_vectors" several times,
35: we have to keep track of the pointers that it has returned */
36: Vec **user_work;
37: PetscInt *mwork_alloc; /* Number of work vectors allocated as part of a work-vector chunck */
38: PetscInt nwork_alloc; /* Number of work vector chunks allocated */
40: PetscInt it; /* Current iteration: inside restart */
41: PetscScalar *nrs; /* temp that holds the coefficients of the Krylov vectors that form the minimum residual solution */
42: Vec sol_temp; /* used to hold temporary solution */
43: } KSP_GMRES;
45: #define HH(a,b) (gmres->hh_origin + (b)*(gmres->max_k+2)+(a))
46: #define HES(a,b) (gmres->hes_origin + (b)*(gmres->max_k+1)+(a))
47: #define CC(a) (gmres->cc_origin + (a))
48: #define SS(a) (gmres->ss_origin + (a))
49: #define GRS(a) (gmres->rs_origin + (a))
51: /* vector names */
52: #define VEC_OFFSET 2
53: #define VEC_TEMP gmres->vecs[0]
54: #define VEC_TEMP_MATOP gmres->vecs[1]
55: #define VEC_VV(i) gmres->vecs[VEC_OFFSET+i]
57: #endif