]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
implement interface to BoomerAMG through PETSc.
authorheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 13 Nov 2010 13:50:55 +0000 (13:50 +0000)
committerheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 13 Nov 2010 13:50:55 +0000 (13:50 +0000)
git-svn-id: https://svn.dealii.org/trunk@22718 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/lac/petsc_precondition.h
deal.II/source/lac/petsc_precondition.cc

index 6d23c6e5b93e13076d9323d9f79c567186d63e2f..6f40b5b1c2ed206506930a0bb99d59ec5e836c34 100644 (file)
@@ -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
 
 
index cceb49b1a69d3c140bd33c5af14f15ce81849f30..15d13c5108d9be913d20d6f0edd05a939d3d48ec 100644 (file)
@@ -16,6 +16,7 @@
 
 #ifdef DEAL_II_USE_PETSC
 
+#  include <base/utilities.h>
 #  include <lac/petsc_matrix_base.h>
 #  include <lac/petsc_vector_base.h>
 #  include <cmath>
@@ -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<char *>(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::

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.