1 /* A. Baker */ 2 /* 3 Private data structure used by the LGMRES method. 4 */ 5 6 #if !defined(__LGMRES) 7 #define __LGMRES 8 9 #include <private/kspimpl.h> /*includes petscksp.h */ 10 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 */ 21 22 PetscScalar *orthogwork; /* holds dot products computed in orthogonalization */ 23 24 /* Work space for computing eigenvalues/singular values */ 25 PetscReal *Dsvd; 26 PetscScalar *Rsvd; 27 28 /* parameters */ 29 PetscReal haptol; /* tolerance used for the "HAPPY BREAK DOWN" */ 30 PetscInt max_k; /* maximum size of the approximation space 31 before restarting */ 32 33 PetscErrorCode (*orthog)(KSP,PetscInt); /* orthogonalization function to use */ 34 KSPGMRESCGSRefinementType cgstype; 35 36 Vec *vecs; /* holds the work vectors */ 37 38 PetscInt q_preallocate; /* 0 = don't pre-allocate space for work vectors */ 39 PetscInt delta_allocate; /* the number of vectors to allocate in each block 40 if not pre-allocated */ 41 PetscInt vv_allocated; /* vv_allocated is the number of allocated lgmres 42 direction vectors */ 43 44 PetscInt vecs_allocated; /* vecs_allocated is the total number of vecs 45 available - used to simplify the dynamic 46 allocation of vectors */ 47 48 Vec **user_work; /* Since we may call the user "obtain_work_vectors" 49 several times, we have to keep track of the pointers 50 that it has returned (so that we may free the 51 storage) */ 52 53 PetscInt *mwork_alloc; /* Number of work vectors allocated as part of 54 a work-vector chunck */ 55 PetscInt nwork_alloc; /* Number of work-vector chunks allocated */ 56 57 58 /* In order to allow the solution to be constructed during the solution 59 process, we need some additional information: */ 60 61 PetscInt it; /* Current iteration */ 62 PetscScalar *nrs; /* temp that holds the coefficients of the 63 Krylov vectors that form the minimum residual 64 solution */ 65 Vec sol_temp; /* used to hold temporary solution */ 66 67 68 /* LGMRES_MOD - make these for the z vectors - new storage for lgmres */ 69 Vec *augvecs; /* holds the error approximation vectors for lgmres. */ 70 Vec **augvecs_user_work; /* same purpose as user_work above, but this one is 71 for our error approx vectors */ 72 /* currently only augvecs_user_work[0] is used, not sure if this will be */ 73 /* extended in the future to use more, or if this is a design bug */ 74 PetscInt aug_vv_allocated; /* aug_vv_allocated is the number of allocated lgmres 75 augmentation vectors */ 76 PetscInt aug_vecs_allocated; /* aug_vecs_allocated is the total number of augmentation vecs 77 available - used to simplify the dynamic 78 allocation of vectors */ 79 PetscScalar *hwork; /* work array to hold Hessenberg product */ 80 81 PetscInt augwork_alloc; /*size of chunk allocated for augmentation vectors */ 82 83 PetscInt aug_dim; /* max number of augmented directions to add */ 84 85 PetscInt aug_ct; /* number of aug. vectors available */ 86 87 PetscInt *aug_order; /*keeps track of order to use aug. vectors*/ 88 89 PetscInt approx_constant; /* = 1 then the approx space at each restart will 90 be size max_k . Therefore, more than (max_k - aug_dim) 91 krylov vectors may be used if less than aug_dim error 92 approximations are available (in the first few restarts, 93 for example) to keep the space a constant size. */ 94 95 PetscInt matvecs; /*keep track of matvecs */ 96 } KSP_LGMRES; 97 98 99 #define HH(a,b) (lgmres->hh_origin + (b)*(lgmres->max_k+2)+(a)) 100 /* HH will be size (max_k+2)*(max_k+1) - think of HH as 101 being stored columnwise (inc. zeros) for access purposes. */ 102 #define HES(a,b) (lgmres->hes_origin + (b)*(lgmres->max_k+1)+(a)) 103 /* HES will be size (max_k + 1) * (max_k + 1) - 104 again, think of HES as being stored columnwise */ 105 #define CC(a) (lgmres->cc_origin + (a)) /* CC will be length (max_k+1) - cosines */ 106 #define SS(a) (lgmres->ss_origin + (a)) /* SS will be length (max_k+1) - sines */ 107 #define GRS(a) (lgmres->rs_origin + (a)) /* GRS will be length (max_k+2) - rt side */ 108 109 /* vector names */ 110 #define VEC_OFFSET 2 111 #define VEC_TEMP lgmres->vecs[0] /* work space */ 112 #define VEC_TEMP_MATOP lgmres->vecs[1] /* work space */ 113 #define VEC_VV(i) lgmres->vecs[VEC_OFFSET+i] /* use to access 114 othog basis vectors */ 115 /*LGMRES_MOD */ 116 #define AUG_OFFSET 1 117 #define AUGVEC(i) lgmres->augvecs[AUG_OFFSET+i] /*error approx vecors */ 118 #define AUG_ORDER(i) lgmres->aug_order[i] /*order in which to augment */ 119 #define A_AUGVEC(i) lgmres->augvecs[AUG_OFFSET+i+lgmres->aug_dim] /*A times error vector */ 120 #define AUG_TEMP lgmres->augvecs[0] /* work vector */ 121 #endif 122 123 124 125