Actual source code: dgmresimpl.h

petsc-3.10.4 2019-02-26
Report Typos and Errors

  4: #define KSPGMRES_NO_MACROS
  5:  #include <../src/ksp/ksp/impls/gmres/gmresimpl.h>
  6:  #include <petscblaslapack.h>

  8: #define KSPDGMRESHEADER \
  9:   /* Data specific to DGMRES */ \
 10:   Vec           *U;     /* Vectors that form the basis of the invariant subspace */ \
 11:   PetscScalar  *T;      /* T=U^T*M^{-1}*A*U */ \
 12:   PetscScalar  *TF;     /* The factors L and U from T = P*L*U */ \
 13:   PetscBLASInt *InvP;   /* Permutation Vector from the LU factorization of T */ \
 14:   PetscInt     neig;    /* number of eigenvalues to extract at each restart */ \
 15:   PetscInt     r;       /* current number of deflated eigenvalues */ \
 16:   PetscInt     max_neig;   /* Maximum number of eigenvalues to deflate */ \
 17:   PetscReal    lambdaN;    /* modulus of the largest eigenvalue of A */ \
 18:   PetscReal    smv;        /* smaller multiple of the remaining allowed number of steps -- used for the adaptive strategy */ \
 19:   PetscBool    force;      /* Force the use of the deflation at the restart */ \
 20:   PetscInt     matvecs;    /* Total number of matrix-vectors */ \
 21:   PetscInt     GreatestEig;   /* Extract the greatest eigenvalues instead */ \
 22:   PetscReal    *wr, *wi, *modul;  /* Real and complex part and modulus of eigenvalues */ \
 23:   PetscScalar  *Q, *Z;  /* Left and right schur/eigenvectors from the QZ algorithm */ \
 24:   PetscInt     *perm;   /* temporary permutation vector */ \
 25:   /* Work spaces */ \
 26:   Vec          *mu;     /* Save the product M^{-1}AU */ \
 27:   PetscScalar  *Sr;     /* Schur vectors to extract */ \
 28:   Vec          *X;      /* Schurs Vectors Sr projected to the entire space */ \
 29:   Vec          *mx;     /* store the product M^{-1}*A*X */ \
 30:   PetscScalar  *umx;    /* store the product U^T*M^{-1}*A*X */ \
 31:   PetscScalar  *xmu;    /* store the product X^T*M^{-1}*A*U */ \
 32:   PetscScalar  *xmx;    /* store the product X^T*M^{-1}*A*X */ \
 33:   PetscScalar  *x1;     /* store the product U^T*x */ \
 34:   PetscScalar  *x2;     /* store the product U^T*x */ \
 35:   PetscScalar  *Sr2;    /* Schur vectors at the improvement step */ \
 36:   PetscScalar  *auau;   /* product of (M*A*U)^T*M*A*U */ \
 37:   PetscScalar  *auu;    /* product of (M*A*U)^T*U */ \
 38:   PetscScalar  *work;   /* work space for LAPACK functions */ \
 39:   PetscBLASInt *iwork;  /* work space for LAPACK functions */ \
 40:   PetscReal    *orth;   /* Coefficients for the orthogonalization */ \
 41:   PetscBool    HasSchur;   /* Indicate if the Schur form had already been computed in this cycle */ \
 42:   PetscBool    improve;    /* 0 = do not improve the eigenvalues; This is an experimental option */

 44: typedef struct {
 45:   KSPGMRESHEADER
 46:   KSPDGMRESHEADER
 47: } KSP_DGMRES;

 49: PETSC_INTERN PetscErrorCode  KSPDGMRESComputeDeflationData(KSP,PetscInt*);

 51: PETSC_EXTERN PetscLogEvent KSP_DGMRESComputeDeflationData;
 52: PETSC_EXTERN PetscLogEvent KSP_DGMRESApplyDeflation;

 54: #define HH(a,b)  (dgmres->hh_origin + (b)*(dgmres->max_k+2)+(a))
 55: #define HES(a,b) (dgmres->hes_origin + (b)*(dgmres->max_k+1)+(a))
 56: #define CC(a)    (dgmres->cc_origin + (a))
 57: #define SS(a)    (dgmres->ss_origin + (a))
 58: #define GRS(a)   (dgmres->rs_origin + (a))

 60: /* vector names */
 61: #define VEC_OFFSET     2
 62: #define VEC_TEMP       dgmres->vecs[0]
 63: #define VEC_TEMP_MATOP dgmres->vecs[1]
 64: #define VEC_VV(i)      dgmres->vecs[VEC_OFFSET+i]

 66: #define EIG_OFFSET             1
 67: #define DGMRES_DEFAULT_EIG     1
 68: #define DGMRES_DEFAULT_MAXEIG 10

 70: #define UU       dgmres->U
 71: #define TT       dgmres->T
 72: #define TTF      dgmres->TF
 73: #define XX       dgmres->X
 74: #define INVP     dgmres->InvP
 75: #define MU       dgmres->mu
 76: #define MX       dgmres->mx
 77: #define UMX      dgmres->umx
 78: #define XMU      dgmres->xmu
 79: #define XMX      dgmres->xmx
 80: #define X1       dgmres->x1
 81: #define X2       dgmres->x2
 82: #define SR       dgmres->Sr
 83: #define SR2      dgmres->Sr2
 84: #define AUAU     dgmres->auau
 85: #define AUU      dgmres->auu
 86: #define MAX_K    dgmres->max_k
 87: #define MAX_NEIG dgmres->max_neig
 88: #define WORK     dgmres->work
 89: #define IWORK    dgmres->iwork
 90: #define ORTH     dgmres->orth
 91: #define SMV      1
 92: #endif