From: Martin Kronbichler Date: Tue, 29 May 2018 14:54:57 +0000 (+0200) Subject: Move all find_active_cell_around_point code to grid_tools_dof_handlers. X-Git-Tag: v9.1.0-rc1~1066^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=06ae3f3f86257ab88c4afcec70c5b9dc68a33268;p=dealii.git Move all find_active_cell_around_point code to grid_tools_dof_handlers. --- diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index 8258bdf239..ecdea29bf6 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -1421,409 +1421,6 @@ namespace GridTools - namespace - { - template class MeshType, int spacedim> - void - find_active_cell_around_point_internal( - const MeshType &mesh, -#ifndef _MSC_VER - std::set::active_cell_iterator> - &searched_cells, - std::set::active_cell_iterator> - &adjacent_cells) -#else - std::set< - typename dealii::internal:: - ActiveCellIterator>::type> - &searched_cells, - std::set< - typename dealii::internal:: - ActiveCellIterator>::type> - &adjacent_cells) -#endif - { -#ifndef _MSC_VER - typedef - typename MeshType::active_cell_iterator cell_iterator; -#else - typedef typename dealii::internal:: - ActiveCellIterator>::type - cell_iterator; -#endif - - // update the searched cells - searched_cells.insert(adjacent_cells.begin(), adjacent_cells.end()); - // now we to collect all neighbors - // of the cells in adjacent_cells we - // have not yet searched. - std::set adjacent_cells_new; - - typename std::set::const_iterator cell = - adjacent_cells.begin(), - endc = - adjacent_cells.end(); - for (; cell != endc; ++cell) - { - std::vector active_neighbors; - get_active_neighbors>(*cell, - active_neighbors); - for (unsigned int i = 0; i < active_neighbors.size(); ++i) - if (searched_cells.find(active_neighbors[i]) == - searched_cells.end()) - adjacent_cells_new.insert(active_neighbors[i]); - } - adjacent_cells.clear(); - adjacent_cells.insert(adjacent_cells_new.begin(), - adjacent_cells_new.end()); - if (adjacent_cells.size() == 0) - { - // we haven't found any other cell that would be a - // neighbor of a previously found cell, but we know - // that we haven't checked all cells yet. that means - // that the domain is disconnected. in that case, - // choose the first previously untouched cell we - // can find - cell_iterator it = mesh.begin_active(); - for (; it != mesh.end(); ++it) - if (searched_cells.find(it) == searched_cells.end()) - { - adjacent_cells.insert(it); - break; - } - } - } - - - - template class MeshType, int spacedim> -#ifndef _MSC_VER - std::pair::active_cell_iterator, - Point> -#else - std::pair< - typename dealii::internal:: - ActiveCellIterator>::type, - Point> -#endif - find_active_cell_around_point_tolerance( - const Mapping & mapping, - const MeshType &mesh, - const Point & p, - const std::vector & marked_vertices, - const double tolerance) - { - typedef typename dealii::internal:: - ActiveCellIterator>::type - active_cell_iterator; - - // The best distance is set to the - // maximum allowable distance from - // the unit cell; we assume a - // max. deviation of the given tolerance - double best_distance = tolerance; - int best_level = -1; - std::pair> best_cell; - - // Find closest vertex and determine - // all adjacent cells - std::vector adjacent_cells_tmp = - find_cells_adjacent_to_vertex( - mesh, find_closest_vertex(mapping, mesh, p, marked_vertices)); - - // Make sure that we have found - // at least one cell adjacent to vertex. - Assert(adjacent_cells_tmp.size() > 0, ExcInternalError()); - - // Copy all the cells into a std::set - std::set adjacent_cells(adjacent_cells_tmp.begin(), - adjacent_cells_tmp.end()); - std::set searched_cells; - - // Determine the maximal number of cells - // in the grid. - // As long as we have not found - // the cell and have not searched - // every cell in the triangulation, - // we keep on looking. - const unsigned int n_active_cells = - mesh.get_triangulation().n_active_cells(); - bool found = false; - unsigned int cells_searched = 0; - while (!found && cells_searched < n_active_cells) - { - typename std::set::const_iterator - cell = adjacent_cells.begin(), - endc = adjacent_cells.end(); - for (; cell != endc; ++cell) - { - try - { - const Point p_cell = - mapping.transform_real_to_unit_cell(*cell, p); - - // calculate the infinity norm of - // the distance vector to the unit cell. - const double dist = - GeometryInfo::distance_to_unit_cell(p_cell); - - // We compare if the point is inside the - // unit cell (or at least not too far - // outside). If it is, it is also checked - // that the cell has a more refined state - if ((dist < best_distance) || - ((dist == best_distance) && - ((*cell)->level() > best_level))) - { - found = true; - best_distance = dist; - best_level = (*cell)->level(); - best_cell = std::make_pair(*cell, p_cell); - } - } - catch ( - typename MappingQGeneric::ExcTransformationFailed - &) - { - // ok, the transformation - // failed presumably - // because the point we - // are looking for lies - // outside the current - // cell. this means that - // the current cell can't - // be the cell around the - // point, so just ignore - // this cell and move on - // to the next - } - } - - // update the number of cells searched - cells_searched += adjacent_cells.size(); - - // if the user provided a custom mask for vertices, - // terminate the search without trying to expand the search - // to all cells of the triangulation, as done below. - if (marked_vertices.size() > 0) - cells_searched = n_active_cells; - - // if we have not found the cell in - // question and have not yet searched every - // cell, we expand our search to - // all the not already searched neighbors of - // the cells in adjacent_cells. This is - // what find_active_cell_around_point_internal - // is for. - if (!found && cells_searched < n_active_cells) - { - find_active_cell_around_point_internal( - mesh, searched_cells, adjacent_cells); - } - } - - AssertThrow(best_cell.first.state() == IteratorState::valid, - ExcPointNotFound(p)); - - return best_cell; - } - } // namespace - - - - template class MeshType, int spacedim> -#ifndef _MSC_VER - typename MeshType::active_cell_iterator -#else - typename dealii::internal:: - ActiveCellIterator>::type -#endif - find_active_cell_around_point(const MeshType &mesh, - const Point & p, - const std::vector & marked_vertices) - { - return find_active_cell_around_point( - StaticMappingQ1::mapping, mesh, p, marked_vertices) - .first; - } - - - - template class MeshType, int spacedim> -#ifndef _MSC_VER - std::pair::active_cell_iterator, Point> -#else - std::pair>::type, - Point> -#endif - find_active_cell_around_point(const Mapping & mapping, - const MeshType &mesh, - const Point & p, - const std::vector & marked_vertices) - { - return find_active_cell_around_point_tolerance( - mapping, mesh, p, marked_vertices, 1e-10); - } - - - - template class MeshType, int spacedim> -#ifndef _MSC_VER - std::vector::active_cell_iterator, - Point>> -#else - std::vector>::type, - Point>> -#endif - find_all_active_cells_around_point(const Mapping & mapping, - const MeshType &mesh, - const Point & p, - const double tolerance, - const std::vector &marked_vertices) - { - // first use the result of the single point function as a guess. In order - // not to make the other find_all_active_cells_around_point more expensive - // and avoid some additional logic there, we first start with one cell as - // given by that other function (that possibly goes through a larger set - // of cells) and later add a list of more cells as appropriate. - std::vector< - std::pair::active_cell_iterator, - Point>> - cells_and_points; - try - { - cells_and_points.push_back(find_active_cell_around_point_tolerance( - mapping, mesh, p, marked_vertices, tolerance)); - } - catch (ExcPointNotFound &) - {} - - if (!cells_and_points.empty()) - { - // check if the given point is on the surface of the unit cell. if yes, - // need to find all neighbors - const Point unit_point = cells_and_points.front().second; - const auto my_cell = cells_and_points.front().first; - Tensor<1, dim> distance_to_center; - unsigned int n_dirs_at_threshold = 0; - unsigned int last_point_at_threshold = numbers::invalid_unsigned_int; - for (unsigned int d = 0; d < dim; ++d) - { - distance_to_center[d] = std::abs(unit_point[d] - 0.5); - if (distance_to_center[d] > 0.5 - tolerance) - { - ++n_dirs_at_threshold; - last_point_at_threshold = d; - } - } - - std::vector::active_cell_iterator> - cells_to_add; - // point is within face -> only need neighbor - if (n_dirs_at_threshold == 1) - { - unsigned int neighbor_index = - 2 * last_point_at_threshold + - (unit_point[last_point_at_threshold] > 0.5 ? 1 : 0); - if (!my_cell->at_boundary(neighbor_index)) - cells_to_add.push_back(my_cell->neighbor(neighbor_index)); - } - // corner point -> use all neighbors - else if (n_dirs_at_threshold == dim) - { - unsigned int local_vertex_index = 0; - for (unsigned int d = 0; d < dim; ++d) - local_vertex_index += (unit_point[d] > 0.5 ? 1 : 0) << d; - std::vector::active_cell_iterator> - cells = find_cells_adjacent_to_vertex( - mesh, my_cell->vertex_index(local_vertex_index)); - for (auto cell : cells) - if (cell != my_cell) - cells_to_add.push_back(cell); - } - // point on line in 3D: We cannot simply take the intersection between - // the two vertices of cels because of hanging nodes. So instead we - // list the vertices around both points and then select the - // appropriate cells according to the result of read_to_unit_cell - // below. - else if (n_dirs_at_threshold == 2) - { - std::pair vertex_indices[3]; - unsigned int count_vertex_indices = 0; - unsigned int free_direction = numbers::invalid_unsigned_int; - for (unsigned int d = 0; d < dim; ++d) - { - if (distance_to_center[d] > 0.5 - tolerance) - { - vertex_indices[count_vertex_indices].first = d; - vertex_indices[count_vertex_indices].second = - unit_point[d] > 0.5 ? 1 : 0; - count_vertex_indices++; - } - else - free_direction = d; - } - - AssertDimension(count_vertex_indices, 2); - Assert(free_direction != numbers::invalid_unsigned_int, - ExcInternalError()); - - const unsigned int first_vertex = - (vertex_indices[0].second << vertex_indices[0].first) + - (vertex_indices[1].second << vertex_indices[1].first); - for (unsigned int d = 0; d < 2; ++d) - { - auto tentative_cells = find_cells_adjacent_to_vertex( - mesh, - my_cell->vertex_index(first_vertex + (d << free_direction))); - for (auto cell : tentative_cells) - { - bool cell_not_yet_present = true; - for (auto other_cell : cells_to_add) - if (cell == other_cell) - { - cell_not_yet_present = false; - break; - } - if (cell_not_yet_present) - cells_to_add.push_back(cell); - } - } - } - - const double original_distance_to_unit_cell = - GeometryInfo::distance_to_unit_cell(unit_point); - for (auto cell : cells_to_add) - { - if (cell != my_cell) - try - { - const Point p_unit = - mapping.transform_real_to_unit_cell(cell, p); - if (GeometryInfo::distance_to_unit_cell(p_unit) < - original_distance_to_unit_cell + tolerance) - cells_and_points.emplace_back(cell, p_unit); - } - catch (typename Mapping::ExcTransformationFailed &) - {} - } - } - std::sort( - cells_and_points.begin(), - cells_and_points.end(), - [](const std::pair::active_cell_iterator, - Point> &a, - const std::pair::active_cell_iterator, - Point> &b) { return a.first < b.first; }); - - return cells_and_points; - } - - - template std::vector>> vertex_to_cell_centers_directions( diff --git a/source/grid/grid_tools.inst.in b/source/grid/grid_tools.inst.in index 3158b88c8d..7abb30f250 100644 --- a/source/grid/grid_tools.inst.in +++ b/source/grid/grid_tools.inst.in @@ -92,17 +92,6 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) deal_II_space_dimension>::active_cell_iterator &, const std::vector &); - template std::vector::active_cell_iterator, - Point>> - find_all_active_cells_around_point( - const Mapping &, - const Triangulation &, - const Point &, - const double, - const std::vector &); - template std::tuple::active_cell_iterator>, @@ -145,34 +134,6 @@ for (deal_II_space_dimension : SPACE_DIMENSIONS) const parallel::distributed::Triangulation &, const Point &p); - - std::pair< - dealii::internal::ActiveCellIterator< - deal_II_space_dimension, - deal_II_space_dimension, - parallel::distributed::Triangulation>::type, - Point> - find_active_cell_around_point( - const Mapping &, - const parallel::distributed::Triangulation &, - const Point &, - const std::vector &); - - std::vector>::type, - Point>> - find_all_active_cells_around_point( - const Mapping &, - const parallel::distributed::Triangulation &, - const Point &, - const double, - const std::vector &); - template unsigned int GridTools::find_closest_vertex( const std::map> &vertices, const Point & p); diff --git a/source/grid/grid_tools_dof_handlers.cc b/source/grid/grid_tools_dof_handlers.cc index d3ddfc3ab6..9e5760e287 100644 --- a/source/grid/grid_tools_dof_handlers.cc +++ b/source/grid/grid_tools_dof_handlers.cc @@ -403,8 +403,144 @@ namespace GridTools } } } + + + + template class MeshType, int spacedim> +#ifndef _MSC_VER + std::pair::active_cell_iterator, + Point> +#else + std::pair< + typename dealii::internal:: + ActiveCellIterator>::type, + Point> +#endif + find_active_cell_around_point_tolerance( + const Mapping & mapping, + const MeshType &mesh, + const Point & p, + const std::vector & marked_vertices, + const double tolerance) + { + typedef typename dealii::internal:: + ActiveCellIterator>::type + active_cell_iterator; + + // The best distance is set to the + // maximum allowable distance from + // the unit cell; we assume a + // max. deviation of the given tolerance + double best_distance = tolerance; + int best_level = -1; + std::pair> best_cell; + + // Find closest vertex and determine + // all adjacent cells + std::vector adjacent_cells_tmp = + find_cells_adjacent_to_vertex( + mesh, find_closest_vertex(mapping, mesh, p, marked_vertices)); + + // Make sure that we have found + // at least one cell adjacent to vertex. + Assert(adjacent_cells_tmp.size() > 0, ExcInternalError()); + + // Copy all the cells into a std::set + std::set adjacent_cells(adjacent_cells_tmp.begin(), + adjacent_cells_tmp.end()); + std::set searched_cells; + + // Determine the maximal number of cells + // in the grid. + // As long as we have not found + // the cell and have not searched + // every cell in the triangulation, + // we keep on looking. + const unsigned int n_active_cells = + mesh.get_triangulation().n_active_cells(); + bool found = false; + unsigned int cells_searched = 0; + while (!found && cells_searched < n_active_cells) + { + typename std::set::const_iterator + cell = adjacent_cells.begin(), + endc = adjacent_cells.end(); + for (; cell != endc; ++cell) + { + try + { + const Point p_cell = + mapping.transform_real_to_unit_cell(*cell, p); + + // calculate the infinity norm of + // the distance vector to the unit cell. + const double dist = + GeometryInfo::distance_to_unit_cell(p_cell); + + // We compare if the point is inside the + // unit cell (or at least not too far + // outside). If it is, it is also checked + // that the cell has a more refined state + if ((dist < best_distance) || + ((dist == best_distance) && + ((*cell)->level() > best_level))) + { + found = true; + best_distance = dist; + best_level = (*cell)->level(); + best_cell = std::make_pair(*cell, p_cell); + } + } + catch ( + typename MappingQGeneric::ExcTransformationFailed + &) + { + // ok, the transformation + // failed presumably + // because the point we + // are looking for lies + // outside the current + // cell. this means that + // the current cell can't + // be the cell around the + // point, so just ignore + // this cell and move on + // to the next + } + } + + // update the number of cells searched + cells_searched += adjacent_cells.size(); + + // if the user provided a custom mask for vertices, + // terminate the search without trying to expand the search + // to all cells of the triangulation, as done below. + if (marked_vertices.size() > 0) + cells_searched = n_active_cells; + + // if we have not found the cell in + // question and have not yet searched every + // cell, we expand our search to + // all the not already searched neighbors of + // the cells in adjacent_cells. This is + // what find_active_cell_around_point_internal + // is for. + if (!found && cells_searched < n_active_cells) + { + find_active_cell_around_point_internal( + mesh, searched_cells, adjacent_cells); + } + } + + AssertThrow(best_cell.first.state() == IteratorState::valid, + ExcPointNotFound(p)); + + return best_cell; + } } // namespace + + template class MeshType, int spacedim> #ifndef _MSC_VER typename MeshType::active_cell_iterator @@ -422,6 +558,7 @@ namespace GridTools } + template class MeshType, int spacedim> #ifndef _MSC_VER std::pair::active_cell_iterator, Point> @@ -435,118 +572,164 @@ namespace GridTools const Point & p, const std::vector & marked_vertices) { - typedef typename dealii::internal:: - ActiveCellIterator>::type - active_cell_iterator; - - // The best distance is set to the - // maximum allowable distance from - // the unit cell; we assume a - // max. deviation of 1e-10 - double best_distance = 1e-10; - int best_level = -1; - std::pair> best_cell; - - // Find closest vertex and determine - // all adjacent cells - std::vector adjacent_cells_tmp = - find_cells_adjacent_to_vertex( - mesh, find_closest_vertex(mapping, mesh, p, marked_vertices)); - - // Make sure that we have found - // at least one cell adjacent to vertex. - Assert(adjacent_cells_tmp.size() > 0, ExcInternalError()); - - // Copy all the cells into a std::set - std::set adjacent_cells(adjacent_cells_tmp.begin(), - adjacent_cells_tmp.end()); - std::set searched_cells; - - // Determine the maximal number of cells - // in the grid. - // As long as we have not found - // the cell and have not searched - // every cell in the triangulation, - // we keep on looking. - const unsigned int n_active_cells = - mesh.get_triangulation().n_active_cells(); - bool found = false; - unsigned int cells_searched = 0; - while (!found && cells_searched < n_active_cells) + return find_active_cell_around_point_tolerance( + mapping, mesh, p, marked_vertices, 1e-10); + } + + + + template class MeshType, int spacedim> +#ifndef _MSC_VER + std::vector::active_cell_iterator, + Point>> +#else + std::vector>::type, + Point>> +#endif + find_all_active_cells_around_point(const Mapping & mapping, + const MeshType &mesh, + const Point & p, + const double tolerance, + const std::vector &marked_vertices) + { + // first use the result of the single point function as a guess. In order + // not to make the other find_all_active_cells_around_point more expensive + // and avoid some additional logic there, we first start with one cell as + // given by that other function (that possibly goes through a larger set + // of cells) and later add a list of more cells as appropriate. + std::vector< + std::pair::active_cell_iterator, + Point>> + cells_and_points; + try { - typename std::set::const_iterator - cell = adjacent_cells.begin(), - endc = adjacent_cells.end(); - for (; cell != endc; ++cell) + cells_and_points.push_back(find_active_cell_around_point_tolerance( + mapping, mesh, p, marked_vertices, tolerance)); + } + catch (ExcPointNotFound &) + {} + + if (!cells_and_points.empty()) + { + // check if the given point is on the surface of the unit cell. if yes, + // need to find all neighbors + const Point unit_point = cells_and_points.front().second; + const auto my_cell = cells_and_points.front().first; + Tensor<1, dim> distance_to_center; + unsigned int n_dirs_at_threshold = 0; + unsigned int last_point_at_threshold = numbers::invalid_unsigned_int; + for (unsigned int d = 0; d < dim; ++d) { - try + distance_to_center[d] = std::abs(unit_point[d] - 0.5); + if (distance_to_center[d] > 0.5 - tolerance) { - const Point p_cell = - mapping.transform_real_to_unit_cell(*cell, p); - - // calculate the infinity norm of - // the distance vector to the unit cell. - const double dist = - GeometryInfo::distance_to_unit_cell(p_cell); - - // We compare if the point is inside the - // unit cell (or at least not too far - // outside). If it is, it is also checked - // that the cell has a more refined state - if ((dist < best_distance) || ((dist == best_distance) && - ((*cell)->level() > best_level))) + ++n_dirs_at_threshold; + last_point_at_threshold = d; + } + } + + std::vector::active_cell_iterator> + cells_to_add; + // point is within face -> only need neighbor + if (n_dirs_at_threshold == 1) + { + unsigned int neighbor_index = + 2 * last_point_at_threshold + + (unit_point[last_point_at_threshold] > 0.5 ? 1 : 0); + if (!my_cell->at_boundary(neighbor_index)) + cells_to_add.push_back(my_cell->neighbor(neighbor_index)); + } + // corner point -> use all neighbors + else if (n_dirs_at_threshold == dim) + { + unsigned int local_vertex_index = 0; + for (unsigned int d = 0; d < dim; ++d) + local_vertex_index += (unit_point[d] > 0.5 ? 1 : 0) << d; + std::vector::active_cell_iterator> + cells = find_cells_adjacent_to_vertex( + mesh, my_cell->vertex_index(local_vertex_index)); + for (auto cell : cells) + if (cell != my_cell) + cells_to_add.push_back(cell); + } + // point on line in 3D: We cannot simply take the intersection between + // the two vertices of cels because of hanging nodes. So instead we + // list the vertices around both points and then select the + // appropriate cells according to the result of read_to_unit_cell + // below. + else if (n_dirs_at_threshold == 2) + { + std::pair vertex_indices[3]; + unsigned int count_vertex_indices = 0; + unsigned int free_direction = numbers::invalid_unsigned_int; + for (unsigned int d = 0; d < dim; ++d) + { + if (distance_to_center[d] > 0.5 - tolerance) { - found = true; - best_distance = dist; - best_level = (*cell)->level(); - best_cell = std::make_pair(*cell, p_cell); + vertex_indices[count_vertex_indices].first = d; + vertex_indices[count_vertex_indices].second = + unit_point[d] > 0.5 ? 1 : 0; + count_vertex_indices++; } + else + free_direction = d; } - catch ( - typename MappingQGeneric::ExcTransformationFailed - &) + + AssertDimension(count_vertex_indices, 2); + Assert(free_direction != numbers::invalid_unsigned_int, + ExcInternalError()); + + const unsigned int first_vertex = + (vertex_indices[0].second << vertex_indices[0].first) + + (vertex_indices[1].second << vertex_indices[1].first); + for (unsigned int d = 0; d < 2; ++d) { - // ok, the transformation - // failed presumably - // because the point we - // are looking for lies - // outside the current - // cell. this means that - // the current cell can't - // be the cell around the - // point, so just ignore - // this cell and move on - // to the next + auto tentative_cells = find_cells_adjacent_to_vertex( + mesh, + my_cell->vertex_index(first_vertex + (d << free_direction))); + for (auto cell : tentative_cells) + { + bool cell_not_yet_present = true; + for (auto other_cell : cells_to_add) + if (cell == other_cell) + { + cell_not_yet_present = false; + break; + } + if (cell_not_yet_present) + cells_to_add.push_back(cell); + } } } - // update the number of cells searched - cells_searched += adjacent_cells.size(); - - // if the user provided a custom mask for vertices, - // terminate the search without trying to expand the search - // to all cells of the triangulation, as done below. - if (marked_vertices.size() > 0) - cells_searched = n_active_cells; - - // if we have not found the cell in - // question and have not yet searched every - // cell, we expand our search to - // all the not already searched neighbors of - // the cells in adjacent_cells. This is - // what find_active_cell_around_point_internal - // is for. - if (!found && cells_searched < n_active_cells) + const double original_distance_to_unit_cell = + GeometryInfo::distance_to_unit_cell(unit_point); + for (auto cell : cells_to_add) { - find_active_cell_around_point_internal( - mesh, searched_cells, adjacent_cells); + if (cell != my_cell) + try + { + const Point p_unit = + mapping.transform_real_to_unit_cell(cell, p); + if (GeometryInfo::distance_to_unit_cell(p_unit) < + original_distance_to_unit_cell + tolerance) + cells_and_points.emplace_back(cell, p_unit); + } + catch (typename Mapping::ExcTransformationFailed &) + {} } } - - AssertThrow(best_cell.first.state() == IteratorState::valid, - ExcPointNotFound(p)); - - return best_cell; + std::sort( + cells_and_points.begin(), + cells_and_points.end(), + [](const std::pair::active_cell_iterator, + Point> &a, + const std::pair::active_cell_iterator, + Point> &b) { return a.first < b.first; }); + + return cells_and_points; } diff --git a/source/grid/grid_tools_dof_handlers.inst.in b/source/grid/grid_tools_dof_handlers.inst.in index f2f0f185b0..bd737f3501 100644 --- a/source/grid/grid_tools_dof_handlers.inst.in +++ b/source/grid/grid_tools_dof_handlers.inst.in @@ -41,6 +41,18 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS; const Point &, const std::vector &); + template std::vector< + std::pair::type, + Point>> + find_all_active_cells_around_point( + const Mapping &, + const X &, + const Point &, + const double, + const std::vector &); + template std::vector< dealii::internal::ActiveCellIterator