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: */
6: #pragma once
8: #include <petsc/private/kspimpl.h>
10: #define KSPGMRESHEADER \
11: /* Hessenberg matrix and orthogonalization information. */ \
12: PetscScalar *hh_origin; /* holds hessenburg matrix that has been multiplied by plane rotations (upper tri) */ \
13: PetscScalar *hes_origin; /* holds the original (unmodified) Hessenberg matrix which may be used to estimate the Singular Values of the matrix */ \
14: PetscScalar *hes_ritz; /* holds the last full Hessenberg matrix to compute (harmonic) Ritz pairs */ \
15: PetscScalar *cc_origin; /* holds cosines for rotation matrices */ \
16: PetscScalar *ss_origin; /* holds sines for rotation matrices */ \
17: PetscScalar *rs_origin; /* holds the right-hand side of the Hessenberg system */ \
18: \
19: PetscScalar *orthogwork; /* holds dot products computed in orthogonalization */ \
20: \
21: /* Work space for computing eigenvalues/singular values */ \
22: PetscReal *Dsvd; \
23: PetscScalar *Rsvd; \
24: \
25: PetscReal haptol; /* tolerance for happy ending */ \
26: PetscInt max_k; /* number of vectors in Krylov space, restart size */ \
27: PetscInt nextra_vecs; /* number of extra vecs needed, e.g. for a pipeline */ \
28: \
29: PetscErrorCode (*orthog)(KSP, PetscInt); \
30: KSPGMRESCGSRefinementType cgstype; \
31: \
32: Vec *vecs; /* the work vectors */ \
33: Vec *vecb; /* holds the last full basis vectors of the Krylov subspace to compute (harmonic) Ritz pairs */ \
34: PetscInt q_preallocate; /* 0=don't preallocate space for work vectors */ \
35: PetscInt delta_allocate; /* number of vectors to preallocaate in each block if not preallocated */ \
36: PetscInt vv_allocated; /* number of allocated gmres direction vectors */ \
37: PetscInt vecs_allocated; /* total number of vecs available */ \
38: /* Since we may call the user "obtain_work_vectors" several times, we have to keep track of the pointers that it has returned */ \
39: Vec **user_work; \
40: PetscInt *mwork_alloc; /* Number of work vectors allocated as part of a work-vector chunk */ \
41: PetscInt nwork_alloc; /* Number of work vector chunks allocated */ \
42: \
43: /* Information for building solution */ \
44: PetscInt it; /* Current iteration: inside restart */ \
45: PetscInt fullcycle; /* Current number of complete cycle */ \
46: PetscScalar *nrs; /* temp that holds the coefficients of the Krylov vectors that form the minimum residual solution */ \
47: Vec sol_temp; /* used to hold temporary solution */ \
48: PetscReal rnorm0; /* residual norm at beginning of the GMRESCycle */ \
49: PetscReal breakdowntol; /* A relative tolerance is used for breakdown check in GMRESCycle */
51: typedef struct {
52: KSPGMRESHEADER
53: } KSP_GMRES;
55: PETSC_INTERN PetscErrorCode KSPView_GMRES(KSP, PetscViewer);
56: PETSC_INTERN PetscErrorCode KSPSetUp_GMRES(KSP);
57: PETSC_INTERN PetscErrorCode KSPSetFromOptions_GMRES(KSP, PetscOptionItems *PetscOptionsObject);
58: PETSC_INTERN PetscErrorCode KSPComputeExtremeSingularValues_GMRES(KSP, PetscReal *, PetscReal *);
59: PETSC_INTERN PetscErrorCode KSPComputeEigenvalues_GMRES(KSP, PetscInt, PetscReal *, PetscReal *, PetscInt *);
60: PETSC_INTERN PetscErrorCode KSPComputeRitz_GMRES(KSP, PetscBool, PetscBool, PetscInt *, Vec[], PetscReal *, PetscReal *);
61: PETSC_INTERN PetscErrorCode KSPReset_GMRES(KSP);
62: PETSC_INTERN PetscErrorCode KSPDestroy_GMRES(KSP);
63: PETSC_INTERN PetscErrorCode KSPGMRESGetNewVectors(KSP, PetscInt);
65: typedef PetscErrorCode (*FCN)(KSP, PetscInt); /* force argument to next function to not be extern C*/
67: PETSC_INTERN PetscErrorCode KSPGMRESSetHapTol_GMRES(KSP, PetscReal);
68: PETSC_INTERN PetscErrorCode KSPGMRESSetPreAllocateVectors_GMRES(KSP);
69: PETSC_INTERN PetscErrorCode KSPGMRESSetRestart_GMRES(KSP, PetscInt);
70: PETSC_INTERN PetscErrorCode KSPGMRESGetRestart_GMRES(KSP, PetscInt *);
71: PETSC_INTERN PetscErrorCode KSPGMRESSetOrthogonalization_GMRES(KSP, FCN);
72: PETSC_INTERN PetscErrorCode KSPGMRESGetOrthogonalization_GMRES(KSP, FCN *);
73: PETSC_INTERN PetscErrorCode KSPGMRESSetCGSRefinementType_GMRES(KSP, KSPGMRESCGSRefinementType);
74: PETSC_INTERN PetscErrorCode KSPGMRESGetCGSRefinementType_GMRES(KSP, KSPGMRESCGSRefinementType *);
76: /* These macros are guarded because they are redefined by derived implementations */
77: #if !defined(KSPGMRES_NO_MACROS)
78: #define HH(a, b) (gmres->hh_origin + (b) * (gmres->max_k + 2) + (a))
79: #define HES(a, b) (gmres->hes_origin + (b) * (gmres->max_k + 1) + (a))
80: #define CC(a) (gmres->cc_origin + (a))
81: #define SS(a) (gmres->ss_origin + (a))
82: #define GRS(a) (gmres->rs_origin + (a))
84: /* vector names */
85: #define VEC_OFFSET 2
86: #define VEC_TEMP gmres->vecs[0]
87: #define VEC_TEMP_MATOP gmres->vecs[1]
88: #define VEC_VV(i) gmres->vecs[VEC_OFFSET + i]
89: #endif