foamToVTK.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 Application
25  foamToVTK
26 
27 Description
28  Legacy VTK file format writer.
29 
30  - Handles volFields, pointFields, surfaceScalarField, surfaceVectorField
31  fields.
32  - Mesh topo changes.
33  - Both ascii and binary.
34  - Single time step writing.
35  - Write subset only.
36  - Automatic decomposition of cells; polygons on boundary undecomposed since
37  handled by vtk.
38 
39 Usage
40  \b foamToVTK [OPTION]
41 
42  Options:
43  - \par -ascii
44  Write VTK data in ASCII format instead of binary.
45 
46  - \par -mesh <name>
47  Use a different mesh name (instead of -region)
48 
49  - \par -fields <fields>
50  Convert selected fields only. For example,
51  \verbatim
52  -fields "( p T U )"
53  \endverbatim
54  The quoting is required to avoid shell expansions and to pass the
55  information as a single argument.
56 
57  - \par -surfaceFields
58  Write surfaceScalarFields (e.g., phi)
59 
60  - \par -cellSet <name>
61  - \par -faceSet <name>
62 
63  - \par -pointSet <name>
64  Restrict conversion to the cellSet, faceSet or pointSet.
65 
66  - \par -nearCellValue
67  Output cell value on patches instead of patch value itself
68 
69  - \par -noInternal
70  Do not generate file for mesh, only for patches
71 
72  - \par -noPointValues
73  No pointFields
74 
75  - \par -noFaceZones
76  No faceZones
77 
78  - \par -noLinks
79  (in parallel) do not link processor files to master
80 
81  - \par poly
82  write polyhedral cells without tet/pyramid decomposition
83 
84  - \par -allPatches
85  Combine all patches into a single file
86 
87  - \par -excludePatches <patchNames>
88  Specify patches (wildcards) to exclude. For example,
89  \verbatim
90  -excludePatches '( inlet_1 inlet_2 "proc.*")'
91  \endverbatim
92  The quoting is required to avoid shell expansions and to pass the
93  information as a single argument. The double quotes denote a regular
94  expression.
95 
96  - \par -useTimeName
97  use the time index in the VTK file name instead of the time index
98 
99 Note
100  mesh subset is handled by vtkMesh. Slight inconsistency in
101  interpolation: on the internal field it interpolates the whole volField
102  to the whole-mesh pointField and then selects only those values it
103  needs for the subMesh (using the fvMeshSubset cellMap(), pointMap()
104  functions). For the patches however it uses the
105  fvMeshSubset.interpolate function to directly interpolate the
106  whole-mesh values onto the subset patch.
107 
108 Note
109  \par new file format:
110  no automatic timestep recognition.
111  However can have .pvd file format which refers to time simulation
112  if XML *.vtu files are available:
113 
114  \verbatim
115  <?xml version="1.0"?>
116  <VTKFile type="Collection"
117  version="0.1"
118  byte_order="LittleEndian"
119  compressor="vtkZLibDataCompressor">
120  <Collection>
121  <DataSet timestep="50" file="pitzDaily_2.vtu"/>
122  <DataSet timestep="100" file="pitzDaily_3.vtu"/>
123  <DataSet timestep="150" file="pitzDaily_4.vtu"/>
124  <DataSet timestep="200" file="pitzDaily_5.vtu"/>
125  <DataSet timestep="250" file="pitzDaily_6.vtu"/>
126  <DataSet timestep="300" file="pitzDaily_7.vtu"/>
127  <DataSet timestep="350" file="pitzDaily_8.vtu"/>
128  <DataSet timestep="400" file="pitzDaily_9.vtu"/>
129  <DataSet timestep="450" file="pitzDaily_10.vtu"/>
130  <DataSet timestep="500" file="pitzDaily_11.vtu"/>
131  </Collection>
132  </VTKFile>
133  \endverbatim
134 
135 \*---------------------------------------------------------------------------*/
136 
137 #include "fvCFD.H"
138 #include "pointMesh.H"
139 #include "volPointInterpolation.H"
140 #include "emptyPolyPatch.H"
141 #include "labelIOField.H"
142 #include "scalarIOField.H"
143 #include "sphericalTensorIOField.H"
144 #include "symmTensorIOField.H"
145 #include "tensorIOField.H"
146 #include "faceZoneMesh.H"
147 #include "Cloud.H"
148 #include "passiveParticle.H"
149 #include "stringListOps.H"
150 
151 #include "vtkMesh.H"
152 #include "readFields.H"
153 #include "writeFuns.H"
154 
155 #include "internalWriter.H"
156 #include "patchWriter.H"
157 #include "lagrangianWriter.H"
158 
159 #include "writeFaceSet.H"
160 #include "writePointSet.H"
161 #include "surfaceMeshWriter.H"
162 #include "writeSurfFields.H"
163 
164 
165 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
166 
167 template<class GeoField>
168 void print(const char* msg, Ostream& os, const PtrList<const GeoField>& flds)
169 {
170  if (flds.size())
171  {
172  os << msg;
173  forAll(flds, i)
174  {
175  os << ' ' << flds[i].name();
176  }
177  os << endl;
178  }
179 }
180 
181 
182 void print(Ostream& os, const wordList& flds)
183 {
184  forAll(flds, i)
185  {
186  os << ' ' << flds[i];
187  }
188  os << endl;
189 }
190 
191 
192 labelList getSelectedPatches
193 (
194  const polyBoundaryMesh& patches,
195  const List<wordRe>& excludePatches
196 )
197 {
198  DynamicList<label> patchIDs(patches.size());
199 
200  Info<< "Combining patches:" << endl;
201 
202  forAll(patches, patchi)
203  {
204  const polyPatch& pp = patches[patchi];
205 
206  if
207  (
208  isType<emptyPolyPatch>(pp)
209  || (Pstream::parRun() && isType<processorPolyPatch>(pp))
210  )
211  {
212  Info<< " discarding empty/processor patch " << patchi
213  << " " << pp.name() << endl;
214  }
215  else if (findStrings(excludePatches, pp.name()))
216  {
217  Info<< " excluding patch " << patchi
218  << " " << pp.name() << endl;
219  }
220  else
221  {
222  patchIDs.append(patchi);
223  Info<< " patch " << patchi << " " << pp.name() << endl;
224  }
225  }
226  return patchIDs.shrink();
227 }
228 
229 
230 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
231 
232 int main(int argc, char *argv[])
233 {
234  argList::addNote
235  (
236  "legacy VTK file format writer"
237  );
238  timeSelector::addOptions();
239 
240  #include "addRegionOption.H"
241  argList::addOption
242  (
243  "fields",
244  "wordList",
245  "only convert the specified fields - eg '(p T U)'"
246  );
247  argList::addOption
248  (
249  "cellSet",
250  "name",
251  "convert a mesh subset corresponding to the specified cellSet"
252  );
253  argList::addOption
254  (
255  "faceSet",
256  "name",
257  "restrict conversion to the specified faceSet"
258  );
259  argList::addOption
260  (
261  "pointSet",
262  "name",
263  "restrict conversion to the specified pointSet"
264  );
265  argList::addBoolOption
266  (
267  "ascii",
268  "write in ASCII format instead of binary"
269  );
270  argList::addBoolOption
271  (
272  "poly",
273  "write polyhedral cells without tet/pyramid decomposition"
274  );
275  argList::addBoolOption
276  (
277  "surfaceFields",
278  "write surfaceScalarFields (e.g., phi)"
279  );
280  argList::addBoolOption
281  (
282  "nearCellValue",
283  "use cell value on patches instead of patch value itself"
284  );
285  argList::addBoolOption
286  (
287  "noInternal",
288  "do not generate file for mesh, only for patches"
289  );
290  argList::addBoolOption
291  (
292  "noPointValues",
293  "no pointFields"
294  );
295  argList::addBoolOption
296  (
297  "allPatches",
298  "combine all patches into a single file"
299  );
300  argList::addOption
301  (
302  "excludePatches",
303  "wordReList",
304  "a list of patches to exclude - eg '( inlet \".*Wall\" )' "
305  );
306  argList::addBoolOption
307  (
308  "noFaceZones",
309  "no faceZones"
310  );
311  argList::addBoolOption
312  (
313  "noLinks",
314  "don't link processor VTK files - parallel only"
315  );
316  argList::addBoolOption
317  (
318  "useTimeName",
319  "use the time name instead of the time index when naming the files"
320  );
321 
322  #include "setRootCase.H"
323  #include "createTime.H"
324 
325  const bool doWriteInternal = !args.optionFound("noInternal");
326  const bool doFaceZones = !args.optionFound("noFaceZones");
327  const bool doLinks = !args.optionFound("noLinks");
328  bool binary = !args.optionFound("ascii");
329  const bool useTimeName = args.optionFound("useTimeName");
330 
331  // Decomposition of polyhedral cells into tets/pyramids cells
332  vtkTopo::decomposePoly = !args.optionFound("poly");
333 
334  if (binary && (sizeof(floatScalar) != 4 || sizeof(label) != 4))
335  {
337  << "Using ASCII rather than binary VTK format because "
338  "floatScalar and/or label are not 4 bytes in size."
339  << nl << endl;
340  binary = false;
341  }
342 
343  const bool nearCellValue = args.optionFound("nearCellValue");
344 
345  if (nearCellValue)
346  {
348  << "Using neighbouring cell value instead of patch value"
349  << nl << endl;
350  }
351 
352  const bool noPointValues = args.optionFound("noPointValues");
353 
354  if (noPointValues)
355  {
357  << "Outputting cell values only" << nl << endl;
358  }
359 
360  const bool allPatches = args.optionFound("allPatches");
361 
362  List<wordRe> excludePatches;
363  if (args.optionFound("excludePatches"))
364  {
365  args.optionLookup("excludePatches")() >> excludePatches;
366 
367  Info<< "Not including patches " << excludePatches << nl << endl;
368  }
369 
370  word cellSetName;
371  word faceSetName;
372  word pointSetName;
373  string vtkName = runTime.caseName();
374 
375  if (args.optionReadIfPresent("cellSet", cellSetName))
376  {
377  vtkName = cellSetName;
378  }
379  else if (Pstream::parRun())
380  {
381  // Strip off leading casename, leaving just processor_DDD ending.
382  string::size_type i = vtkName.rfind("processor");
383 
384  if (i != string::npos)
385  {
386  vtkName = vtkName.substr(i);
387  }
388  }
389  args.optionReadIfPresent("faceSet", faceSetName);
390  args.optionReadIfPresent("pointSet", pointSetName);
391 
392 
393 
394  instantList timeDirs = timeSelector::select0(runTime, args);
395 
396  #include "createNamedMesh.H"
397 
398  // VTK/ directory in the case
399  fileName fvPath(runTime.path()/"VTK");
400 
401  // Directory of mesh (region0 gets filtered out)
402  fileName regionPrefix = "";
403 
404  if (regionName != polyMesh::defaultRegion)
405  {
406  fvPath = fvPath/regionName;
407  regionPrefix = regionName;
408  }
409 
410  if (isDir(fvPath))
411  {
412  if
413  (
414  args.optionFound("time")
415  || args.optionFound("latestTime")
416  || cellSetName.size()
417  || faceSetName.size()
418  || pointSetName.size()
419  || regionName != polyMesh::defaultRegion
420  )
421  {
422  Info<< "Keeping old VTK files in " << fvPath << nl << endl;
423  }
424  else
425  {
426  Info<< "Deleting old VTK files in " << fvPath << nl << endl;
427 
428  rmDir(fvPath);
429  }
430  }
431 
432  mkDir(fvPath);
433 
434 
435  // Mesh wrapper; does subsetting and decomposition
436  vtkMesh vMesh(mesh, cellSetName);
437 
438 
439  // Scan for all possible lagrangian clouds
440  HashSet<fileName> allCloudDirs;
441  forAll(timeDirs, timeI)
442  {
443  runTime.setTime(timeDirs[timeI], timeI);
444  fileNameList cloudDirs
445  (
446  readDir
447  (
448  runTime.timePath()/regionPrefix/cloud::prefix,
449  fileType::directory
450  )
451  );
452  forAll(cloudDirs, i)
453  {
454  IOobjectList sprayObjs
455  (
456  mesh,
457  runTime.timeName(),
458  cloud::prefix/cloudDirs[i]
459  );
460 
461  IOobject* positionsPtr = sprayObjs.lookup(word("positions"));
462 
463  if (positionsPtr)
464  {
465  if (allCloudDirs.insert(cloudDirs[i]))
466  {
467  Info<< "At time: " << runTime.timeName()
468  << " detected cloud directory : " << cloudDirs[i]
469  << endl;
470  }
471  }
472  }
473  }
474 
475 
476  forAll(timeDirs, timeI)
477  {
478  runTime.setTime(timeDirs[timeI], timeI);
479 
480  Info<< "Time: " << runTime.timeName() << endl;
481 
482  word timeDesc =
483  useTimeName ? runTime.timeName() : Foam::name(runTime.timeIndex());
484 
485  // Check for new polyMesh/ and update mesh, fvMeshSubset and cell
486  // decomposition.
487  polyMesh::readUpdateState meshState = vMesh.readUpdate();
488 
489  const fvMesh& mesh = vMesh.mesh();
490 
491  if
492  (
493  meshState == polyMesh::TOPO_CHANGE
494  || meshState == polyMesh::TOPO_PATCH_CHANGE
495  )
496  {
497  Info<< " Read new mesh" << nl << endl;
498  }
499 
500 
501  // If faceSet: write faceSet only (as polydata)
502  if (faceSetName.size())
503  {
504  // Load the faceSet
505  faceSet set(mesh, faceSetName);
506 
507  // Filename as if patch with same name.
508  mkDir(fvPath/set.name());
509 
510  fileName patchFileName
511  (
512  fvPath/set.name()/set.name()
513  + "_"
514  + timeDesc
515  + ".vtk"
516  );
517 
518  Info<< " FaceSet : " << patchFileName << endl;
519 
520  writeFaceSet(binary, vMesh, set, patchFileName);
521 
522  continue;
523  }
524  // If pointSet: write pointSet only (as polydata)
525  if (pointSetName.size())
526  {
527  // Load the pointSet
528  pointSet set(mesh, pointSetName);
529 
530  // Filename as if patch with same name.
531  mkDir(fvPath/set.name());
532 
533  fileName patchFileName
534  (
535  fvPath/set.name()/set.name()
536  + "_"
537  + timeDesc
538  + ".vtk"
539  );
540 
541  Info<< " pointSet : " << patchFileName << endl;
542 
543  writePointSet(binary, vMesh, set, patchFileName);
544 
545  continue;
546  }
547 
548 
549  // Search for list of objects for this time
550  IOobjectList objects(mesh, runTime.timeName());
551 
552  HashSet<word> selectedFields;
553  bool specifiedFields = args.optionReadIfPresent
554  (
555  "fields",
556  selectedFields
557  );
558 
559  // Construct the vol fields (on the original mesh if subsetted)
560 
561  PtrList<const volScalarField> vsf;
562  PtrList<const volVectorField> vvf;
563  PtrList<const volSphericalTensorField> vSpheretf;
564  PtrList<const volSymmTensorField> vSymmtf;
565  PtrList<const volTensorField> vtf;
566 
567  if (!specifiedFields || selectedFields.size())
568  {
569  readFields(vMesh, vMesh.baseMesh(), objects, selectedFields, vsf);
570  print(" volScalarFields :", Info, vsf);
571 
572  readFields(vMesh, vMesh.baseMesh(), objects, selectedFields, vvf);
573  print(" volVectorFields :", Info, vvf);
574 
575  readFields
576  (
577  vMesh,
578  vMesh.baseMesh(),
579  objects,
580  selectedFields,
581  vSpheretf
582  );
583  print(" volSphericalTensorFields :", Info, vSpheretf);
584 
585  readFields
586  (
587  vMesh,
588  vMesh.baseMesh(),
589  objects,
590  selectedFields,
591  vSymmtf
592  );
593  print(" volSymmTensorFields :", Info, vSymmtf);
594 
595  readFields(vMesh, vMesh.baseMesh(), objects, selectedFields, vtf);
596  print(" volTensorFields :", Info, vtf);
597  }
598 
599  label nVolFields =
600  vsf.size()
601  + vvf.size()
602  + vSpheretf.size()
603  + vSymmtf.size()
604  + vtf.size();
605 
606 
607  // Construct pointMesh only if necessary since constructs edge
608  // addressing (expensive on polyhedral meshes)
609  if (noPointValues)
610  {
611  Info<< " pointScalarFields : switched off"
612  << " (\"-noPointValues\" (at your option)\n";
613  Info<< " pointVectorFields : switched off"
614  << " (\"-noPointValues\" (at your option)\n";
615  }
616 
617  PtrList<const pointScalarField> psf;
618  PtrList<const pointVectorField> pvf;
619  PtrList<const pointSphericalTensorField> pSpheretf;
620  PtrList<const pointSymmTensorField> pSymmtf;
621  PtrList<const pointTensorField> ptf;
622 
623  if (!noPointValues && !(specifiedFields && selectedFields.empty()))
624  {
625  readFields
626  (
627  vMesh,
628  pointMesh::New(vMesh.baseMesh()),
629  objects,
630  selectedFields,
631  psf
632  );
633  print(" pointScalarFields :", Info, psf);
634 
635  readFields
636  (
637  vMesh,
638  pointMesh::New(vMesh.baseMesh()),
639  objects,
640  selectedFields,
641  pvf
642  );
643  print(" pointVectorFields :", Info, pvf);
644 
645  readFields
646  (
647  vMesh,
648  pointMesh::New(vMesh.baseMesh()),
649  objects,
650  selectedFields,
651  pSpheretf
652  );
653  print(" pointSphericalTensorFields :", Info, pSpheretf);
654 
655  readFields
656  (
657  vMesh,
658  pointMesh::New(vMesh.baseMesh()),
659  objects,
660  selectedFields,
661  pSymmtf
662  );
663  print(" pointSymmTensorFields :", Info, pSymmtf);
664 
665  readFields
666  (
667  vMesh,
668  pointMesh::New(vMesh.baseMesh()),
669  objects,
670  selectedFields,
671  ptf
672  );
673  print(" pointTensorFields :", Info, ptf);
674  }
675  Info<< endl;
676 
677  label nPointFields =
678  psf.size()
679  + pvf.size()
680  + pSpheretf.size()
681  + pSymmtf.size()
682  + ptf.size();
683 
684  if (doWriteInternal)
685  {
686  // Create file and write header
687  fileName vtkFileName
688  (
689  fvPath/vtkName
690  + "_"
691  + timeDesc
692  + ".vtk"
693  );
694 
695  Info<< " Internal : " << vtkFileName << endl;
696 
697  // Write mesh
698  internalWriter writer(vMesh, binary, vtkFileName);
699 
700  // VolFields + cellID
701  writeFuns::writeCellDataHeader
702  (
703  writer.os(),
704  vMesh.nFieldCells(),
705  1+nVolFields
706  );
707 
708  // Write cellID field
709  writer.writeCellIDs();
710 
711  // Write volFields
712  writer.write(vsf);
713  writer.write(vvf);
714  writer.write(vSpheretf);
715  writer.write(vSymmtf);
716  writer.write(vtf);
717 
718  if (!noPointValues)
719  {
720  writeFuns::writePointDataHeader
721  (
722  writer.os(),
723  vMesh.nFieldPoints(),
724  nVolFields+nPointFields
725  );
726 
727  // pointFields
728  writer.write(psf);
729  writer.write(pvf);
730  writer.write(pSpheretf);
731  writer.write(pSymmtf);
732  writer.write(ptf);
733 
734  // Interpolated volFields
735  volPointInterpolation pInterp(mesh);
736  writer.write(pInterp, vsf);
737  writer.write(pInterp, vvf);
738  writer.write(pInterp, vSpheretf);
739  writer.write(pInterp, vSymmtf);
740  writer.write(pInterp, vtf);
741  }
742  }
743 
744  //---------------------------------------------------------------------
745  //
746  // Write surface fields
747  //
748  //---------------------------------------------------------------------
749 
750  if (args.optionFound("surfaceFields"))
751  {
752  PtrList<const surfaceScalarField> ssf;
753  readFields
754  (
755  vMesh,
756  vMesh.baseMesh(),
757  objects,
758  selectedFields,
759  ssf
760  );
761  print(" surfScalarFields :", Info, ssf);
762 
763  PtrList<const surfaceVectorField> svf;
764  readFields
765  (
766  vMesh,
767  vMesh.baseMesh(),
768  objects,
769  selectedFields,
770  svf
771  );
772  print(" surfVectorFields :", Info, svf);
773 
774  if (ssf.size() + svf.size() > 0)
775  {
776  // Rework the scalar fields into vectorfields.
777  label sz = svf.size();
778 
779  svf.setSize(sz+ssf.size());
780 
781  surfaceVectorField n(mesh.Sf()/mesh.magSf());
782 
783  forAll(ssf, i)
784  {
785  surfaceVectorField* ssfiPtr = (ssf[i]*n).ptr();
786  ssfiPtr->rename(ssf[i].name());
787  svf.set(sz+i, ssfiPtr);
788  }
789  ssf.clear();
790 
791  mkDir(fvPath / "surfaceFields");
792 
793  fileName surfFileName
794  (
795  fvPath
796  /"surfaceFields"
797  /"surfaceFields"
798  + "_"
799  + timeDesc
800  + ".vtk"
801  );
802 
804  (
805  binary,
806  vMesh,
807  surfFileName,
808  svf
809  );
810  }
811  }
812 
813 
814  //---------------------------------------------------------------------
815  //
816  // Write patches (POLYDATA file, one for each patch)
817  //
818  //---------------------------------------------------------------------
819 
820  const polyBoundaryMesh& patches = mesh.boundaryMesh();
821 
822  if (allPatches)
823  {
824  mkDir(fvPath/"allPatches");
825 
826  fileName patchFileName;
827 
828  if (vMesh.useSubMesh())
829  {
830  patchFileName =
831  fvPath/"allPatches"/cellSetName
832  + "_"
833  + timeDesc
834  + ".vtk";
835  }
836  else
837  {
838  patchFileName =
839  fvPath/"allPatches"/"allPatches"
840  + "_"
841  + timeDesc
842  + ".vtk";
843  }
844 
845  Info<< " Combined patches : " << patchFileName << endl;
846 
847  patchWriter writer
848  (
849  vMesh,
850  binary,
851  nearCellValue,
852  patchFileName,
853  getSelectedPatches(patches, excludePatches)
854  );
855 
856  // VolFields + patchID
857  writeFuns::writeCellDataHeader
858  (
859  writer.os(),
860  writer.nFaces(),
861  1+nVolFields
862  );
863 
864  // Write patchID field
865  writer.writePatchIDs();
866 
867  // Write volFields
868  writer.write(vsf);
869  writer.write(vvf);
870  writer.write(vSpheretf);
871  writer.write(vSymmtf);
872  writer.write(vtf);
873 
874  if (!noPointValues)
875  {
876  writeFuns::writePointDataHeader
877  (
878  writer.os(),
879  writer.nPoints(),
880  nPointFields
881  );
882 
883  // Write pointFields
884  writer.write(psf);
885  writer.write(pvf);
886  writer.write(pSpheretf);
887  writer.write(pSymmtf);
888  writer.write(ptf);
889 
890  // no interpolated volFields since I cannot be bothered to
891  // create the patchInterpolation for all subpatches.
892  }
893  }
894  else
895  {
896  forAll(patches, patchi)
897  {
898  const polyPatch& pp = patches[patchi];
899 
900  if (!findStrings(excludePatches, pp.name()))
901  {
902  mkDir(fvPath/pp.name());
903 
904  fileName patchFileName;
905 
906  if (vMesh.useSubMesh())
907  {
908  patchFileName =
909  fvPath/pp.name()/cellSetName
910  + "_"
911  + timeDesc
912  + ".vtk";
913  }
914  else
915  {
916  patchFileName =
917  fvPath/pp.name()/pp.name()
918  + "_"
919  + timeDesc
920  + ".vtk";
921  }
922 
923  Info<< " Patch : " << patchFileName << endl;
924 
925  patchWriter writer
926  (
927  vMesh,
928  binary,
929  nearCellValue,
930  patchFileName,
931  labelList(1, patchi)
932  );
933 
934  if (!isA<emptyPolyPatch>(pp))
935  {
936  // VolFields + patchID
937  writeFuns::writeCellDataHeader
938  (
939  writer.os(),
940  writer.nFaces(),
941  1+nVolFields
942  );
943 
944  // Write patchID field
945  writer.writePatchIDs();
946 
947  // Write volFields
948  writer.write(vsf);
949  writer.write(vvf);
950  writer.write(vSpheretf);
951  writer.write(vSymmtf);
952  writer.write(vtf);
953 
954  if (!noPointValues)
955  {
956  writeFuns::writePointDataHeader
957  (
958  writer.os(),
959  writer.nPoints(),
960  nVolFields
961  + nPointFields
962  );
963 
964  // Write pointFields
965  writer.write(psf);
966  writer.write(pvf);
967  writer.write(pSpheretf);
968  writer.write(pSymmtf);
969  writer.write(ptf);
970 
971  PrimitivePatchInterpolation<primitivePatch> pInter
972  (
973  pp
974  );
975 
976  // Write interpolated volFields
977  writer.write(pInter, vsf);
978  writer.write(pInter, vvf);
979  writer.write(pInter, vSpheretf);
980  writer.write(pInter, vSymmtf);
981  writer.write(pInter, vtf);
982  }
983  }
984  }
985  }
986  }
987 
988  //---------------------------------------------------------------------
989  //
990  // Write faceZones (POLYDATA file, one for each zone)
991  //
992  //---------------------------------------------------------------------
993 
994  if (doFaceZones)
995  {
996  PtrList<const surfaceScalarField> ssf;
997  readFields
998  (
999  vMesh,
1000  vMesh.baseMesh(),
1001  objects,
1002  selectedFields,
1003  ssf
1004  );
1005  print(" surfScalarFields :", Info, ssf);
1006 
1007  PtrList<const surfaceVectorField> svf;
1008  readFields
1009  (
1010  vMesh,
1011  vMesh.baseMesh(),
1012  objects,
1013  selectedFields,
1014  svf
1015  );
1016  print(" surfVectorFields :", Info, svf);
1017 
1018  const faceZoneMesh& zones = mesh.faceZones();
1019 
1020  forAll(zones, zoneI)
1021  {
1022  const faceZone& fz = zones[zoneI];
1023 
1024  mkDir(fvPath/fz.name());
1025 
1026  fileName patchFileName;
1027 
1028  if (vMesh.useSubMesh())
1029  {
1030  patchFileName =
1031  fvPath/fz.name()/cellSetName
1032  + "_"
1033  + timeDesc
1034  + ".vtk";
1035  }
1036  else
1037  {
1038  patchFileName =
1039  fvPath/fz.name()/fz.name()
1040  + "_"
1041  + timeDesc
1042  + ".vtk";
1043  }
1044 
1045  Info<< " FaceZone : " << patchFileName << endl;
1046 
1048  (
1049  IndirectList<face>(mesh.faces(), fz),
1050  mesh.points()
1051  );
1052 
1053  surfaceMeshWriter writer
1054  (
1055  binary,
1056  pp,
1057  fz.name(),
1058  patchFileName
1059  );
1060 
1061  // Number of fields
1062  writeFuns::writeCellDataHeader
1063  (
1064  writer.os(),
1065  pp.size(),
1066  ssf.size()+svf.size()
1067  );
1068 
1069  writer.write(ssf);
1070  writer.write(svf);
1071  }
1072  }
1073 
1074 
1075 
1076  //---------------------------------------------------------------------
1077  //
1078  // Write lagrangian data
1079  //
1080  //---------------------------------------------------------------------
1081 
1082  forAllConstIter(HashSet<fileName>, allCloudDirs, iter)
1083  {
1084  const fileName& cloudName = iter.key();
1085 
1086  // Always create the cloud directory.
1087  mkDir(fvPath/cloud::prefix/cloudName);
1088 
1089  fileName lagrFileName
1090  (
1091  fvPath/cloud::prefix/cloudName/cloudName
1092  + "_" + timeDesc + ".vtk"
1093  );
1094 
1095  Info<< " Lagrangian: " << lagrFileName << endl;
1096 
1097 
1098  IOobjectList sprayObjs
1099  (
1100  mesh,
1101  runTime.timeName(),
1102  cloud::prefix/cloudName
1103  );
1104 
1105  IOobject* positionsPtr = sprayObjs.lookup(word("positions"));
1106 
1107  if (positionsPtr)
1108  {
1109  wordList labelNames(sprayObjs.names(labelIOField::typeName));
1110  Info<< " labels :";
1111  print(Info, labelNames);
1112 
1113  wordList scalarNames(sprayObjs.names(scalarIOField::typeName));
1114  Info<< " scalars :";
1115  print(Info, scalarNames);
1116 
1117  wordList vectorNames(sprayObjs.names(vectorIOField::typeName));
1118  Info<< " vectors :";
1119  print(Info, vectorNames);
1120 
1121  wordList sphereNames
1122  (
1123  sprayObjs.names
1124  (
1125  sphericalTensorIOField::typeName
1126  )
1127  );
1128  Info<< " spherical tensors :";
1129  print(Info, sphereNames);
1130 
1131  wordList symmNames
1132  (
1133  sprayObjs.names
1134  (
1135  symmTensorIOField::typeName
1136  )
1137  );
1138  Info<< " symm tensors :";
1139  print(Info, symmNames);
1140 
1141  wordList tensorNames(sprayObjs.names(tensorIOField::typeName));
1142  Info<< " tensors :";
1143  print(Info, tensorNames);
1144 
1145  lagrangianWriter writer
1146  (
1147  vMesh,
1148  binary,
1149  lagrFileName,
1150  cloudName,
1151  false
1152  );
1153 
1154  // Write number of fields
1155  writer.writeParcelHeader
1156  (
1157  labelNames.size()
1158  + scalarNames.size()
1159  + vectorNames.size()
1160  + sphereNames.size()
1161  + symmNames.size()
1162  + tensorNames.size()
1163  );
1164 
1165  // Fields
1166  writer.writeIOField<label>(labelNames);
1167  writer.writeIOField<scalar>(scalarNames);
1168  writer.writeIOField<vector>(vectorNames);
1169  writer.writeIOField<sphericalTensor>(sphereNames);
1170  writer.writeIOField<symmTensor>(symmNames);
1171  writer.writeIOField<tensor>(tensorNames);
1172  }
1173  else
1174  {
1175  lagrangianWriter writer
1176  (
1177  vMesh,
1178  binary,
1179  lagrFileName,
1180  cloudName,
1181  true
1182  );
1183 
1184  // Write number of fields
1185  writer.writeParcelHeader(0);
1186  }
1187  }
1188  }
1189 
1190 
1191  //---------------------------------------------------------------------
1192  //
1193  // Link parallel outputs back to undecomposed case for ease of loading
1194  //
1195  //---------------------------------------------------------------------
1196 
1197  if (Pstream::parRun() && doLinks)
1198  {
1199  mkDir(runTime.path()/".."/"VTK");
1200  chDir(runTime.path()/".."/"VTK");
1201 
1202  Info<< "Linking all processor files to " << runTime.path()/".."/"VTK"
1203  << endl;
1204 
1205  // Get list of vtk files
1206  fileName procVTK
1207  (
1208  fileName("..")
1209  /"processor" + name(Pstream::myProcNo())
1210  /"VTK"
1211  );
1212 
1213  fileNameList dirs(readDir(procVTK, fileType::directory));
1214  label sz = dirs.size();
1215  dirs.setSize(sz+1);
1216  dirs[sz] = ".";
1217 
1218  forAll(dirs, i)
1219  {
1220  fileNameList subFiles(readDir(procVTK/dirs[i], fileType::file));
1221 
1222  forAll(subFiles, j)
1223  {
1224  fileName procFile(procVTK/dirs[i]/subFiles[j]);
1225 
1226  if (exists(procFile))
1227  {
1228  string cmd
1229  (
1230  "ln -s "
1231  + procFile
1232  + " "
1233  + "processor"
1234  + name(Pstream::myProcNo())
1235  + "_"
1236  + procFile.name()
1237  );
1238  if (system(cmd.c_str()) == -1)
1239  {
1241  << "Could not execute command " << cmd << endl;
1242  }
1243  }
1244  }
1245  }
1246  }
1247 
1248  Info<< "End\n" << endl;
1249 
1250  return 0;
1251 }
1252 
1253 
1254 // ************************************************************************* //
List< instant > instantList
List of instants.
Definition: instantList.H:42
bool exists(const fileName &, const bool checkVariants=true, const bool followLink=true)
Does the name exist (as directory or file) in the file system?
Definition: POSIX.C:520
ZoneMesh< faceZone, polyMesh > faceZoneMesh
A ZoneMesh with the type faceZone.
#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
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const HashSet< word > &selectedFields, LIFOStack< regIOobject *> &storedObjects)
Read the selected GeometricFields of the specified type.
Definition: ReadFields.C:244
Foam::word regionName
void writePointSet(const bool binary, const primitiveMesh &mesh, const topoSet &set, const fileName &fileName)
Write pointSet to vtk polydata file.
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
engineTime & runTime
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:108
bool chDir(const fileName &dir)
Change the current directory to the one given and return true,.
Definition: POSIX.C:284
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
Operations on lists of strings.
PrimitivePatch< IndirectList< face >, const pointField & > indirectPrimitivePatch
Foam::indirectPrimitivePatch.
bool optionReadIfPresent(const word &opt, T &) const
Read a value from the named option if present.
Definition: argListI.H:198
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
bool findStrings(const wordReListMatcher &matcher, const std::string &str)
Return true if string matches one of the regular expressions.
Definition: stringListOps.H:52
void writeSurfFields(const bool binary, const vtkMesh &vMesh, const fileName &fileName, const UPtrList< const surfaceVectorField > &surfVectorFields)
dynamicFvMesh & mesh
SymmTensor< scalar > symmTensor
SymmTensor of scalars.
Definition: symmTensor.H:48
bool isDir(const fileName &, const bool followLink=true)
Does the name exist as a directory in the file system?
Definition: POSIX.C:539
float floatScalar
Float precision floating point scalar type.
Definition: floatScalar.H:52
Write faceSet to vtk polydata file. Only one data which is original faceID.
void writeFaceSet(const bool binary, const vtkMesh &vMesh, const faceSet &set, const fileName &fileName)
List< label > labelList
A List of labels.
Definition: labelList.H:56
fileNameList readDir(const fileName &, const fileType=fileType::file, const bool filterVariants=true, const bool followLink=true)
Read a directory and return the entries as a string list.
Definition: POSIX.C:662
graph_traits< Graph >::vertices_size_type size_type
Definition: SloanRenumber.C:73
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
bool rmDir(const fileName &)
Remove a directory and its contents.
Definition: POSIX.C:1051
Write a patch with its data.
static const char nl
Definition: Ostream.H:265
objects
virtual void rename(const word &newName)
Rename.
Definition: regIOobject.C:410
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:290
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
List< word > wordList
A List of words.
Definition: fileName.H:54
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
label n
const word cloudName(propsDict.lookup("cloudName"))
SphericalTensor< scalar > sphericalTensor
SphericalTensor of scalars.
Foam::argList args(argc, argv)
List< fileName > fileNameList
A List of fileNames.
Definition: fileNameList.H:50
Foam::faceZoneMesh.
Tensor< scalar > tensor
Tensor of scalars.
Definition: tensor.H:51
int system(const std::string &command)
Execute the specified command.
Definition: POSIX.C:1234
IStringStream optionLookup(const word &opt) const
Return an IStringStream from the named option.
Definition: argListI.H:114