From: Daniel Arndt Date: Wed, 7 Dec 2022 22:15:41 +0000 (+0000) Subject: Revert changing due to boost X-Git-Tag: v9.5.0-rc1~697^2~19 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=8d43cf50ebba4335a767ad9f71d027d21b9785d3;p=dealii.git Revert changing due to boost --- 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, unsigned int>; + std::function 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, 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, 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::is_inside_unit_cell( + reference_locations[0])) + { + current_cell = cell; + found_cell = true; + break; + } + } } if (!found_cell)