FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 <finiteVolume/fvMesh.H>
28 #include <finiteVolume/volFields.H>
29 #include <OpenFOAM/OFstream.H>
30 #include <OpenFOAM/IOmanip.H>
31 #include "itoa.H"
32 #include "ensightWriteBinary.H"
33 
34 using namespace Foam;
35 
36 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
37 
38 void writeData(const scalarField& sf, OFstream& ensightFile)
39 {
40  forAll(sf, i)
41  {
42  if (mag( sf[i] ) >= scalar(floatScalarVSMALL))
43  {
44  ensightFile << setw(12) << sf[i] << nl;
45  }
46  else
47  {
48  ensightFile << setw(12) << scalar(0) << nl;
49  }
50  }
51 }
52 
53 
54 template<class Type>
55 scalarField map
56 (
57  const Field<Type>& vf,
58  const labelList& map,
59  const label cmpt
60 )
61 {
62  scalarField mf(map.size());
63 
64  forAll(map, i)
65  {
66  mf[i] = component(vf[map[i]], cmpt);
67  }
68 
69  return mf;
70 }
71 
72 
73 template<class Type>
74 scalarField map
75 (
76  const Field<Type>& vf,
77  const labelList& map1,
78  const labelList& map2,
79  const label cmpt
80 )
81 {
82  scalarField mf(map1.size() + map2.size());
83 
84  forAll(map1, i)
85  {
86  mf[i] = component(vf[map1[i]], cmpt);
87  }
88 
89  label offset = map1.size();
90 
91  forAll(map2, i)
92  {
93  mf[i + offset] = component(vf[map2[i]], cmpt);
94  }
95 
96  return mf;
97 }
98 
99 
100 template<class Type>
101 void writeAllData
102 (
103  const char* key,
104  const Field<Type>& vf,
105  const labelList& prims,
106  const label nPrims,
107  OFstream& ensightFile
108 )
109 {
110  if (nPrims)
111  {
112  if (Pstream::master())
113  {
114  ensightFile << key << nl;
115 
116  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
117  {
118  writeData(map(vf, prims, cmpt), ensightFile);
119 
120  for (int slave=1; slave<Pstream::nProcs(); slave++)
121  {
122  IPstream fromSlave(Pstream::scheduled, slave);
123  scalarField data(fromSlave);
124  writeData(data, ensightFile);
125  }
126  }
127  }
128  else
129  {
130  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
131  {
133  toMaster<< map(vf, prims, cmpt);
134  }
135  }
136  }
137 }
138 
139 
140 template<class Type>
141 void writeAllDataBinary
142 (
143  const char* key,
144  const Field<Type>& vf,
145  const labelList& prims,
146  const label nPrims,
147  std::ofstream& ensightFile
148 )
149 {
150  if (nPrims)
151  {
152  if (Pstream::master())
153  {
154  writeEnsDataBinary(key,ensightFile);
155 
156  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
157  {
158  writeEnsDataBinary(map(vf, prims, cmpt), ensightFile);
159 
160  for (int slave=1; slave<Pstream::nProcs(); slave++)
161  {
162  IPstream fromSlave(Pstream::scheduled, slave);
163  scalarField data(fromSlave);
164  writeEnsDataBinary(data, ensightFile);
165  }
166  }
167  }
168  else
169  {
170  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
171  {
173  toMaster<< map(vf, prims, cmpt);
174  }
175  }
176  }
177 }
178 
179 
180 template<class Type>
181 void writeAllFaceData
182 (
183  const char* key,
184  const labelList& prims,
185  const label nPrims,
186  const Field<Type>& pf,
187  const labelList& patchProcessors,
188  OFstream& ensightFile
189 )
190 {
191  if (nPrims)
192  {
193  if (Pstream::master())
194  {
195  ensightFile << key << nl;
196 
197  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
198  {
199  writeData(map(pf, prims, cmpt), ensightFile);
200 
201  forAll (patchProcessors, i)
202  {
203  if (patchProcessors[i] != 0)
204  {
205  label slave = patchProcessors[i];
206  IPstream fromSlave(Pstream::scheduled, slave);
207  scalarField pf(fromSlave);
208 
209  writeData(pf, ensightFile);
210  }
211  }
212  }
213  }
214  else
215  {
216  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
217  {
219  toMaster<< map(pf, prims, cmpt);
220  }
221  }
222  }
223 }
224 
225 
226 template<class Type>
227 void writeAllFaceDataBinary
228 (
229  const char* key,
230  const labelList& prims,
231  const label nPrims,
232  const Field<Type>& pf,
233  const labelList& patchProcessors,
234  std::ofstream& ensightFile
235 )
236 {
237  if (nPrims)
238  {
239  if (Pstream::master())
240  {
241  writeEnsDataBinary(key,ensightFile);
242 
243  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
244  {
245  writeEnsDataBinary(map(pf, prims, cmpt), ensightFile);
246 
247  forAll (patchProcessors, i)
248  {
249  if (patchProcessors[i] != 0)
250  {
251  label slave = patchProcessors[i];
252  IPstream fromSlave(Pstream::scheduled, slave);
253  scalarField pf(fromSlave);
254 
255  writeEnsDataBinary(pf, ensightFile);
256  }
257  }
258  }
259  }
260  else
261  {
262  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
263  {
265  toMaster<< map(pf, prims, cmpt);
266  }
267  }
268  }
269 }
270 
271 
272 template<class Type>
273 bool writePatchField
274 (
275  const Foam::Field<Type>& pf,
276  const Foam::label patchi,
277  const Foam::label ensightPatchI,
278  const Foam::faceSets& boundaryFaceSet,
280  const Foam::labelList& patchProcessors,
281  Foam::OFstream& ensightFile
282 )
283 {
284  if (nfp.nTris || nfp.nQuads || nfp.nPolys)
285  {
286  if (Pstream::master())
287  {
288  ensightFile
289  << "part" << nl
290  << setw(10) << ensightPatchI << nl;
291  }
292 
293  writeAllFaceData
294  (
295  "tria3",
296  boundaryFaceSet.tris,
297  nfp.nTris,
298  pf,
299  patchProcessors,
300  ensightFile
301  );
302 
303  writeAllFaceData
304  (
305  "quad4",
306  boundaryFaceSet.quads,
307  nfp.nQuads,
308  pf,
309  patchProcessors,
310  ensightFile
311  );
312 
313  writeAllFaceData
314  (
315  "nsided",
316  boundaryFaceSet.polys,
317  nfp.nPolys,
318  pf,
319  patchProcessors,
320  ensightFile
321  );
322 
323  return true;
324  }
325  else
326  {
327  return false;
328  }
329 }
330 
331 
332 template<class Type>
333 bool writePatchFieldBinary
334 (
335  const Foam::Field<Type>& pf,
336  const Foam::label patchi,
337  const Foam::label ensightPatchI,
338  const Foam::faceSets& boundaryFaceSet,
340  const Foam::labelList& patchProcessors,
341  std::ofstream& ensightFile
342 )
343 {
344  if (nfp.nTris || nfp.nQuads || nfp.nPolys)
345  {
346  if (Pstream::master())
347  {
348  writeEnsDataBinary("part",ensightFile);
349  writeEnsDataBinary(ensightPatchI,ensightFile);
350  }
351 
352  writeAllFaceDataBinary
353  (
354  "tria3",
355  boundaryFaceSet.tris,
356  nfp.nTris,
357  pf,
358  patchProcessors,
359  ensightFile
360  );
361 
362  writeAllFaceDataBinary
363  (
364  "quad4",
365  boundaryFaceSet.quads,
366  nfp.nQuads,
367  pf,
368  patchProcessors,
369  ensightFile
370  );
371 
372  writeAllFaceDataBinary
373  (
374  "nsided",
375  boundaryFaceSet.polys,
376  nfp.nPolys,
377  pf,
378  patchProcessors,
379  ensightFile
380  );
381 
382  return true;
383  }
384  else
385  {
386  return false;
387  }
388 }
389 
390 
391 template<class Type>
392 void writePatchField
393 (
394  const Foam::word& fieldName,
395  const Foam::Field<Type>& pf,
396  const Foam::word& patchName,
397  const Foam::ensightMesh& eMesh,
398  const Foam::fileName& postProcPath,
399  const Foam::word& prepend,
400  const Foam::label timeIndex,
401  Foam::Ostream& ensightCaseFile
402 )
403 {
404  const Time& runTime = eMesh.mesh().time();
405 
406  const List<faceSets>& boundaryFaceSets = eMesh.boundaryFaceSets();
407  const wordList& allPatchNames = eMesh.allPatchNames();
408  const List<labelList>& allPatchProcs = eMesh.allPatchProcs();
410  nPatchPrims = eMesh.nPatchPrims();
411 
412  label ensightPatchI = eMesh.patchPartOffset();
413 
414  label patchi = -1;
415 
416  forAll(allPatchNames, i)
417  {
418  if (allPatchNames[i] == patchName)
419  {
420  patchi = i;
421  break;
422  }
423  ensightPatchI++;
424  }
425 
426 
427  const labelList& patchProcessors = allPatchProcs[patchi];
428 
429  word pfName = patchName + '.' + fieldName;
430 
431  word timeFile = prepend + itoa(timeIndex);
432 
433  OFstream *ensightFilePtr = NULL;
434  if (Pstream::master())
435  {
436  if (timeIndex == 0)
437  {
438  ensightCaseFile.setf(ios_base::left);
439 
440  ensightCaseFile
442  << " per element: 1 "
443  << setw(15) << pfName
444  << (' ' + prepend + "***." + pfName).c_str()
445  << nl;
446  }
447 
448  // set the filename of the ensight file
449  fileName ensightFileName(timeFile + "." + pfName);
450  ensightFilePtr = new OFstream
451  (
452  postProcPath/ensightFileName,
453  runTime.writeFormat(),
454  runTime.writeVersion(),
455  runTime.writeCompression()
456  );
457  }
458 
459  OFstream& ensightFile = *ensightFilePtr;
460 
461  if (Pstream::master())
462  {
464  }
465 
466  if (patchi >= 0)
467  {
469  (
470  pf,
471  patchi,
472  ensightPatchI,
473  boundaryFaceSets[patchi],
474  nPatchPrims.find(patchName)(),
475  patchProcessors,
476  ensightFile
477  );
478  }
479  else
480  {
481  faceSets nullFaceSets;
482 
484  (
485  Field<Type>(),
486  -1,
487  ensightPatchI,
488  nullFaceSets,
489  nPatchPrims.find(patchName)(),
490  patchProcessors,
491  ensightFile
492  );
493  }
494 
495  if (Pstream::master())
496  {
497  delete ensightFilePtr;
498  }
499 }
500 
501 
502 template<class Type>
503 void ensightFieldAscii
504 (
506  const Foam::ensightMesh& eMesh,
507  const Foam::fileName& postProcPath,
508  const Foam::word& prepend,
509  const Foam::label timeIndex,
510  Foam::Ostream& ensightCaseFile
511 )
512 {
513  Info<< "Converting field " << fieldObject.name() << endl;
514 
515  word timeFile = prepend + itoa(timeIndex);
516 
517  const fvMesh& mesh = eMesh.mesh();
518  const Time& runTime = mesh.time();
519 
520  const cellSets& meshCellSets = eMesh.meshCellSets();
521  const List<faceSets>& boundaryFaceSets = eMesh.boundaryFaceSets();
522  const wordList& allPatchNames = eMesh.allPatchNames();
523  const List<labelList>& allPatchProcs = eMesh.allPatchProcs();
524  const wordHashSet& patchNames = eMesh.patchNames();
526  nPatchPrims = eMesh.nPatchPrims();
527 
528  const labelList& tets = meshCellSets.tets;
529  const labelList& pyrs = meshCellSets.pyrs;
530  const labelList& prisms = meshCellSets.prisms;
531  const labelList& wedges = meshCellSets.wedges;
532  const labelList& hexes = meshCellSets.hexes;
533  const labelList& polys = meshCellSets.polys;
534 
535  OFstream *ensightFilePtr = NULL;
536  if (Pstream::master())
537  {
538  // set the filename of the ensight file
539  fileName ensightFileName(timeFile + "." + fieldObject.name());
540  ensightFilePtr = new OFstream
541  (
542  postProcPath/ensightFileName,
543  runTime.writeFormat(),
544  runTime.writeVersion(),
545  runTime.writeCompression()
546  );
547  }
548 
549  OFstream& ensightFile = *ensightFilePtr;
550 
551  GeometricField<Type, fvPatchField, volMesh> vf(fieldObject, mesh);
552 
553  if (patchNames.empty())
554  {
555  if (Pstream::master())
556  {
557  if (timeIndex == 0)
558  {
559  ensightCaseFile.setf(ios_base::left);
560 
561  ensightCaseFile
563  << " per element: 1 "
564  << setw(15) << vf.name()
565  << (' ' + prepend + "***." + vf.name()).c_str()
566  << nl;
567  }
568 
569  ensightFile
570  << pTraits<Type>::typeName << nl
571  << "part" << nl
572  << setw(10) << 1 << nl;
573 
574  ensightFile.setf(ios_base::scientific, ios_base::floatfield);
575  ensightFile.precision(5);
576  }
577 
578  if (meshCellSets.nHexesWedges)
579  {
580  if (Pstream::master())
581  {
582  ensightFile << "hexa8" << nl;
583 
584  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
585  {
586  writeData
587  (
588  map(vf, hexes, wedges, cmpt),
589  ensightFile
590  );
591 
592  for (int slave=1; slave<Pstream::nProcs(); slave++)
593  {
594  IPstream fromSlave(Pstream::scheduled, slave);
595  scalarField data(fromSlave);
596  writeData(data, ensightFile);
597  }
598  }
599  }
600  else
601  {
602  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
603  {
605  toMaster<< map(vf, hexes, wedges, cmpt);
606  }
607  }
608  }
609 
610  writeAllData("penta6", vf, prisms, meshCellSets.nPrisms, ensightFile);
611  writeAllData("pyramid5", vf, pyrs, meshCellSets.nPyrs, ensightFile);
612  writeAllData("tetra4", vf, tets, meshCellSets.nTets, ensightFile);
613  writeAllData("nfaced", vf, polys, meshCellSets.nPolys, ensightFile);
614  }
615 
616  label ensightPatchI = eMesh.patchPartOffset();
617 
618  forAll(allPatchNames, patchi)
619  {
620  const word& patchName = allPatchNames[patchi];
621  const labelList& patchProcessors = allPatchProcs[patchi];
622 
623  if (patchNames.empty() || patchNames.found(patchName))
624  {
625  if (mesh.boundary()[patchi].size())
626  {
627  if
628  (
630  (
631  vf.boundaryField()[patchi],
632  patchi,
633  ensightPatchI,
634  boundaryFaceSets[patchi],
635  nPatchPrims.find(patchName)(),
636  patchProcessors,
637  ensightFile
638  )
639  )
640  {
641  ensightPatchI++;
642  }
643 
644  }
645  else if (Pstream::master())
646  {
647  faceSets nullFaceSet;
648 
649  if
650  (
652  (
653  Field<Type>(),
654  -1,
655  ensightPatchI,
656  nullFaceSet,
657  nPatchPrims.find(patchName)(),
658  patchProcessors,
659  ensightFile
660  )
661  )
662  {
663  ensightPatchI++;
664  }
665  }
666  }
667  }
668 
669  if (Pstream::master())
670  {
671  delete ensightFilePtr;
672  }
673 }
674 
675 
676 template<class Type>
677 void ensightFieldBinary
678 (
679  const Foam::IOobject& fieldObject,
680  const Foam::ensightMesh& eMesh,
681  const Foam::fileName& postProcPath,
682  const Foam::word& prepend,
683  const Foam::label timeIndex,
684  Foam::Ostream& ensightCaseFile
685 )
686 {
687  Info<< "Converting field (binary) " << fieldObject.name() << endl;
688 
689  word timeFile = prepend + itoa(timeIndex);
690 
691  const fvMesh& mesh = eMesh.mesh();
692  //const Time& runTime = mesh.time();
693 
694  const cellSets& meshCellSets = eMesh.meshCellSets();
695  const List<faceSets>& boundaryFaceSets = eMesh.boundaryFaceSets();
696  const wordList& allPatchNames = eMesh.allPatchNames();
697  const List<labelList>& allPatchProcs = eMesh.allPatchProcs();
698  const wordHashSet& patchNames = eMesh.patchNames();
700  nPatchPrims = eMesh.nPatchPrims();
701 
702  const labelList& tets = meshCellSets.tets;
703  const labelList& pyrs = meshCellSets.pyrs;
704  const labelList& prisms = meshCellSets.prisms;
705  const labelList& wedges = meshCellSets.wedges;
706  const labelList& hexes = meshCellSets.hexes;
707  const labelList& polys = meshCellSets.polys;
708 
709  std::ofstream *ensightFilePtr = NULL;
710  if (Pstream::master())
711  {
712  // set the filename of the ensight file
713  fileName ensightFileName(timeFile + "." + fieldObject.name());
714  ensightFilePtr = new std::ofstream
715  (
716  (postProcPath/ensightFileName).c_str(),
717  ios_base::out | ios_base::binary | ios_base::trunc
718  );
719  // Check on file opened?
720  }
721 
722  std::ofstream& ensightFile = *ensightFilePtr;
723 
724  GeometricField<Type, fvPatchField, volMesh> vf(fieldObject, mesh);
725 
726  if (patchNames.empty())
727  {
728  if (Pstream::master())
729  {
730  if (timeIndex == 0)
731  {
732  ensightCaseFile.setf(ios_base::left);
733 
734  ensightCaseFile
736  << " per element: 1 "
737  << setw(15) << vf.name()
738  << (' ' + prepend + "***." + vf.name()).c_str()
739  << nl;
740  }
741 
743  writeEnsDataBinary("part",ensightFile);
744  writeEnsDataBinary(1,ensightFile);
745  }
746 
747  if (meshCellSets.nHexesWedges)
748  {
749  if (Pstream::master())
750  {
751  writeEnsDataBinary("hexa8",ensightFile);
752 
753  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
754  {
756  (
757  map(vf, hexes, wedges, cmpt),
758  ensightFile
759  );
760 
761  for (int slave=1; slave<Pstream::nProcs(); slave++)
762  {
763  IPstream fromSlave(Pstream::scheduled, slave);
764  scalarField data(fromSlave);
765  writeEnsDataBinary(data, ensightFile);
766  }
767  }
768  }
769  else
770  {
771  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
772  {
774  toMaster<< map(vf, hexes, wedges, cmpt);
775  }
776  }
777  }
778 
779  writeAllDataBinary
780  (
781  "penta6",
782  vf,
783  prisms,
784  meshCellSets.nPrisms,
785  ensightFile
786  );
787 
788  writeAllDataBinary
789  (
790  "pyramid5",
791  vf,
792  pyrs,
793  meshCellSets.nPyrs,
794  ensightFile
795  );
796 
797  writeAllDataBinary
798  (
799  "tetra4",
800  vf,
801  tets,
802  meshCellSets.nTets,
803  ensightFile
804  );
805 
806  writeAllDataBinary
807  (
808  "nfaced",
809  vf,
810  polys,
811  meshCellSets.nPolys,
812  ensightFile
813  );
814  }
815 
816  label ensightPatchI = eMesh.patchPartOffset();
817 
818  forAll(allPatchNames, patchi)
819  {
820  const word& patchName = allPatchNames[patchi];
821  const labelList& patchProcessors = allPatchProcs[patchi];
822 
823  if (patchNames.empty() || patchNames.found(patchName))
824  {
825  if (mesh.boundary()[patchi].size())
826  {
827  if
828  (
829  writePatchFieldBinary
830  (
831  vf.boundaryField()[patchi],
832  patchi,
833  ensightPatchI,
834  boundaryFaceSets[patchi],
835  nPatchPrims.find(patchName)(),
836  patchProcessors,
837  ensightFile
838  )
839  )
840  {
841  ensightPatchI++;
842  }
843 
844  }
845  else if (Pstream::master())
846  {
847  faceSets nullFaceSet;
848 
849  if
850  (
851  writePatchFieldBinary
852  (
853  Field<Type>(),
854  -1,
855  ensightPatchI,
856  nullFaceSet,
857  nPatchPrims.find(patchName)(),
858  patchProcessors,
859  ensightFile
860  )
861  )
862  {
863  ensightPatchI++;
864  }
865  }
866  }
867  }
868 
869  if (Pstream::master())
870  {
871  ensightFile.close();
872  }
873 }
874 
875 
876 template<class Type>
877 void ensightField
878 (
879  const Foam::IOobject& fieldObject,
880  const Foam::ensightMesh& eMesh,
881  const Foam::fileName& postProcPath,
882  const Foam::word& prepend,
883  const Foam::label timeIndex,
884  const bool binary,
885  Foam::Ostream& ensightCaseFile
886 )
887 {
888  if (binary)
889  {
890  ensightFieldBinary<Type>
891  (
892  fieldObject,
893  eMesh,
894  postProcPath,
895  prepend,
896  timeIndex,
897  ensightCaseFile
898  );
899  }
900  else
901  {
902  ensightFieldAscii<Type>
903  (
904  fieldObject,
905  eMesh,
906  postProcPath,
907  prepend,
908  timeIndex,
909  ensightCaseFile
910  );
911  }
912 }
913 
914 
915 // ************************ vim: set sw=4 sts=4 et: ************************ //