]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add some documentation.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 21 Aug 2008 20:48:49 +0000 (20:48 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 21 Aug 2008 20:48:49 +0000 (20:48 +0000)
git-svn-id: https://svn.dealii.org/trunk@16643 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 9d1df99a2de3553da3f51af555370cd8e49c7b89..b0aff0cec055fb45a751db06e5befeff276838c6 100755 (executable)
@@ -37,46 +37,49 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace TrilinosWrappers
 {
-                                  // This implements an AMG
-                                  // preconditioner based on the
-                                  // Trilinos ML implementation.
-                                  // What this class does is twofold.
-                                  // When the constructor of the class
-                                  // is invoked, a ML preconditioner
-                                  // object is created based on the
-                                  // DoFHandler and matrix
-                                  // that we want the preconditioner to
-                                  // be based on. A call of
-                                  // the respective <code>vmult</code>
-                                  // function does call the respective
-                                  // operation in the Trilinos package,
-                                  // where it is called 
-                                  // <code>ApplyInverse</code>.
-                                  // There are a few pecularities in
-                                  // the constructor. Since the
-                                  // Trilinos objects we want to use are
-                                  // heavily dependent on Epetra objects,
-                                  // the fundamental construction
-                                  // routines for vectors and 
-                                  // matrices in Trilinos, we do a 
-                                  // copy of our deal.II preconditioner
-                                  // matrix to a Epetra matrix. This 
-                                  // is of course not optimal, but for
-                                  // the time being there is no direct
-                                  // support for our data interface.
-                                  // When doing this time-consuming 
-                                  // operation, we can still profit 
-                                  // from the fact that some of the
-                                  // entries in the preconditioner matrix
-                                  // are zero and hence can be 
-                                  // neglected.
+  
+/**  
+ * This class implements an algebraic multigrid (AMG) preconditioner
+ * based on the Trilinos ML implementation.  What this class does is
+ * twofold.  When the initialize() function is invoked, a ML
+ * preconditioner object is created based on the matrix
+ * that we want the preconditioner to be based on. A call of the
+ * respective <code>vmult</code> function does call the respective
+ * operation in the Trilinos package, where it is called
+ * <code>ApplyInverse</code>. Use of this class is explained in the
+ * @ref step_31 "step-31" tutorial program.
+ *
+ * There are a few pecularities in initialize(). Since the Trilinos
+ * objects we want to use are heavily dependent on Epetra objects, the
+ * fundamental construction routines for vectors and matrices in
+ * Trilinos, we do a copy of our deal.II preconditioner matrix to a
+ * Epetra matrix. This is of course not optimal, but for the time
+ * being there is no direct support for our data interface.  When
+ * doing this time-consuming operation, we can still profit from the
+ * fact that some of the entries in the preconditioner matrix are zero
+ * and hence can be neglected.
+ *
+ * @author Martin Kronbichler, 2008
+ */
   class PreconditionAMG : public Subscriptor
   {
     public:
+                                      /**
+                                       * Constructor.
+                                       */
       PreconditionAMG ();
+
+                                      /**
+                                       * Destructor.
+                                       */
       ~PreconditionAMG ();
-    
+
+                                      /**
+                                       * Let Trilinos compute a
+                                       * multilevel hierarchy for the
+                                       * solution of a linear system
+                                       * with the given matrix.
+                                       */
       void initialize (const dealii::SparseMatrix<double> &matrix,
                       const std::vector<double>  &null_space,
                       const unsigned int          null_space_dimension,
@@ -85,15 +88,30 @@ namespace TrilinosWrappers
                       const bool                  output_details,
                       const double                drop_tolerance = 1e-13);
 
+                                      /**
+                                       * Multiply the source vector
+                                       * with the preconditioner
+                                       * represented by this object,
+                                       * and return it in the
+                                       * destination vector.
+                                       */
       void vmult (dealii::Vector<double>       &dst,
                  const dealii::Vector<double> &src) const;
 
     private:
-    
+
+                                      /**
+                                       * Pointer to the Trilinos AMG object.
+                                       */
       boost::shared_ptr<ML_Epetra::MultiLevelPreconditioner> multigrid_operator;
-    
+
       Epetra_SerialComm  communicator;
       boost::shared_ptr<Epetra_Map>       Map;
+
+                                      /**
+                                       * A copy of the deal.II matrix
+                                       * into Trilinos format.
+                                       */
       boost::shared_ptr<Epetra_CrsMatrix> Matrix;
   };
 }
index c025380f31bc02a359b1b7a8291ac0bdbce3c320..26a1361fa2e54a866b903726735e1fa2b5eb5529 100755 (executable)
@@ -37,19 +37,21 @@ namespace TrilinosWrappers
   {}
 
 
+  
   PreconditionAMG::~PreconditionAMG ()
   {}
 
   
-  void PreconditionAMG::initialize (
-    const dealii::SparseMatrix<double> &matrix,
-    const std::vector<double>  &null_space,
-    const unsigned int          null_space_dimension,
-    const bool                  elliptic,
-    const bool                  higher_order_elements,
-    const bool                  output_details,
-    const double                drop_tolerance
-  )
+  
+  void
+  PreconditionAMG::
+  initialize (const dealii::SparseMatrix<double> &matrix,
+             const std::vector<double>  &null_space,
+             const unsigned int          null_space_dimension,
+             const bool                  elliptic,
+             const bool                  higher_order_elements,
+             const bool                  output_details,
+             const double                drop_tolerance)
   {
     Assert (drop_tolerance >= 0,
            ExcMessage ("Drop tolerance must be a non-negative number."));

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.