18 #ifdef DEAL_II_WITH_PETSC
28 # include <petscconf.h>
61 PetscErrorCode ierr = PCDestroy(&
pc);
71 PetscErrorCode ierr = PCApply(
pc, src, dst);
80 PetscErrorCode ierr = PCApplyTranspose(
pc, src, dst);
89 PetscErrorCode ierr = PCSetUp(
pc);
96 static MPI_Comm
comm = PETSC_COMM_SELF;
97 MPI_Comm pcomm = PetscObjectComm(
reinterpret_cast<PetscObject
>(
pc));
98 if (pcomm != MPI_COMM_NULL)
111 PetscErrorCode ierr = PetscObjectGetComm(
112 reinterpret_cast<PetscObject
>(
static_cast<const Mat &
>(
matrix)), &
comm);
117 # if DEAL_II_PETSC_VERSION_LT(3, 5, 0)
129 PetscErrorCode ierr = PCCreate(
comm, &
pc);
154 PetscErrorCode ierr = PCCreate(
comm, &
pc);
176 PetscErrorCode ierr = PCSetType(
pc,
const_cast<char *
>(PCJACOBI));
179 ierr = PCSetFromOptions(
pc);
205 const MPI_Comm &
comm,
211 PetscErrorCode ierr = PCCreate(
comm, &
pc);
232 PetscErrorCode ierr = PCSetType(
pc,
const_cast<char *
>(PCBJACOBI));
235 ierr = PCSetFromOptions(
pc);
286 PetscErrorCode ierr = PCSetType(
pc,
const_cast<char *
>(PCSOR));
293 ierr = PCSetFromOptions(
pc);
330 PetscErrorCode ierr = PCSetType(
pc,
const_cast<char *
>(PCSOR));
338 ierr = PCSORSetSymmetric(
pc, SOR_SYMMETRIC_SWEEP);
341 ierr = PCSetFromOptions(
pc);
378 PetscErrorCode ierr = PCSetType(
pc,
const_cast<char *
>(PCICC));
385 ierr = PCSetFromOptions(
pc);
422 PetscErrorCode ierr = PCSetType(
pc,
const_cast<char *
>(PCILU));
429 ierr = PCSetFromOptions(
pc);
437 const bool symmetric_operator,
438 const double strong_threshold,
439 const double max_row_sum,
440 const unsigned int aggressive_coarsening_num_levels,
441 const bool output_details,
445 const unsigned int n_sweeps_coarse,
447 const unsigned int max_iter,
449 : symmetric_operator(symmetric_operator)
450 , strong_threshold(strong_threshold)
451 , max_row_sum(max_row_sum)
452 , aggressive_coarsening_num_levels(aggressive_coarsening_num_levels)
453 , output_details(output_details)
454 , relaxation_type_up(relaxation_type_up)
455 , relaxation_type_down(relaxation_type_down)
456 , relaxation_type_coarse(relaxation_type_coarse)
457 , n_sweeps_coarse(n_sweeps_coarse)
465 # ifdef DEAL_II_PETSC_WITH_HYPRE
476 std::string string_type;
478 switch (relaxation_type)
481 string_type =
"Jacobi";
485 string_type =
"sequential-Gauss-Seidel";
489 string_type =
"seqboundary-Gauss-Seidel";
492 string_type =
"SOR/Jacobi";
496 string_type =
"backward-SOR/Jacobi";
500 string_type =
"symmetric-SOR/Jacobi";
504 string_type =
" l1scaled-SOR/Jacobi";
508 string_type =
"Gaussian-elimination";
512 string_type =
"l1-Gauss-Seidel";
516 string_type =
"backward-l1-Gauss-Seidel";
522 string_type =
"Chebyshev";
525 string_type =
"FCF-Jacobi";
529 string_type =
"l1scaled-Jacobi";
532 string_type =
"None";
551 const MPI_Comm &
comm,
557 PetscErrorCode ierr = PCCreate(
comm, &
pc);
560 # ifdef DEAL_II_PETSC_WITH_HYPRE
565 ExcMessage(
"Your PETSc installation does not include a copy of "
566 "the hypre package necessary for this preconditioner."));
585 # ifdef DEAL_II_PETSC_WITH_HYPRE
586 PetscErrorCode ierr = PCSetType(
pc,
const_cast<char *
>(PCHYPRE));
589 ierr = PCHYPRESetType(
pc,
"boomeramg");
633 auto relaxation_type_is_symmetric =
650 ExcMessage(
"Use a symmetric smoother for relaxation_type_up"));
655 ExcMessage(
"Use a symmetric smoother for relaxation_type_down"));
660 ExcMessage(
"Use a symmetric smoother for relaxation_type_coarse"));
681 ierr = PCSetFromOptions(
pc);
685 ExcMessage(
"Your PETSc installation does not include a copy of "
686 "the hypre package necessary for this preconditioner."));
696 # ifdef DEAL_II_PETSC_WITH_HYPRE
706 (void)additional_data_;
708 ExcMessage(
"Your PETSc installation does not include a copy of "
709 "the hypre package necessary for this preconditioner."));
718 const unsigned int n_levels,
719 const double threshold,
721 const bool output_details)
724 , threshold(threshold)
726 , output_details(output_details)
754 # ifdef DEAL_II_PETSC_WITH_HYPRE
757 PetscErrorCode ierr = PCSetType(
pc,
const_cast<char *
>(PCHYPRE));
760 ierr = PCHYPRESetType(
pc,
"parasails");
771 "ParaSails parameter symmetric can only be equal to 0, 1, 2!"));
773 std::stringstream ssStream;
779 ssStream <<
"nonsymmetric";
791 ssStream <<
"nonsymmetric,SPD";
799 "ParaSails parameter symmetric can only be equal to 0, 1, 2!"));
815 ierr = PCSetFromOptions(
pc);
821 ExcMessage(
"Your PETSc installation does not include a copy of "
822 "the hypre package necessary for this preconditioner."));
853 PetscErrorCode ierr = PCSetType(
pc,
const_cast<char *
>(PCNONE));
856 ierr = PCSetFromOptions(
pc);
864 const double zero_pivot,
865 const double damping)
867 , zero_pivot(zero_pivot)
897 PetscErrorCode ierr = PCSetType(
pc,
const_cast<char *
>(PCLU));
910 ierr = PCSetFromOptions(
pc);
918 const bool use_vertices,
919 const bool use_edges,
920 const bool use_faces,
923 : use_vertices(use_vertices)
924 , use_edges(use_edges)
925 , use_faces(use_faces)
947 PetscErrorCode ierr = PCCreate(
comm, &
pc);
969 # if DEAL_II_PETSC_VERSION_GTE(3, 10, 0)
970 PetscErrorCode ierr = PCSetType(pc,
const_cast<char *
>(PCBDDC));
976 MatType current_type;
980 ierr = PCGetOperators(pc, &
A, &P);
982 ierr = PCGetUseAmat(pc, &flg);
985 ierr = MatGetType(flg ?
A : P, ¤t_type);
988 strcmp(current_type, MATIS) == 0,
990 "Matrix must be of IS type. For this, the variant of reinit that includes the active dofs must be used."));
994 std::stringstream ssStream;
996 if (additional_data.use_vertices)
1000 if (additional_data.use_edges)
1004 if (additional_data.use_faces)
1008 if (additional_data.symmetric)
1012 if (additional_data.coords.size() > 0)
1016 std::vector<PetscReal> coords_petsc(additional_data.coords.size() *
1018 for (
unsigned int i = 0, j = 0; i < additional_data.coords.size(); ++i)
1020 for (j = 0; j < dim; ++j)
1021 coords_petsc[dim * i + j] = additional_data.coords[i][j];
1024 ierr = PCSetCoordinates(pc,
1026 additional_data.coords.size(),
1027 coords_petsc.data());
1033 ierr = PCSetCoordinates(pc, 0, 0,
nullptr);
1038 ierr = PCSetFromOptions(pc);
1042 false,
ExcMessage(
"BDDC preconditioner requires PETSc 3.10.0 or newer"));
1055 additional_data = additional_data_;
1057 create_pc_with_mat(matrix_);
1076 PetscErrorCode ierr;
1079 ierr = PCDestroy(&
pc);
1084 ierr = PCSetType(
pc, PCSHELL);
1086 ierr = PCShellSetContext(
pc,
static_cast<void *
>(
this));
1092 ierr = PCShellSetName(
pc,
"deal.II user solve");
1100 PetscErrorCode ierr;
1101 # if DEAL_II_PETSC_VERSION_LT(3, 5, 0)
1102 ierr = PCSetOperators(
pc,
matrix,
matrix, DIFFERENT_NONZERO_PATTERN);
1114 PetscFunctionBeginUser;
1115 PetscErrorCode ierr = PCShellGetContext(ppc, &ctx);
1121 "std::function vmult"));
1125 user->vmult(dst, src);
1126 PetscFunctionReturn(0);
1134 PetscFunctionBeginUser;
1135 PetscErrorCode ierr = PCShellGetContext(ppc, &ctx);
1141 "std::function vmultT"));
1145 user->vmult(dst, src);
1146 PetscFunctionReturn(0);
AdditionalData additional_data
void Tvmult(VectorBase &dst, const VectorBase &src) const
const PC & get_pc() const
const MPI_Comm & get_mpi_communicator() const
virtual ~PreconditionBase()
void vmult(VectorBase &dst, const VectorBase &src) const
void create_pc_with_comm(const MPI_Comm &)
void create_pc_with_mat(const MatrixBase &)
PreconditionBlockJacobi()
AdditionalData additional_data
AdditionalData additional_data
AdditionalData additional_data
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
AdditionalData additional_data
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
AdditionalData additional_data
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
AdditionalData additional_data
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
AdditionalData additional_data
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
AdditionalData additional_data
AdditionalData additional_data
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
AdditionalData additional_data
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
static int pcapply_transpose(PC pc, Vec src, Vec dst)
static int pcapply(PC pc, Vec src, Vec dst)
PreconditionShell()=default
void initialize(const MPI_Comm &comm)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcInvalidState()
static ::ExceptionBase & ExcFunctionNotProvided(std::string arg1)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
@ matrix
Contents is actually a matrix.
@ symmetric
Matrix is symmetric.
void set_option_value(const std::string &name, const std::string &value)
AdditionalData(const bool use_vertices=true, const bool use_edges=false, const bool use_faces=false, const bool symmetric=false, const std::vector< Point< dim >> coords={})
RelaxationType relaxation_type_up
unsigned int n_sweeps_coarse
AdditionalData(const bool symmetric_operator=false, const double strong_threshold=0.25, const double max_row_sum=0.9, const unsigned int aggressive_coarsening_num_levels=0, const bool output_details=false, const RelaxationType relaxation_type_up=RelaxationType::SORJacobi, const RelaxationType relaxation_type_down=RelaxationType::SORJacobi, const RelaxationType relaxation_type_coarse=RelaxationType::GaussianElimination, const unsigned int n_sweeps_coarse=1, const double tol=0.0, const unsigned int max_iter=1, const bool w_cycle=false)
unsigned int aggressive_coarsening_num_levels
RelaxationType relaxation_type_down
RelaxationType relaxation_type_coarse
AdditionalData(const unsigned int levels=0)
AdditionalData(const unsigned int levels=0)
AdditionalData(const double pivoting=1.e-6, const double zero_pivot=1.e-12, const double damping=0.0)
AdditionalData(const unsigned int symmetric=1, const unsigned int n_levels=1, const double threshold=0.1, const double filter=0.05, const bool output_details=false)
AdditionalData(const double omega=1)
AdditionalData(const double omega=1)