]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add an initial version of some preconditioners based on Trilinos' Ifpack functionality.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 3 Sep 2008 01:05:15 +0000 (01:05 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 3 Sep 2008 01:05:15 +0000 (01:05 +0000)
git-svn-id: https://svn.dealii.org/trunk@16716 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/step-31/step-31.cc
deal.II/lac/include/lac/trilinos_precondition.h [new file with mode: 0755]
deal.II/lac/source/trilinos_precondition.cc [new file with mode: 0755]

index cfd21ca5e15299cacd1cee6ee3a14b021618d027..57d80394d08f9d19222e3a1355209b4e96bc2735 100644 (file)
@@ -30,6 +30,7 @@
 #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>
@@ -1905,9 +1906,11 @@ void BoussinesqFlowProblem<dim>::solve ()
                                  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);
diff --git a/deal.II/lac/include/lac/trilinos_precondition.h b/deal.II/lac/include/lac/trilinos_precondition.h
new file mode 100755 (executable)
index 0000000..43bcbb1
--- /dev/null
@@ -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 <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     ---------------------------*/
diff --git a/deal.II/lac/source/trilinos_precondition.cc b/deal.II/lac/source/trilinos_precondition.cc
new file mode 100755 (executable)
index 0000000..c980901
--- /dev/null
@@ -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 <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

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.