From: bangerth Date: Thu, 21 Aug 2008 20:41:29 +0000 (+0000) Subject: Rename the preconditioner to make it conform with our usual naming scheme. Use it... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=e636e5adf5a640a443892d034f361086e59e0ed4;p=dealii-svn.git Rename the preconditioner to make it conform with our usual naming scheme. Use it from step-31. git-svn-id: https://svn.dealii.org/trunk@16642 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/examples/step-31/step-31.cc b/deal.II/examples/step-31/step-31.cc index 1f7fb24e13..2e8fbf7d86 100644 --- a/deal.II/examples/step-31/step-31.cc +++ b/deal.II/examples/step-31/step-31.cc @@ -29,7 +29,7 @@ #include #include #include -#include +#include #include #include @@ -60,13 +60,6 @@ #include // This is Trilinos -#include -#include -#include -#include -#include -#include -#include // Next, we import all deal.II // names into global namespace @@ -75,211 +68,6 @@ using namespace dealii; // @sect3{Defining the AMG preconditioner} - // 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. -class PreconditionerTrilinosAmg : public Subscriptor -{ - public: - PreconditionerTrilinosAmg (); - - void initialize (const SparseMatrix &matrix, - const std::vector &null_space, - const unsigned int null_space_dimension, - const bool higher_order_elements, - const bool elliptic, - const bool output_details, - const double drop_tolerance = 1e-13); - - void vmult (Vector &dst, - const Vector &src) const; - - private: - - boost::shared_ptr multigrid_operator; - - Epetra_SerialComm communicator; - boost::shared_ptr Map; - boost::shared_ptr Matrix; -}; - - -PreconditionerTrilinosAmg::PreconditionerTrilinosAmg () -{} - -void PreconditionerTrilinosAmg::initialize ( - const 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.")); - - const unsigned int n_rows = matrix.m(); - const SparsityPattern *sparsity_pattern = &(matrix.get_sparsity_pattern()); - - // Init Epetra Matrix, avoid - // storing the nonzero elements. - { - Map.reset (new Epetra_Map(n_rows, 0, communicator)); - - std::vector row_lengths (n_rows); - for (SparseMatrix::const_iterator p = matrix.begin(); - p != matrix.end(); ++p) - if (std::abs(p->value()) > drop_tolerance) - ++row_lengths[p->row()]; - - Matrix.reset (new Epetra_CrsMatrix(Copy, *Map, &row_lengths[0], true)); - - const unsigned int max_nonzero_entries - = *std::max_element (row_lengths.begin(), row_lengths.end()); - - std::vector values(max_nonzero_entries, 0); - std::vector row_indices(max_nonzero_entries); - - for (unsigned int row=0; row::const_iterator p = matrix.begin(row); - p != matrix.end(row); ++p) - if (std::abs(p->value()) > drop_tolerance) - { - row_indices[index] = p->column(); - values[index] = p->value(); - ++index; - } - - Assert (index == static_cast(row_lengths[row]), - ExcMessage("Filtering out zeros could not " - "be successfully finished!")); - - Matrix->InsertGlobalValues(row, row_lengths[row], - &values[0], &row_indices[0]); - } - - Matrix->FillComplete(); - } - - // Build the AMG preconditioner. - Teuchos::ParameterList parameter_list; - - // The implementation is able - // to distinguish between - // matrices from elliptic problems - // and convection dominated - // problems. We use the standard - // options for elliptic problems, - // except that we use a - // Chebyshev smoother instead - // of a symmetric Gauss-Seidel - // smoother. For most elliptic - // problems, Chebyshev is better - // than Gauss-Seidel (SSOR). - if (elliptic) - { - ML_Epetra::SetDefaults("SA",parameter_list); - parameter_list.set("smoother: type", "Chebyshev"); - parameter_list.set("smoother: sweeps", 4); - } - else - { - ML_Epetra::SetDefaults("NSSA",parameter_list); - parameter_list.set("aggregation: type", "Uncoupled"); - parameter_list.set("aggregation: block scaling", true); - } - - if (output_details) - parameter_list.set("ML output", 10); - else - parameter_list.set("ML output", 0); - - if (higher_order_elements) - parameter_list.set("aggregation: type", "MIS"); - - Assert (n_rows * null_space_dimension == null_space.size(), - ExcDimensionMismatch(n_rows * null_space_dimension, - null_space.size())); - - if (null_space_dimension > 1) - { - parameter_list.set("null space: type", "pre-computed"); - parameter_list.set("null space: dimension", int(null_space_dimension)); - parameter_list.set("null space: vectors", (double *)&null_space[0]); - } - - multigrid_operator = boost::shared_ptr - (new ML_Epetra::MultiLevelPreconditioner(*Matrix, parameter_list, true)); - - if (output_details) - multigrid_operator->PrintUnused(0); -} - - // For the implementation of the - // vmult function we - // note that invoking a call of - // the Trilinos preconditioner - // requires us to use Epetra vectors - // as well. Luckily, it is sufficient - // to provide a view, i.e., feed - // Trilinos with a pointer to the - // data, so we avoid copying the - // content of the vectors during - // the iteration. In the declaration - // of the right hand side, we need - // to cast the source vector (that - // is const in all deal.II - // calls) to non-constant value, as - // this is the way Trilinos wants to - // have them. -void PreconditionerTrilinosAmg::vmult (Vector &dst, - const Vector &src) const -{ - Epetra_Vector LHS (View, *Map, dst.begin()); - Epetra_Vector RHS (View, *Map, const_cast(src.begin())); - - const int res = multigrid_operator->ApplyInverse (RHS, LHS); - - Assert (res == 0, - ExcMessage ("Trilinos AMG MultiLevel preconditioner returned " - "with an error!")); -} - // @sect3{Equation data} @@ -768,7 +556,7 @@ class BoussinesqFlowProblem double old_time_step; unsigned int timestep_number; - boost::shared_ptr Amg_preconditioner; + boost::shared_ptr Amg_preconditioner; boost::shared_ptr > Mp_preconditioner; bool rebuild_stokes_matrix; @@ -1386,8 +1174,8 @@ BoussinesqFlowProblem::build_stokes_preconditioner () // boost library. assemble_stokes_preconditioner (); - Amg_preconditioner = boost::shared_ptr - (new PreconditionerTrilinosAmg()); + Amg_preconditioner = boost::shared_ptr + (new TrilinosWrappers::PreconditionAMG()); const unsigned int n_u = stokes_preconditioner_matrix.block(0,0).m(); std::vector null_space (dim * n_u, 0.); @@ -2089,7 +1877,7 @@ void BoussinesqFlowProblem::solve () /*BlockSchurPreconditioner::type, SparseILU > preconditioner (stokes_matrix, mp_inverse, *A_preconditioner);*/ - LinearSolvers::BlockSchurPreconditioner > preconditioner (stokes_matrix, mp_inverse, *Amg_preconditioner); diff --git a/deal.II/lac/include/lac/trilinos_preconditioner_amg.h b/deal.II/lac/include/lac/trilinos_precondition_amg.h similarity index 96% rename from deal.II/lac/include/lac/trilinos_preconditioner_amg.h rename to deal.II/lac/include/lac/trilinos_precondition_amg.h index 9af75d80b2..9d1df99a2d 100755 --- a/deal.II/lac/include/lac/trilinos_preconditioner_amg.h +++ b/deal.II/lac/include/lac/trilinos_precondition_amg.h @@ -71,11 +71,11 @@ namespace TrilinosWrappers // entries in the preconditioner matrix // are zero and hence can be // neglected. - class PreconditionerTrilinosAmg : public Subscriptor + class PreconditionAMG : public Subscriptor { public: - PreconditionerTrilinosAmg (); - ~PreconditionerTrilinosAmg (); + PreconditionAMG (); + ~PreconditionAMG (); void initialize (const dealii::SparseMatrix &matrix, const std::vector &null_space, diff --git a/deal.II/lac/source/trilinos_preconditioner_amg.cc b/deal.II/lac/source/trilinos_precondition_amg.cc similarity index 95% rename from deal.II/lac/source/trilinos_preconditioner_amg.cc rename to deal.II/lac/source/trilinos_precondition_amg.cc index 72cb844b5a..c025380f31 100755 --- a/deal.II/lac/source/trilinos_preconditioner_amg.cc +++ b/deal.II/lac/source/trilinos_precondition_amg.cc @@ -12,7 +12,7 @@ //--------------------------------------------------------------------------- -#include +#include #include #include @@ -33,10 +33,15 @@ DEAL_II_NAMESPACE_OPEN namespace TrilinosWrappers { - PreconditionerTrilinosAmg::PreconditionerTrilinosAmg () + PreconditionAMG::PreconditionAMG () {} - void PreconditionerTrilinosAmg::initialize ( + + PreconditionAMG::~PreconditionAMG () + {} + + + void PreconditionAMG::initialize ( const dealii::SparseMatrix &matrix, const std::vector &null_space, const unsigned int null_space_dimension, @@ -165,8 +170,8 @@ namespace TrilinosWrappers // calls) to non-constant value, as // this is the way Trilinos wants to // have them. - void PreconditionerTrilinosAmg::vmult (dealii::Vector &dst, - const dealii::Vector &src) const + void PreconditionAMG::vmult (dealii::Vector &dst, + const dealii::Vector &src) const { Epetra_Vector LHS (View, *Map, dst.begin()); Epetra_Vector RHS (View, *Map, const_cast(src.begin()));