applications/utilities/postProcessing/dataConversion/foamToVTK/foamToVTK.C

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