]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Also instantiate for value_type float
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Fri, 14 Sep 2018 15:34:33 +0000 (17:34 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Wed, 20 Feb 2019 12:46:19 +0000 (13:46 +0100)
cmake/configure/configure_2_trilinos.cmake
include/deal.II/lac/trilinos_sparse_matrix.h
source/lac/trilinos_sparse_matrix.cc

index 661c739b7c7f14f1099984052b7b23b1048cebe9..c159fe255270ddf74655ca2562b943591bd32ce5 100644 (file)
@@ -228,7 +228,7 @@ MACRO(FEATURE_TRILINOS_CONFIGURE_EXTERNAL)
   IF (TRILINOS_WITH_MPI)
     SET(DEAL_II_EXPAND_EPETRA_VECTOR "LinearAlgebra::EpetraWrappers::Vector")
     SET(DEAL_II_EXPAND_TPETRA_VECTOR_DOUBLE "LinearAlgebra::TpetraWrappers::Vector<double>")
-    #SET(DEAL_II_EXPAND_TPETRA_VECTOR_FLOAT "LinearAlgebra::TpetraWrappers::Vector<float>")
+    SET(DEAL_II_EXPAND_TPETRA_VECTOR_FLOAT "LinearAlgebra::TpetraWrappers::Vector<float>")
   ENDIF()
   IF(${DEAL_II_TRILINOS_WITH_SACADO})
     # Note: Only CMake 3.0 and greater support line continuation with the "\" character
index b158a224d493619427b1647fc8e397acd29f6c02..2ae42c5ef1576dfa1936ba1252a79d96750aa79f 100644 (file)
@@ -1325,13 +1325,6 @@ namespace TrilinosWrappers
         const Number *   values,
         const bool       elide_zero_values = false);
 
-    void
-    set(const size_type       row,
-        const size_type       n_cols,
-        const size_type *     col_indices,
-        const TrilinosScalar *values,
-        const bool            elide_zero_values = false);
-
     /**
      * Add @p value to the element (<i>i,j</i>).
      *
@@ -1588,7 +1581,13 @@ namespace TrilinosWrappers
      * distributed. Otherwise, an exception will be thrown.
      */
     template <typename VectorType>
-    void
+    typename std::enable_if<std::is_same<typename VectorType::value_type,
+                                         TrilinosScalar>::value>::type
+    vmult(VectorType &dst, const VectorType &src) const;
+
+    template <typename VectorType>
+    typename std::enable_if<!std::is_same<typename VectorType::value_type,
+                                          TrilinosScalar>::value>::type
     vmult(VectorType &dst, const VectorType &src) const;
 
     /**
@@ -1602,7 +1601,13 @@ namespace TrilinosWrappers
      * see the discussion about @p VectorType in vmult().
      */
     template <typename VectorType>
-    void
+    typename std::enable_if<std::is_same<typename VectorType::value_type,
+                                         TrilinosScalar>::value>::type
+    Tvmult(VectorType &dst, const VectorType &src) const;
+
+    template <typename VectorType>
+    typename std::enable_if<!std::is_same<typename VectorType::value_type,
+                                          TrilinosScalar>::value>::type
     Tvmult(VectorType &dst, const VectorType &src) const;
 
     /**
@@ -2992,6 +2997,32 @@ namespace TrilinosWrappers
   // Inline the set() and add() functions, since they will be called
   // frequently, and the compiler can optimize away some unnecessary loops
   // when the sizes are given at compile time.
+  template <>
+  void
+  SparseMatrix::set<TrilinosScalar>(const size_type       row,
+                                    const size_type       n_cols,
+                                    const size_type *     col_indices,
+                                    const TrilinosScalar *values,
+                                    const bool            elide_zero_values);
+
+
+
+  template <typename Number>
+  void
+  SparseMatrix::set(const size_type  row,
+                    const size_type  n_cols,
+                    const size_type *col_indices,
+                    const Number *   values,
+                    const bool       elide_zero_values)
+  {
+    std::vector<TrilinosScalar> trilinos_values(n_cols);
+    std::copy(values, values + n_cols, trilinos_values.begin());
+    this->set(
+      row, n_cols, col_indices, trilinos_values.data(), elide_zero_values);
+  }
+
+
+
   inline void
   SparseMatrix::set(const size_type      i,
                     const size_type      j,
@@ -3298,6 +3329,13 @@ namespace TrilinosWrappers
     } // namespace LinearOperatorImplementation
   }   // namespace internal
 
+  template <>
+  void
+  SparseMatrix::set<TrilinosScalar>(const size_type       row,
+                                    const size_type       n_cols,
+                                    const size_type *     col_indices,
+                                    const TrilinosScalar *values,
+                                    const bool            elide_zero_values);
 #    endif // DOXYGEN
 
 } /* namespace TrilinosWrappers */
index 0996f76053ee05c4c842e24636f2a0bb3d33913a..250a6d0519183b2490b48ee89836b2ac5d40595a 100644 (file)
@@ -1540,28 +1540,13 @@ namespace TrilinosWrappers
 
 
 
-  template <typename Number>
+  template <>
   void
-  SparseMatrix::set(const size_type  row,
-                    const size_type  n_cols,
-                    const size_type *col_indices,
-                    const Number *   values,
-                    const bool       elide_zero_values)
-  {
-    std::vector<TrilinosScalar> trilinos_values(n_cols);
-    std::copy(values, values + n_cols, trilinos_values.begin());
-    this->set(
-      row, n_cols, col_indices, trilinos_values.data(), elide_zero_values);
-  }
-
-
-
-  void
-  SparseMatrix::set(const size_type       row,
-                    const size_type       n_cols,
-                    const size_type *     col_indices,
-                    const TrilinosScalar *values,
-                    const bool            elide_zero_values)
+  SparseMatrix::set<TrilinosScalar>(const size_type       row,
+                                    const size_type       n_cols,
+                                    const size_type *     col_indices,
+                                    const TrilinosScalar *values,
+                                    const bool            elide_zero_values)
   {
     AssertIndexRange(row, this->m());
 
@@ -2124,7 +2109,8 @@ namespace TrilinosWrappers
 
 
   template <typename VectorType>
-  void
+  typename std::enable_if<
+    std::is_same<typename VectorType::value_type, TrilinosScalar>::value>::type
   SparseMatrix::vmult(VectorType &dst, const VectorType &src) const
   {
     Assert(&src != &dst, ExcSourceEqualsDestination());
@@ -2157,7 +2143,18 @@ namespace TrilinosWrappers
 
 
   template <typename VectorType>
-  void
+  typename std::enable_if<
+    !std::is_same<typename VectorType::value_type, TrilinosScalar>::value>::type
+  SparseMatrix::vmult(VectorType & /*dst*/, const VectorType & /*src*/) const
+  {
+    AssertThrow(false, ExcNotImplemented());
+  }
+
+
+
+  template <typename VectorType>
+  typename std::enable_if<
+    std::is_same<typename VectorType::value_type, TrilinosScalar>::value>::type
   SparseMatrix::Tvmult(VectorType &dst, const VectorType &src) const
   {
     Assert(&src != &dst, ExcSourceEqualsDestination());
@@ -2186,6 +2183,16 @@ namespace TrilinosWrappers
 
 
 
+  template <typename VectorType>
+  typename std::enable_if<
+    !std::is_same<typename VectorType::value_type, TrilinosScalar>::value>::type
+  SparseMatrix::Tvmult(VectorType & /*dst*/, const VectorType & /*src*/) const
+  {
+    AssertThrow(false, ExcNotImplemented());
+  }
+
+
+
   template <typename VectorType>
   void
   SparseMatrix::vmult_add(VectorType &dst, const VectorType &src) const
@@ -3508,10 +3515,10 @@ namespace TrilinosWrappers
     dealii::LinearAlgebra::TpetraWrappers::Vector<double> &,
     const dealii::LinearAlgebra::TpetraWrappers::Vector<double> &) const;
 
-  /*template void
-     SparseMatrix::vmult(
-       dealii::LinearAlgebra::TpetraWrappers::Vector<float> &,
-       const dealii::LinearAlgebra::TpetraWrappers::Vector<float> &) const;*/
+  template void
+  SparseMatrix::vmult(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float> &) const;
 
   template void
   SparseMatrix::vmult(
@@ -3537,6 +3544,11 @@ namespace TrilinosWrappers
     dealii::LinearAlgebra::TpetraWrappers::Vector<double> &,
     const dealii::LinearAlgebra::TpetraWrappers::Vector<double> &) const;
 
+  template void
+  SparseMatrix::Tvmult(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float> &) const;
+
   template void
   SparseMatrix::Tvmult(
     dealii::LinearAlgebra::EpetraWrappers::Vector &,
@@ -3561,6 +3573,11 @@ namespace TrilinosWrappers
     dealii::LinearAlgebra::TpetraWrappers::Vector<double> &,
     const dealii::LinearAlgebra::TpetraWrappers::Vector<double> &) const;
 
+  template void
+  SparseMatrix::vmult_add(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float> &) const;
+
   template void
   SparseMatrix::vmult_add(
     dealii::LinearAlgebra::EpetraWrappers::Vector &,
@@ -3585,6 +3602,11 @@ namespace TrilinosWrappers
     dealii::LinearAlgebra::TpetraWrappers::Vector<double> &,
     const dealii::LinearAlgebra::TpetraWrappers::Vector<double> &) const;
 
+  template void
+  SparseMatrix::Tvmult_add(
+    dealii::LinearAlgebra::TpetraWrappers::Vector<float> &,
+    const dealii::LinearAlgebra::TpetraWrappers::Vector<float> &) const;
+
   template void
   SparseMatrix::Tvmult_add(
     dealii::LinearAlgebra::EpetraWrappers::Vector &,

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.