FEFieldFunction: Detect case when some points were not found

diff --git a/include/deal.II/numerics/fe_field_function.templates.h b/include/deal.II/numerics/fe_field_function.templates.h
index 90dc040011fd2bf3ecc74e4c00d0c841c839a6ad..8c3da383f76e077a430db23e20ecc48338cb3815 100644 (file)
--- a/include/deal.II/numerics/fe_field_function.templates.h
+++ b/include/deal.II/numerics/fe_field_function.templates.h
@@ -523,7 +523,9 @@ namespace Functions
   {
     // Calling the GridTools routine and preparing output
     auto cell_qpoint_map =
-      GridTools::compute_point_locations(cache, points, cell_hint.get());
+      GridTools::compute_point_locations_try_all(cache,
+                                                 points,
+                                                 cell_hint.get());
     const auto &tria_cells = std::get<0>(cell_qpoint_map);
     cells.resize(tria_cells.size());
     unsigned int i = 0;
@@ -531,6 +533,10 @@ namespace Functions
       cells[i++] = typename DoFHandlerType::cell_iterator(*c, dh);
     qpoints = std::get<1>(cell_qpoint_map);
     maps    = std::get<2>(cell_qpoint_map);
+
+    // Ensure that we found all points
+    AssertThrow(std::get<3>(cell_qpoint_map).empty(),
+                VectorTools::ExcPointNotAvailableHere());
     return cells.size();
   }