Speed up the function MatrixFree::n_active_entries_per_cell_batch

diff --git a/include/deal.II/matrix_free/matrix_free.h b/include/deal.II/matrix_free/matrix_free.h
index 498d0e828e9402462a8976deba994317a20ff8ec..903d551b0407bc4aeaa841387eef4715a6da2a9b 100644 (file)
--- a/include/deal.II/matrix_free/matrix_free.h
+++ b/include/deal.II/matrix_free/matrix_free.h
@@ -2626,16 +2626,14 @@ inline unsigned int
 MatrixFree<dim, Number, VectorizedArrayType>::n_active_entries_per_cell_batch(
   const unsigned int cell_batch_index) const
 {
+  Assert(!dof_info.empty(), ExcNotInitialized());
   AssertIndexRange(cell_batch_index, task_info.cell_partition_data.back());
-  unsigned int n_lanes = VectorizedArrayType::size();
-  while (n_lanes > 1 &&
-         cell_level_index[cell_batch_index * VectorizedArrayType::size() +
-                          n_lanes - 1] ==
-           cell_level_index[cell_batch_index * VectorizedArrayType::size() +
-                            n_lanes - 2])
-    --n_lanes;
-  AssertIndexRange(n_lanes - 1, VectorizedArrayType::size());
-  return n_lanes;
+  const std::vector<unsigned char> &n_lanes_filled =
+    dof_info[0].n_vectorization_lanes_filled
+      [internal::MatrixFreeFunctions::DoFInfo::dof_access_cell];
+  AssertIndexRange(cell_batch_index, n_lanes_filled.size());
+
+  return n_lanes_filled[cell_batch_index];
 }
 
 
@@ -2646,13 +2644,12 @@ MatrixFree<dim, Number, VectorizedArrayType>::n_active_entries_per_face_batch(
   const unsigned int face_batch_index) const
 {
   AssertIndexRange(face_batch_index, face_info.faces.size());
-  unsigned int n_lanes = VectorizedArrayType::size();
-  while (n_lanes > 1 &&
-         face_info.faces[face_batch_index].cells_interior[n_lanes - 1] ==
-           numbers::invalid_unsigned_int)
-    --n_lanes;
-  AssertIndexRange(n_lanes - 1, VectorizedArrayType::size());
-  return n_lanes;
+  Assert(!dof_info.empty(), ExcNotInitialized());
+  const std::vector<unsigned char> &n_lanes_filled =
+    dof_info[0].n_vectorization_lanes_filled
+      [internal::MatrixFreeFunctions::DoFInfo::dof_access_face_interior];
+  AssertIndexRange(face_batch_index, n_lanes_filled.size());
+  return n_lanes_filled[face_batch_index];
 }
 
 

diff --git a/include/deal.II/matrix_free/matrix_free.templates.h b/include/deal.II/matrix_free/matrix_free.templates.h
index e6adbab2b07fc04a95c350e7884123e5b85c3fa2..5cf7a73d1114094616a7cff0ab786200db035db2 100644 (file)
--- a/include/deal.II/matrix_free/matrix_free.templates.h
+++ b/include/deal.II/matrix_free/matrix_free.templates.h
@@ -514,11 +514,24 @@ MatrixFree<dim, Number, VectorizedArrayType>::internal_reinit(
           dof_info[i].dofs_per_cell.push_back(
             dof_handler[i]->get_fe(0).n_dofs_per_cell());
 
+          const unsigned int n_regular_cells = cell_level_index.size();
           // if indices are not initialized, the cell_level_index might not be
           // divisible by the vectorization length. But it must be for
           // mapping_info...
           while (cell_level_index.size() % VectorizedArrayType::size() != 0)
             cell_level_index.push_back(cell_level_index.back());
+
+          // adjust lengths for vectorization
+          dof_info[i]
+            .n_vectorization_lanes_filled
+              [internal::MatrixFreeFunctions::DoFInfo::dof_access_cell]
+            .resize(cell_level_index.size() / VectorizedArrayType::size(),
+                    VectorizedArrayType::size());
+          if (n_regular_cells < cell_level_index.size())
+            dof_info[i].n_vectorization_lanes_filled
+              [internal::MatrixFreeFunctions::DoFInfo::dof_access_cell]
+              [n_regular_cells / VectorizedArrayType::size()] =
+              n_regular_cells % VectorizedArrayType::size();
         }
     }