Actual source code: telescope.h

  1: #pragma once

  3: /* Telescope */
  4: typedef enum {
  5:   TELESCOPE_DEFAULT = 0,
  6:   TELESCOPE_DMDA,
  7:   TELESCOPE_DMPLEX,
  8:   TELESCOPE_COARSEDM
  9: } PCTelescopeType;

 11: typedef struct _PC_Telescope *PC_Telescope;
 12: struct _PC_Telescope {
 13:   PetscSubcomm     psubcomm;
 14:   PetscSubcommType subcommtype;
 15:   MPI_Comm         subcomm;
 16:   PetscInt         redfactor; /* factor to reduce comm size by */
 17:   KSP              ksp;
 18:   IS               isin;
 19:   VecScatter       scatter;
 20:   Vec              xred, yred, xtmp;
 21:   Mat              Bred;
 22:   PetscBool        ignore_dm, ignore_kspcomputeoperators, use_coarse_dm;
 23:   PCTelescopeType  sr_type;
 24:   void            *dm_ctx;
 25:   PetscErrorCode (*pctelescope_setup_type)(PC, PC_Telescope);
 26:   PetscErrorCode (*pctelescope_matcreate_type)(PC, PC_Telescope, MatReuse, Mat *);
 27:   PetscErrorCode (*pctelescope_matnullspacecreate_type)(PC, PC_Telescope, Mat);
 28:   PetscErrorCode (*pctelescope_reset_type)(PC);
 29: };

 31: /* DMDA */
 32: typedef struct {
 33:   DM        dmrepart;
 34:   Mat       permutation;
 35:   Vec       xp;
 36:   PetscInt  Mp_re, Np_re, Pp_re;
 37:   PetscInt *range_i_re, *range_j_re, *range_k_re;
 38:   PetscInt *start_i_re, *start_j_re, *start_k_re;
 39: } PC_Telescope_DMDACtx;

 41: static inline PetscBool PetscSubcomm_isActiveRank(PetscSubcomm scomm)
 42: {
 43:   if (scomm->color == 0) return PETSC_TRUE;
 44:   else return PETSC_FALSE;
 45: }

 47: static inline PetscBool PCTelescope_isActiveRank(PC_Telescope sred)
 48: {
 49:   if (sred->psubcomm) return PetscSubcomm_isActiveRank(sred->psubcomm);
 50:   else {
 51:     if (sred->subcomm != MPI_COMM_NULL) return PETSC_TRUE;
 52:     else return PETSC_FALSE;
 53:   }
 54: }

 56: PetscErrorCode PCTelescopeSetUp_dmda(PC, PC_Telescope);
 57: PetscErrorCode PCTelescopeMatCreate_dmda(PC, PC_Telescope, MatReuse, Mat *);
 58: PetscErrorCode PCTelescopeMatNullSpaceCreate_dmda(PC, PC_Telescope, Mat);
 59: PetscErrorCode PCApply_Telescope_dmda(PC, Vec, Vec);
 60: PetscErrorCode PCApplyRichardson_Telescope_dmda(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *);
 61: PetscErrorCode PCReset_Telescope_dmda(PC);
 62: PetscErrorCode PCTelescopeSetUp_CoarseDM(PC, PC_Telescope);
 63: PetscErrorCode PCApply_Telescope_CoarseDM(PC, Vec, Vec);
 64: PetscErrorCode PCTelescopeMatNullSpaceCreate_CoarseDM(PC, PC_Telescope, Mat);
 65: PetscErrorCode PCReset_Telescope_CoarseDM(PC);
 66: PetscErrorCode PCApplyRichardson_Telescope_CoarseDM(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *);
 67: PetscErrorCode DMView_DA_Short(DM, PetscViewer);