From 8d43cf50ebba4335a767ad9f71d027d21b9785d3 Mon Sep 17 00:00:00 2001
From: Daniel Arndt <arndtd@ornl.gov>
Date: Wed, 7 Dec 2022 22:15:41 +0000
Subject: [PATCH] Revert changing due to boost

---
 source/grid/grid_tools.cc            | 25 +++++++++++++++++++++-
 source/particles/particle_handler.cc | 32 +++++++++++++++++++++++++++-
 2 files changed, 55 insertions(+), 2 deletions(-)

diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc
index 9d5da20c39..dc2aa6242a 100644
--- a/source/grid/grid_tools.cc
+++ b/source/grid/grid_tools.cc
@@ -2963,7 +2963,30 @@ namespace GridTools
           {
             if (!used_vertices_rtree.empty())
               {
-                Assert(false, ExcInternalError());
+                // If we have an rtree at our disposal, use it.
+                using ValueType = std::pair<Point<spacedim>, unsigned int>;
+                std::function<bool(const ValueType &)> marked;
+                if (marked_vertices.size() == mesh.n_vertices())
+                  marked = [&marked_vertices](const ValueType &value) -> bool {
+                    return marked_vertices[value.second];
+                  };
+                else
+                  marked = [](const ValueType &) -> bool { return true; };
+
+                std::vector<std::pair<Point<spacedim>, unsigned int>> res;
+                used_vertices_rtree.query(
+                  boost::geometry::index::nearest(p, 1) &&
+                    boost::geometry::index::satisfies(marked),
+                  std::back_inserter(res));
+
+                // Searching for a point which is located outside the
+                // triangulation results in res.size() = 0
+                Assert(res.size() < 2,
+                       dealii::ExcMessage("There can not be multiple results"));
+
+                if (res.size() > 0)
+                  if (any_cell_marked(vertex_to_cells[res[0].second]))
+                    closest_vertex_index = res[0].second;
               }
             else
               {
diff --git a/source/particles/particle_handler.cc b/source/particles/particle_handler.cc
index 20f260cc6b..e55fe4da7f 100644
--- a/source/particles/particle_handler.cc
+++ b/source/particles/particle_handler.cc
@@ -1395,7 +1395,37 @@ namespace Particles
 
           if (!found_cell)
             {
-              Assert(false, ExcInternalError());
+              // The particle is not in a neighbor of the old cell.
+              // Look for the new cell in the whole local domain.
+              // This case is rare.
+              std::vector<std::pair<Point<spacedim>, unsigned int>>
+                closest_vertex_in_domain;
+              triangulation_cache->get_used_vertices_rtree().query(
+                boost::geometry::index::nearest(out_particle->get_location(),
+                                                1),
+                std::back_inserter(closest_vertex_in_domain));
+
+              // We should have one and only one result
+              AssertDimension(closest_vertex_in_domain.size(), 1);
+              const unsigned int closest_vertex_index_in_domain =
+                closest_vertex_in_domain[0].second;
+
+              // Search all of the cells adjacent to the closest vertex of the
+              // domain. Most likely we will find the particle in them.
+              for (const auto &cell :
+                   vertex_to_cells[closest_vertex_index_in_domain])
+                {
+                  mapping->transform_points_real_to_unit_cell(
+                    cell, real_locations, reference_locations);
+
+                  if (GeometryInfo<dim>::is_inside_unit_cell(
+                        reference_locations[0]))
+                    {
+                      current_cell = cell;
+                      found_cell   = true;
+                      break;
+                    }
+                }
             }
 
           if (!found_cell)
-- 
2.39.5