OGS
MeshLib::VtkMappedMeshSource Class Referencefinal

Detailed Description

Adapter which maps a MeshLib::Mesh to a vtkUnstructuredGridAlgorithm. Allows for zero-copy access of the mesh from the visualization side.

Definition at line 44 of file VtkMappedMeshSource.h.

#include <VtkMappedMeshSource.h>

Inheritance diagram for MeshLib::VtkMappedMeshSource:
[legend]
Collaboration diagram for MeshLib::VtkMappedMeshSource:
[legend]

Public Member Functions

 vtkTypeMacro (VtkMappedMeshSource, vtkUnstructuredGridAlgorithm)
 
void PrintSelf (std::ostream &os, vtkIndent indent) override
 
 VtkMappedMeshSource (const VtkMappedMeshSource &)=delete
 
void operator= (const VtkMappedMeshSource &)=delete
 
void SetMesh (const MeshLib::Mesh *mesh)
 Sets the mesh. Calling is mandatory. More...
 
const MeshLib::MeshGetMesh () const
 Returns the mesh. More...
 

Static Public Member Functions

static VtkMappedMeshSourceNew ()
 

Protected Member Functions

 VtkMappedMeshSource ()
 
int ProcessRequest (vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
 
int RequestData (vtkInformation *, vtkInformationVector **, vtkInformationVector *) override
 
int RequestInformation (vtkInformation *, vtkInformationVector **, vtkInformationVector *) override
 

Private Member Functions

template<typename T >
void addProperty (MeshLib::PropertyVector< T > const &property) const
 Adds a zero-copy vtk array wrapper. More...
 

Private Attributes

const MeshLib::Mesh_mesh {}
 
int NumberOfDimensions {0}
 
int NumberOfNodes {0}
 
vtkNew< vtkPoints > Points
 
vtkNew< vtkPointData > PointData
 
vtkNew< vtkCellData > CellData
 
vtkNew< vtkFieldData > FieldData
 

Constructor & Destructor Documentation

◆ VtkMappedMeshSource() [1/2]

MeshLib::VtkMappedMeshSource::VtkMappedMeshSource ( const VtkMappedMeshSource )
delete

◆ VtkMappedMeshSource() [2/2]

MeshLib::VtkMappedMeshSource::VtkMappedMeshSource ( )
protected

Definition at line 40 of file VtkMappedMeshSource.cpp.

42 {
43  this->SetNumberOfInputPorts(0);
44 }

Member Function Documentation

◆ addProperty()

template<typename T >
void MeshLib::VtkMappedMeshSource::addProperty ( MeshLib::PropertyVector< T > const &  property) const
private

Adds a zero-copy vtk array wrapper.

Definition at line 231 of file VtkMappedMeshSource.cpp.

233 {
234  vtkNew<vtkAOSDataArrayTemplate<T>> dataArray;
235  const bool hasArrayOwnership = false;
236  dataArray->SetArray(const_cast<T*>(property.data()),
237  static_cast<vtkIdType>(property.size()),
238  static_cast<int>(!hasArrayOwnership));
239  dataArray->SetNumberOfComponents(property.getNumberOfGlobalComponents());
240  dataArray->SetName(property.getPropertyName().c_str());
241 
243  {
244  this->PointData->AddArray(dataArray.GetPointer());
245  }
246  else if (property.getMeshItemType() == MeshLib::MeshItemType::Cell)
247  {
248  this->CellData->AddArray(dataArray.GetPointer());
249  }
250  else if (property.getMeshItemType() ==
252  {
253  this->FieldData->AddArray(dataArray.GetPointer());
254  }
255 }
std::string const & getPropertyName() const
MeshItemType getMeshItemType() const
int getNumberOfGlobalComponents() const
std::size_t size() const
vtkNew< vtkPointData > PointData
vtkNew< vtkFieldData > FieldData
vtkNew< vtkCellData > CellData

References MeshLib::Cell, CellData, FieldData, MeshLib::PropertyVectorBase::getMeshItemType(), MeshLib::PropertyVectorBase::getNumberOfGlobalComponents(), MeshLib::PropertyVectorBase::getPropertyName(), MeshLib::IntegrationPoint, MeshLib::Node, PointData, and MeshLib::PropertyVector< PROP_VAL_TYPE >::size().

Referenced by RequestData().

◆ GetMesh()

const MeshLib::Mesh* MeshLib::VtkMappedMeshSource::GetMesh ( ) const
inline

Returns the mesh.

Definition at line 62 of file VtkMappedMeshSource.h.

62 { return _mesh; }

References _mesh.

Referenced by VtkVisPipeline::addPipelineItem(), MeshItem::getMesh(), and MainWindow::showMeshQualitySelectionDialog().

◆ New()

static VtkMappedMeshSource* MeshLib::VtkMappedMeshSource::New ( )
static

Referenced by MeshItem::MeshItem(), and RequestData().

◆ operator=()

void MeshLib::VtkMappedMeshSource::operator= ( const VtkMappedMeshSource )
delete

◆ PrintSelf()

void MeshLib::VtkMappedMeshSource::PrintSelf ( std::ostream &  os,
vtkIndent  indent 
)
override

◆ ProcessRequest()

int MeshLib::VtkMappedMeshSource::ProcessRequest ( vtkInformation *  request,
vtkInformationVector **  inputVector,
vtkInformationVector *  outputVector 
)
overrideprotected

Definition at line 46 of file VtkMappedMeshSource.cpp.

49 {
50  if (request->Has(vtkDemandDrivenPipeline::REQUEST_DATA()))
51  {
52  return this->RequestData(request, inputVector, outputVector);
53  }
54 
55  if (request->Has(vtkDemandDrivenPipeline::REQUEST_INFORMATION()))
56  {
57  return this->RequestInformation(request, inputVector, outputVector);
58  }
59 
60  return this->Superclass::ProcessRequest(request, inputVector, outputVector);
61 }
int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override
int RequestInformation(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override

References RequestData(), and RequestInformation().

◆ RequestData()

int MeshLib::VtkMappedMeshSource::RequestData ( vtkInformation *  ,
vtkInformationVector **  ,
vtkInformationVector *  outputVector 
)
overrideprotected

Definition at line 63 of file VtkMappedMeshSource.cpp.

66 {
67  vtkSmartPointer<vtkInformation> outInfo =
68  outputVector->GetInformationObject(0);
69  vtkSmartPointer<vtkUnstructuredGrid> output =
70  vtkUnstructuredGrid::SafeDownCast(
71  outInfo->Get(vtkDataObject::DATA_OBJECT()));
72 
73  if (outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_PIECE_NUMBER()) >
74  0)
75  {
76  return 1;
77  }
78 
79  // Points
80  this->Points->Reset();
81 
82  vtkSmartPointer<VtkMeshNodalCoordinatesTemplate<double>> nodeCoords =
83  vtkSmartPointer<VtkMeshNodalCoordinatesTemplate<double>>::New();
84  nodeCoords->SetNodes(_mesh->getNodes());
85  this->Points->SetData(nodeCoords);
86  output->SetPoints(this->Points.GetPointer());
87 
88  // Cells
89  auto elems = _mesh->getElements();
90  output->Allocate(elems.size());
91  for (auto& cell : elems)
92  {
93  auto cellType = OGSToVtkCellType(cell->getCellType());
94 
95  const MeshLib::Element* const elem = cell;
96  const unsigned numNodes(elem->getNumberOfNodes());
97  const MeshLib::Node* const* nodes = cell->getNodes();
98  vtkSmartPointer<vtkIdList> ptIds = vtkSmartPointer<vtkIdList>::New();
99  ptIds->SetNumberOfIds(numNodes);
100 
101  for (unsigned i = 0; i < numNodes; ++i)
102  {
103  ptIds->SetId(i, nodes[i]->getID());
104  }
105 
106  if (cellType == VTK_WEDGE)
107  {
108  for (unsigned i = 0; i < 3; ++i)
109  {
110  const auto prism_swap_id = ptIds->GetId(i);
111  ptIds->SetId(i, ptIds->GetId(i + 3));
112  ptIds->SetId(i + 3, prism_swap_id);
113  }
114  }
115  else if (cellType == VTK_QUADRATIC_WEDGE)
116  {
117  std::array<vtkIdType, 15> ogs_nodeIds{};
118  for (unsigned i = 0; i < 15; ++i)
119  {
120  ogs_nodeIds[i] = ptIds->GetId(i);
121  }
122  for (unsigned i = 0; i < 3; ++i)
123  {
124  ptIds->SetId(i, ogs_nodeIds[i + 3]);
125  ptIds->SetId(i + 3, ogs_nodeIds[i]);
126  }
127  for (unsigned i = 0; i < 3; ++i)
128  {
129  ptIds->SetId(6 + i, ogs_nodeIds[8 - i]);
130  }
131  for (unsigned i = 0; i < 3; ++i)
132  {
133  ptIds->SetId(9 + i, ogs_nodeIds[14 - i]);
134  }
135  ptIds->SetId(12, ogs_nodeIds[9]);
136  ptIds->SetId(13, ogs_nodeIds[11]);
137  ptIds->SetId(14, ogs_nodeIds[10]);
138  }
139 
140  output->InsertNextCell(cellType, ptIds);
141  }
142 
143  // Arrays
144  for (auto [name, property] : _mesh->getProperties())
145  {
146  if (auto p = dynamic_cast<PropertyVector<double>*>(property))
147  {
148  addProperty(*p);
149  }
150  else if (auto p = dynamic_cast<PropertyVector<float>*>(property))
151  {
152  addProperty(*p);
153  }
154  else if (auto p = dynamic_cast<PropertyVector<int>*>(property))
155  {
156  addProperty(*p);
157  }
158  else if (auto p = dynamic_cast<PropertyVector<unsigned>*>(property))
159  {
160  addProperty(*p);
161  }
162  else if (auto p = dynamic_cast<PropertyVector<long>*>(property))
163  {
164  addProperty(*p);
165  }
166  else if (auto p = dynamic_cast<PropertyVector<long long>*>(property))
167  {
168  addProperty(*p);
169  }
170  else if (auto p =
171  dynamic_cast<PropertyVector<unsigned long>*>(property))
172  {
173  addProperty(*p);
174  }
175  else if (auto p = dynamic_cast<PropertyVector<unsigned long long>*>(
176  property))
177  {
178  addProperty(*p);
179  }
180  else if (auto p = dynamic_cast<PropertyVector<std::size_t>*>(property))
181  {
182  addProperty(*p);
183  }
184  else if (auto p = dynamic_cast<PropertyVector<char>*>(property))
185  {
186  addProperty(*p);
187  }
188  else if (auto p =
189  dynamic_cast<PropertyVector<unsigned char>*>(property))
190  {
191  addProperty(*p);
192  }
193  else
194  {
195  OGS_FATAL(
196  "Mesh property '{:s}' of unhandled data type '{:s}'. Please "
197  "check the data type of the mesh properties. The available "
198  "data types are:"
199  "\n\t double,"
200  "\n\t float,"
201  "\n\t int,"
202  "\n\t unsigned,"
203  "\n\t long,"
204  "\n\t long long,"
205  "\n\t unsigned long long,"
206  "\n\t char,",
207  "\n\t unsigned char.",
208  property->getPropertyName(),
209  typeid(*property).name());
210  }
211  }
212 
213  output->GetPointData()->ShallowCopy(this->PointData.GetPointer());
214  output->GetCellData()->ShallowCopy(this->CellData.GetPointer());
215  output->GetFieldData()->ShallowCopy(this->FieldData.GetPointer());
216  return 1;
217 }
#define OGS_FATAL(...)
Definition: Error.h:26
int OGSToVtkCellType(MeshLib::CellType ogs)
Definition: VtkOGSEnum.cpp:17
virtual unsigned getNumberOfNodes() const =0
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition: Mesh.h:95
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:98
Properties & getProperties()
Definition: Mesh.h:123
void addProperty(MeshLib::PropertyVector< T > const &property) const
Adds a zero-copy vtk array wrapper.
static VtkMappedMeshSource * New()

References _mesh, addProperty(), CellData, FieldData, MeshLib::Mesh::getElements(), MeshLib::Mesh::getNodes(), MeshLib::Element::getNumberOfNodes(), MeshLib::Mesh::getProperties(), MaterialPropertyLib::name, New(), OGS_FATAL, OGSToVtkCellType(), PointData, and Points.

Referenced by ProcessRequest().

◆ RequestInformation()

int MeshLib::VtkMappedMeshSource::RequestInformation ( vtkInformation *  ,
vtkInformationVector **  ,
vtkInformationVector *   
)
overrideprotected

Definition at line 219 of file VtkMappedMeshSource.cpp.

223 {
224  this->NumberOfDimensions = 3;
226 
227  return 1;
228 }
std::size_t getNumberOfNodes() const
Get the number of nodes.
Definition: Mesh.h:89

References _mesh, MeshLib::Mesh::getNumberOfNodes(), NumberOfDimensions, and NumberOfNodes.

Referenced by ProcessRequest().

◆ SetMesh()

void MeshLib::VtkMappedMeshSource::SetMesh ( const MeshLib::Mesh mesh)
inline

Sets the mesh. Calling is mandatory.

Definition at line 55 of file VtkMappedMeshSource.h.

56  {
57  this->_mesh = mesh;
58  this->Modified();
59  }

References _mesh.

◆ vtkTypeMacro()

MeshLib::VtkMappedMeshSource::vtkTypeMacro ( VtkMappedMeshSource  ,
vtkUnstructuredGridAlgorithm   
)

Member Data Documentation

◆ _mesh

const MeshLib::Mesh* MeshLib::VtkMappedMeshSource::_mesh {}
private

◆ CellData

vtkNew<vtkCellData> MeshLib::VtkMappedMeshSource::CellData
private

Definition at line 89 of file VtkMappedMeshSource.h.

Referenced by addProperty(), and RequestData().

◆ FieldData

vtkNew<vtkFieldData> MeshLib::VtkMappedMeshSource::FieldData
private

Definition at line 90 of file VtkMappedMeshSource.h.

Referenced by addProperty(), and RequestData().

◆ NumberOfDimensions

int MeshLib::VtkMappedMeshSource::NumberOfDimensions {0}
private

Definition at line 84 of file VtkMappedMeshSource.h.

Referenced by RequestInformation().

◆ NumberOfNodes

int MeshLib::VtkMappedMeshSource::NumberOfNodes {0}
private

Definition at line 85 of file VtkMappedMeshSource.h.

Referenced by RequestInformation().

◆ PointData

vtkNew<vtkPointData> MeshLib::VtkMappedMeshSource::PointData
private

Definition at line 88 of file VtkMappedMeshSource.h.

Referenced by addProperty(), and RequestData().

◆ Points

vtkNew<vtkPoints> MeshLib::VtkMappedMeshSource::Points
private

Definition at line 87 of file VtkMappedMeshSource.h.

Referenced by RequestData().


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