ensightField.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 | Copyright (C) 2011-2016 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::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  {
124  toMaster<< vf.component(cmpt);
125  }
126  }
127  }
128 }
129 
130 
131 template<class Type>
132 bool writePatchField
133 (
134  const Field<Type>& pf,
135  const label patchi,
136  const label ensightPatchi,
137  const faceSets& boundaryFaceSet,
138  const ensightMesh::nFacePrimitives& nfp,
139  ensightStream& ensightFile
140 )
141 {
142  if (nfp.nTris || nfp.nQuads || nfp.nPolys)
143  {
144  if (Pstream::master())
145  {
146  ensightFile.writePartHeader(ensightPatchi);
147  }
148 
149  writeField
150  (
151  "tria3",
152  Field<Type>(pf, boundaryFaceSet.tris),
153  ensightFile
154  );
155 
156  writeField
157  (
158  "quad4",
159  Field<Type>(pf, boundaryFaceSet.quads),
160  ensightFile
161  );
162 
163  writeField
164  (
165  "nsided",
166  Field<Type>(pf, boundaryFaceSet.polys),
167  ensightFile
168  );
169 
170  return true;
171  }
172  else
173  {
174  return false;
175  }
176 }
177 
178 
179 template<class Type>
180 void writePatchField
181 (
182  const word& fieldName,
183  const Field<Type>& pf,
184  const word& patchName,
185  const ensightMesh& eMesh,
186  const fileName& postProcPath,
187  const word& prepend,
188  const label timeIndex,
189  const bool binary,
190  Ostream& ensightCaseFile
191 )
192 {
193  const Time& runTime = eMesh.mesh().time();
194 
195  const List<faceSets>& boundaryFaceSets = eMesh.boundaryFaceSets();
196  const wordList& allPatchNames = eMesh.allPatchNames();
198  nPatchPrims = eMesh.nPatchPrims();
199 
200  label ensightPatchi = eMesh.patchPartOffset();
201 
202  label patchi = -1;
203 
204  forAll(allPatchNames, i)
205  {
206  if (allPatchNames[i] == patchName)
207  {
208  patchi = i;
209  break;
210  }
211  ensightPatchi++;
212  }
213 
214 
215  word pfName = patchName + '.' + fieldName;
216 
217  word timeFile = prepend + itoa(timeIndex);
218 
219  ensightStream* ensightFilePtr = NULL;
220  if (Pstream::master())
221  {
222  if (timeIndex == 0)
223  {
224  ensightCaseFile.setf(ios_base::left);
225 
226  ensightCaseFile
228  << " per element: 1 "
229  << setw(15) << pfName
230  << (' ' + prepend + "****." + pfName).c_str()
231  << nl;
232  }
233 
234  // set the filename of the ensight file
235  fileName ensightFileName(timeFile + "." + pfName);
236 
237  if (binary)
238  {
239  ensightFilePtr = new ensightBinaryStream
240  (
241  postProcPath/ensightFileName,
242  runTime
243  );
244  }
245  else
246  {
247  ensightFilePtr = new ensightAsciiStream
248  (
249  postProcPath/ensightFileName,
250  runTime
251  );
252  }
253  }
254 
255  ensightStream& ensightFile = *ensightFilePtr;
256 
257  if (Pstream::master())
258  {
260  }
261 
262  if (patchi >= 0)
263  {
265  (
266  pf,
267  patchi,
268  ensightPatchi,
269  boundaryFaceSets[patchi],
270  nPatchPrims.find(patchName)(),
271  ensightFile
272  );
273  }
274  else
275  {
276  faceSets nullFaceSets;
277 
279  (
280  Field<Type>(),
281  -1,
282  ensightPatchi,
283  nullFaceSets,
284  nPatchPrims.find(patchName)(),
285  ensightFile
286  );
287  }
288 
289  if (Pstream::master())
290  {
291  delete ensightFilePtr;
292  }
293 }
294 
295 
296 template<class Type>
297 void ensightField
298 (
300  const ensightMesh& eMesh,
301  const fileName& postProcPath,
302  const word& prepend,
303  const label timeIndex,
304  const bool binary,
305  Ostream& ensightCaseFile
306 )
307 {
308  Info<< "Converting field " << vf.name() << endl;
309 
310  word timeFile = prepend + itoa(timeIndex);
311 
312  const fvMesh& mesh = eMesh.mesh();
313  const Time& runTime = mesh.time();
314 
315  const cellSets& meshCellSets = eMesh.meshCellSets();
316  const List<faceSets>& boundaryFaceSets = eMesh.boundaryFaceSets();
317  const wordList& allPatchNames = eMesh.allPatchNames();
318  const wordHashSet& patchNames = eMesh.patchNames();
320  nPatchPrims = eMesh.nPatchPrims();
321  const List<faceSets>& faceZoneFaceSets = eMesh.faceZoneFaceSets();
322  const wordHashSet& faceZoneNames = eMesh.faceZoneNames();
324  nFaceZonePrims = eMesh.nFaceZonePrims();
325 
326  const labelList& tets = meshCellSets.tets;
327  const labelList& pyrs = meshCellSets.pyrs;
328  const labelList& prisms = meshCellSets.prisms;
329  const labelList& wedges = meshCellSets.wedges;
330  const labelList& hexes = meshCellSets.hexes;
331  const labelList& polys = meshCellSets.polys;
332 
333  ensightStream* ensightFilePtr = NULL;
334  if (Pstream::master())
335  {
336  // set the filename of the ensight file
337  fileName ensightFileName(timeFile + "." + vf.name());
338 
339  if (binary)
340  {
341  ensightFilePtr = new ensightBinaryStream
342  (
343  postProcPath/ensightFileName,
344  runTime
345  );
346  }
347  else
348  {
349  ensightFilePtr = new ensightAsciiStream
350  (
351  postProcPath/ensightFileName,
352  runTime
353  );
354  }
355  }
356 
357  ensightStream& ensightFile = *ensightFilePtr;
358 
359  if (Pstream::master())
360  {
361  if (timeIndex == 0)
362  {
363  ensightCaseFile.setf(ios_base::left);
364 
365  ensightCaseFile
367  << " per element: 1 "
368  << setw(15) << vf.name()
369  << (' ' + prepend + "****." + vf.name()).c_str()
370  << nl;
371  }
372 
374  }
375 
376  if (patchNames.empty())
377  {
378  eMesh.barrier();
379 
380  if (Pstream::master())
381  {
382  ensightFile.writePartHeader(1);
383  }
384 
385  writeField
386  (
387  "hexa8",
388  map(vf, hexes, wedges),
389  ensightFile
390  );
391 
392  writeField
393  (
394  "penta6",
395  Field<Type>(vf, prisms),
396  ensightFile
397  );
398 
399  writeField
400  (
401  "pyramid5",
402  Field<Type>(vf, pyrs),
403  ensightFile
404  );
405 
406  writeField
407  (
408  "tetra4",
409  Field<Type>(vf, tets),
410  ensightFile
411  );
412 
413  writeField
414  (
415  "nfaced",
416  Field<Type>(vf, polys),
417  ensightFile
418  );
419  }
420 
421  label ensightPatchi = eMesh.patchPartOffset();
422 
423  forAll(allPatchNames, patchi)
424  {
425  const word& patchName = allPatchNames[patchi];
426 
427  eMesh.barrier();
428 
429  if (patchNames.empty() || patchNames.found(patchName))
430  {
431  if
432  (
434  (
435  vf.boundaryField()[patchi],
436  patchi,
437  ensightPatchi,
438  boundaryFaceSets[patchi],
439  nPatchPrims.find(patchName)(),
440  ensightFile
441  )
442  )
443  {
444  ensightPatchi++;
445  }
446  }
447  }
448 
449  // write faceZones, if requested
450  if (faceZoneNames.size())
451  {
452  // Interpolates cell values to faces - needed only when exporting
453  // faceZones...
455  (
457  );
458 
459  forAllConstIter(wordHashSet, faceZoneNames, iter)
460  {
461  const word& faceZoneName = iter.key();
462 
463  eMesh.barrier();
464 
465  label zoneID = mesh.faceZones().findZoneID(faceZoneName);
466 
467  const faceZone& fz = mesh.faceZones()[zoneID];
468 
469  // Prepare data to write
470  label nIncluded = 0;
471  forAll(fz, i)
472  {
473  if (eMesh.faceToBeIncluded(fz[i]))
474  {
475  ++nIncluded;
476  }
477  }
478 
479  Field<Type> values(nIncluded);
480 
481  // Loop on the faceZone and store the needed field values
482  label j = 0;
483  forAll(fz, i)
484  {
485  label facei = fz[i];
486  if (mesh.isInternalFace(facei))
487  {
488  values[j] = sf[facei];
489  ++j;
490  }
491  else
492  {
493  if (eMesh.faceToBeIncluded(facei))
494  {
495  label patchi = mesh.boundaryMesh().whichPatch(facei);
496  const polyPatch& pp = mesh.boundaryMesh()[patchi];
497  label patchFacei = pp.whichFace(facei);
498  Type value = sf.boundaryField()[patchi][patchFacei];
499  values[j] = value;
500  ++j;
501  }
502  }
503  }
504 
505  if
506  (
508  (
509  values,
510  zoneID,
511  ensightPatchi,
512  faceZoneFaceSets[zoneID],
513  nFaceZonePrims.find(faceZoneName)(),
514  ensightFile
515  )
516  )
517  {
518  ensightPatchi++;
519  }
520  }
521  }
522  if (Pstream::master())
523  {
524  delete ensightFilePtr;
525  }
526 }
527 
528 
529 template<class Type>
530 void ensightPointField
531 (
533  const ensightMesh& eMesh,
534  const fileName& postProcPath,
535  const word& prepend,
536  const label timeIndex,
537  const bool binary,
538  Ostream& ensightCaseFile
539 )
540 {
541  Info<< "Converting field " << pf.name() << endl;
542 
543  word timeFile = prepend + itoa(timeIndex);
544 
545  const fvMesh& mesh = eMesh.mesh();
546  const wordList& allPatchNames = eMesh.allPatchNames();
547  const wordHashSet& patchNames = eMesh.patchNames();
548  const wordHashSet& faceZoneNames = eMesh.faceZoneNames();
549 
550 
551  ensightStream* ensightFilePtr = NULL;
552  if (Pstream::master())
553  {
554  // set the filename of the ensight file
555  fileName ensightFileName(timeFile + "." + pf.name());
556 
557  if (binary)
558  {
559  ensightFilePtr = new ensightBinaryStream
560  (
561  postProcPath/ensightFileName,
562  mesh.time()
563  );
564  }
565  else
566  {
567  ensightFilePtr = new ensightAsciiStream
568  (
569  postProcPath/ensightFileName,
570  mesh.time()
571  );
572  }
573  }
574 
575  ensightStream& ensightFile = *ensightFilePtr;
576 
577  if (Pstream::master())
578  {
579  if (timeIndex == 0)
580  {
581  ensightCaseFile.setf(ios_base::left);
582 
583  ensightCaseFile
585  << " per node: 1 "
586  << setw(15) << pf.name()
587  << (' ' + prepend + "****." + pf.name()).c_str()
588  << nl;
589  }
590 
591  ensightFile.write(ensightPTraits<Type>::typeName);
592  }
593 
594  if (eMesh.patchNames().empty())
595  {
596  eMesh.barrier();
597 
598  if (Pstream::master())
599  {
600  ensightFile.writePartHeader(1);
601  }
602 
603  writeField
604  (
605  "coordinates",
607  ensightFile
608  );
609  }
610 
611 
612  label ensightPatchi = eMesh.patchPartOffset();
613 
614  forAll(allPatchNames, patchi)
615  {
616  const word& patchName = allPatchNames[patchi];
617 
618  eMesh.barrier();
619 
620  if (patchNames.empty() || patchNames.found(patchName))
621  {
622  const fvPatch& p = mesh.boundary()[patchi];
623  if
624  (
626  > 0
627  )
628  {
629  // Renumber the patch points/faces into unique points
630  labelList pointToGlobal;
631  labelList uniqueMeshPointLabels;
632  autoPtr<globalIndex> globalPointsPtr =
633  mesh.globalData().mergePoints
634  (
635  p.patch().meshPoints(),
636  p.patch().meshPointMap(),
637  pointToGlobal,
638  uniqueMeshPointLabels
639  );
640 
641  if (Pstream::master())
642  {
643  ensightFile.writePartHeader(ensightPatchi);
644  }
645 
646  writeField
647  (
648  "coordinates",
649  Field<Type>(pf.primitiveField(), uniqueMeshPointLabels),
650  ensightFile
651  );
652 
653  ensightPatchi++;
654  }
655  }
656  }
657 
658  // write faceZones, if requested
659  if (faceZoneNames.size())
660  {
661  forAllConstIter(wordHashSet, faceZoneNames, iter)
662  {
663  const word& faceZoneName = iter.key();
664 
665  eMesh.barrier();
666 
667  label zoneID = mesh.faceZones().findZoneID(faceZoneName);
668 
669  const faceZone& fz = mesh.faceZones()[zoneID];
670 
671  if (returnReduce(fz().nPoints(), sumOp<label>()) > 0)
672  {
673  // Renumber the faceZone points/faces into unique points
674  labelList pointToGlobal;
675  labelList uniqueMeshPointLabels;
676  autoPtr<globalIndex> globalPointsPtr =
677  mesh.globalData().mergePoints
678  (
679  fz().meshPoints(),
680  fz().meshPointMap(),
681  pointToGlobal,
682  uniqueMeshPointLabels
683  );
684 
685  if (Pstream::master())
686  {
687  ensightFile.writePartHeader(ensightPatchi);
688  }
689 
690  writeField
691  (
692  "coordinates",
694  (
695  pf.primitiveField(),
696  uniqueMeshPointLabels
697  ),
698  ensightFile
699  );
700 
701  ensightPatchi++;
702  }
703  }
704  }
705 
706  if (Pstream::master())
707  {
708  delete ensightFilePtr;
709  }
710 }
711 
712 
713 template<class Type>
714 void ensightField
715 (
717  const ensightMesh& eMesh,
718  const fileName& postProcPath,
719  const word& prepend,
720  const label timeIndex,
721  const bool binary,
722  const bool nodeValues,
723  Ostream& ensightCaseFile
724 )
725 {
726  if (nodeValues)
727  {
729  (
731  );
732  pfld.ref().rename(vf.name());
733 
734  ensightPointField<Type>
735  (
736  pfld,
737  eMesh,
738  postProcPath,
739  prepend,
740  timeIndex,
741  binary,
742  ensightCaseFile
743  );
744  }
745  else
746  {
747  ensightField<Type>
748  (
749  vf,
750  eMesh,
751  postProcPath,
752  prepend,
753  timeIndex,
754  binary,
755  ensightCaseFile
756  );
757  }
758 }
759 
760 
761 // ************************************************************************* //
A HashTable with keys but without contents.
Definition: HashSet.H:59
labelList polys
Definition: cellSets.H:59
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const List< faceSets > & faceZoneFaceSets() const
Definition: ensightMesh.H:311
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
uint8_t direction
Definition: direction.H:46
static int masterNo()
Process index of the master.
Definition: UPstream.H:405
const cellSets & meshCellSets() const
Definition: ensightMesh.H:286
autoPtr< globalIndex > mergePoints(labelList &pointToGlobal, labelList &uniquePoints) const
Helper for merging (collocated!) mesh point data.
A class for handling file names.
Definition: fileName.H:69
label whichFace(const label l) const
Return label of face in patch from global face label.
Definition: polyPatch.H:377
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)
labelList hexes
Definition: cellSets.H:58
ios_base::fmtflags setf(const ios_base::fmtflags f)
Set flags of stream.
Definition: IOstream.H:496
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > linearInterpolate(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: linear.H:107
word itoa(const label n)
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
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.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
labelList quads
Definition: faceSets.H:53
bool faceToBeIncluded(const label facei) const
When exporting faceZones, check if a given face has to be included.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
virtual void write(const char *)=0
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:411
const labelList & uniquePointMap() const
Local points that are unique.
Definition: ensightMesh.H:348
const Map< label > & meshPointMap() const
Mesh point map. Given the global point index find its.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Generic GeometricField class.
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
labelList pyrs
Definition: cellSets.H:55
Abstract base class for writing Ensight data.
Definition: ensightStream.H:50
labelList polys
Definition: faceSets.H:54
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:143
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
const HashTable< nFacePrimitives > & nFaceZonePrims() const
Definition: ensightMesh.H:321
static const volPointInterpolation & New(const fvMesh &mesh)
Input inter-processor communications stream.
Definition: IPstream.H:50
labelList tets
Definition: cellSets.H:54
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:138
label findZoneID(const word &zoneName) const
Find zone index given a name.
Definition: ZoneMesh.C:341
Conversion of OpenFOAM pTraits into the Ensight equivalent.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
label patchPartOffset() const
The ensight part id for the first patch.
Definition: ensightMesh.H:327
Pre-declare SubField and related Field type.
Definition: Field.H:57
const List< faceSets > & boundaryFaceSets() const
Definition: ensightMesh.H:291
A class for handling words, derived from string.
Definition: word.H:59
label nPoints
const HashTable< nFacePrimitives > & nPatchPrims() const
Definition: ensightMesh.H:306
virtual label size() const
Return size.
Definition: fvPatch.H:161
labelList tris
Definition: faceSets.H:52
const wordHashSet & faceZoneNames() const
Definition: ensightMesh.H:316
bool hasSubMesh() const
Have subMesh?
An STL-conforming hash table.
Definition: HashTable.H:61
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
bool empty() const
Return true if the hash table is empty.
Definition: HashTableI.H:72
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:109
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1140
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
Istream and Ostream manipulators taking arguments.
const fvMesh & mesh() const
Definition: ensightMesh.H:281
labelList prisms
Definition: cellSets.H:56
static const char nl
Definition: Ostream.H:262
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Definition: Field.C:650
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 label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:399
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 Mesh & mesh() const
Return mesh.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const wordHashSet & patchNames() const
Definition: ensightMesh.H:301
label timeIndex
Definition: getTimeIndex.H:4
const wordList & allPatchNames() const
Definition: ensightMesh.H:296
messageStream Info
const faceZoneMesh & faceZones() const
Return face zone mesh.
Definition: polyMesh.H:463
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:53
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:54
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
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.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:545
labelList wedges
Definition: cellSets.H:57
const word & name() const
Return name.
Definition: IOobject.H:260
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243