27 #include "vtkPVFoamReader.h" 32 #include "vtkOpenFOAMPoints.H" 36 #include "vtkCellArray.h" 37 #include "vtkIdTypeArray.h" 38 #include "vtkUnstructuredGrid.h" 42 vtkUnstructuredGrid* Foam::vtkPVFoam::volumeVTKMesh
45 polyDecomp& decompInfo
59 Info<<
"<beg> Foam::vtkPVFoam::volumeVTKMesh" <<
endl;
72 const labelList& owner = mesh.faceOwner();
74 labelList& superCells = decompInfo.superCells();
75 labelList& addPointCellLabels = decompInfo.addPointCellLabels();
79 if (!reader_->GetUseVTKPolyhedron())
83 Info<<
"... scanning for polyhedra" <<
endl;
88 const cellModel& model = cellShapes[celli].model();
100 const cell& cFaces = mesh.cells()[celli];
104 const face& f = mesh.faces()[cFaces[cFacei]];
108 f.nTrianglesQuads(mesh.points(), nTris, nQuads);
110 nAddCells += nQuads + nTris;
121 addPointCellLabels.setSize(nAddPoints);
128 Info<<
" mesh nCells = " << mesh.nCells() <<
nl 129 <<
" nPoints = " << mesh.nPoints() <<
nl 130 <<
" nAddCells = " << nAddCells <<
nl 131 <<
" nAddPoints = " << nAddPoints <<
endl;
134 superCells.setSize(mesh.nCells() + nAddCells);
138 Info<<
"... converting points" <<
endl;
143 vtkpoints->Allocate(mesh.nPoints() + nAddPoints);
155 Info<<
"... converting cells" <<
endl;
158 vtkmesh->Allocate(mesh.nCells() + nAddCells);
161 label addPointi = 0, addCelli = 0;
165 vtkIdType nodeIds[8];
169 DynamicList<vtkIdType> faceStream(256);
173 const cellShape& cellShape = cellShapes[celli];
174 const cellModel& cellModel = cellShape.model();
176 superCells[addCelli++] = celli;
178 if (cellModel == tet)
180 for (
int j = 0; j < 4; j++)
182 nodeIds[j] = cellShape[j];
184 vtkmesh->InsertNextCell
191 else if (cellModel == pyr)
193 for (
int j = 0; j < 5; j++)
195 nodeIds[j] = cellShape[j];
197 vtkmesh->InsertNextCell
204 else if (cellModel == prism)
208 nodeIds[0] = cellShape[0];
209 nodeIds[1] = cellShape[2];
210 nodeIds[2] = cellShape[1];
211 nodeIds[3] = cellShape[3];
212 nodeIds[4] = cellShape[5];
213 nodeIds[5] = cellShape[4];
215 vtkmesh->InsertNextCell
222 else if (cellModel == tetWedge && !reader_->GetUseVTKPolyhedron())
226 nodeIds[0] = cellShape[0];
227 nodeIds[1] = cellShape[2];
228 nodeIds[2] = cellShape[1];
229 nodeIds[3] = cellShape[3];
230 nodeIds[4] = cellShape[4];
231 nodeIds[5] = cellShape[3];
233 vtkmesh->InsertNextCell
240 else if (cellModel == wedge)
244 nodeIds[0] = cellShape[0];
245 nodeIds[1] = cellShape[1];
246 nodeIds[2] = cellShape[2];
247 nodeIds[3] = cellShape[2];
248 nodeIds[4] = cellShape[3];
249 nodeIds[5] = cellShape[4];
250 nodeIds[6] = cellShape[5];
251 nodeIds[7] = cellShape[6];
253 vtkmesh->InsertNextCell
260 else if (cellModel == hex)
262 for (
int j = 0; j < 8; j++)
264 nodeIds[j] = cellShape[j];
266 vtkmesh->InsertNextCell
273 else if (reader_->GetUseVTKPolyhedron())
276 const labelList& cFaces = mesh.cells()[celli];
278 #ifdef HAS_VTK_POLYHEDRON 279 vtkIdType nFaces = cFaces.
size();
280 vtkIdType nLabels = nFaces;
285 const face& f = mesh.faces()[cFaces[cFacei]];
293 faceStream.reserve(nLabels + nFaces);
297 const face& f = mesh.faces()[cFaces[cFacei]];
298 const bool isOwner = (owner[cFaces[cFacei]] == celli);
299 const label nFacePoints = f.size();
302 faceStream.append(nFacePoints);
308 faceStream.append(f[fp]);
317 faceStream.append(f[fp]);
322 vtkmesh->InsertNextCell(VTK_POLYHEDRON, nFaces, faceStream.data());
328 HashSet<vtkIdType, Hash<label>> hashUniqId(2*256);
332 const face& f = mesh.faces()[cFaces[cFacei]];
336 hashUniqId.insert(f[fp]);
341 faceStream = hashUniqId.sortedToc();
343 vtkmesh->InsertNextCell
345 VTK_CONVEX_POINT_SET,
346 vtkIdType(faceStream.size()),
356 addPointCellLabels[addPointi] = celli;
359 const label newVertexLabel = mesh.nPoints() + addPointi;
363 bool substituteCell =
true;
365 const labelList& cFaces = mesh.cells()[celli];
368 const face& f = mesh.faces()[cFaces[cFacei]];
369 const bool isOwner = (owner[cFaces[cFacei]] == celli);
374 f.nTrianglesQuads(mesh.points(), nTris, nQuads);
381 f.trianglesQuads(mesh.points(), trii, quadi, triFcs, quadFcs);
387 substituteCell =
false;
391 superCells[addCelli++] = celli;
394 const face& quad = quadFcs[quadI];
405 nodeIds[0] = quad[3];
406 nodeIds[1] = quad[2];
407 nodeIds[2] = quad[1];
408 nodeIds[3] = quad[0];
412 nodeIds[0] = quad[0];
413 nodeIds[1] = quad[1];
414 nodeIds[2] = quad[2];
415 nodeIds[3] = quad[3];
417 nodeIds[4] = newVertexLabel;
418 vtkmesh->InsertNextCell
430 substituteCell =
false;
434 superCells[addCelli++] = celli;
437 const face& tri = triFcs[triI];
452 nodeIds[3] = newVertexLabel;
454 vtkmesh->InsertNextCell
467 vtkmesh->SetPoints(vtkpoints);
472 Info<<
"<end> Foam::vtkPVFoam::volumeVTKMesh" <<
endl;
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
static void printMemory()
Simple memory used debugging information.
static const cellModel * lookup(const word &)
Look up a model by name and return a pointer to the model or NULL.
void size(const label)
Override size to be inconsistent with allocated storage.
Ostream & endl(Ostream &os)
Add newline and flush stream.
void vtkInsertNextOpenFOAMPoint(vtkPoints *points, const Foam::point &p)
#define forAllReverse(list, i)
Reverse loop across all elements in list.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
List< cellShape > cellShapeList
List of cellShapes and PtrList of List of cellShape.
List< label > labelList
A List of labels.