From 6225139ff3ef6273575d57656bad05827d62df62 Mon Sep 17 00:00:00 2001
From: Martin Kronbichler <kronbichler@lnm.mw.tum.de>
Date: Wed, 13 Mar 2019 17:16:04 +0100
Subject: [PATCH] Fix MatrixFree::loop() for the case where VectorType is not a
 vector

---
 include/deal.II/matrix_free/matrix_free.h    | 23 ++++++++++++++++----
 source/matrix_free/task_info.cc              |  1 +
 tests/matrix_free/matrix_free_type_traits.cc | 11 ++++++++++
 3 files changed, 31 insertions(+), 4 deletions(-)

diff --git a/include/deal.II/matrix_free/matrix_free.h b/include/deal.II/matrix_free/matrix_free.h
index 4877be670b..65a2f533f5 100644
--- a/include/deal.II/matrix_free/matrix_free.h
+++ b/include/deal.II/matrix_free/matrix_free.h
@@ -1436,7 +1436,6 @@ public:
    * Queries whether or not the geometry-related information for the cells has
    * been set.
    */
-
   bool
   mapping_initialized() const;
 
@@ -3321,13 +3320,15 @@ namespace internal
 
 
     /**
-     * Zero out vector region for vector that do _not_ support
-     * exchange on a subset of DoFs <==> begin() + ind == local_element(ind)
+     * Zero out vector region for vector that do _not_ support exchange on a
+     * subset of DoFs <==> begin() + ind == local_element(ind) but are still a
+     * vector type
      */
     template <
       typename VectorType,
       typename std::enable_if<!has_exchange_on_subset<VectorType>::value,
-                              VectorType>::type * = nullptr>
+                              VectorType>::type * = nullptr,
+      typename VectorType::value_type *           = nullptr>
     void
     zero_vector_region(const unsigned int range_index, VectorType &vec) const
     {
@@ -3341,6 +3342,20 @@ namespace internal
 
 
 
+    /**
+     * Zero out vector region for non-vector types, i.e., classes that do not
+     * have VectorType::value_type
+     */
+    void
+    zero_vector_region(const unsigned int, ...) const
+    {
+      Assert(false,
+             ExcNotImplemented("Zeroing is only implemented for vector types "
+                               "which provide VectorType::value_type"));
+    }
+
+
+
     const dealii::MatrixFree<dim, Number> &matrix_free;
     const typename dealii::MatrixFree<dim, Number>::DataAccessOnFaces
          vector_face_access;
diff --git a/source/matrix_free/task_info.cc b/source/matrix_free/task_info.cc
index d60db786c8..203fcad6ea 100644
--- a/source/matrix_free/task_info.cc
+++ b/source/matrix_free/task_info.cc
@@ -631,6 +631,7 @@ namespace internal
       n_blocks             = 0;
       scheme               = none;
       partition_row_index.clear();
+      partition_row_index.resize(2);
       cell_partition_data.clear();
       face_partition_data.clear();
       boundary_partition_data.clear();
diff --git a/tests/matrix_free/matrix_free_type_traits.cc b/tests/matrix_free/matrix_free_type_traits.cc
index fb0f597c0c..4ea0771d0b 100644
--- a/tests/matrix_free/matrix_free_type_traits.cc
+++ b/tests/matrix_free/matrix_free_type_traits.cc
@@ -90,5 +90,16 @@ main()
           << "unsigned int = "
           << internal::is_serial_or_dummy<unsigned int>::value << std::endl;
 
+  // check that MatrixFree::cell_loop can run for non-vector types
+  MatrixFree<2> matrix_free;
+  int           dummy = 0;
+  matrix_free.cell_loop(
+    std::function<void(const MatrixFree<2> &,
+                       int &,
+                       const int &,
+                       const std::pair<unsigned int, unsigned int> &)>(),
+    dummy,
+    dummy);
+
   deallog << "OK" << std::endl;
 }
-- 
2.39.5