Actual source code: fgmresp.h
1: /*$Id: fgmresp.h,v 1.7 2001/08/07 03:03:55 balay Exp $*/
2: /*
3: Private data structure used by the FGMRES method.
4: */
9: #include src/sles/ksp/kspimpl.h
11: typedef struct {
12: /* Hessenberg matrix and orthogonalization information. */
13: PetscScalar *hh_origin; /* holds hessenburg matrix that has been
14: multiplied by plane rotations (upper tri) */
15: PetscScalar *hes_origin; /* holds the original (unmodified) hessenberg matrix
16: which may be used to estimate the Singular Values
17: of the matrix */
18: PetscScalar *cc_origin; /* holds cosines for rotation matrices */
19: PetscScalar *ss_origin; /* holds sines for rotation matrices */
20: PetscScalar *rs_origin; /* holds the right-hand-side of the Hessenberg system */
22: /* Work space for computing eigenvalues/singular values */
23: PetscReal *Dsvd;
24: PetscScalar *Rsvd;
25:
26: /* parameters */
27: PetscReal haptol; /* tolerance used for the "HAPPY BREAK DOWN" */
28: int max_k; /* maximum number of Krylov directions to find
29: before restarting */
31: int (*orthog)(KSP,int); /* orthogonalization function to use */
32:
33: Vec *vecs; /* holds the work vectors */
34:
35: int q_preallocate; /* 0 = don't pre-allocate space for work vectors */
36: int delta_allocate; /* the number of vectors to allocate in each block
37: if not pre-allocated */
38: int vv_allocated; /* vv_allocated is the number of allocated fgmres
39: direction vectors */
40:
41: int vecs_allocated; /* vecs_allocated is the total number of vecs
42: available - used to simplify the dynamic
43: allocation of vectors */
44:
45: Vec **user_work; /* Since we may call the user "obtain_work_vectors"
46: several times, we have to keep track of the pointers
47: that it has returned (so that we may free the
48: storage) */
50: int *mwork_alloc; /* Number of work vectors allocated as part of
51: a work-vector chunck */
52: int nwork_alloc; /* Number of work-vector chunks allocated */
55: /* In order to allow the solution to be constructed during the solution
56: process, we need some additional information: */
58: int it; /* Current iteration */
59: PetscScalar *nrs; /* temp that holds the coefficients of the
60: Krylov vectors that form the minimum residual
61: solution */
62: Vec sol_temp; /* used to hold temporary solution */
65: /* new storage for fgmres */
66: Vec *prevecs; /* holds the preconditioned basis vectors for fgmres.
67: We will allocate these at the same time as vecs
68: above (and in the same "chunks". */
69: Vec **prevecs_user_work; /* same purpose as user_work above, but this one is
70: for our preconditioned vectors */
72: /* we need a function for interacting with the pcfamily */
73:
74: int (*modifypc)(KSP,int,int,PetscReal,void*); /* function to modify the preconditioner*/
75: int (*modifydestroy)(void*);
76: void *modifyctx;
77: } KSP_FGMRES;
80: #define HH(a,b) (fgmres->hh_origin + (b)*(fgmres->max_k+2)+(a))
81: /* HH will be size (max_k+2)*(max_k+1) - think of HH as
82: being stored columnwise for access purposes. */
83: #define HES(a,b) (fgmres->hes_origin + (b)*(fgmres->max_k+1)+(a))
84: /* HES will be size (max_k + 1) * (max_k + 1) -
85: again, think of HES as being stored columnwise */
86: #define CC(a) (fgmres->cc_origin + (a)) /* CC will be length (max_k+1) - cosines */
87: #define SS(a) (fgmres->ss_origin + (a)) /* SS will be length (max_k+1) - sines */
88: #define RS(a) (fgmres->rs_origin + (a)) /* RS will be length (max_k+2) - rt side */
90: /* vector names */
91: #define VEC_OFFSET 2
92: #define VEC_SOLN ksp->vec_sol /* solution */
93: #define VEC_RHS ksp->vec_rhs /* right-hand side */
94: #define VEC_TEMP fgmres->vecs[0] /* work space */
95: #define VEC_TEMP_MATOP fgmres->vecs[1] /* work space */
96: #define VEC_VV(i) fgmres->vecs[VEC_OFFSET+i] /* use to access
97: othog basis vectors */
98: #define PREVEC(i) fgmres->prevecs[i] /* use to access
99: preconditioned basis */
101: #endif