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
72 const labelList& owner = mesh.faceOwner();
74 labelList& superCells = decompInfo.superCells();
75 labelList& addPointCellLabels = decompInfo.addPointCellLabels();
79 if (!reader_->GetUseVTKPolyhedron())
84 <<
"... scanning for polyhedra" <<
endl;
89 const cellModel& model = cellShapes[celli].model();
101 const cell& cFaces = mesh.cells()[celli];
105 const face& f = mesh.faces()[cFaces[cFacei]];
109 f.nTrianglesQuads(mesh.points(), nTris, nQuads);
111 nAddCells += nQuads + nTris;
122 addPointCellLabels.setSize(nAddPoints);
129 Info<<
" mesh nCells = " << mesh.nCells() <<
nl 130 <<
" nPoints = " << mesh.nPoints() <<
nl 131 <<
" nAddCells = " << nAddCells <<
nl 132 <<
" nAddPoints = " << nAddPoints <<
endl;
135 superCells.setSize(mesh.nCells() + nAddCells);
139 Info<<
" ... converting points" <<
endl;
144 vtkpoints->Allocate(mesh.nPoints() + nAddPoints);
156 Info<<
" ... converting cells" <<
endl;
159 vtkmesh->Allocate(mesh.nCells() + nAddCells);
162 label addPointi = 0, addCelli = 0;
166 vtkIdType nodeIds[8];
170 DynamicList<vtkIdType> faceStream(256);
174 const cellShape& cellShape = cellShapes[celli];
175 const cellModel& cellModel = cellShape.model();
177 superCells[addCelli++] = celli;
179 if (cellModel == tet)
181 for (
int j = 0; j < 4; j++)
183 nodeIds[j] = cellShape[j];
185 vtkmesh->InsertNextCell
192 else if (cellModel == pyr)
194 for (
int j = 0; j < 5; j++)
196 nodeIds[j] = cellShape[j];
198 vtkmesh->InsertNextCell
205 else if (cellModel == prism)
209 nodeIds[0] = cellShape[0];
210 nodeIds[1] = cellShape[2];
211 nodeIds[2] = cellShape[1];
212 nodeIds[3] = cellShape[3];
213 nodeIds[4] = cellShape[5];
214 nodeIds[5] = cellShape[4];
216 vtkmesh->InsertNextCell
223 else if (cellModel == tetWedge && !reader_->GetUseVTKPolyhedron())
227 nodeIds[0] = cellShape[0];
228 nodeIds[1] = cellShape[2];
229 nodeIds[2] = cellShape[1];
230 nodeIds[3] = cellShape[3];
231 nodeIds[4] = cellShape[4];
232 nodeIds[5] = cellShape[3];
234 vtkmesh->InsertNextCell
241 else if (cellModel == wedge)
245 nodeIds[0] = cellShape[0];
246 nodeIds[1] = cellShape[1];
247 nodeIds[2] = cellShape[2];
248 nodeIds[3] = cellShape[2];
249 nodeIds[4] = cellShape[3];
250 nodeIds[5] = cellShape[4];
251 nodeIds[6] = cellShape[5];
252 nodeIds[7] = cellShape[6];
254 vtkmesh->InsertNextCell
261 else if (cellModel == hex)
263 for (
int j = 0; j < 8; j++)
265 nodeIds[j] = cellShape[j];
267 vtkmesh->InsertNextCell
274 else if (reader_->GetUseVTKPolyhedron())
277 const labelList& cFaces = mesh.cells()[celli];
279 #ifdef HAS_VTK_POLYHEDRON 280 vtkIdType nFaces = cFaces.
size();
281 vtkIdType nLabels = nFaces;
286 const face& f = mesh.faces()[cFaces[cFacei]];
294 faceStream.reserve(nLabels + nFaces);
298 const face& f = mesh.faces()[cFaces[cFacei]];
299 const bool isOwner = (owner[cFaces[cFacei]] == celli);
300 const label nFacePoints = f.size();
303 faceStream.append(nFacePoints);
309 faceStream.append(f[fp]);
318 faceStream.append(f[fp]);
323 vtkmesh->InsertNextCell(VTK_POLYHEDRON, nFaces, faceStream.data());
329 HashSet<vtkIdType, Hash<label>> hashUniqId(2*256);
333 const face& f = mesh.faces()[cFaces[cFacei]];
337 hashUniqId.insert(f[fp]);
342 faceStream = hashUniqId.sortedToc();
344 vtkmesh->InsertNextCell
346 VTK_CONVEX_POINT_SET,
347 vtkIdType(faceStream.size()),
357 addPointCellLabels[addPointi] = celli;
360 const label newVertexLabel = mesh.nPoints() + addPointi;
364 bool substituteCell =
true;
366 const labelList& cFaces = mesh.cells()[celli];
369 const face& f = mesh.faces()[cFaces[cFacei]];
370 const bool isOwner = (owner[cFaces[cFacei]] == celli);
375 f.nTrianglesQuads(mesh.points(), nTris, nQuads);
382 f.trianglesQuads(mesh.points(), trii, quadi, triFcs, quadFcs);
388 substituteCell =
false;
392 superCells[addCelli++] = celli;
395 const face& quad = quadFcs[quadI];
406 nodeIds[0] = quad[3];
407 nodeIds[1] = quad[2];
408 nodeIds[2] = quad[1];
409 nodeIds[3] = quad[0];
413 nodeIds[0] = quad[0];
414 nodeIds[1] = quad[1];
415 nodeIds[2] = quad[2];
416 nodeIds[3] = quad[3];
418 nodeIds[4] = newVertexLabel;
419 vtkmesh->InsertNextCell
431 substituteCell =
false;
435 superCells[addCelli++] = celli;
438 const face& tri = triFcs[triI];
453 nodeIds[3] = newVertexLabel;
455 vtkmesh->InsertNextCell
468 vtkmesh->SetPoints(vtkpoints);
#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 const cellModel * lookup(const word &)
Look up a model by name and return a pointer to the model or nullptr.
static void printMemory()
Simple memory used debugging information.
void vtkInsertNextOpenFOAMPoint(vtkPoints *points, const Foam::point &p)
void size(const label)
Override size to be inconsistent with allocated storage.
Ostream & endl(Ostream &os)
Add newline and flush stream.
#define forAllReverse(list, i)
Reverse loop across all elements in list.
List< cellShape > cellShapeList
List of cellShapes and PtrList of List of cellShape.
Pre-declare SubField and related Field type.
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
List< label > labelList
A List of labels.
#define InfoInFunction
Report an information message using Foam::Info.