]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implemented a TrilinosWrapper function to a Chebyshev preconditioner.
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 13 Nov 2008 15:12:28 +0000 (15:12 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 13 Nov 2008 15:12:28 +0000 (15:12 +0000)
git-svn-id: https://svn.dealii.org/trunk@17567 0785d39b-7218-0410-832d-ea1e28bc413d

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

index ca32752a89cdd4605393987b6e8bb4ff62c2a9b9..7ad20e145aabc99ee07db0b5d0a0fcf16b277735 100755 (executable)
@@ -36,6 +36,7 @@
 
 // forward declarations
 class Ifpack_Preconditioner;
+class Ifpack_Chebyshev;
 namespace ML_Epetra
 {
   class MultiLevelPreconditioner;
@@ -901,7 +902,7 @@ namespace TrilinosWrappers
  * @ingroup Preconditioners
  * @author Martin Kronbichler, 2008
  */
-  class PreconditionBlockDirect : public PreconditionBase
+  class PreconditionBlockwiseDirect : public PreconditionBase
   {
     public:
                                        /**
@@ -948,6 +949,107 @@ namespace TrilinosWrappers
   
 
 
+
+
+
+/**
+ * A wrapper class for a Chebyshev preconditioner for Trilinos matrices.
+ *
+ * The AdditionalData data structure allows to set preconditioner
+ * options.
+ *
+ * @ingroup TrilinosWrappers
+ * @ingroup Preconditioners
+ * @author Martin Kronbichler, 2008
+ */
+  class PreconditionChebyshev : public PreconditionBase
+  {
+    public:
+                                       /**
+                                        * Standardized data struct to
+                                        * pipe additional parameters
+                                        * to the preconditioner.
+                                        */      
+      struct AdditionalData
+      {
+                                       /**
+                                       * Constructor.
+                                       */
+       AdditionalData (const unsigned int degree           = 1,
+                       const double       max_eigenvalue   = 10.,
+                       const double       eigenvalue_ratio = 30.,
+                       const double       min_eigenvalue   = 1.,
+                       const double       min_diagonal     = 1e-12,
+                       const bool         nonzero_starting = false);
+       
+                                       /**
+                                       * This determines the degree of the
+                                       * Chebyshev polynomial.
+                                       */
+       unsigned int degree;
+
+                                       /**
+                                       * This sets the maximum eigenvalue
+                                       * of the matrix, which needs to be
+                                       * set properly for appropriate
+                                       * performance of the Chebyshev
+                                       * preconditioner.
+                                       */
+       double max_eigenvalue;
+
+                                       /**
+                                       * This sets the ratio between the
+                                       * maximum and the minimum
+                                       * eigenvalue.
+                                       */
+       double eigenvalue_ratio;
+
+                                       /**
+                                       * This sets the minimum eigenvalue,
+                                       * which is an optional parameter
+                                       * only used internally for checking
+                                       * whether we use an identity matrix.
+                                       */
+       double min_eigenvalue;
+
+                                       /**
+                                       * This sets a threshold below which
+                                       * the diagonal element will not be
+                                       * inverted in the Chebyshev
+                                       * algorithm.
+                                       */
+       double min_diagonal;
+
+                                       /**
+                                       * This flag let the preconditioner
+                                       * start at a nonzero value when
+                                       * generating the Chebyshev
+                                       * polynomial.
+                                       */
+       double nonzero_starting;
+      };
+
+                                       /**
+                                        * Initialize function. Takes
+                                        * the matrix which is used to
+                                        * form the preconditioner, and
+                                        * additional flags if there
+                                        * are any.
+                                        */
+      void initialize (const SparseMatrix   &matrix,
+                      const AdditionalData &additional_data = AdditionalData());
+
+    private:
+                                      /**
+                                       * This is a pointer to the
+                                       * Ifpack data contained in
+                                       * this preconditioner.
+                                       */
+      Teuchos::RCP<Ifpack_Chebyshev> ifpack;
+  };
+  
+
+
 /**
  * This class implements an algebraic multigrid (AMG) preconditioner
  * based on the Trilinos ML implementation, which is a black-box
index 5545321c48ea0887eb20e7966afcf04f47dc9d35..6c4a54a4d12f57789b7ffbbfff064cb3475c4c84 100755 (executable)
@@ -19,6 +19,7 @@
 #ifdef DEAL_II_USE_TRILINOS
 
 #include <Ifpack.h>
+#include <Ifpack_Chebyshev.h>
 #include <Teuchos_ParameterList.hpp>
 #include <Epetra_Vector.h>
 #include <Epetra_MultiVector.h>
@@ -372,7 +373,7 @@ namespace TrilinosWrappers
 
 /* ---------------------- PreconditionBlockDirect --------------------- */
 
-  PreconditionBlockDirect::AdditionalData::
+  PreconditionBlockwiseDirect::AdditionalData::
   AdditionalData (const unsigned int overlap)
                   :
                   overlap  (overlap)
@@ -381,8 +382,8 @@ namespace TrilinosWrappers
 
 
   void
-  PreconditionBlockDirect::initialize (const SparseMatrix   &matrix,
-                                      const AdditionalData &additional_data)
+  PreconditionBlockwiseDirect::initialize (const SparseMatrix   &matrix,
+                                          const AdditionalData &additional_data)
   {
     preconditioner.release();
     ifpack.release();
@@ -411,6 +412,67 @@ namespace TrilinosWrappers
 
 
 
+/* ---------------------- PreconditionBlockDirect --------------------- */
+
+  PreconditionChebyshev::AdditionalData::
+  AdditionalData (const unsigned int degree,
+                 const double       max_eigenvalue,
+                 const double       eigenvalue_ratio,
+                 const double       min_eigenvalue,
+                 const double       min_diagonal,
+                 const bool         nonzero_starting)
+                  :
+                  degree  (degree),
+                 max_eigenvalue (max_eigenvalue),
+                 eigenvalue_ratio (eigenvalue_ratio),
+                 min_eigenvalue (min_eigenvalue),
+                 min_diagonal (min_diagonal),
+                 nonzero_starting (nonzero_starting)
+  {}
+
+
+
+  void
+  PreconditionChebyshev::initialize (const SparseMatrix   &matrix,
+                                    const AdditionalData &additional_data)
+  {
+    preconditioner.release();
+    ifpack.release();
+
+    ifpack = Teuchos::rcp (new Ifpack_Chebyshev(&*matrix.matrix));
+    Assert (&*ifpack != 0, ExcMessage ("Trilinos could not create this "
+                                      "preconditioner"));
+
+    int ierr;
+
+    Teuchos::ParameterList parameter_list;
+    parameter_list.set ("chebyshev: ratio eigenvalue",
+                       additional_data.eigenvalue_ratio);
+    parameter_list.set ("chebyshev: min eigenvalue",
+                       additional_data.min_eigenvalue);
+    parameter_list.set ("chebyshev: max eigenvalue",
+                       additional_data.max_eigenvalue);
+    parameter_list.set ("chebyshev: degree",
+                       (int)additional_data.degree);
+    parameter_list.set ("chebyshev: min diagonal value",
+                       additional_data.min_diagonal);
+    parameter_list.set ("chebyshev: zero starting solution",
+                       !additional_data.nonzero_starting);
+
+    ierr = ifpack->SetParameters(parameter_list);
+    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+
+    ierr = ifpack->Initialize();
+    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+
+    ierr = ifpack->Compute();
+    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+
+    preconditioner = Teuchos::rcp (&*ifpack, false);
+  }
+
+
+
 /* -------------------------- PreconditionAMG -------------------------- */
 
   PreconditionAMG::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.