]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Move PreconditionBase::vmult in h file.
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 19 May 2009 06:11:08 +0000 (06:11 +0000)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 19 May 2009 06:11:08 +0000 (06:11 +0000)
git-svn-id: https://svn.dealii.org/trunk@18867 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/trilinos_precondition.h
deal.II/lac/source/trilinos_precondition.cc

index 26ab08b814aaa7befc46f6b6383f241487327b5b..0e5ee7303f3766add25d3ee145da2c7fe67cdba7 100755 (executable)
@@ -16,9 +16,9 @@
 
 #include <base/config.h>
 #include <base/subscriptor.h>
-
 #include <base/std_cxx1x/shared_ptr.h>
 
+#include <lac/trilinos_vector_base.h>
 
 #ifdef DEAL_II_USE_TRILINOS
 
@@ -31,6 +31,7 @@
 
 #include <Teuchos_RCP.hpp>
 #include <Epetra_Operator.h>
+#include <Epetra_Vector.h>
 
 // forward declarations
 class Ifpack_Preconditioner;
@@ -53,8 +54,7 @@ template <typename number> class Vector;
 
 namespace TrilinosWrappers
 {
-
-  class VectorBase;
+  // forward declarations
   class SparseMatrix;
   class BlockSparseMatrix;
   class SolverBase;
@@ -113,8 +113,9 @@ namespace TrilinosWrappers
                                        * in the Trilinos wrapper
                                        * class.
                                        */
-      void vmult (dealii::Vector<double>       &dst,
-                 const dealii::Vector<double> &src) const;
+      template <typename number>
+      void vmult (dealii::Vector<number>       &dst,
+                 const dealii::Vector<number> &src) const;
 
                                        /**
                                        * Exception.
@@ -1465,8 +1466,73 @@ namespace TrilinosWrappers
       std_cxx1x::shared_ptr<SparseMatrix> Matrix;
   };
 
+
+
+// -------------------------- inline and template functions ----------------------
+
+
+#ifndef DOXYGEN
+
+  inline
+  void
+  PreconditionBase::vmult (VectorBase       &dst,
+                          const VectorBase &src) const
+  {
+    Assert (dst.vector_partitioner().SameAs(preconditioner->OperatorRangeMap()),
+           ExcNonMatchingMaps("dst"));
+    Assert (src.vector_partitioner().SameAs(preconditioner->OperatorDomainMap()),
+           ExcNonMatchingMaps("src"));
+    
+    const int ierr = preconditioner->ApplyInverse (src.trilinos_vector(), 
+                                                  dst.trilinos_vector());
+    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+  }
+
+
+                                       // For the implementation of
+                                       // the <code>vmult</code>
+                                       // function with deal.II data
+                                       // structures we note that
+                                       // invoking a call of the
+                                       // Trilinos preconditioner
+                                       // requires us to use Epetra
+                                       // vectors as well. We do this
+                                       // by providing a view, i.e.,
+                                       // feed Trilinos with a
+                                       // pointer to the data, so we
+                                       // avoid copying the content
+                                       // of the vectors during the
+                                       // iteration (this function is
+                                       // only useful when used in
+                                       // serial anyway). In the
+                                       // declaration of the right
+                                       // hand side, we need to cast
+                                       // the source vector (that is
+                                       // <code>const</code> in all
+                                       // deal.II calls) to
+                                       // non-constant value, as this
+                                       // is the way Trilinos wants
+                                       // to have them.
+  template <typename number>
+  inline
+  void PreconditionBase::vmult (dealii::Vector<number>       &dst,
+                               const dealii::Vector<number> &src) const
+  {
+    
+    Epetra_Vector LHS (View, preconditioner->OperatorDomainMap(),
+                      dst.begin());
+    Epetra_Vector RHS (View, preconditioner->OperatorRangeMap(),
+                      const_cast<double*>(src.begin()));
+  
+    const int ierr = preconditioner->ApplyInverse (RHS, LHS);
+    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
+  }
+
+#endif
+
 }
 
+
 /*@}*/
 
 
index 659a444cec1817fc2f760d0591d01976081c298b..646cb9f5da614c427811d5a31d4eb6dd97afb792 100755 (executable)
@@ -15,7 +15,6 @@
 #include <lac/vector.h>
 #include <lac/sparse_matrix.h>
 #include <lac/trilinos_precondition.h>
-#include <lac/trilinos_vector_base.h>
 #include <lac/trilinos_sparse_matrix.h>
 
 #ifdef DEAL_II_USE_TRILINOS
@@ -23,9 +22,7 @@
 #include <Ifpack.h>
 #include <Ifpack_Chebyshev.h>
 #include <Teuchos_ParameterList.hpp>
-#include <Epetra_Vector.h>
 #include <Epetra_MultiVector.h>
-#include <Epetra_Vector.h>
 #include <ml_include.h>
 #include <ml_MultiLevelPreconditioner.h>
 
@@ -63,60 +60,6 @@ namespace TrilinosWrappers
 
 
 
-  void
-  PreconditionBase::vmult (VectorBase       &dst,
-                          const VectorBase &src) const
-  {
-    Assert (dst.vector_partitioner().SameAs(preconditioner->OperatorRangeMap()),
-           ExcNonMatchingMaps("dst"));
-    Assert (src.vector_partitioner().SameAs(preconditioner->OperatorDomainMap()),
-           ExcNonMatchingMaps("src"));
-    
-    const int ierr = preconditioner->ApplyInverse (src.trilinos_vector(), 
-                                                  dst.trilinos_vector());
-    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
-  }
-
-
-                                       // For the implementation of
-                                       // the <code>vmult</code>
-                                       // function with deal.II data
-                                       // structures we note that
-                                       // invoking a call of the
-                                       // Trilinos preconditioner
-                                       // requires us to use Epetra
-                                       // vectors as well. We do this
-                                       // by providing a view, i.e.,
-                                       // feed Trilinos with a
-                                       // pointer to the data, so we
-                                       // avoid copying the content
-                                       // of the vectors during the
-                                       // iteration (this function is
-                                       // only useful when used in
-                                       // serial anyway). In the
-                                       // declaration of the right
-                                       // hand side, we need to cast
-                                       // the source vector (that is
-                                       // <code>const</code> in all
-                                       // deal.II calls) to
-                                       // non-constant value, as this
-                                       // is the way Trilinos wants
-                                       // to have them.
-  void PreconditionBase::vmult (dealii::Vector<double>       &dst,
-                               const dealii::Vector<double> &src) const
-  {
-    
-    Epetra_Vector LHS (View, preconditioner->OperatorDomainMap(),
-                      dst.begin());
-    Epetra_Vector RHS (View, preconditioner->OperatorRangeMap(),
-                      const_cast<double*>(src.begin()));
-  
-    const int ierr = preconditioner->ApplyInverse (RHS, LHS);
-    AssertThrow (ierr == 0, ExcTrilinosError(ierr));
-  }
-
-
-
 /* -------------------------- PreconditionJacobi -------------------------- */
 
   PreconditionJacobi::AdditionalData::
@@ -699,8 +642,7 @@ namespace TrilinosWrappers
   }
   
   
-  void PreconditionAMG::
-  reinit ()
+  void PreconditionAMG::reinit ()
   {
     multilevel_operator->ReComputePreconditioner();
   }

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.