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. The beginning of this
  4:  data structure MUST be identical to the beginning of KSP_GMRES since they
  5:  share several functions!
  6: */


 11:  #include src/ksp/ksp/kspimpl.h

 13: typedef struct {
 14:     /* Hessenberg matrix and orthogonalization information. */
 15:     PetscScalar *hh_origin;   /* holds hessenburg matrix that has been
 16:                             multiplied by plane rotations (upper tri) */
 17:     PetscScalar *hes_origin;  /* holds the original (unmodified) hessenberg matrix 
 18:                             which may be used to estimate the Singular Values
 19:                             of the matrix */
 20:     PetscScalar *cc_origin;   /* holds cosines for rotation matrices */
 21:     PetscScalar *ss_origin;   /* holds sines for rotation matrices */
 22:     PetscScalar *rs_origin;   /* holds the right-hand-side of the Hessenberg system */

 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:     int    max_k;             /* maximum number of Krylov directions to find 
 31:                                  before restarting */

 33:     int   (*orthog)(KSP,int); /* orthogonalization function to use */
 34:     KSPGMRESCGSRefinementType cgstype;
 35: 
 36:     Vec   *vecs;              /* holds the work vectors */
 37: 
 38:     int    q_preallocate;     /* 0 = don't pre-allocate space for work vectors */
 39:     int    delta_allocate;    /* the number of vectors to allocate in each block 
 40:                                  if not pre-allocated */
 41:     int    vv_allocated;      /* vv_allocated is the number of allocated fgmres 
 42:                                  direction vectors */
 43: 
 44:     int    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) */

 53:     int    *mwork_alloc;      /* Number of work vectors allocated as part of
 54:                                  a work-vector chunck */
 55:     int    nwork_alloc;       /* Number of work-vector chunks allocated */


 58:     /* In order to allow the solution to be constructed during the solution
 59:        process, we need some additional information: */

 61:     int    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 */


 68:     /* new storage for fgmres */
 69:     Vec  *prevecs;            /* holds the preconditioned basis vectors for fgmres.  
 70:                                 We will allocate these at the same time as vecs 
 71:                                 above (and in the same "chunks". */
 72:     Vec  **prevecs_user_work; /* same purpose as user_work above, but this one is
 73:                                 for our preconditioned vectors */

 75:     /* we need a function for interacting with the pcfamily */
 76: 
 77:     int    (*modifypc)(KSP,int,int,PetscReal,void*);  /* function to modify the preconditioner*/
 78:     int    (*modifydestroy)(void*);
 79:     void   *modifyctx;
 80: } KSP_FGMRES;

 82: #define HH(a,b)  (fgmres->hh_origin + (b)*(fgmres->max_k+2)+(a)) 
 83:                  /* HH will be size (max_k+2)*(max_k+1)  -  think of HH as 
 84:                     being stored columnwise for access purposes. */
 85: #define HES(a,b) (fgmres->hes_origin + (b)*(fgmres->max_k+1)+(a)) 
 86:                   /* HES will be size (max_k + 1) * (max_k + 1) - 
 87:                      again, think of HES as being stored columnwise */
 88: #define CC(a)    (fgmres->cc_origin + (a)) /* CC will be length (max_k+1) - cosines */
 89: #define SS(a)    (fgmres->ss_origin + (a)) /* SS will be length (max_k+1) - sines */
 90: #define RS(a)    (fgmres->rs_origin + (a)) /* RS will be length (max_k+2) - rt side */

 92: /* vector names */
 93: #define VEC_OFFSET     2
 94: #define VEC_SOLN       ksp->vec_sol                  /* solution */ 
 95: #define VEC_RHS        ksp->vec_rhs                  /* right-hand side */
 96: #define VEC_TEMP       fgmres->vecs[0]               /* work space */  
 97: #define VEC_TEMP_MATOP fgmres->vecs[1]               /* work space */
 98: #define VEC_VV(i)      fgmres->vecs[VEC_OFFSET+i]    /* use to access
 99:                                                         othog basis vectors */
100: #define PREVEC(i)      fgmres->prevecs[i]            /* use to access 
101:                                                         preconditioned basis */

103: #endif