FreeFOAM The Cross-Platform CFD Toolkit
foamToVTK.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 Application
25  foamToVTK
26 
27 Description
28  Legacy VTK file format writer.
29 
30  - handles volScalar, volVector, pointScalar, pointVector, surfaceScalar
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 
41  - foamToVTK [OPTION]
42 
43 
44  @param -ascii \n
45  Write VTK data in ASCII format instead of binary.
46 
47  @param -mesh <name>\n
48  Use a different mesh name (instead of -region)
49 
50  @param -fields <fields>\n
51  Convert selected fields only. For example,
52  @verbatim
53  -fields "( p T U )"
54  @endverbatim
55  The quoting is required to avoid shell expansions and to pass the
56  information as a single argument.
57 
58  @param -surfaceFields \n
59  Write surfaceScalarFields (e.g., phi)
60 
61  @param -cellSet <name>\n
62  @param -faceSet <name>\n
63  @param -pointSet <name>\n.
64  Restrict conversion to the cellSet, faceSet or pointSet.
65 
66  @param -nearCellValue \n
67  Output cell value on patches instead of patch value itself.
68 
69  @param -noInternal \n
70  Do not generate file for mesh, only for patches.
71 
72  @param -noPointValues \n
73  No pointFields.
74 
75  @param -noFaceZones \n
76  No faceZones.
77 
78  @param -noLinks \n
79  (in parallel) do not link processor files to master.
80 
81  @param -allPatches \n
82  Combine all patches into a single file.
83 
84  @param -excludePatches <patchNames>\n
85  Specify patches to exclude. For example,
86  @verbatim
87  -excludePatches "( inlet_1 inlet_2 )"
88  @endverbatim
89  The quoting is required to avoid shell expansions and to pass the
90  information as a single argument.
91 
92  @param -useTimeName \n
93  use the time index in the VTK file name instead of the time index
94 
95  @param -case <dir> \n
96  Case directory.
97 
98  @param -noZero \n
99  Ignore time step 0.
100 
101  @param -constant \n
102  Include the constant directory.
103 
104  @param -latestTime \n
105  Only apply to the latest time step.
106 
107  @param -time <time>\n
108  Apply only to specified time.
109 
110  @param -parallel \n
111  Run in parallel.
112 
113  @param -help \n
114  Display help message.
115 
116  @param -doc \n
117  Display Doxygen API documentation page for this application.
118 
119  @param -srcDoc \n
120  Display Doxygen source documentation page for this application.
121 
122 Note
123  mesh subset is handled by vtkMesh. Slight inconsistency in
124  interpolation: on the internal field it interpolates the whole volfield
125  to the whole-mesh pointField and then selects only those values it
126  needs for the subMesh (using the fvMeshSubset cellMap(), pointMap()
127  functions). For the patches however it uses the
128  fvMeshSubset.interpolate function to directly interpolate the
129  whole-mesh values onto the subset patch.
130 
131 \*---------------------------------------------------------------------------*/
132 
133 #include <finiteVolume/fvCFD.H>
134 #include <OpenFOAM/pointMesh.H>
136 #include <OpenFOAM/emptyPolyPatch.H>
137 #include <OpenFOAM/labelIOField.H>
138 #include <OpenFOAM/scalarIOField.H>
141 #include <OpenFOAM/tensorIOField.H>
142 #include <OpenFOAM/faceZoneMesh.H>
143 #include <lagrangian/Cloud.H>
145 
146 #include "vtkMesh.H"
147 #include "readFields.H"
148 #include "writeFuns.H"
149 
150 #include "internalWriter.H"
151 #include "patchWriter.H"
152 #include "lagrangianWriter.H"
153 
154 #include "writeFaceSet.H"
155 #include "writePointSet.H"
156 #include "writePatchGeom.H"
157 #include "writeSurfFields.H"
158 
159 
160 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
161 
162 
163 static const label VTK_TETRA = 10;
164 static const label VTK_PYRAMID = 14;
165 static const label VTK_WEDGE = 13;
166 static const label VTK_HEXAHEDRON = 12;
167 
168 
169 template<class GeoField>
170 void print(const char* msg, Ostream& os, const PtrList<GeoField>& flds)
171 {
172  if (flds.size())
173  {
174  os << msg;
175  forAll(flds, i)
176  {
177  os<< ' ' << flds[i].name();
178  }
179  os << endl;
180  }
181 }
182 
183 
184 void print(Ostream& os, const wordList& flds)
185 {
186  forAll(flds, i)
187  {
188  os<< ' ' << flds[i];
189  }
190  os << endl;
191 }
192 
193 
194 labelList getSelectedPatches
195 (
196  const polyBoundaryMesh& patches,
197  const HashSet<word>& excludePatches
198 )
199 {
200  DynamicList<label> patchIDs(patches.size());
201 
202  Info<< "Combining patches:" << endl;
203 
204  forAll(patches, patchI)
205  {
206  const polyPatch& pp = patches[patchI];
207 
208  if
209  (
210  isA<emptyPolyPatch>(pp)
211  || (Pstream::parRun() && isA<processorPolyPatch>(pp))
212  )
213  {
214  Info<< " discarding empty/processor patch " << patchI
215  << " " << pp.name() << endl;
216  }
217  else if (!excludePatches.found(pp.name()))
218  {
219  patchIDs.append(patchI);
220  Info<< " patch " << patchI << " " << pp.name() << endl;
221  }
222  }
223  return patchIDs.shrink();
224 }
225 
226 
227 
228 
229 // Main program:
230 
231 int main(int argc, char *argv[])
232 {
233  timeSelector::addOptions();
234 
235 # include <OpenFOAM/addRegionOption.H>
236 
237  argList::validOptions.insert("fields", "fields");
238  argList::validOptions.insert("cellSet", "cellSet name");
239  argList::validOptions.insert("faceSet", "faceSet name");
240  argList::validOptions.insert("pointSet", "pointSet name");
241  argList::validOptions.insert("ascii","");
242  argList::validOptions.insert("surfaceFields","");
243  argList::validOptions.insert("nearCellValue","");
244  argList::validOptions.insert("noInternal","");
245  argList::validOptions.insert("noPointValues","");
246  argList::validOptions.insert("allPatches","");
247  argList::validOptions.insert("excludePatches","patches to exclude");
248  argList::validOptions.insert("noFaceZones","");
249  argList::validOptions.insert("noLinks","");
250  argList::validOptions.insert("useTimeName","");
251 
252 # include <OpenFOAM/setRootCase.H>
253 # include <OpenFOAM/createTime.H>
254 
255  bool doWriteInternal = !args.optionFound("noInternal");
256  bool doFaceZones = !args.optionFound("noFaceZones");
257  bool doLinks = !args.optionFound("noLinks");
258  bool binary = !args.optionFound("ascii");
259  bool useTimeName = args.optionFound("useTimeName");
260 
261  if (binary && (sizeof(floatScalar) != 4 || sizeof(label) != 4))
262  {
264  << "floatScalar and/or label are not 4 bytes in size" << nl
265  << "Hence cannot use binary VTK format. Please use -ascii"
266  << exit(FatalError);
267  }
268 
269  bool nearCellValue = args.optionFound("nearCellValue");
270 
271  if (nearCellValue)
272  {
274  << "Using neighbouring cell value instead of patch value"
275  << nl << endl;
276  }
277 
278  bool noPointValues = args.optionFound("noPointValues");
279 
280  if (noPointValues)
281  {
283  << "Outputting cell values only" << nl << endl;
284  }
285 
286  bool allPatches = args.optionFound("allPatches");
287 
288  HashSet<word> excludePatches;
289  if (args.optionFound("excludePatches"))
290  {
291  args.optionLookup("excludePatches")() >> excludePatches;
292 
293  Info<< "Not including patches " << excludePatches << nl << endl;
294  }
295 
296  word cellSetName;
297  string vtkName;
298 
299  if (args.optionFound("cellSet"))
300  {
301  cellSetName = args.option("cellSet");
302  vtkName = cellSetName;
303  }
304  else if (Pstream::parRun())
305  {
306  // Strip off leading casename, leaving just processor_DDD ending.
307  vtkName = runTime.caseName();
308 
309  string::size_type i = vtkName.rfind("processor");
310 
311  if (i != string::npos)
312  {
313  vtkName = vtkName.substr(i);
314  }
315  }
316  else
317  {
318  vtkName = runTime.caseName();
319  }
320 
321 
322  instantList timeDirs = timeSelector::select0(runTime, args);
323 
324 # include <OpenFOAM/createNamedMesh.H>
325 
326  // VTK/ directory in the case
327  fileName fvPath(runTime.path()/"VTK");
328  // Directory of mesh (region0 gets filtered out)
329  fileName regionPrefix = "";
330 
331  if (regionName != polyMesh::defaultRegion)
332  {
333  fvPath = fvPath/regionName;
334  regionPrefix = regionName;
335  }
336 
337  if (isDir(fvPath))
338  {
339  if
340  (
341  args.optionFound("time")
342  || args.optionFound("latestTime")
343  || cellSetName.size()
344  || regionName != polyMesh::defaultRegion
345  )
346  {
347  Info<< "Keeping old VTK files in " << fvPath << nl << endl;
348  }
349  else
350  {
351  Info<< "Deleting old VTK files in " << fvPath << nl << endl;
352 
353  rmDir(fvPath);
354  }
355  }
356 
357  mkDir(fvPath);
358 
359 
360  // mesh wrapper; does subsetting and decomposition
361  vtkMesh vMesh(mesh, cellSetName);
362 
363 
364  // Scan for all possible lagrangian clouds
365  HashSet<fileName> allCloudDirs;
366  forAll(timeDirs, timeI)
367  {
368  runTime.setTime(timeDirs[timeI], timeI);
369  fileNameList cloudDirs
370  (
371  readDir
372  (
373  runTime.timePath()/regionPrefix/cloud::prefix,
374  fileName::DIRECTORY
375  )
376  );
377  forAll(cloudDirs, i)
378  {
379  IOobjectList sprayObjs
380  (
381  mesh,
382  runTime.timeName(),
383  cloud::prefix/cloudDirs[i]
384  );
385 
386  IOobject* positionsPtr = sprayObjs.lookup("positions");
387 
388  if (positionsPtr)
389  {
390  if (allCloudDirs.insert(cloudDirs[i]))
391  {
392  Info<< "At time: " << runTime.timeName()
393  << " detected cloud directory : " << cloudDirs[i]
394  << endl;
395  }
396  }
397  }
398  }
399 
400 
401  forAll(timeDirs, timeI)
402  {
403  runTime.setTime(timeDirs[timeI], timeI);
404 
405  Info<< "Time: " << runTime.timeName() << endl;
406 
407  word timeDesc = "";
408  if (useTimeName)
409  {
410  timeDesc = runTime.timeName();
411  }
412  else
413  {
414  timeDesc = name(runTime.timeIndex());
415  }
416 
417  // Check for new polyMesh/ and update mesh, fvMeshSubset and cell
418  // decomposition.
419  polyMesh::readUpdateState meshState = vMesh.readUpdate();
420 
421  const fvMesh& mesh = vMesh.mesh();
422 
423  if
424  (
425  meshState == polyMesh::TOPO_CHANGE
426  || meshState == polyMesh::TOPO_PATCH_CHANGE
427  )
428  {
429  Info<< " Read new mesh" << nl << endl;
430  }
431 
432 
433  // If faceSet: write faceSet only (as polydata)
434  if (args.optionFound("faceSet"))
435  {
436  // Load the faceSet
437  faceSet set(mesh, args.option("faceSet"));
438 
439  // Filename as if patch with same name.
440  mkDir(fvPath/set.name());
441 
442  fileName patchFileName
443  (
444  fvPath/set.name()/set.name()
445  + "_"
446  + timeDesc
447  + ".vtk"
448  );
449 
450  Info<< " FaceSet : " << patchFileName << endl;
451 
452  writeFaceSet(binary, vMesh, set, patchFileName);
453 
454  continue;
455  }
456  // If pointSet: write pointSet only (as polydata)
457  if (args.optionFound("pointSet"))
458  {
459  // Load the pointSet
460  pointSet set(mesh, args.option("pointSet"));
461 
462  // Filename as if patch with same name.
463  mkDir(fvPath/set.name());
464 
465  fileName patchFileName
466  (
467  fvPath/set.name()/set.name()
468  + "_"
469  + timeDesc
470  + ".vtk"
471  );
472 
473  Info<< " pointSet : " << patchFileName << endl;
474 
475  writePointSet(binary, vMesh, set, patchFileName);
476 
477  continue;
478  }
479 
480 
481  // Search for list of objects for this time
482  IOobjectList objects(mesh, runTime.timeName());
483 
484  HashSet<word> selectedFields;
485  if (args.optionFound("fields"))
486  {
487  args.optionLookup("fields")() >> selectedFields;
488  }
489 
490  // Construct the vol fields (on the original mesh if subsetted)
491 
493  readFields(vMesh, vMesh.baseMesh(), objects, selectedFields, vsf);
494  print(" volScalarFields :", Info, vsf);
495 
497  readFields(vMesh, vMesh.baseMesh(), objects, selectedFields, vvf);
498  print(" volVectorFields :", Info, vvf);
499 
501  readFields(vMesh, vMesh.baseMesh(), objects, selectedFields, vSpheretf);
502  print(" volSphericalTensorFields :", Info, vSpheretf);
503 
505  readFields(vMesh, vMesh.baseMesh(), objects, selectedFields, vSymmtf);
506  print(" volSymmTensorFields :", Info, vSymmtf);
507 
509  readFields(vMesh, vMesh.baseMesh(), objects, selectedFields, vtf);
510  print(" volTensorFields :", Info, vtf);
511 
512  label nVolFields =
513  vsf.size()
514  + vvf.size()
515  + vSpheretf.size()
516  + vSymmtf.size()
517  + vtf.size();
518 
519 
520  // Construct pointMesh only if nessecary since constructs edge
521  // addressing (expensive on polyhedral meshes)
522  if (noPointValues)
523  {
524  Info<< " pointScalarFields : switched off"
525  << " (\"-noPointValues\" option)\n";
526  Info<< " pointVectorFields : switched off"
527  << " (\"-noPointValues\" option)\n";
528  }
529 
535 
536  if (!noPointValues)
537  {
538  readFields
539  (
540  vMesh,
541  pointMesh::New(vMesh.baseMesh()),
542  objects,
543  selectedFields,
544  psf
545  );
546  print(" pointScalarFields :", Info, psf);
547 
548  readFields
549  (
550  vMesh,
551  pointMesh::New(vMesh.baseMesh()),
552  objects,
553  selectedFields,
554  pvf
555  );
556  print(" pointVectorFields :", Info, pvf);
557 
558  readFields
559  (
560  vMesh,
561  pointMesh::New(vMesh.baseMesh()),
562  objects,
563  selectedFields,
564  pSpheretf
565  );
566  print(" pointSphericalTensorFields :", Info, pSpheretf);
567 
568  readFields
569  (
570  vMesh,
571  pointMesh::New(vMesh.baseMesh()),
572  objects,
573  selectedFields,
574  pSymmtf
575  );
576  print(" pointSymmTensorFields :", Info, pSymmtf);
577 
578  readFields
579  (
580  vMesh,
581  pointMesh::New(vMesh.baseMesh()),
582  objects,
583  selectedFields,
584  ptf
585  );
586  print(" pointTensorFields :", Info, ptf);
587  }
588  Info<< endl;
589 
590  label nPointFields =
591  psf.size()
592  + pvf.size()
593  + pSpheretf.size()
594  + pSymmtf.size()
595  + ptf.size();
596 
597  if (doWriteInternal)
598  {
599  //
600  // Create file and write header
601  //
602  fileName vtkFileName
603  (
604  fvPath/vtkName
605  + "_"
606  + timeDesc
607  + ".vtk"
608  );
609 
610  Info<< " Internal : " << vtkFileName << endl;
611 
612  // Write mesh
613  internalWriter writer(vMesh, binary, vtkFileName);
614 
615  // VolFields + cellID
616  writeFuns::writeCellDataHeader
617  (
618  writer.os(),
619  vMesh.nFieldCells(),
620  1+nVolFields
621  );
622 
623  // Write cellID field
624  writer.writeCellIDs();
625 
626  // Write volFields
627  writer.write(vsf);
628  writer.write(vvf);
629  writer.write(vSpheretf);
630  writer.write(vSymmtf);
631  writer.write(vtf);
632 
633  if (!noPointValues)
634  {
635  writeFuns::writePointDataHeader
636  (
637  writer.os(),
638  vMesh.nFieldPoints(),
639  nVolFields+nPointFields
640  );
641 
642  // pointFields
643  writer.write(psf);
644  writer.write(pvf);
645  writer.write(pSpheretf);
646  writer.write(pSymmtf);
647  writer.write(ptf);
648 
649  // Interpolated volFields
650  volPointInterpolation pInterp(mesh);
651  writer.write(pInterp, vsf);
652  writer.write(pInterp, vvf);
653  writer.write(pInterp, vSpheretf);
654  writer.write(pInterp, vSymmtf);
655  writer.write(pInterp, vtf);
656  }
657  }
658 
659  //---------------------------------------------------------------------
660  //
661  // Write surface fields
662  //
663  //---------------------------------------------------------------------
664 
665  if (args.optionFound("surfaceFields"))
666  {
668  readFields
669  (
670  vMesh,
671  vMesh.baseMesh(),
672  objects,
673  selectedFields,
674  ssf
675  );
676  print(" surfScalarFields :", Info, ssf);
677 
679  readFields
680  (
681  vMesh,
682  vMesh.baseMesh(),
683  objects,
684  selectedFields,
685  svf
686  );
687  print(" surfVectorFields :", Info, svf);
688 
689  if (ssf.size() + svf.size() > 0)
690  {
691  // Rework the scalar fields into vectorfields.
692  label sz = svf.size();
693 
694  svf.setSize(sz+ssf.size());
695 
696  surfaceVectorField n(mesh.Sf()/mesh.magSf());
697 
698  forAll(ssf, i)
699  {
700  svf.set(sz+i, ssf[i]*n);
701  svf[sz+i].rename(ssf[i].name());
702  }
703  ssf.clear();
704 
705  mkDir(fvPath / "surfaceFields");
706 
707  fileName surfFileName
708  (
709  fvPath
710  /"surfaceFields"
711  /"surfaceFields"
712  + "_"
713  + timeDesc
714  + ".vtk"
715  );
716 
718  (
719  binary,
720  vMesh,
721  surfFileName,
722  svf
723  );
724  }
725  }
726 
727 
728  //---------------------------------------------------------------------
729  //
730  // Write patches (POLYDATA file, one for each patch)
731  //
732  //---------------------------------------------------------------------
733 
734  const polyBoundaryMesh& patches = mesh.boundaryMesh();
735 
736  if (allPatches)
737  {
738  mkDir(fvPath/"allPatches");
739 
740  fileName patchFileName;
741 
742  if (vMesh.useSubMesh())
743  {
744  patchFileName =
745  fvPath/"allPatches"/cellSetName
746  + "_"
747  + timeDesc
748  + ".vtk";
749  }
750  else
751  {
752  patchFileName =
753  fvPath/"allPatches"/"allPatches"
754  + "_"
755  + timeDesc
756  + ".vtk";
757  }
758 
759  Info<< " Combined patches : " << patchFileName << endl;
760 
762  (
763  vMesh,
764  binary,
765  nearCellValue,
766  patchFileName,
767  getSelectedPatches(patches, excludePatches)
768  );
769 
770  // VolFields + patchID
771  writeFuns::writeCellDataHeader
772  (
773  writer.os(),
774  writer.nFaces(),
775  1+nVolFields
776  );
777 
778  // Write patchID field
779  writer.writePatchIDs();
780 
781  // Write volFields
782  writer.write(vsf);
783  writer.write(vvf);
784  writer.write(vSpheretf);
785  writer.write(vSymmtf);
786  writer.write(vtf);
787 
788  if (!noPointValues)
789  {
790  writeFuns::writePointDataHeader
791  (
792  writer.os(),
793  writer.nPoints(),
794  nPointFields
795  );
796 
797  // Write pointFields
798  writer.write(psf);
799  writer.write(pvf);
800  writer.write(pSpheretf);
801  writer.write(pSymmtf);
802  writer.write(ptf);
803 
804  // no interpolated volFields since I cannot be bothered to
805  // create the patchInterpolation for all subpatches.
806  }
807  }
808  else
809  {
810  forAll(patches, patchI)
811  {
812  const polyPatch& pp = patches[patchI];
813 
814  if (!excludePatches.found(pp.name()))
815  {
816  mkDir(fvPath/pp.name());
817 
818  fileName patchFileName;
819 
820  if (vMesh.useSubMesh())
821  {
822  patchFileName =
823  fvPath/pp.name()/cellSetName
824  + "_"
825  + timeDesc
826  + ".vtk";
827  }
828  else
829  {
830  patchFileName =
831  fvPath/pp.name()/pp.name()
832  + "_"
833  + timeDesc
834  + ".vtk";
835  }
836 
837  Info<< " Patch : " << patchFileName << endl;
838 
840  (
841  vMesh,
842  binary,
843  nearCellValue,
844  patchFileName,
845  labelList(1, patchI)
846  );
847 
848  if (!isA<emptyPolyPatch>(pp))
849  {
850  // VolFields + patchID
851  writeFuns::writeCellDataHeader
852  (
853  writer.os(),
854  writer.nFaces(),
855  1+nVolFields
856  );
857 
858  // Write patchID field
859  writer.writePatchIDs();
860 
861  // Write volFields
862  writer.write(vsf);
863  writer.write(vvf);
864  writer.write(vSpheretf);
865  writer.write(vSymmtf);
866  writer.write(vtf);
867 
868  if (!noPointValues)
869  {
870  writeFuns::writePointDataHeader
871  (
872  writer.os(),
873  writer.nPoints(),
874  nVolFields
875  + nPointFields
876  );
877 
878  // Write pointFields
879  writer.write(psf);
880  writer.write(pvf);
881  writer.write(pSpheretf);
882  writer.write(pSymmtf);
883  writer.write(ptf);
884 
886  (
887  pp
888  );
889 
890  // Write interpolated volFields
891  writer.write(pInter, vsf);
892  writer.write(pInter, vvf);
893  writer.write(pInter, vSpheretf);
894  writer.write(pInter, vSymmtf);
895  writer.write(pInter, vtf);
896  }
897  }
898  }
899  }
900  }
901 
902  //---------------------------------------------------------------------
903  //
904  // Write faceZones (POLYDATA file, one for each zone)
905  //
906  //---------------------------------------------------------------------
907 
908  if (doFaceZones)
909  {
910  const faceZoneMesh& zones = mesh.faceZones();
911 
912  forAll(zones, zoneI)
913  {
914  const faceZone& pp = zones[zoneI];
915 
916  mkDir(fvPath/pp.name());
917 
918  fileName patchFileName;
919 
920  if (vMesh.useSubMesh())
921  {
922  patchFileName =
923  fvPath/pp.name()/cellSetName
924  + "_"
925  + timeDesc
926  + ".vtk";
927  }
928  else
929  {
930  patchFileName =
931  fvPath/pp.name()/pp.name()
932  + "_"
933  + timeDesc
934  + ".vtk";
935  }
936 
937  Info<< " FaceZone : " << patchFileName << endl;
938 
939  std::ofstream str(patchFileName.c_str());
940 
941  writeFuns::writeHeader(str, binary, pp.name());
942  str << "DATASET POLYDATA" << std::endl;
943 
945  (
946  binary,
947  pp().localFaces(),
948  pp().localPoints(),
949  str
950  );
951  }
952  }
953 
954 
955 
956  //---------------------------------------------------------------------
957  //
958  // Write lagrangian data
959  //
960  //---------------------------------------------------------------------
961 
962  forAllConstIter(HashSet<fileName>, allCloudDirs, iter)
963  {
964  const fileName& cloudName = iter.key();
965 
966  // Always create the cloud directory.
967  mkDir(fvPath/cloud::prefix/cloudName);
968 
969  fileName lagrFileName
970  (
971  fvPath/cloud::prefix/cloudName/cloudName
972  + "_" + timeDesc + ".vtk"
973  );
974 
975  Info<< " Lagrangian: " << lagrFileName << endl;
976 
977 
978  IOobjectList sprayObjs
979  (
980  mesh,
981  runTime.timeName(),
982  cloud::prefix/cloudName
983  );
984 
985  IOobject* positionsPtr = sprayObjs.lookup("positions");
986 
987  if (positionsPtr)
988  {
989  wordList labelNames(sprayObjs.names(labelIOField::typeName));
990  Info<< " labels :";
991  print(Info, labelNames);
992 
993  wordList scalarNames(sprayObjs.names(scalarIOField::typeName));
994  Info<< " scalars :";
995  print(Info, scalarNames);
996 
997  wordList vectorNames(sprayObjs.names(vectorIOField::typeName));
998  Info<< " vectors :";
999  print(Info, vectorNames);
1000 
1001  wordList sphereNames
1002  (
1003  sprayObjs.names
1004  (
1005  sphericalTensorIOField::typeName
1006  )
1007  );
1008  Info<< " spherical tensors :";
1009  print(Info, sphereNames);
1010 
1011  wordList symmNames
1012  (
1013  sprayObjs.names
1014  (
1015  symmTensorIOField::typeName
1016  )
1017  );
1018  Info<< " symm tensors :";
1019  print(Info, symmNames);
1020 
1021  wordList tensorNames(sprayObjs.names(tensorIOField::typeName));
1022  Info<< " tensors :";
1023  print(Info, tensorNames);
1024 
1026  (
1027  vMesh,
1028  binary,
1029  lagrFileName,
1030  cloudName,
1031  false
1032  );
1033 
1034  // Write number of fields
1035  writer.writeParcelHeader
1036  (
1037  labelNames.size()
1038  + scalarNames.size()
1039  + vectorNames.size()
1040  + sphereNames.size()
1041  + symmNames.size()
1042  + tensorNames.size()
1043  );
1044 
1045  // Fields
1046  writer.writeIOField<label>(labelNames);
1047  writer.writeIOField<scalar>(scalarNames);
1048  writer.writeIOField<vector>(vectorNames);
1049  writer.writeIOField<sphericalTensor>(sphereNames);
1050  writer.writeIOField<symmTensor>(symmNames);
1051  writer.writeIOField<tensor>(tensorNames);
1052  }
1053  else
1054  {
1056  (
1057  vMesh,
1058  binary,
1059  lagrFileName,
1060  cloudName,
1061  true
1062  );
1063 
1064  // Write number of fields
1065  writer.writeParcelHeader(0);
1066  }
1067  }
1068  }
1069 
1070 
1071  //---------------------------------------------------------------------
1072  //
1073  // Link parallel outputs back to undecomposed case for ease of loading
1074  //
1075  //---------------------------------------------------------------------
1076 
1077  if (Pstream::parRun() && doLinks)
1078  {
1079  mkDir(runTime.path()/".."/"VTK");
1080  chDir(runTime.path()/".."/"VTK");
1081 
1082  Info<< "Linking all processor files to " << runTime.path()/".."/"VTK"
1083  << endl;
1084 
1085  // Get list of vtk files
1086  fileName procVTK
1087  (
1088  fileName("..")
1089  /"processor" + name(Pstream::myProcNo())
1090  /"VTK"
1091  );
1092 
1093  fileNameList dirs(readDir(procVTK, fileName::DIRECTORY));
1094  label sz = dirs.size();
1095  dirs.setSize(sz+1);
1096  dirs[sz] = ".";
1097 
1098  forAll(dirs, i)
1099  {
1100  fileNameList subFiles(readDir(procVTK/dirs[i], fileName::FILE));
1101 
1102  forAll(subFiles, j)
1103  {
1104  fileName procFile(procVTK/dirs[i]/subFiles[j]);
1105 
1106  if (exists(procFile))
1107  {
1108  string cmd
1109  (
1110  "ln -s "
1111  + procFile
1112  + " "
1113  + "processor"
1114  + name(Pstream::myProcNo())
1115  + "_"
1116  + procFile.name()
1117  );
1118  if (system(cmd.c_str()) == -1)
1119  {
1121  << "Could not execute command " << cmd << endl;
1122  }
1123  }
1124  }
1125  }
1126  }
1127 
1128  Info<< "End\n" << endl;
1129 
1130  return 0;
1131 }
1132 
1133 
1134 // ************************ vim: set sw=4 sts=4 et: ************************ //