From 403f0779417b474ae7527262c8ec601a3c4e41d0 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Wed, 25 Nov 2020 16:14:12 +0100 Subject: [PATCH] MatrixFree: use active_fe_indices within FaceComparator --- .../deal.II/matrix_free/face_setup_internal.h | 61 +++++++++++++++++-- .../matrix_free/matrix_free.templates.h | 12 ++-- 2 files changed, 64 insertions(+), 9 deletions(-) diff --git a/include/deal.II/matrix_free/face_setup_internal.h b/include/deal.II/matrix_free/face_setup_internal.h index bf154269a2..a2f0c5462e 100644 --- a/include/deal.II/matrix_free/face_setup_internal.h +++ b/include/deal.II/matrix_free/face_setup_internal.h @@ -1023,8 +1023,11 @@ namespace internal * to batch similar faces together for vectorization. */ inline bool - compare_faces_for_vectorization(const FaceToCellTopology<1> &face1, - const FaceToCellTopology<1> &face2) + compare_faces_for_vectorization( + const FaceToCellTopology<1> & face1, + const FaceToCellTopology<1> & face2, + const std::vector &active_fe_indices, + const unsigned int length) { if (face1.interior_face_no != face2.interior_face_no) return false; @@ -1034,6 +1037,19 @@ namespace internal return false; if (face1.face_orientation != face2.face_orientation) return false; + + if (active_fe_indices.size() > 0) + { + if (active_fe_indices[face1.cells_interior[0] / length] != + active_fe_indices[face2.cells_interior[0] / length]) + return false; + + if (face2.cells_exterior[0] != numbers::invalid_unsigned_int) + if (active_fe_indices[face1.cells_exterior[0] / length] != + active_fe_indices[face2.cells_exterior[0] / length]) + return false; + } + return true; } @@ -1048,10 +1064,37 @@ namespace internal template struct FaceComparator { + FaceComparator(const std::vector &active_fe_indices) + : active_fe_indices(active_fe_indices) + {} + bool operator()(const FaceToCellTopology &face1, const FaceToCellTopology &face2) const { + // check if active fe indices match + if (active_fe_indices.size() > 0) + { + // ... for interior faces + if (active_fe_indices[face1.cells_interior[0] / length] < + active_fe_indices[face2.cells_interior[0] / length]) + return true; + else if (active_fe_indices[face1.cells_interior[0] / length] > + active_fe_indices[face2.cells_interior[0] / length]) + return false; + + // ... for exterior faces + if (face2.cells_exterior[0] != numbers::invalid_unsigned_int) + { + if (active_fe_indices[face1.cells_exterior[0] / length] < + active_fe_indices[face2.cells_exterior[0] / length]) + return true; + else if (active_fe_indices[face1.cells_exterior[0] / length] > + active_fe_indices[face2.cells_exterior[0] / length]) + return false; + } + } + for (unsigned int i = 0; i < length; ++i) if (face1.cells_interior[i] < face2.cells_interior[i]) return true; @@ -1079,6 +1122,9 @@ namespace internal return false; } + + private: + const std::vector &active_fe_indices; }; @@ -1089,7 +1135,8 @@ namespace internal const std::vector> &faces_in, const std::vector & hard_vectorization_boundary, std::vector & face_partition_data, - std::vector> &faces_out) + std::vector> &faces_out, + const std::vector & active_fe_indices) { FaceToCellTopology macro_face; std::vector> faces_type; @@ -1119,7 +1166,9 @@ namespace internal { // Compare current face with first face of type type if (compare_faces_for_vectorization(faces_in[face], - faces_in[face_type[0]])) + faces_in[face_type[0]], + active_fe_indices, + vectorization_width)) { face_type.push_back(face); goto face_found; @@ -1131,9 +1180,11 @@ namespace internal } // insert new faces in sorted list to get good data locality + FaceComparator face_comparator( + active_fe_indices); std::set, FaceComparator> - new_faces; + new_faces(face_comparator); for (const auto &face_type : faces_type) { macro_face.interior_face_no = diff --git a/include/deal.II/matrix_free/matrix_free.templates.h b/include/deal.II/matrix_free/matrix_free.templates.h index 01809a8ea9..c0dce20168 100644 --- a/include/deal.II/matrix_free/matrix_free.templates.h +++ b/include/deal.II/matrix_free/matrix_free.templates.h @@ -1403,7 +1403,8 @@ MatrixFree::initialize_indices( face_setup.inner_faces, hard_vectorization_boundary, task_info.face_partition_data, - face_info.faces); + face_info.faces, + dof_info[0].cell_active_fe_index); // on boundary faces, we must also respect the vectorization boundary of // the inner faces because we might have dependencies on ghosts of @@ -1412,7 +1413,8 @@ MatrixFree::initialize_indices( face_setup.boundary_faces, hard_vectorization_boundary, task_info.boundary_partition_data, - face_info.faces); + face_info.faces, + dof_info[0].cell_active_fe_index); // for the other ghosted faces, there are no scheduling restrictions hard_vectorization_boundary.clear(); @@ -1422,7 +1424,8 @@ MatrixFree::initialize_indices( face_setup.inner_ghost_faces, hard_vectorization_boundary, task_info.ghost_face_partition_data, - face_info.faces); + face_info.faces, + dof_info[0].cell_active_fe_index); hard_vectorization_boundary.clear(); hard_vectorization_boundary.resize( task_info.refinement_edge_face_partition_data.size(), false); @@ -1430,7 +1433,8 @@ MatrixFree::initialize_indices( face_setup.refinement_edge_faces, hard_vectorization_boundary, task_info.refinement_edge_face_partition_data, - face_info.faces); + face_info.faces, + dof_info[0].cell_active_fe_index); cell_level_index.resize( cell_level_index.size() + -- 2.39.5