From: heister Date: Sat, 13 Nov 2010 13:50:55 +0000 (+0000) Subject: implement interface to BoomerAMG through PETSc. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=779321368f41d419b67b3ec00adccaa4454369a2;p=dealii-svn.git implement interface to BoomerAMG through PETSc. git-svn-id: https://svn.dealii.org/trunk@22718 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/lac/petsc_precondition.h b/deal.II/include/deal.II/lac/petsc_precondition.h index 6d23c6e5b9..6f40b5b1c2 100644 --- a/deal.II/include/deal.II/lac/petsc_precondition.h +++ b/deal.II/include/deal.II/lac/petsc_precondition.h @@ -594,9 +594,143 @@ namespace PETScWrappers virtual void set_preconditioner_type (PC &pc) const; }; + + + +/** + * A class that implements the interface to use the BoomerAMG algebraic + * multigrid preconditioner from the HYPRE suite. Note that PETSc has to be + * configured with HYPRE (e.g. with --download-hypre=1). + * + * The preconditioner does support parallel distributed computations. See @ref + * step-40 for an example. + * + * @ingroup PETScWrappers + * @author Timo Heister, 2010 + */ + class PreconditionBoomerAMG : public PreconditionerBase + { + public: + /** + * Standardized data struct to + * pipe additional flags to the + * preconditioner. + */ + struct AdditionalData + { + /** + * Constructor. Note that BoomerAMG + * offers a lot more options to set + * than what is exposed here. + */ + 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 + ); + + /** + * Set this flag to true if you + * have a symmetric system matrix + * and you want to use a solver + * which asumes a symmetric + * preconditioner like CG. The + * relaxation is done with + * SSOR/Jacobi when set to true and + * with SOR/Jacobi otherwise. + */ + bool symmetric_operator; + + /** + * Threshold of when nodes are + * considered strongly + * connected. See + * HYPRE_BoomerAMGSetStrongThreshold(). Recommended + * values are 0.25 for 2d and 0.5 + * for 3d problems, but it is + * problem dependent. + */ + double strong_threshold; + + /** + * If set to a value smaller than + * 1.0 then diagonally dominant + * parts of the matrix are treated + * as having no strongly connected + * nodes. If the row sum weighted + * by the diagonal entry is bigger + * than the given value, it is + * considered diagonally + * dominant. This feature is turned + * of by setting the value to + * 1.0. This is the default as some + * matrices can result in having + * only diagonally dominant entries + * and thus no multigrid levels are + * constructed. The default in + * BoomerAMG for this is 0.9. When + * you try this, check for a + * reasonable number of levels + * created. + */ + double max_row_sum; + + /** + * Number of levels of aggressive + * coarsening. Increasing this + * value reduces the construction + * time and memory requirements but + * may decrease effectiveness.*/ + unsigned int aggressive_coarsening_num_levels; + + /** + * Setting this flag to true + * produces debug output from + * HYPRE, when the preconditioner + * is constructed. + */ + bool output_details; + + + + + }; + + + + /** + * Constructor. Take the matrix which + * is used to form the preconditioner, + * and additional flags if there are + * any. + */ + PreconditionBoomerAMG (const MatrixBase &matrix, + const AdditionalData &additional_data = AdditionalData()); + + protected: + /** + * Store a copy of the flags for this + * particular preconditioner. + */ + const AdditionalData additional_data; + + /** + * Function that takes a Krylov + * Subspace Preconditioner context + * object, and sets the type of + * preconditioner that is appropriate + * for the present class. + */ + virtual void set_preconditioner_type (PC &pc) const; + }; + + } + DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/source/lac/petsc_precondition.cc b/deal.II/source/lac/petsc_precondition.cc index cceb49b1a6..15d13c5108 100644 --- a/deal.II/source/lac/petsc_precondition.cc +++ b/deal.II/source/lac/petsc_precondition.cc @@ -16,6 +16,7 @@ #ifdef DEAL_II_USE_PETSC +# include # include # include # include @@ -286,6 +287,78 @@ namespace PETScWrappers } +/* ----------------- PreconditionBoomerAMG -------------------- */ + + PreconditionBoomerAMG::AdditionalData:: + AdditionalData(const bool symmetric_operator, + const double strong_threshold, + const double max_row_sum, + const unsigned int aggressive_coarsening_num_levels, + const bool output_details + ) + : + symmetric_operator(symmetric_operator), + strong_threshold(strong_threshold), + max_row_sum(max_row_sum), + aggressive_coarsening_num_levels(aggressive_coarsening_num_levels), + output_details(output_details) + {} + + + PreconditionBoomerAMG::PreconditionBoomerAMG (const MatrixBase &matrix, + const AdditionalData &additional_data) + : + PreconditionerBase (matrix), + additional_data (additional_data) + {} + + + void + PreconditionBoomerAMG::set_preconditioner_type (PC &pc) const + { + // set the right type for the + // preconditioner + int ierr; + ierr = PCSetType (pc, const_cast(PCHYPRE)); + AssertThrow (ierr == 0, ExcPETScError(ierr)); + + ierr = PCHYPRESetType(pc, "boomeramg"); + AssertThrow (ierr == 0, ExcPETScError(ierr)); + + if (additional_data.output_details) + PetscOptionsSetValue("-pc_hypre_boomeramg_print_statistics","1"); + + PetscOptionsSetValue("-pc_hypre_boomeramg_agg_nl", + Utilities::int_to_string( + additional_data.aggressive_coarsening_num_levels + ).c_str()); + + std::stringstream ssStream; + ssStream << additional_data.max_row_sum; + PetscOptionsSetValue("-pc_hypre_boomeramg_max_row_sum", ssStream.str().c_str()); + + ssStream.str(""); // empty the stringstream + ssStream << additional_data.strong_threshold; + PetscOptionsSetValue("-pc_hypre_boomeramg_strong_threshold", ssStream.str().c_str()); + + if (additional_data.symmetric_operator) + { + PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_up", "symmetric-SOR/Jacobi"); + PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_down", "symmetric-SOR/Jacobi"); + PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_coarse", "Gaussian-elimination"); + } + else + { + PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_up", "SOR/Jacobi"); + PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_down", "SOR/Jacobi"); + PetscOptionsSetValue("-pc_hypre_boomeramg_relax_type_coarse", "Gaussian-elimination"); + } + + ierr = PCSetFromOptions (pc); + AssertThrow (ierr == 0, ExcPETScError(ierr)); + } + + /* ----------------- PreconditionLU -------------------- */ PreconditionLU::AdditionalData::