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 -------------------------------------------------------------------------------
10 License
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 
30 inline 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 
64 inline 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 
99 inline vtkSmartPointer<vtkUnstructuredGrid>
101 (
102  const fvMesh& mesh,
103  const bool decompPoly
104 )
105 {
107  (
108  #ifdef VTK_CELL_ARRAY_V2
109  vtk::vtuSizing::contentType::INTERNAL2
110  #else
111  vtk::vtuSizing::contentType::INTERNAL1
112  #endif
113  );
114 
115  vtk::vtuSizing sizing(mesh, decompPoly);
116 
118 
120 
121  UList<uint8_t> cellTypesUL
122  (
124  );
125 
128 
129  auto cellOffsets = vtkSmartPointer<vtkIdTypeArray>::New();
130  auto faceLocations = vtkSmartPointer<vtkIdTypeArray>::New();
131 
132  const auto nConnect
133  (
134  sizing.sizeOf(output, vtk::vtuSizing::slotType::CELLS)
135  );
136 
137  UList<vtkIdType> cellOffsetsUL
138  (
140  (
141  cellOffsets,
142  sizing.sizeOf(output, vtk::vtuSizing::slotType::CELLS_OFFSETS)
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,
172  sizing.sizeOf(output, vtk::vtuSizing::slotType::FACES)
173  )
174  );
175 
176  UList<vtkIdType> faceLocationsUL
177  (
179  (
180  faceLocations,
181  sizing.sizeOf(output, vtk::vtuSizing::slotType::FACES_OFFSETS)
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 // ************************************************************************* //
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::output
static Ostream & output(Ostream &os, const IntRange< T > &range)
Definition: IntRanges.C:66
Foam::vtk::vtuSizing::contentType
contentType
Types of content that the storage may represent.
Definition: foamVtuSizing.H:208
Foam::vtk::Tools::asUList
UList< uint8_t > asUList(vtkUnsignedCharArray *array, const label size)
Wrap vtkUnsignedCharArray as a UList.
Definition: foamVtkToolsI.H:33
Foam::vtk::vtuAdaptor::points
vtkSmartPointer< vtkPoints > points(const fvMesh &mesh) const
The vtk points for the mesh (and decomposition)
Definition: foamVtkVtuAdaptorI.H:32
Foam::Field< vector >
Foam::vtk::vtuSizing::nFieldCells
label nFieldCells() const noexcept
Number of field cells = nCells + nAddCells.
Definition: foamVtuSizingI.H:99
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:85
Foam::vtk::vtuSizing
Sizing descriptions and routines for transcribing an OpenFOAM volume mesh into a VTK unstructured gri...
Definition: foamVtuSizing.H:201
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
Foam::vtk::vtuSizing::populateInternal
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.
Definition: foamVtuSizing.C:757
cellTypes
const labelList & cellTypes
Definition: setCellMask.H:33
Foam::Vector< scalar >
Foam::UList< label >
points
const pointField & points
Definition: gmvOutputHeader.H:1
cells
const cellShapeList & cells
Definition: gmvOutputHeader.H:3
Foam::UList::size
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
Foam::foamVtkMeshMaps
Bookkeeping for mesh subsetting and/or polyhedral cell decomposition. Although the main use case is f...
Definition: foamVtkMeshMaps.H:58
Foam::vtk::vtuSizing::sizeOf
label sizeOf(const enum contentType output, const enum slotType slot) const
Return the required size for the storage slot.
Definition: foamVtuSizing.C:471
Foam::vtk::vtuAdaptor::internal
vtkSmartPointer< vtkUnstructuredGrid > internal(const fvMesh &mesh, const bool decompPoly=false)
Internal mesh as vtkUnstructuredGrid.
Definition: foamVtkVtuAdaptorI.H:101