Actual source code: petschpddm.h

  1: #ifndef PETSCHPDDM_H
  2: #define PETSCHPDDM_H

  4: #include <petsc/private/kspimpl.h>

  6: #define PETSC_KSPHPDDM_DEFAULT_PRECISION \
  7:   (PetscDefined(USE_REAL_SINGLE) ? KSP_HPDDM_PRECISION_SINGLE : (PetscDefined(USE_REAL_DOUBLE) ? KSP_HPDDM_PRECISION_DOUBLE : (PetscDefined(USE_REAL___FLOAT128) ? KSP_HPDDM_PRECISION_QUADRUPLE : KSP_HPDDM_PRECISION_HALF)))
  8: #define PETSC_PCHPDDM_MAXLEVELS 9
  9: PETSC_EXTERN PetscLogEvent  PC_HPDDM_PtAP;
 10: PETSC_EXTERN PetscLogEvent  PC_HPDDM_PtBP;
 11: PETSC_EXTERN PetscLogEvent  PC_HPDDM_Next;
 12: PETSC_INTERN PetscErrorCode HPDDMLoadDL_Private(PetscBool *);

 14: namespace HPDDM
 15: {
 16: template <class>
 17: class Schwarz; /* forward definitions of two needed HPDDM classes */
 18: class PETScOperator;
 19: } // namespace HPDDM

 21: struct PC_HPDDM_Level {
 22:   VecScatter                   scatter;   /* scattering from PETSc nonoverlapping numbering to HPDDM overlapping */
 23:   Vec                         *v[2];      /* working vectors */
 24:   Mat                          V[3];      /* working matrices */
 25:   KSP                          ksp;       /* KSP coupling the action of pc and P */
 26:   PC                           pc;        /* inner fine-level PC, acting like a multigrid smoother */
 27:   HPDDM::Schwarz<PetscScalar> *P;         /* coarse-level HPDDM solver */
 28:   Vec                          D;         /* partition of unity */
 29:   PetscReal                    threshold; /* threshold for selecting local deflation vectors */
 30:   PetscInt                     nu;        /* number of local deflation vectors */
 31:   const struct PC_HPDDM       *parent;    /* parent PC */
 32: };

 34: struct PC_HPDDM {
 35:   PC_HPDDM_Level            **levels;                                       /* array of shells */
 36:   Mat                         aux;                                          /* local auxiliary matrix defined at the finest level on PETSC_COMM_SELF */
 37:   Mat                         B;                                            /* right-hand side matrix defined at the finest level on PETSC_COMM_SELF */
 38:   Vec                         normal;                                       /* temporary Vec when preconditioning the normal equations with KSPLSQR */
 39:   IS                          is;                                           /* global numbering of the auxiliary matrix */
 40:   PetscInt                    N;                                            /* number of levels */
 41:   PCHPDDMCoarseCorrectionType correction;                                   /* type of coarse correction */
 42:   PetscBool                   Neumann;                                      /* aux is the local Neumann matrix? */
 43:   PetscBool                   log_separate;                                 /* separate events for each level? */
 44:   PetscBool                   share;                                        /* shared subdomain KSP between SLEPc and PETSc? */
 45:   PetscBool                   deflation;                                    /* aux is the local deflation space? */
 46:   PetscErrorCode (*setup)(Mat, PetscReal, Vec, Vec, PetscReal, IS, void *); /* setup function for the auxiliary matrix */
 47:   void *setup_ctx;                                                          /* context for setup */
 48: };

 50: struct KSP_HPDDM {
 51:   HPDDM::PETScOperator *op;
 52:   PetscReal             rcntl[1];
 53:   int                   icntl[2];
 54:   unsigned short        scntl[2];
 55:   char                  cntl[5];
 56:   KSPHPDDMPrecision     precision;
 57: };

 59: PETSC_INTERN const char HPDDMCitation[];
 60: PETSC_INTERN PetscBool  HPDDMCite;

 62: #include <HPDDM.hpp>

 64: #endif /* PETSCHPDDM_H */