ensightOutputVolFieldTemplates.C
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-2021 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
29#include "ensightMesh.H"
30
31#include "fvMesh.H"
32#include "linear.H"
34#include "interpolation.H"
35#include "processorFvPatch.H"
36#include "DynamicField.H"
37
38// * * * * * * * * * * * * * * * * Detail * * * * * * * * * * * * * * * * * //
39
40template<class Type>
42(
45 const ensightMesh& ensMesh
46)
47{
48 bool parallel = Pstream::parRun();
49
50 const fvMesh& mesh = vf.mesh();
51 const polyBoundaryMesh& bmesh = mesh.boundaryMesh();
52
53 const Map<ensightCells>& cellZoneParts = ensMesh.cellZoneParts();
54 const Map<ensightFaces>& faceZoneParts = ensMesh.faceZoneParts();
55 const Map<ensightFaces>& boundaryParts = ensMesh.boundaryParts();
56
57
58 // Write internalMesh and cellZones - sorted by index
59 for (const label zoneId : cellZoneParts.sortedToc())
60 {
61 const ensightCells& part = cellZoneParts[zoneId];
62
63 ensightOutput::writeField(os, vf, part, parallel);
64 }
65
66
67 // Write patches - sorted by index
68 for (const label patchId : boundaryParts.sortedToc())
69 {
70 const ensightFaces& part = boundaryParts[patchId];
71
72 if (patchId < 0 || patchId >= bmesh.size())
73 {
74 // Future handling of combined patches?
75 continue;
76 }
77
78 const label patchStart = bmesh[patchId].start();
79
80 // Either use a flat boundary field for all patches,
81 // or patch-local face ids
82
83 // Operate on a copy
84 ensightFaces localPart(part);
85
86 // Change from global faceIds to patch-local
87 localPart.decrFaceIds(patchStart);
88
89 ensightOutput::writeField
90 (
91 os,
93 localPart,
94 parallel
95 );
96 }
97
98
99 // No face zones data
100 if (faceZoneParts.empty())
101 {
102 return true;
103 }
104
105
106 // Flat boundary field
107 // similar to volPointInterpolation::flatBoundaryField()
108
109 Field<Type> flat(mesh.nBoundaryFaces(), Zero);
110
111 const fvBoundaryMesh& bm = mesh.boundary();
112 forAll(vf.boundaryField(), patchi)
113 {
114 const polyPatch& pp = bm[patchi].patch();
115 const auto& bf = vf.boundaryField()[patchi];
116
117 if (isA<processorFvPatch>(bm[patchi]))
118 {
119 // Use average value for processor faces
120 // own cell value = patchInternalField
121 // nei cell value = evaluated boundary values
123 (
124 flat,
125 bf.size(),
126 pp.offset()
127 ) = (0.5 * (bf.patchInternalField() + bf));
128 }
129 else if (!isA<emptyFvPatch>(bm[patchi]))
130 {
132 (
133 flat,
134 bf.size(),
135 pp.offset()
136 ) = bf;
137 }
138 }
139
140
141 // Interpolate cell values to faces
142 // - only needed when exporting faceZones...
143
146
147 const auto& sfld = tsfld();
148
149
150 // Local output buffer
151 label maxLen = 0;
152
153 forAllConstIters(faceZoneParts, iter)
154 {
155 maxLen = max(maxLen, iter.val().size());
156 }
157
158 DynamicField<Type> values(maxLen);
159
160
161 //
162 // Write requested faceZones - sorted by index
163 //
164 for (const label zoneId : faceZoneParts.sortedToc())
165 {
166 const ensightFaces& part = faceZoneParts[zoneId];
167
168 // Loop over face ids to store the needed field values
169 // - internal faces use linear interpolation
170 // - boundary faces use the corresponding patch value
171
172 // Local copy of the field
173 values.resize(part.size());
174 values = Zero;
175
176 auto valIter = values.begin();
177
178 for (const label faceId : part.faceIds())
179 {
180 *valIter =
181 (
182 mesh.isInternalFace(faceId)
183 ? sfld[faceId]
184 : flat[faceId - mesh.nInternalFaces()]
185 );
186
187 ++valIter;
188 }
189
190 // The field is already in the proper element order
191 // - just need its corresponding sub-fields
192 ensightOutput::Detail::writeFaceSubField(os, values, part, parallel);
193 }
194
195 return true;
196}
197
198
199template<class Type>
201(
204 const ensightMesh& ensMesh
205)
206{
207 bool parallel = Pstream::parRun();
208
209 const polyMesh& mesh = ensMesh.mesh();
210
211 const Map<ensightCells>& cellZoneParts = ensMesh.cellZoneParts();
212 const Map<ensightFaces>& faceZoneParts = ensMesh.faceZoneParts();
213 const Map<ensightFaces>& boundaryParts = ensMesh.boundaryParts();
214
215 //
216 // Write internalMesh and cellZones - sorted by index
217 //
218 for (const label zoneId : cellZoneParts.sortedToc())
219 {
220 const ensightCells& part = cellZoneParts[zoneId];
221
222 if (Pstream::master())
223 {
224 os.beginPart(part.index());
225 }
226
227 labelList uniquePointLabels;
228 part.uniqueMeshPoints(mesh, uniquePointLabels, parallel);
229
230 ensightOutput::Detail::writeFieldComponents
231 (
232 os,
233 ensightFile::coordinates,
234 UIndirectList<Type>(pf.internalField(), uniquePointLabels),
235 parallel
236 );
237 }
238
239
240 //
241 // Write patches - sorted by index
242 //
243 for (const label patchId : boundaryParts.sortedToc())
244 {
245 const ensightFaces& part = boundaryParts[patchId];
246
247 if (Pstream::master())
248 {
249 os.beginPart(part.index());
250 }
251
252 labelList uniquePointLabels;
253 part.uniqueMeshPoints(mesh, uniquePointLabels, parallel);
254
255 const auto& bfld = pf.boundaryField()[patchId];
256
257 // Only valuePointPatchField is actually derived from Field
258 const auto* vpp = isA<Field<Type>>(bfld);
259 if (vpp)
260 {
261 // Require patch local indices, not mesh point labels.
262 // But need to use polyPatch meshPointMap() to recover the
263 // local indices since the ensight output will have jumbled
264 // the face output order
265
266 const polyPatch& pp = mesh.boundaryMesh()[patchId];
267
268 for (label& pointi : uniquePointLabels)
269 {
270 pointi = pp.meshPointMap()[pointi];
271 }
272
273 ensightOutput::Detail::writeFieldComponents
274 (
275 os,
276 ensightFile::coordinates,
277 UIndirectList<Type>(*vpp, uniquePointLabels),
278 parallel
279 );
280 }
281 else
282 {
283 ensightOutput::Detail::writeFieldComponents
284 (
285 os,
286 ensightFile::coordinates,
287 UIndirectList<Type>(pf.internalField(), uniquePointLabels),
288 parallel
289 );
290 }
291 }
292
293 //
294 // Write requested faceZones - sorted by index
295 //
296 for (const label zoneId : faceZoneParts.sortedToc())
297 {
298 const ensightFaces& part = faceZoneParts[zoneId];
299
300 if (Pstream::master())
301 {
302 os.beginPart(part.index());
303 }
304
305 // CAVEAT - does not properly handle valuePointPatchField,
306 // uses internalField only
307
308 {
309 labelList uniquePointLabels;
310 part.uniqueMeshPoints(mesh, uniquePointLabels, parallel);
311
312 ensightOutput::Detail::writeFieldComponents
313 (
314 os,
315 ensightFile::coordinates,
316 UIndirectList<Type>(pf.internalField(), uniquePointLabels),
317 parallel
318 );
319 }
320 }
321
322 return true;
323}
324
325
326// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
327
328template<class Type>
330(
333 const ensightMesh& ensMesh,
334 const bool nodeValues
335)
336{
337 if (nodeValues)
338 {
340 (
341 volPointInterpolation::New(vf.mesh()).interpolate(vf)
342 );
343 pfld.ref().checkOut();
344 pfld.ref().rename(vf.name());
345
346 return ensightOutput::writePointField<Type>(os, pfld, ensMesh);
347 }
348
349 return ensightOutput::writeVolField<Type>(os, vf, ensMesh);
350}
351
352
353// ************************************************************************* //
Y[inertIndex] max(0.0)
const Mesh & mesh() const
Return mesh.
Dynamically sized Field.
Definition: DynamicField.H:65
Generic templated field type.
Definition: Field.H:82
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const Internal & internalField() const
Return a const-reference to the dimensioned internal field.
List< Key > sortedToc() const
The table of contents (the keys) in sorted order.
Definition: HashTable.C:137
bool empty() const noexcept
True if the hash table is empty.
Definition: HashTableI.H:59
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:65
A HashTable to objects of type <T> with a label key.
Definition: Map.H:60
const Map< label > & meshPointMap() const
Mesh point map.
A List obtained as a section of another List.
Definition: SubList.H:70
A List with indirect addressing. Like IndirectList but does not store addressing.
Definition: IndirectList.H:79
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
label size() const noexcept
The number of elements in the list.
Definition: UPtrListI.H:106
Sorting/classification of cells (3D) into corresponding ensight element types.
Definition: ensightCells.H:58
label uniqueMeshPoints(const polyMesh &mesh, labelList &uniqueMeshPointLabels, bool parallel) const
Globally unique mesh points. Required when writing point fields.
Sorting/classification of faces (2D) into corresponding ensight types.
Definition: ensightFaces.H:75
void decrFaceIds(const label off)
Decrease face ids by specified offset value.
label uniqueMeshPoints(const polyMesh &mesh, labelList &uniqueMeshPointLabels, bool parallel) const
const labelList & faceIds() const noexcept
Processor-local face ids of all elements.
Definition: ensightFacesI.H:79
label size(const elemType etype) const
Processor-local size of the specified element type.
Definition: ensightFacesI.H:67
Ensight output with specialized write() for strings, integers and floats. Correctly handles binary wr...
Definition: ensightFile.H:55
Encapsulation of volume meshes for writing in ensight format. It manages cellZones,...
Definition: ensightMesh.H:83
const Map< ensightFaces > & boundaryParts() const noexcept
Face elements per selected patch, lookup by patch index.
Definition: ensightMesh.H:184
const Map< ensightCells > & cellZoneParts() const noexcept
Face elements per selected patch, lookup by patch index.
Definition: ensightMesh.H:170
const polyMesh & mesh() const noexcept
Reference to the underlying polyMesh.
Definition: ensightMesh.H:159
const Map< ensightFaces > & faceZoneParts() const noexcept
Face elements per faceZone, lookup by zone index.
Definition: ensightMesh.H:177
label index() const noexcept
The index in a list (0-based)
Definition: ensightPart.H:124
Foam::fvBoundaryMesh.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
label start() const noexcept
The start label of boundary faces in the polyMesh face list.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:75
label offset() const
The offset where this patch starts in the boundary face list.
Definition: polyPatch.C:309
A class for managing temporary objects.
Definition: tmp.H:65
T & ref() const
Definition: tmpI.H:227
dynamicFvMesh & mesh
A collection of functions for writing volField content in ensight format.
OBJstream os(runTime.globalPath()/outputName)
label patchId(-1)
label faceId(-1)
bool writeVolField(ensightFile &os, const GeometricField< Type, fvPatchField, volMesh > &vf, const ensightMesh &ensMesh)
Write volume field component-wise.
bool writePointField(ensightFile &os, const GeometricField< Type, pointPatchField, pointMesh > &pf, const ensightMesh &ensMesh)
Write point field component-wise.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > linearInterpolate(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: linear.H:112
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
Definition: stdFoam.H:278