foamVtkVtuAdaptorI.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2016-2020 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
29
30inline vtkSmartPointer<vtkPoints>
32(
33 const fvMesh& mesh
34) const
35{
36 // Convert OpenFOAM mesh vertices to VTK
37 auto vtkpoints = vtkSmartPointer<vtkPoints>::New();
38
39 // Normal points
40 const pointField& pts = mesh.points();
41
42 // Additional cell centres
43 const labelUList& addPoints = this->additionalIds();
44
45 vtkpoints->SetNumberOfPoints(pts.size() + addPoints.size());
46
47 // Normal points
48 vtkIdType pointId = 0;
49 for (const point& p : pts)
50 {
51 vtkpoints->SetPoint(pointId++, p.v_);
52 }
53
54 // Cell centres
55 for (const label meshCelli : addPoints)
56 {
57 vtkpoints->SetPoint(pointId++, mesh.cellCentres()[meshCelli].v_);
58 }
59
60 return vtkpoints;
61}
62
63
64inline vtkSmartPointer<vtkPoints>
66(
67 const fvMesh& mesh,
68 const labelUList& pointMap
69) const
70{
71 // Convert OpenFOAM mesh vertices to VTK
72 auto vtkpoints = vtkSmartPointer<vtkPoints>::New();
73
74 // Normal points
75 const pointField& pts = mesh.points();
76
77 // Additional cell centres
78 const labelUList& addPoints = this->additionalIds();
79
80 vtkpoints->SetNumberOfPoints(pointMap.size() + addPoints.size());
81
82 // Normal points
83 vtkIdType pointId = 0;
84 for (const label meshPointi : pointMap)
85 {
86 vtkpoints->SetPoint(pointId++, pts[meshPointi].v_);
87 }
88
89 // Cell centres
90 for (const label meshCelli : addPoints)
91 {
92 vtkpoints->SetPoint(pointId++, mesh.cellCentres()[meshCelli].v_);
93 }
94
95 return vtkpoints;
96}
97
98
99inline vtkSmartPointer<vtkUnstructuredGrid>
101(
102 const fvMesh& mesh,
103 const bool decompPoly
104)
105{
107 (
108 #ifdef VTK_CELL_ARRAY_V2
110 #else
112 #endif
113 );
114
115 vtk::vtuSizing sizing(mesh, decompPoly);
116
117 auto vtkmesh = vtkSmartPointer<vtkUnstructuredGrid>::New();
118
119 auto cellTypes = vtkSmartPointer<vtkUnsignedCharArray>::New();
120
121 UList<uint8_t> cellTypesUL
122 (
124 );
125
126 auto cells = vtkSmartPointer<vtkCellArray>::New();
127 auto faces = vtkSmartPointer<vtkIdTypeArray>::New();
128
129 auto cellOffsets = vtkSmartPointer<vtkIdTypeArray>::New();
130 auto faceLocations = vtkSmartPointer<vtkIdTypeArray>::New();
131
132 const auto nConnect
133 (
135 );
136
137 UList<vtkIdType> cellOffsetsUL
138 (
140 (
141 cellOffsets,
143 )
144 );
145
146 #ifdef VTK_CELL_ARRAY_V2
147
148 auto cellConnect = vtkSmartPointer<vtkIdTypeArray>::New();
149
150 UList<vtkIdType> cellsUL
151 (
152 vtk::Tools::asUList(cellConnect, nConnect)
153 );
154
155 #else
156
157 cells->GetData()->SetNumberOfTuples(sizing.nFieldCells());
158
159 UList<vtkIdType> cellsUL
160 (
161 cells->WritePointer(sizing.nFieldCells(), nConnect),
162 nConnect
163 );
164
165 #endif
166
167 UList<vtkIdType> facesUL
168 (
170 (
171 faces,
173 )
174 );
175
176 UList<vtkIdType> faceLocationsUL
177 (
179 (
180 faceLocations,
182 )
183 );
184
185
186 sizing.populateInternal
187 (
188 mesh,
189 cellTypesUL,
190 cellsUL, cellOffsetsUL,
191 facesUL, faceLocationsUL,
192 static_cast<foamVtkMeshMaps&>(*this),
193 output
194 );
195
196
197 // Convert OpenFOAM mesh vertices to VTK
198 // - can only do this *after* populating the decompInfo with cell-ids
199 // for any additional points (ie, mesh cell-centres)
200 vtkmesh->SetPoints(this->points(mesh));
201
202 #ifdef VTK_CELL_ARRAY_V2
203
204 // Move into a vtkCellArray
205 cells->SetData(cellOffsets, cellConnect);
206
207 if (facesUL.size())
208 {
209 vtkmesh->SetCells(cellTypes, cells, faceLocations, faces);
210 }
211 else
212 {
213 vtkmesh->SetCells(cellTypes, cells, nullptr, nullptr);
214 }
215 #else
216
217 if (facesUL.size())
218 {
219 vtkmesh->SetCells(cellTypes, cellOffsets, cells, faceLocations, faces);
220 }
221 else
222 {
223 vtkmesh->SetCells(cellTypes, cellOffsets, cells, nullptr, nullptr);
224 }
225
226 #endif
227
228 return vtkmesh;
229}
230
231
232// ************************************************************************* //
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
Bookkeeping for mesh subsetting and/or polyhedral cell decomposition. Although the main use case is f...
const labelList & additionalIds() const noexcept
Any additional (user) labels.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1083
const vectorField & cellCentres() const
Sizing descriptions and routines for transcribing an OpenFOAM volume mesh into a VTK unstructured gri...
label sizeOf(const enum contentType output, const enum slotType slot) const
Return the required size for the storage slot.
@ FACES_OFFSETS
Faces end-offsets (XML) or locations (INTERNAL1)
@ CELLS
Cell connectivity (ALL)
@ FACES
Face-stream (XML, INTERNAL)
contentType
Types of content that the storage may represent.
@ INTERNAL2
Internal vtkUnstructuredGrid content, VTK_CELL_ARRAY_V2.
@ INTERNAL1
Internal vtkUnstructuredGrid content.
void populateInternal(const polyMesh &mesh, UList< uint8_t > &cellTypes, UList< int > &connectivity, UList< int > &offsets, UList< int > &faces, UList< int > &facesOffsets, foamVtkMeshMaps &maps, const enum contentType output) const
Populate lists for Internal VTK format.
label nFieldCells() const noexcept
Number of field cells = nCells + nAddCells.
volScalarField & p
dynamicFvMesh & mesh
const pointField & points
const cellShapeList & cells
const wordList internal
Standard dimensioned field types (scalar, vector, tensor, etc)
UList< uint8_t > asUList(vtkUnsignedCharArray *array, const label size)
Wrap vtkUnsignedCharArray as a UList.
Definition: foamVtkToolsI.H:33
static Ostream & output(Ostream &os, const IntRange< T > &range)
Definition: IntRanges.C:66
const labelList & cellTypes
Definition: setCellMask.H:33