From: bangerth Date: Wed, 3 Sep 2008 01:05:15 +0000 (+0000) Subject: Add an initial version of some preconditioners based on Trilinos' Ifpack functionality. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=50c0ed7c6951cfcc948981d8b132111bdae16f79;p=dealii-svn.git Add an initial version of some preconditioners based on Trilinos' Ifpack functionality. git-svn-id: https://svn.dealii.org/trunk@16716 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 cfd21ca5e1..57d80394d0 100644 --- a/deal.II/examples/step-31/step-31.cc +++ b/deal.II/examples/step-31/step-31.cc @@ -30,6 +30,7 @@ #include #include #include +#include #include #include @@ -1905,9 +1906,11 @@ void BoussinesqFlowProblem::solve () 1e-8*temperature_rhs.l2_norm()); SolverCG cg (solver_control); -//TODO: Need to do something about this! + TrilinosWrappers::PreconditionSSOR preconditioner (temperature_matrix, + 1.2); cg.solve (temperature_matrix, temperature_solution, - temperature_rhs, PreconditionIdentity()); + temperature_rhs, + preconditioner); // produce a consistent temperature field temperature_constraints.distribute (temperature_solution); diff --git a/deal.II/lac/include/lac/trilinos_precondition.h b/deal.II/lac/include/lac/trilinos_precondition.h new file mode 100755 index 0000000000..43bcbb138a --- /dev/null +++ b/deal.II/lac/include/lac/trilinos_precondition.h @@ -0,0 +1,90 @@ +//--------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2008 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//--------------------------------------------------------------------------- +#ifndef __deal2__trilinos_precondition_h +#define __deal2__trilinos_precondition_h + + +#include +#include + + +#ifdef DEAL_II_USE_TRILINOS + +// forward declarations +class Ifpack_Preconditioner; + + +DEAL_II_NAMESPACE_OPEN + +namespace TrilinosWrappers +{ + class Vector; + class SparseMatrix; + + + class PreconditionSSOR + { + public: + public: + /** + * Standardized data struct to + * pipe additional flags to the + * preconditioner. + */ + struct AdditionalData + { + /** + * Constructor. By default, + * set the damping parameter + * to one. + */ + AdditionalData (const double omega = 1); + + /** + * Relaxation parameter. + */ + double omega; + }; + + /** + * Constructor. Take the matrix which + * is used to form the preconditioner, + * and additional flags if there are + * any. + */ + PreconditionSSOR (const SparseMatrix &matrix, + const AdditionalData &additional_data = AdditionalData()); + + /** + * Apply the preconditioner. + */ + void vmult (Vector &dst, + const Vector &src) const; + + private: + boost::shared_ptr preconditioner; + }; + + +} + + + +DEAL_II_NAMESPACE_CLOSE + +#endif // DEAL_II_USE_TRILINOS + +/*---------------------------- trilinos_precondition_base.h ---------------------------*/ + +#endif +/*---------------------------- trilinos_precondition_base.h ---------------------------*/ diff --git a/deal.II/lac/source/trilinos_precondition.cc b/deal.II/lac/source/trilinos_precondition.cc new file mode 100755 index 0000000000..c980901f24 --- /dev/null +++ b/deal.II/lac/source/trilinos_precondition.cc @@ -0,0 +1,73 @@ +//--------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2008 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//--------------------------------------------------------------------------- + + +#include +#include +#include + +#include + + + +#ifdef DEAL_II_USE_TRILINOS + + +DEAL_II_NAMESPACE_OPEN + +namespace TrilinosWrappers +{ + + PreconditionSSOR::AdditionalData:: + AdditionalData (const double omega) + : + omega (omega) + {} + + + PreconditionSSOR::PreconditionSSOR (const SparseMatrix &matrix, + const AdditionalData &additional_data) + : + preconditioner (Ifpack().Create ("point relaxation", + &*matrix.matrix, 0)) + { + Assert (preconditioner != 0, ExcMessage ("Trilinos could not create this preconditioner")); + + int ierr; + + Teuchos::ParameterList parameter_list; + parameter_list.set ("relaxation: type", "symmetric Gauss-Seidel"); + parameter_list.set ("relaxation: damping factor", additional_data.omega); + + ierr = preconditioner->SetParameters(parameter_list); + AssertThrow (ierr == 0, ExcTrilinosError(ierr)); + + ierr = preconditioner->Initialize(); + AssertThrow (ierr == 0, ExcTrilinosError(ierr)); + + ierr = preconditioner->Compute(); + AssertThrow (ierr == 0, ExcTrilinosError(ierr)); + } + + + void + PreconditionSSOR::vmult (Vector &dst, + const Vector &src) const + { + preconditioner->ApplyInverse (*src.vector, *dst.vector); + } +} + +DEAL_II_NAMESPACE_CLOSE + +#endif // DEAL_II_USE_TRILINOS