From: Denis Davydov <davydden@gmail.com>
Date: Thu, 21 Feb 2019 18:09:04 +0000 (+0100)
Subject: FEEvaluation: switch between internal::vector_access() via SFINAE
X-Git-Tag: v9.1.0-rc1~323^2~7
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=bf6a24cc20ed8f3cdcf276abfc7b8a6d08b68417;p=dealii.git

FEEvaluation: switch between internal::vector_access() via SFINAE
---

diff --git a/include/deal.II/matrix_free/fe_evaluation.h b/include/deal.II/matrix_free/fe_evaluation.h
index 8c3cc1d3bd..dc97327d07 100644
--- a/include/deal.II/matrix_free/fe_evaluation.h
+++ b/include/deal.II/matrix_free/fe_evaluation.h
@@ -3439,7 +3439,10 @@ namespace internal
 
 
   // access to generic vectors that have operator ().
-  template <typename VectorType>
+  // FIXME: this is wrong for Trilinos/Petsc MPI vectors
+  template <typename VectorType,
+            typename std::enable_if<!has_local_element<VectorType>::value,
+                                    VectorType>::type * = nullptr>
   inline typename VectorType::value_type &
   vector_access(VectorType &vec, const unsigned int entry)
   {
@@ -3451,19 +3454,20 @@ namespace internal
   // access to distributed MPI vectors that have a local_element(uint)
   // method to access data in local index space, which is what we use in
   // DoFInfo and hence in read_dof_values etc.
-  template <typename Number>
-  inline Number &
-  vector_access(LinearAlgebra::distributed::Vector<Number> &vec,
-                const unsigned int                          entry)
+  template <typename VectorType,
+            typename std::enable_if<has_local_element<VectorType>::value,
+                                    VectorType>::type * = nullptr>
+  inline typename VectorType::value_type &
+  vector_access(VectorType &vec, const unsigned int entry)
   {
     return vec.local_element(entry);
   }
 
 
 
-  // this is to make sure that the parallel partitioning in the
-  // LinearAlgebra::distributed::Vector is really the same as stored in
-  // MatrixFree
+  // this is to make sure that the parallel partitioning in VectorType
+  // is really the same as stored in MatrixFree
+  // FIXME: this is incorrect for PETSc/Trilinos MPI vectors
   template <typename VectorType>
   inline void
   check_vector_compatibility(
@@ -3476,6 +3480,10 @@ namespace internal
     AssertDimension(vec.size(), dof_info.vector_partitioner->size());
   }
 
+
+  // this is to make sure that the parallel partitioning in the
+  // LinearAlgebra::distributed::Vector is really the same as stored in
+  // MatrixFree
   template <typename Number>
   inline void
   check_vector_compatibility(
diff --git a/tests/matrix_free/fe_evaluation_type_traits.cc b/tests/matrix_free/fe_evaluation_type_traits.cc
index 68c837d5e1..458a625e8a 100644
--- a/tests/matrix_free/fe_evaluation_type_traits.cc
+++ b/tests/matrix_free/fe_evaluation_type_traits.cc
@@ -24,12 +24,80 @@
 
 #include "../tests.h"
 
+// dummy class we use to check typetraits and internal function.
+// this one mimics LA::d::Vec
+template <typename Number>
+class Dummy
+{
+public:
+  using value_type = Number;
+
+  Number
+  local_element(const unsigned int) const
+  {
+    deallog << "Dummy::local_element() const" << std::endl;
+    return Number();
+  }
+
+  Number &
+  local_element(const unsigned int)
+  {
+    deallog << "Dummy::local_element()" << std::endl;
+    return dummy;
+  }
+
+  Number
+  operator()(const unsigned int) const
+  {
+    deallog << "Dummy::operator() const" << std::endl;
+    return Number();
+  }
+
+  Number &
+  operator()(const unsigned int)
+  {
+    deallog << "Dummy::operator()" << std::endl;
+    return dummy;
+  }
+
+private:
+  Number dummy;
+};
+
+
+template <typename Number>
+class Dummy2
+{
+public:
+  using value_type = Number;
+
+  Number
+  operator()(const unsigned int) const
+  {
+    deallog << "Dummy2::operator() const" << std::endl;
+    return Number();
+  }
+
+  Number &
+  operator()(const unsigned int)
+  {
+    deallog << "Dummy2::operator()" << std::endl;
+    return dummy;
+  }
+
+private:
+  Number dummy;
+};
+
 
 int
 main()
 {
   initlog();
 
+  Dummy<double>  dummy;
+  Dummy2<double> dummy2;
+
   deallog << "has_local_element:" << std::endl
           << "LinearAlgebra::distributed::Vector = "
           << internal::has_local_element<
@@ -37,7 +105,18 @@ main()
           << std::endl
           << "TrilinosWrappers::MPI::Vector = "
           << internal::has_local_element<TrilinosWrappers::MPI::Vector>::value
+          << std::endl
+          << "Dummy = " << internal::has_local_element<Dummy<double>>::value
+          << std::endl
+          << "Dummy2 = " << internal::has_local_element<Dummy2<double>>::value
           << std::endl;
 
+  // now check internal::vector_access wrapper
+  // we expect local_element() to be called
+  deallog << "internal::vector_access:" << std::endl;
+  internal::vector_access(dummy, 0);
+  internal::vector_access(dummy2, 0);
+
+
   deallog << "OK" << std::endl;
 }
diff --git a/tests/matrix_free/fe_evaluation_type_traits.with_trilinos=true.output b/tests/matrix_free/fe_evaluation_type_traits.with_trilinos=true.output
index b62ca9b86d..78f0765434 100644
--- a/tests/matrix_free/fe_evaluation_type_traits.with_trilinos=true.output
+++ b/tests/matrix_free/fe_evaluation_type_traits.with_trilinos=true.output
@@ -2,4 +2,9 @@
 DEAL::has_local_element:
 DEAL::LinearAlgebra::distributed::Vector = 1
 DEAL::TrilinosWrappers::MPI::Vector = 0
+DEAL::Dummy = 1
+DEAL::Dummy2 = 0
+DEAL::internal::vector_access:
+DEAL::Dummy::local_element()
+DEAL::Dummy2::operator()
 DEAL::OK