ensightField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "ensightField.H"
27 #include "fvMesh.H"
28 #include "volFields.H"
29 #include "OFstream.H"
30 #include "IOmanip.H"
31 #include "itoa.H"
32 #include "volPointInterpolation.H"
33 #include "ensightBinaryStream.H"
34 #include "ensightAsciiStream.H"
35 #include "globalIndex.H"
36 #include "ensightPTraits.H"
37 
38 using namespace Foam;
39 
40 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
41 
42 template<class Type>
45 (
46  const fvMeshSubset& meshSubsetter,
48 )
49 {
50  if (meshSubsetter.hasSubMesh())
51  {
53  (
54  meshSubsetter.interpolate(vf)
55  );
56  tfld.ref().checkOut();
57  tfld.ref().rename(vf.name());
58  return tfld;
59  }
60  else
61  {
62  return vf;
63  }
64 }
65 
66 
67 template<class Type>
68 Field<Type> map
69 (
70  const Field<Type>& vf,
71  const labelList& map1,
72  const labelList& map2
73 )
74 {
75  Field<Type> mf(map1.size() + map2.size());
76 
77  forAll(map1, i)
78  {
79  mf[i] = vf[map1[i]];
80  }
81 
82  label offset = map1.size();
83 
84  forAll(map2, i)
85  {
86  mf[i + offset] = vf[map2[i]];
87  }
88 
89  return mf;
90 }
91 
92 
93 template<class Type>
94 void writeField
95 (
96  const char* key,
97  const Field<Type>& vf,
98  ensightStream& ensightFile
99 )
100 {
101  if (returnReduce(vf.size(), sumOp<label>()) > 0)
102  {
103  if (Pstream::master())
104  {
105  ensightFile.write(key);
106 
107  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
108  {
109  ensightFile.write(vf.component(cmpt));
110 
111  for (int slave=1; slave<Pstream::nProcs(); slave++)
112  {
113  IPstream fromSlave(Pstream::commsTypes::scheduled, slave);
114  scalarField slaveData(fromSlave);
115  ensightFile.write(slaveData);
116  }
117  }
118  }
119  else
120  {
121  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
122  {
123  OPstream toMaster
124  (
127  );
128  toMaster<< vf.component(cmpt);
129  }
130  }
131  }
132 }
133 
134 
135 template<class Type>
136 bool writePatchField
137 (
138  const Field<Type>& pf,
139  const label patchi,
140  const label ensightPatchi,
141  const faceSets& boundaryFaceSet,
142  const ensightMesh::nFacePrimitives& nfp,
143  ensightStream& ensightFile
144 )
145 {
146  if (nfp.nTris || nfp.nQuads || nfp.nPolys)
147  {
148  if (Pstream::master())
149  {
150  ensightFile.writePartHeader(ensightPatchi);
151  }
152 
153  writeField
154  (
155  "tria3",
156  Field<Type>(pf, boundaryFaceSet.tris),
157  ensightFile
158  );
159 
160  writeField
161  (
162  "quad4",
163  Field<Type>(pf, boundaryFaceSet.quads),
164  ensightFile
165  );
166 
167  writeField
168  (
169  "nsided",
170  Field<Type>(pf, boundaryFaceSet.polys),
171  ensightFile
172  );
173 
174  return true;
175  }
176  else
177  {
178  return false;
179  }
180 }
181 
182 
183 template<class Type>
184 void writePatchField
185 (
186  const word& fieldName,
187  const Field<Type>& pf,
188  const word& patchName,
189  const ensightMesh& eMesh,
190  const fileName& postProcPath,
191  const word& prepend,
192  const label timeIndex,
193  const bool binary,
194  Ostream& ensightCaseFile
195 )
196 {
197  const Time& runTime = eMesh.mesh().time();
198 
199  const List<faceSets>& boundaryFaceSets = eMesh.boundaryFaceSets();
200  const wordList& allPatchNames = eMesh.allPatchNames();
202  nPatchPrims = eMesh.nPatchPrims();
203 
204  label ensightPatchi = eMesh.patchPartOffset();
205 
206  label patchi = -1;
207 
208  forAll(allPatchNames, i)
209  {
210  if (allPatchNames[i] == patchName)
211  {
212  patchi = i;
213  break;
214  }
215  ensightPatchi++;
216  }
217 
218 
219  word pfName = patchName + '.' + fieldName;
220 
221  word timeFile = prepend + itoa(timeIndex);
222 
223  ensightStream* ensightFilePtr = nullptr;
224  if (Pstream::master())
225  {
226  if (timeIndex == 0)
227  {
228  ensightCaseFile.setf(ios_base::left);
229 
230  ensightCaseFile
232  << " per element: 1 "
233  << setw(15) << pfName
234  << (' ' + prepend + "****." + pfName).c_str()
235  << nl;
236  }
237 
238  // set the filename of the ensight file
239  fileName ensightFileName(timeFile + "." + pfName);
240 
241  if (binary)
242  {
243  ensightFilePtr = new ensightBinaryStream
244  (
245  postProcPath/ensightFileName,
246  runTime
247  );
248  }
249  else
250  {
251  ensightFilePtr = new ensightAsciiStream
252  (
253  postProcPath/ensightFileName,
254  runTime
255  );
256  }
257  }
258 
259  ensightStream& ensightFile = *ensightFilePtr;
260 
261  if (Pstream::master())
262  {
264  }
265 
266  if (patchi >= 0)
267  {
269  (
270  pf,
271  patchi,
272  ensightPatchi,
273  boundaryFaceSets[patchi],
274  nPatchPrims.find(patchName)(),
275  ensightFile
276  );
277  }
278  else
279  {
280  faceSets nullFaceSets;
281 
283  (
284  Field<Type>(),
285  -1,
286  ensightPatchi,
287  nullFaceSets,
288  nPatchPrims.find(patchName)(),
289  ensightFile
290  );
291  }
292 
293  if (Pstream::master())
294  {
295  delete ensightFilePtr;
296  }
297 }
298 
299 
300 template<class Type>
301 void ensightField
302 (
304  const ensightMesh& eMesh,
305  const fileName& postProcPath,
306  const word& prepend,
307  const label timeIndex,
308  const bool binary,
309  Ostream& ensightCaseFile
310 )
311 {
312  Info<< "Converting field " << vf.name() << endl;
313 
314  word timeFile = prepend + itoa(timeIndex);
315 
316  const fvMesh& mesh = eMesh.mesh();
317  const Time& runTime = mesh.time();
318 
319  const cellSets& meshCellSets = eMesh.meshCellSets();
320  const List<faceSets>& boundaryFaceSets = eMesh.boundaryFaceSets();
321  const wordList& allPatchNames = eMesh.allPatchNames();
322  const wordHashSet& patchNames = eMesh.patchNames();
324  nPatchPrims = eMesh.nPatchPrims();
325  const List<faceSets>& faceZoneFaceSets = eMesh.faceZoneFaceSets();
326  const wordHashSet& faceZoneNames = eMesh.faceZoneNames();
328  nFaceZonePrims = eMesh.nFaceZonePrims();
329 
330  const labelList& tets = meshCellSets.tets;
331  const labelList& pyrs = meshCellSets.pyrs;
332  const labelList& prisms = meshCellSets.prisms;
333  const labelList& wedges = meshCellSets.wedges;
334  const labelList& hexes = meshCellSets.hexes;
335  const labelList& polys = meshCellSets.polys;
336 
337  ensightStream* ensightFilePtr = nullptr;
338  if (Pstream::master())
339  {
340  // set the filename of the ensight file
341  fileName ensightFileName(timeFile + "." + vf.name());
342 
343  if (binary)
344  {
345  ensightFilePtr = new ensightBinaryStream
346  (
347  postProcPath/ensightFileName,
348  runTime
349  );
350  }
351  else
352  {
353  ensightFilePtr = new ensightAsciiStream
354  (
355  postProcPath/ensightFileName,
356  runTime
357  );
358  }
359  }
360 
361  ensightStream& ensightFile = *ensightFilePtr;
362 
363  if (Pstream::master())
364  {
365  if (timeIndex == 0)
366  {
367  ensightCaseFile.setf(ios_base::left);
368 
369  ensightCaseFile
371  << " per element: 1 "
372  << setw(15) << vf.name()
373  << (' ' + prepend + "****." + vf.name()).c_str()
374  << nl;
375  }
376 
378  }
379 
380  if (patchNames.empty())
381  {
382  eMesh.barrier();
383 
384  if (Pstream::master())
385  {
386  ensightFile.writePartHeader(1);
387  }
388 
389  writeField
390  (
391  "hexa8",
392  map(vf, hexes, wedges),
393  ensightFile
394  );
395 
396  writeField
397  (
398  "penta6",
399  Field<Type>(vf, prisms),
400  ensightFile
401  );
402 
403  writeField
404  (
405  "pyramid5",
406  Field<Type>(vf, pyrs),
407  ensightFile
408  );
409 
410  writeField
411  (
412  "tetra4",
413  Field<Type>(vf, tets),
414  ensightFile
415  );
416 
417  writeField
418  (
419  "nfaced",
420  Field<Type>(vf, polys),
421  ensightFile
422  );
423  }
424 
425  label ensightPatchi = eMesh.patchPartOffset();
426 
427  forAll(allPatchNames, patchi)
428  {
429  const word& patchName = allPatchNames[patchi];
430 
431  eMesh.barrier();
432 
433  if (patchNames.empty() || patchNames.found(patchName))
434  {
435  if
436  (
438  (
439  vf.boundaryField()[patchi],
440  patchi,
441  ensightPatchi,
442  boundaryFaceSets[patchi],
443  nPatchPrims.find(patchName)(),
444  ensightFile
445  )
446  )
447  {
448  ensightPatchi++;
449  }
450  }
451  }
452 
453  // write faceZones, if requested
454  if (faceZoneNames.size())
455  {
456  // Interpolates cell values to faces - needed only when exporting
457  // faceZones...
459  (
461  );
462 
463  forAllConstIter(wordHashSet, faceZoneNames, iter)
464  {
465  const word& faceZoneName = iter.key();
466 
467  eMesh.barrier();
468 
469  label zoneID = mesh.faceZones().findZoneID(faceZoneName);
470 
471  const faceZone& fz = mesh.faceZones()[zoneID];
472 
473  // Prepare data to write
474  label nIncluded = 0;
475  forAll(fz, i)
476  {
477  if (eMesh.faceToBeIncluded(fz[i]))
478  {
479  ++nIncluded;
480  }
481  }
482 
483  Field<Type> values(nIncluded);
484 
485  // Loop on the faceZone and store the needed field values
486  label j = 0;
487  forAll(fz, i)
488  {
489  label facei = fz[i];
490  if (mesh.isInternalFace(facei))
491  {
492  values[j] = sf[facei];
493  ++j;
494  }
495  else
496  {
497  if (eMesh.faceToBeIncluded(facei))
498  {
499  label patchi = mesh.boundaryMesh().whichPatch(facei);
500  const polyPatch& pp = mesh.boundaryMesh()[patchi];
501  label patchFacei = pp.whichFace(facei);
502  Type value = sf.boundaryField()[patchi][patchFacei];
503  values[j] = value;
504  ++j;
505  }
506  }
507  }
508 
509  if
510  (
512  (
513  values,
514  zoneID,
515  ensightPatchi,
516  faceZoneFaceSets[zoneID],
517  nFaceZonePrims.find(faceZoneName)(),
518  ensightFile
519  )
520  )
521  {
522  ensightPatchi++;
523  }
524  }
525  }
526  if (Pstream::master())
527  {
528  delete ensightFilePtr;
529  }
530 }
531 
532 
533 template<class Type>
534 void ensightPointField
535 (
537  const ensightMesh& eMesh,
538  const fileName& postProcPath,
539  const word& prepend,
540  const label timeIndex,
541  const bool binary,
542  Ostream& ensightCaseFile
543 )
544 {
545  Info<< "Converting field " << pf.name() << endl;
546 
547  word timeFile = prepend + itoa(timeIndex);
548 
549  const fvMesh& mesh = eMesh.mesh();
550  const wordList& allPatchNames = eMesh.allPatchNames();
551  const wordHashSet& patchNames = eMesh.patchNames();
552  const wordHashSet& faceZoneNames = eMesh.faceZoneNames();
553 
554 
555  ensightStream* ensightFilePtr = nullptr;
556  if (Pstream::master())
557  {
558  // set the filename of the ensight file
559  fileName ensightFileName(timeFile + "." + pf.name());
560 
561  if (binary)
562  {
563  ensightFilePtr = new ensightBinaryStream
564  (
565  postProcPath/ensightFileName,
566  mesh.time()
567  );
568  }
569  else
570  {
571  ensightFilePtr = new ensightAsciiStream
572  (
573  postProcPath/ensightFileName,
574  mesh.time()
575  );
576  }
577  }
578 
579  ensightStream& ensightFile = *ensightFilePtr;
580 
581  if (Pstream::master())
582  {
583  if (timeIndex == 0)
584  {
585  ensightCaseFile.setf(ios_base::left);
586 
587  ensightCaseFile
589  << " per node: 1 "
590  << setw(15) << pf.name()
591  << (' ' + prepend + "****." + pf.name()).c_str()
592  << nl;
593  }
594 
595  ensightFile.write(ensightPTraits<Type>::typeName);
596  }
597 
598  if (eMesh.patchNames().empty())
599  {
600  eMesh.barrier();
601 
602  if (Pstream::master())
603  {
604  ensightFile.writePartHeader(1);
605  }
606 
607  writeField
608  (
609  "coordinates",
611  ensightFile
612  );
613  }
614 
615 
616  label ensightPatchi = eMesh.patchPartOffset();
617 
618  forAll(allPatchNames, patchi)
619  {
620  const word& patchName = allPatchNames[patchi];
621 
622  eMesh.barrier();
623 
624  if (patchNames.empty() || patchNames.found(patchName))
625  {
626  const fvPatch& p = mesh.boundary()[patchi];
627  if
628  (
630  > 0
631  )
632  {
633  // Renumber the patch points/faces into unique points
634  labelList pointToGlobal;
635  labelList uniqueMeshPointLabels;
636  autoPtr<globalIndex> globalPointsPtr =
637  mesh.globalData().mergePoints
638  (
639  p.patch().meshPoints(),
640  p.patch().meshPointMap(),
641  pointToGlobal,
642  uniqueMeshPointLabels
643  );
644 
645  if (Pstream::master())
646  {
647  ensightFile.writePartHeader(ensightPatchi);
648  }
649 
650  writeField
651  (
652  "coordinates",
653  Field<Type>(pf.primitiveField(), uniqueMeshPointLabels),
654  ensightFile
655  );
656 
657  ensightPatchi++;
658  }
659  }
660  }
661 
662  // write faceZones, if requested
663  if (faceZoneNames.size())
664  {
665  forAllConstIter(wordHashSet, faceZoneNames, iter)
666  {
667  const word& faceZoneName = iter.key();
668 
669  eMesh.barrier();
670 
671  label zoneID = mesh.faceZones().findZoneID(faceZoneName);
672 
673  const faceZone& fz = mesh.faceZones()[zoneID];
674 
675  if (returnReduce(fz().nPoints(), sumOp<label>()) > 0)
676  {
677  // Renumber the faceZone points/faces into unique points
678  labelList pointToGlobal;
679  labelList uniqueMeshPointLabels;
680  autoPtr<globalIndex> globalPointsPtr =
681  mesh.globalData().mergePoints
682  (
683  fz().meshPoints(),
684  fz().meshPointMap(),
685  pointToGlobal,
686  uniqueMeshPointLabels
687  );
688 
689  if (Pstream::master())
690  {
691  ensightFile.writePartHeader(ensightPatchi);
692  }
693 
694  writeField
695  (
696  "coordinates",
698  (
699  pf.primitiveField(),
700  uniqueMeshPointLabels
701  ),
702  ensightFile
703  );
704 
705  ensightPatchi++;
706  }
707  }
708  }
709 
710  if (Pstream::master())
711  {
712  delete ensightFilePtr;
713  }
714 }
715 
716 
717 template<class Type>
718 void ensightField
719 (
721  const ensightMesh& eMesh,
722  const fileName& postProcPath,
723  const word& prepend,
724  const label timeIndex,
725  const bool binary,
726  const bool nodeValues,
727  Ostream& ensightCaseFile
728 )
729 {
730  if (nodeValues)
731  {
733  (
735  );
736  pfld.ref().rename(vf.name());
737 
738  ensightPointField<Type>
739  (
740  pfld,
741  eMesh,
742  postProcPath,
743  prepend,
744  timeIndex,
745  binary,
746  ensightCaseFile
747  );
748  }
749  else
750  {
751  ensightField<Type>
752  (
753  vf,
754  eMesh,
755  postProcPath,
756  prepend,
757  timeIndex,
758  binary,
759  ensightCaseFile
760  );
761  }
762 }
763 
764 
765 // ************************************************************************* //
const HashTable< nFacePrimitives > & nPatchPrims() const
Definition: ensightMesh.H:303
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
A HashTable with keys but without contents.
Definition: HashSet.H:59
labelList polys
Definition: cellSets.H:59
label findZoneID(const word &zoneName) const
Find zone index given a name.
Definition: MeshZones.C:341
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
static int masterNo()
Process index of the master.
Definition: UPstream.H:417
const word & name() const
Return name.
Definition: IOobject.H:303
A class for handling file names.
Definition: fileName.H:79
void ensightField(const Foam::GeometricField< Type, Foam::fvPatchField, Foam::volMesh > &vf, const Foam::ensightMesh &eMesh, const Foam::fileName &postProcPath, const Foam::word &prepend, const Foam::label timeIndex, const bool binary, const bool nodeValues, Foam::Ostream &ensightCaseFile)
const wordHashSet & faceZoneNames() const
Definition: ensightMesh.H:313
labelList hexes
Definition: cellSets.H:58
ios_base::fmtflags setf(const ios_base::fmtflags f)
Set flags of stream.
Definition: IOstream.H:493
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > linearInterpolate(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: linear.H:108
const Boundary & boundaryField() const
Return const-reference to the boundary field.
word itoa(const label n)
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
label patchPartOffset() const
The ensight part id for the first patch.
Definition: ensightMesh.H:324
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
Foam::tmp< Foam::GeometricField< Type, Foam::fvPatchField, Foam::volMesh > > volField(const Foam::fvMeshSubset &, const Foam::GeometricField< Type, Foam::fvPatchField, Foam::volMesh > &vf)
Wrapper to get hold of the field or the subsetted field.
uint8_t direction
Definition: direction.H:45
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
bool empty() const
Return true if the hash table is empty.
Definition: HashTableI.H:72
labelList quads
Definition: faceSets.H:53
const labelList & uniquePointMap() const
Local points that are unique.
Definition: ensightMesh.H:345
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:62
virtual void write(const char *)=0
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Generic GeometricField class.
autoPtr< globalIndex > mergePoints(labelList &pointToGlobal, labelList &uniquePoints) const
Helper for merging (collocated!) mesh point data.
labelList pyrs
Definition: cellSets.H:55
Abstract base class for writing Ensight data.
Definition: ensightStream.H:50
labelList polys
Definition: faceSets.H:54
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
Input inter-processor communications stream.
Definition: IPstream.H:50
const cellSets & meshCellSets() const
Definition: ensightMesh.H:283
labelList tets
Definition: cellSets.H:54
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:142
const Map< label > & meshPointMap() const
Mesh point map. Given the global point index find its.
Conversion of OpenFOAM pTraits into the Ensight equivalent.
const List< faceSets > & boundaryFaceSets() const
Definition: ensightMesh.H:288
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:113
Pre-declare SubField and related Field type.
Definition: Field.H:56
A class for handling words, derived from string.
Definition: word.H:59
label nPoints
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:138
const wordHashSet & patchNames() const
Definition: ensightMesh.H:298
const fvMesh & mesh() const
Definition: ensightMesh.H:278
labelList tris
Definition: faceSets.H:52
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1394
An STL-conforming hash table.
Definition: HashTable.H:61
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Definition: Field.C:441
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
virtual label size() const
Return size.
Definition: fvPatch.H:156
const HashTable< nFacePrimitives > & nFaceZonePrims() const
Definition: ensightMesh.H:318
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
Istream and Ostream manipulators taking arguments.
labelList prisms
Definition: cellSets.H:56
static const char nl
Definition: Ostream.H:260
const Mesh & mesh() const
Return mesh.
const wordList & allPatchNames() const
Definition: ensightMesh.H:293
static void barrier()
Helper to cause barrier. Necessary on Quadrics.
Post-processing mesh subset tool. Given the original mesh and the list of selected cells...
Definition: fvMeshSubset.H:73
Output inter-processor communications stream.
Definition: OPstream.H:50
volScalarField sf(fieldObject, mesh)
static tmp< GeometricField< Type, fvPatchField, volMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &, const fvMesh &sMesh, const labelList &patchMap, const labelList &cellMap, const labelList &faceMap)
Map volume field.
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
void writePatchField(const Foam::word &fieldName, const Foam::Field< Type > &pf, const Foam::word &patchName, const Foam::ensightMesh &eMesh, const Foam::fileName &postProcPath, const Foam::word &prepend, const Foam::label timeIndex, Foam::Ostream &ensightCaseFile)
label patchi
const meshFaceZones & faceZones() const
Return face zones.
Definition: polyMesh.H:476
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
label timeIndex
Definition: getTimeIndex.H:4
const List< faceSets > & faceZoneFaceSets() const
Definition: ensightMesh.H:308
messageStream Info
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Omanip< int > setw(const int i)
Definition: IOmanip.H:199
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:64
virtual void writePartHeader(const label)=0
A class for managing temporary objects.
Definition: PtrList.H:53
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
bool faceToBeIncluded(const label facei) const
When exporting faceZones, check if a given face has to be included.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
labelList wedges
Definition: cellSets.H:57
Namespace for OpenFOAM.
bool hasSubMesh() const
Have subMesh?
label whichFace(const label l) const
Return label of face in patch from global face label.
Definition: polyPatch.H:389
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540