Loading [MathJax]/extensions/tex2jax.js
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 43 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 82 of file FindCellsForPoints.h.

84 {
85 auto const cell_ids = findCellIdsImpl(coords, tolerance);
86
87 vtkNew<vtkExtractCells> extract_cells;
88 extract_cells->SetInputData(bulk_mesh_);
89 extract_cells->SetCellList(cell_ids);
90 extract_cells->Update();
91 auto const filtered_bulk_mesh = extract_cells->GetOutput();
92
93 return filtered_bulk_mesh;
94 }
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 71 of file FindCellsForPoints.h.

73 {
74 using namespace ranges;
75 auto const cell_ids = findCellIdsImpl(coords, tolerance);
76 return views::iota(0, cell_ids->GetNumberOfIds()) |
77 views::transform([&](vtkIdType i)
78 { return cell_ids->GetId(i); }) |
79 to_vector;
80 }

References findCellIdsImpl().

◆ findCellIdsImpl()

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

Definition at line 97 of file FindCellsForPoints.h.

99 {
100 std::array<double, 6> bbox = {
101 coords[0] - tolerance, coords[0] + tolerance,
102 coords[1] - tolerance, coords[1] + tolerance,
103 coords[2] - tolerance, coords[2] + tolerance};
104
105 vtkNew<vtkIdList> cell_ids;
106 locator_->FindCellsWithinBounds(bbox.data(), cell_ids);
107
108 return cell_ids;
109 }
vtkSmartPointer< vtkCellTreeLocator > locator_
constexpr ranges::views::view_closure coords
Definition Mesh.h:234

References locator_.

Referenced by find(), and findCellIds().

◆ initialize()

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

Definition at line 46 of file FindCellsForPoints.h.

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

References bulk_mesh_, and locator_.

Referenced by ApplicationUtils::computeNaturalCoords().

Member Data Documentation

◆ bulk_mesh_

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

Definition at line 111 of file FindCellsForPoints.h.

Referenced by find(), and initialize().

◆ locator_

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

Definition at line 112 of file FindCellsForPoints.h.

Referenced by findCellIdsImpl(), and initialize().


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