#include <lac/sparse_direct.h>
#include <lac/sparse_ilu.h>
#include <lac/trilinos_sparse_matrix.h>
+#include <lac/trilinos_precondition.h>
#include <lac/trilinos_precondition_amg.h>
#include <grid/tria.h>
1e-8*temperature_rhs.l2_norm());
SolverCG<TrilinosWrappers::Vector> 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);
--- /dev/null
+//---------------------------------------------------------------------------
+// $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 <base/config.h>
+#include <boost/shared_ptr.hpp>
+
+
+#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<Ifpack_Preconditioner> preconditioner;
+ };
+
+
+}
+
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif // DEAL_II_USE_TRILINOS
+
+/*---------------------------- trilinos_precondition_base.h ---------------------------*/
+
+#endif
+/*---------------------------- trilinos_precondition_base.h ---------------------------*/
--- /dev/null
+//---------------------------------------------------------------------------
+// $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 <lac/trilinos_precondition.h>
+#include <lac/trilinos_vector.h>
+#include <lac/trilinos_sparse_matrix.h>
+
+#include <Ifpack.h>
+
+
+
+#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