From: bangerth Date: Thu, 21 Aug 2008 20:48:49 +0000 (+0000) Subject: Add some documentation. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=6c6f9c5f425aab91bfa93d6051986a8a0042dd8d;p=dealii-svn.git Add some documentation. git-svn-id: https://svn.dealii.org/trunk@16643 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/trilinos_precondition_amg.h b/deal.II/lac/include/lac/trilinos_precondition_amg.h index 9d1df99a2d..b0aff0cec0 100755 --- a/deal.II/lac/include/lac/trilinos_precondition_amg.h +++ b/deal.II/lac/include/lac/trilinos_precondition_amg.h @@ -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 vmult - // function does call the respective - // operation in the Trilinos package, - // where it is called - // ApplyInverse. - - // 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 vmult function does call the respective + * operation in the Trilinos package, where it is called + * ApplyInverse. 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 &matrix, const std::vector &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 &dst, const dealii::Vector &src) const; private: - + + /** + * Pointer to the Trilinos AMG object. + */ boost::shared_ptr multigrid_operator; - + Epetra_SerialComm communicator; boost::shared_ptr Map; + + /** + * A copy of the deal.II matrix + * into Trilinos format. + */ boost::shared_ptr Matrix; }; } diff --git a/deal.II/lac/source/trilinos_precondition_amg.cc b/deal.II/lac/source/trilinos_precondition_amg.cc index c025380f31..26a1361fa2 100755 --- a/deal.II/lac/source/trilinos_precondition_amg.cc +++ b/deal.II/lac/source/trilinos_precondition_amg.cc @@ -37,19 +37,21 @@ namespace TrilinosWrappers {} + PreconditionAMG::~PreconditionAMG () {} - void PreconditionAMG::initialize ( - const dealii::SparseMatrix &matrix, - const std::vector &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 &matrix, + const std::vector &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."));