1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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.
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
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/>.
24 \*---------------------------------------------------------------------------*/
26 #include "ensightField.H"
28 #include "volFields.H"
32 #include "volPointInterpolation.H"
33 #include "ensightBinaryStream.H"
34 #include "ensightAsciiStream.H"
35 #include "globalIndex.H"
39 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
44 const Field<Type>& vf,
45 const labelList& map1,
49 Field<Type> mf(map1.size() + map2.size());
56 label offset = map1.size();
60 mf[i + offset] = vf[map2[i]];
71 const Field<Type>& vf,
72 ensightStream& ensightFile
75 if (returnReduce(vf.size(), sumOp<label>()) > 0)
77 if (Pstream::master())
79 ensightFile.write(key);
81 for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
83 ensightFile.write(vf.component(cmpt));
85 for (int slave=1; slave<Pstream::nProcs(); slave++)
87 IPstream fromSlave(Pstream::scheduled, slave);
88 scalarField slaveData(fromSlave);
89 ensightFile.write(slaveData);
95 for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
97 OPstream toMaster(Pstream::scheduled, Pstream::masterNo());
98 toMaster<< vf.component(cmpt);
108 const Field<Type>& pf,
110 const label ensightPatchI,
111 const faceSets& boundaryFaceSet,
112 const ensightMesh::nFacePrimitives& nfp,
113 ensightStream& ensightFile
116 if (nfp.nTris || nfp.nQuads || nfp.nPolys)
118 if (Pstream::master())
120 ensightFile.writePartHeader(ensightPatchI);
126 Field<Type>(pf, boundaryFaceSet.tris),
133 Field<Type>(pf, boundaryFaceSet.quads),
140 Field<Type>(pf, boundaryFaceSet.polys),
156 const word& fieldName,
157 const Field<Type>& pf,
158 const word& patchName,
159 const ensightMesh& eMesh,
160 const fileName& postProcPath,
162 const label timeIndex,
164 Ostream& ensightCaseFile
167 const Time& runTime = eMesh.mesh().time();
169 const List<faceSets>& boundaryFaceSets = eMesh.boundaryFaceSets();
170 const wordList& allPatchNames = eMesh.allPatchNames();
171 const HashTable<ensightMesh::nFacePrimitives>&
172 nPatchPrims = eMesh.nPatchPrims();
174 label ensightPatchI = eMesh.patchPartOffset();
178 forAll(allPatchNames, i)
180 if (allPatchNames[i] == patchName)
189 word pfName = patchName + '.' + fieldName;
191 word timeFile = prepend + itoa(timeIndex);
193 ensightStream* ensightFilePtr = NULL;
194 if (Pstream::master())
198 ensightCaseFile.setf(ios_base::left);
201 << pTraits<Type>::typeName
202 << " per element: 1 "
203 << setw(15) << pfName
204 << (' ' + prepend + "***." + pfName).c_str()
208 // set the filename of the ensight file
209 fileName ensightFileName(timeFile + "." + pfName);
213 ensightFilePtr = new ensightBinaryStream
215 postProcPath/ensightFileName,
221 ensightFilePtr = new ensightAsciiStream
223 postProcPath/ensightFileName,
229 ensightStream& ensightFile = *ensightFilePtr;
231 if (Pstream::master())
233 ensightFile.write(pTraits<Type>::typeName);
243 boundaryFaceSets[patchi],
244 nPatchPrims.find(patchName)(),
250 faceSets nullFaceSets;
258 nPatchPrims.find(patchName)(),
263 if (Pstream::master())
265 delete ensightFilePtr;
273 const GeometricField<Type, fvPatchField, volMesh>& vf,
274 const ensightMesh& eMesh,
275 const fileName& postProcPath,
277 const label timeIndex,
279 Ostream& ensightCaseFile
282 Info<< "Converting field " << vf.name() << endl;
284 word timeFile = prepend + itoa(timeIndex);
286 const fvMesh& mesh = eMesh.mesh();
287 const Time& runTime = mesh.time();
289 const cellSets& meshCellSets = eMesh.meshCellSets();
290 const List<faceSets>& boundaryFaceSets = eMesh.boundaryFaceSets();
291 const wordList& allPatchNames = eMesh.allPatchNames();
292 const wordHashSet& patchNames = eMesh.patchNames();
293 const HashTable<ensightMesh::nFacePrimitives>&
294 nPatchPrims = eMesh.nPatchPrims();
295 const List<faceSets>& faceZoneFaceSets = eMesh.faceZoneFaceSets();
296 const wordHashSet& faceZoneNames = eMesh.faceZoneNames();
297 const HashTable<ensightMesh::nFacePrimitives>&
298 nFaceZonePrims = eMesh.nFaceZonePrims();
300 const labelList& tets = meshCellSets.tets;
301 const labelList& pyrs = meshCellSets.pyrs;
302 const labelList& prisms = meshCellSets.prisms;
303 const labelList& wedges = meshCellSets.wedges;
304 const labelList& hexes = meshCellSets.hexes;
305 const labelList& polys = meshCellSets.polys;
307 ensightStream* ensightFilePtr = NULL;
308 if (Pstream::master())
310 // set the filename of the ensight file
311 fileName ensightFileName(timeFile + "." + vf.name());
315 ensightFilePtr = new ensightBinaryStream
317 postProcPath/ensightFileName,
323 ensightFilePtr = new ensightAsciiStream
325 postProcPath/ensightFileName,
331 ensightStream& ensightFile = *ensightFilePtr;
333 if (patchNames.empty())
337 if (Pstream::master())
341 ensightCaseFile.setf(ios_base::left);
344 << pTraits<Type>::typeName
345 << " per element: 1 "
346 << setw(15) << vf.name()
347 << (' ' + prepend + "***." + vf.name()).c_str()
351 ensightFile.write(pTraits<Type>::typeName);
352 ensightFile.writePartHeader(1);
358 map(vf, hexes, wedges),
365 Field<Type>(vf, prisms),
372 Field<Type>(vf, pyrs),
379 Field<Type>(vf, tets),
386 Field<Type>(vf, polys),
391 label ensightPatchI = eMesh.patchPartOffset();
393 forAll(allPatchNames, patchi)
395 const word& patchName = allPatchNames[patchi];
399 if (patchNames.empty() || patchNames.found(patchName))
405 vf.boundaryField()[patchi],
408 boundaryFaceSets[patchi],
409 nPatchPrims.find(patchName)(),
419 // write faceZones, if requested
420 if (faceZoneNames.size())
422 // Interpolates cell values to faces - needed only when exporting
424 GeometricField<Type, fvsPatchField, surfaceMesh> sf
426 linearInterpolate(vf)
429 forAllConstIter(wordHashSet, faceZoneNames, iter)
431 const word& faceZoneName = iter.key();
435 label zoneID = mesh.faceZones().findZoneID(faceZoneName);
437 const faceZone& fz = mesh.faceZones()[zoneID];
439 // Prepare data to write
443 if (eMesh.faceToBeIncluded(fz[i]))
449 Field<Type> values(nIncluded);
451 // Loop on the faceZone and store the needed field values
456 if (mesh.isInternalFace(faceI))
458 values[j] = sf[faceI];
463 if (eMesh.faceToBeIncluded(faceI))
465 label patchI = mesh.boundaryMesh().whichPatch(faceI);
466 const polyPatch& pp = mesh.boundaryMesh()[patchI];
467 label patchFaceI = pp.whichFace(faceI);
468 Type value = sf.boundaryField()[patchI][patchFaceI];
482 faceZoneFaceSets[zoneID],
483 nFaceZonePrims.find(faceZoneName)(),
492 if (Pstream::master())
494 delete ensightFilePtr;
500 void ensightPointField
502 const GeometricField<Type, pointPatchField, pointMesh>& pf,
503 const ensightMesh& eMesh,
504 const fileName& postProcPath,
506 const label timeIndex,
508 Ostream& ensightCaseFile
511 Info<< "Converting field " << pf.name() << endl;
513 word timeFile = prepend + itoa(timeIndex);
515 const fvMesh& mesh = eMesh.mesh();
516 const wordList& allPatchNames = eMesh.allPatchNames();
517 const wordHashSet& patchNames = eMesh.patchNames();
518 const wordHashSet& faceZoneNames = eMesh.faceZoneNames();
521 ensightStream* ensightFilePtr = NULL;
522 if (Pstream::master())
524 // set the filename of the ensight file
525 fileName ensightFileName(timeFile + "." + pf.name());
529 ensightFilePtr = new ensightBinaryStream
531 postProcPath/ensightFileName,
537 ensightFilePtr = new ensightAsciiStream
539 postProcPath/ensightFileName,
545 ensightStream& ensightFile = *ensightFilePtr;
547 if (eMesh.patchNames().empty())
551 if (Pstream::master())
555 ensightCaseFile.setf(ios_base::left);
558 << pTraits<Type>::typeName
560 << setw(15) << pf.name()
561 << (' ' + prepend + "***." + pf.name()).c_str()
565 ensightFile.write(pTraits<Type>::typeName);
566 ensightFile.writePartHeader(1);
572 Field<Type>(pf.internalField(), eMesh.uniquePointMap()),
578 label ensightPatchI = eMesh.patchPartOffset();
580 forAll(allPatchNames, patchi)
582 const word& patchName = allPatchNames[patchi];
586 if (patchNames.empty() || patchNames.found(patchName))
588 const fvPatch& p = mesh.boundary()[patchi];
591 returnReduce(p.size(), sumOp<label>())
595 // Renumber the patch points/faces into unique points
596 labelList pointToGlobal;
597 labelList uniqueMeshPointLabels;
598 autoPtr<globalIndex> globalPointsPtr =
599 mesh.globalData().mergePoints
601 p.patch().meshPoints(),
602 p.patch().meshPointMap(),
604 uniqueMeshPointLabels
607 if (Pstream::master())
609 ensightFile.writePartHeader(ensightPatchI);
615 Field<Type>(pf.internalField(), uniqueMeshPointLabels),
624 // write faceZones, if requested
625 if (faceZoneNames.size())
627 forAllConstIter(wordHashSet, faceZoneNames, iter)
629 const word& faceZoneName = iter.key();
633 label zoneID = mesh.faceZones().findZoneID(faceZoneName);
635 const faceZone& fz = mesh.faceZones()[zoneID];
637 if (returnReduce(fz().nPoints(), sumOp<label>()) > 0)
639 // Renumber the faceZone points/faces into unique points
640 labelList pointToGlobal;
641 labelList uniqueMeshPointLabels;
642 autoPtr<globalIndex> globalPointsPtr =
643 mesh.globalData().mergePoints
648 uniqueMeshPointLabels
651 if (Pstream::master())
653 ensightFile.writePartHeader(ensightPatchI);
662 uniqueMeshPointLabels
672 if (Pstream::master())
674 delete ensightFilePtr;
682 const IOobject& fieldObject,
683 const ensightMesh& eMesh,
684 const fileName& postProcPath,
686 const label timeIndex,
688 const bool nodeValues,
689 Ostream& ensightCaseFile
693 GeometricField<Type, fvPatchField, volMesh> vf(fieldObject, eMesh.mesh());
697 tmp<GeometricField<Type, pointPatchField, pointMesh> > pfld
699 volPointInterpolation::New(eMesh.mesh()).interpolate(vf)
701 pfld().rename(vf.name());
703 ensightPointField<Type>
730 // ************************************************************************* //