OGS
ApplicationUtils::FindCellsForPoints Class Reference

Detailed Description

An algorithm finding mesh elements within a certain tolerance from a given point.

The algorithm will find a single element if for points lying within a mesh element (unless the tolerance is very high), and will find multiple elements for points on element faces/edges/corners.

Definition at line 37 of file FindCellsForPoints.h.

#include <FindCellsForPoints.h>

Public Member Functions

void initialize (vtkUnstructuredGrid *const bulk_mesh)
std::vector< vtkIdType > findCellIds (Eigen::Vector3d const &coords, double const tolerance)
vtkSmartPointer< vtkUnstructuredGrid > find (Eigen::Vector3d const &coords, double const tolerance)

Private Member Functions

vtkSmartPointer< vtkIdList > findCellIdsImpl (Eigen::Vector3d const &coords, double const tolerance)

Private Attributes

vtkSmartPointer< vtkDataSet > bulk_mesh_
vtkSmartPointer< vtkCellTreeLocator > locator_

Member Function Documentation

◆ find()

vtkSmartPointer< vtkUnstructuredGrid > ApplicationUtils::FindCellsForPoints::find ( Eigen::Vector3d const & coords,
double const tolerance )
inline

Definition at line 76 of file FindCellsForPoints.h.

78 {
79 auto const cell_ids = findCellIdsImpl(coords, tolerance);
80
81 vtkNew<vtkExtractCells> extract_cells;
82 extract_cells->SetInputData(bulk_mesh_);
83 extract_cells->SetCellList(cell_ids);
84 extract_cells->Update();
85 auto const filtered_bulk_mesh = extract_cells->GetOutput();
86
87 return filtered_bulk_mesh;
88 }
vtkSmartPointer< vtkDataSet > bulk_mesh_
vtkSmartPointer< vtkIdList > findCellIdsImpl(Eigen::Vector3d const &coords, double const tolerance)

References bulk_mesh_, and findCellIdsImpl().

Referenced by ApplicationUtils::computeNaturalCoords().

◆ findCellIds()

std::vector< vtkIdType > ApplicationUtils::FindCellsForPoints::findCellIds ( Eigen::Vector3d const & coords,
double const tolerance )
inline

Definition at line 65 of file FindCellsForPoints.h.

67 {
68 using namespace ranges;
69 auto const cell_ids = findCellIdsImpl(coords, tolerance);
70 return views::iota(0, cell_ids->GetNumberOfIds()) |
71 views::transform([&](vtkIdType i)
72 { return cell_ids->GetId(i); }) |
73 to_vector;
74 }

References findCellIdsImpl().

◆ findCellIdsImpl()

vtkSmartPointer< vtkIdList > ApplicationUtils::FindCellsForPoints::findCellIdsImpl ( Eigen::Vector3d const & coords,
double const tolerance )
inlineprivate

Definition at line 91 of file FindCellsForPoints.h.

93 {
94 std::array<double, 6> bbox = {
95 coords[0] - tolerance, coords[0] + tolerance,
96 coords[1] - tolerance, coords[1] + tolerance,
97 coords[2] - tolerance, coords[2] + tolerance};
98
99 vtkNew<vtkIdList> cell_ids;
100 locator_->FindCellsWithinBounds(bbox.data(), cell_ids);
101
102 return cell_ids;
103 }
vtkSmartPointer< vtkCellTreeLocator > locator_
constexpr ranges::views::view_closure coords
Definition Mesh.h:223

References locator_.

Referenced by find(), and findCellIds().

◆ initialize()

void ApplicationUtils::FindCellsForPoints::initialize ( vtkUnstructuredGrid *const bulk_mesh)
inline

Definition at line 40 of file FindCellsForPoints.h.

41 {
42#if VTK_VERSION_STRIPPED < 940
43 vtkNew<vtkIdFilter> id_filter;
44 id_filter->SetInputData(bulk_mesh);
45 id_filter->SetCellIdsArrayName("bulk_element_ids");
46 id_filter->SetPointIds(false);
47 id_filter->SetCellIds(true);
48 id_filter->Update();
49 bulk_mesh_ = id_filter->GetOutput();
50#else
51 vtkNew<vtkGenerateIds> generate_ids;
52 generate_ids->SetInputData(bulk_mesh);
53 generate_ids->SetCellIdsArrayName("bulk_element_ids");
54 generate_ids->PointIdsOff();
55 generate_ids->CellIdsOn();
56 generate_ids->Update();
58 static_cast<vtkUnstructuredGrid*>(generate_ids->GetOutput());
59#endif
60 locator_ = vtkSmartPointer<vtkCellTreeLocator>::New();
61 locator_->SetDataSet(bulk_mesh_);
62 locator_->BuildLocator();
63 }

References bulk_mesh_, and locator_.

Referenced by ApplicationUtils::computeNaturalCoords().

Member Data Documentation

◆ bulk_mesh_

vtkSmartPointer<vtkDataSet> ApplicationUtils::FindCellsForPoints::bulk_mesh_
private

Definition at line 105 of file FindCellsForPoints.h.

Referenced by find(), and initialize().

◆ locator_

vtkSmartPointer<vtkCellTreeLocator> ApplicationUtils::FindCellsForPoints::locator_
private

Definition at line 106 of file FindCellsForPoints.h.

Referenced by findCellIdsImpl(), and initialize().


The documentation for this class was generated from the following file: