applications/utilities/mesh/manipulation/setSet/setSet.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 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 Description
  25     Manipulate a cell/face/point/ set or zone interactively.
  26
  27 \*---------------------------------------------------------------------------*/
  28
  29 #include "argList.H"
  30 #include "Time.H"
  31 #include "polyMesh.H"
  32 #include "globalMeshData.H"
  33 #include "IStringStream.H"
  34 #include "cellSet.H"
  35 #include "faceSet.H"
  36 #include "pointSet.H"
  37 #include "topoSetSource.H"
  38 #include "OFstream.H"
  39 #include "IFstream.H"
  40 #include "demandDrivenData.H"
  41 #include "writePatch.H"
  42 #include "writePointSet.H"
  43 #include "IOobjectList.H"
  44 #include "cellZoneSet.H"
  45 #include "faceZoneSet.H"
  46 #include "pointZoneSet.H"
  47 #include "timeSelector.H"
  48
  49 #include <stdio.h>
  50
  51
  52 #ifdef HAS_READLINE
  53 # include <readline/readline.h>
  54 # include <readline/history.h>
  55 #endif
  56
  57 using namespace Foam;
  58
  59 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  60
  61
  62 #ifdef HAS_READLINE
  63 static const char* historyFile = ".setSet";
  64 #endif
  65
  66
  67 // Write set to VTK readable files
  68 void writeVTK
  69 (
  70     const polyMesh& mesh,
  71     const topoSet& currentSet,
  72     const fileName& vtkName
  73 )
  74 {
  75     if (isA<faceSet>(currentSet))
  76     {
  77         // Faces of set with OpenFOAM faceID as value
  78
  79         faceList setFaces(currentSet.size());
  80         labelList faceValues(currentSet.size());
  81         label setFaceI = 0;
  82
  83         forAllConstIter(topoSet, currentSet, iter)
  84         {
  85             setFaces[setFaceI] = mesh.faces()[iter.key()];
  86             faceValues[setFaceI] = iter.key();
  87             setFaceI++;
  88         }
  89
  90         primitiveFacePatch fp(setFaces, mesh.points());
  91
  92         writePatch
  93         (
  94             true,
  95             currentSet.name(),
  96             fp,
  97             "faceID",
  98             faceValues,
  99             mesh.time().path()/vtkName
 100         );
 101     }
 102     else if (isA<cellSet>(currentSet))
 103     {
 104         // External faces of cellset with OpenFOAM cellID as value
 105
 106         Map<label> cellFaces(currentSet.size());
 107
 108         forAllConstIter(cellSet, currentSet, iter)
 109         {
 110             label cellI = iter.key();
 111
 112             const cell& cFaces = mesh.cells()[cellI];
 113
 114             forAll(cFaces, i)
 115             {
 116                 label faceI = cFaces[i];
 117
 118                 if (mesh.isInternalFace(faceI))
 119                 {
 120                     label otherCellI = mesh.faceOwner()[faceI];
 121
 122                     if (otherCellI == cellI)
 123                     {
 124                         otherCellI = mesh.faceNeighbour()[faceI];
 125                     }
 126
 127                     if (!currentSet.found(otherCellI))
 128                     {
 129                         cellFaces.insert(faceI, cellI);
 130                     }
 131                 }
 132                 else
 133                 {
 134                     cellFaces.insert(faceI, cellI);
 135                 }
 136             }
 137         }
 138
 139         faceList setFaces(cellFaces.size());
 140         labelList faceValues(cellFaces.size());
 141         label setFaceI = 0;
 142
 143         forAllConstIter(Map<label>, cellFaces, iter)
 144         {
 145             setFaces[setFaceI] = mesh.faces()[iter.key()];
 146             faceValues[setFaceI] = iter();              // Cell ID
 147             setFaceI++;
 148         }
 149
 150         primitiveFacePatch fp(setFaces, mesh.points());
 151
 152         writePatch
 153         (
 154             true,
 155             currentSet.name(),
 156             fp,
 157             "cellID",
 158             faceValues,
 159             mesh.time().path()/vtkName
 160         );
 161     }
 162     else if (isA<pointSet>(currentSet))
 163     {
 164         writePointSet
 165         (
 166             true,
 167             mesh,
 168             currentSet,
 169             mesh.time().path()/vtkName
 170         );
 171     }
 172     else
 173     {
 174         WarningIn
 175         (
 176             "void writeVTK"
 177             "(const polyMesh& mesh, const topoSet& currentSet,"
 178             "const fileName& vtkName)"
 179         )   << "Don't know how to handle set of type " << currentSet.type()
 180             << endl;
 181     }
 182 }
 183
 184
 185 void printHelp(Ostream& os)
 186 {
 187     os  << "Please type 'help', 'list', 'quit', 'time ddd'"
 188         << " or a set command after prompt." << endl
 189         << "'list' will show all current cell/face/point sets." << endl
 190         << "'time ddd' will change the current time." << endl
 191         << endl
 192         << "A set command should be of the following form" << endl
 193         << endl
 194         << "    cellSet|faceSet|pointSet <setName> <action> <source>"
 195         << endl
 196         << endl
 197         << "The <action> is one of" << endl
 198         << "    list            - prints the contents of the set" << endl
 199         << "    clear           - clears the set" << endl
 200         << "    invert          - inverts the set" << endl
 201         << "    remove          - remove the set" << endl
 202         << "    new <source>    - sets to set to the source set" << endl
 203         << "    add <source>    - adds all elements from the source set" << endl
 204         << "    delete <source> - deletes      ,," << endl
 205         << "    subset <source> - combines current set with the source set"
 206         << endl
 207         << endl
 208         << "The sources come in various forms. Type a wrong source"
 209         << " to see all the types available." << endl
 210         << endl
 211         << "Example: pick up all cells connected by point or face to patch"
 212         << " movingWall" << endl
 213         << endl
 214         << "Pick up all faces of patch:" << endl
 215         << "    faceSet f0 new patchToFace movingWall" << endl
 216         << "Add faces 0,1,2:" << endl
 217         << "    faceSet f0 add labelToFace (0 1 2)" << endl
 218         << "Pick up all points used by faces in faceSet f0:" << endl
 219         << "    pointSet p0 new faceToPoint f0 all" << endl
 220         << "Pick up cell which has any face in f0:" << endl
 221         << "    cellSet c0 new faceToCell f0 any" << endl
 222         << "Add cells which have any point in p0:" << endl
 223         << "    cellSet c0 add pointToCell p0 any" << endl
 224         << "List set:" << endl
 225         << "    cellSet c0 list" << endl
 226         << endl
 227         << "Zones can be set using zoneSets from corresponding sets:" << endl
 228         << "    cellZoneSet c0Zone new setToCellZone c0" << endl
 229         << "    faceZoneSet f0Zone new setToFaceZone f0" << endl
 230         << endl
 231         << "or if orientation is important:" << endl
 232         << "    faceZoneSet f0Zone new setsToFaceZone f0 c0" << endl
 233         << endl
 234         << "ZoneSets can be manipulated using the general actions:" << endl
 235         << "    list            - prints the contents of the set" << endl
 236         << "    clear           - clears the set" << endl
 237         << "    invert          - inverts the set (undefined orientation)"
 238         << endl
 239         << "    remove          - remove the set" << endl
 240         << endl;
 241 }
 242
 243
 244 void printAllSets(const polyMesh& mesh, Ostream& os)
 245 {
 246     IOobjectList objects
 247     (
 248         mesh,
 249         mesh.time().findInstance
 250         (
 251             polyMesh::meshSubDir/"sets",
 252             word::null,
 253             IOobject::READ_IF_PRESENT,
 254             mesh.facesInstance()
 255         ),
 256         polyMesh::meshSubDir/"sets"
 257     );
 258     IOobjectList cellSets(objects.lookupClass(cellSet::typeName));
 259     if (cellSets.size())
 260     {
 261         os  << "cellSets:" << endl;
 262         forAllConstIter(IOobjectList, cellSets, iter)
 263         {
 264             cellSet set(*iter());
 265             os  << '\t' << set.name() << "\tsize:" << set.size() << endl;
 266         }
 267     }
 268     IOobjectList faceSets(objects.lookupClass(faceSet::typeName));
 269     if (faceSets.size())
 270     {
 271         os  << "faceSets:" << endl;
 272         forAllConstIter(IOobjectList, faceSets, iter)
 273         {
 274             faceSet set(*iter());
 275             os  << '\t' << set.name() << "\tsize:" << set.size() << endl;
 276         }
 277     }
 278     IOobjectList pointSets(objects.lookupClass(pointSet::typeName));
 279     if (pointSets.size())
 280     {
 281         os  << "pointSets:" << endl;
 282         forAllConstIter(IOobjectList, pointSets, iter)
 283         {
 284             pointSet set(*iter());
 285             os  << '\t' << set.name() << "\tsize:" << set.size() << endl;
 286         }
 287     }
 288
 289     const cellZoneMesh& cellZones = mesh.cellZones();
 290     if (cellZones.size())
 291     {
 292         os  << "cellZones:" << endl;
 293         forAll(cellZones, i)
 294         {
 295             const cellZone& zone = cellZones[i];
 296             os  << '\t' << zone.name() << "\tsize:" << zone.size() << endl;
 297         }
 298     }
 299     const faceZoneMesh& faceZones = mesh.faceZones();
 300     if (faceZones.size())
 301     {
 302         os  << "faceZones:" << endl;
 303         forAll(faceZones, i)
 304         {
 305             const faceZone& zone = faceZones[i];
 306             os  << '\t' << zone.name() << "\tsize:" << zone.size() << endl;
 307         }
 308     }
 309     const pointZoneMesh& pointZones = mesh.pointZones();
 310     if (pointZones.size())
 311     {
 312         os  << "pointZones:" << endl;
 313         forAll(pointZones, i)
 314         {
 315             const pointZone& zone = pointZones[i];
 316             os  << '\t' << zone.name() << "\tsize:" << zone.size() << endl;
 317         }
 318     }
 319
 320     os  << endl;
 321 }
 322
 323
 324 template<class ZoneType>
 325 void removeZone
 326 (
 327     ZoneMesh<ZoneType, polyMesh>& zones,
 328     const word& setName
 329 )
 330 {
 331     label zoneID = zones.findZoneID(setName);
 332
 333     if (zoneID != -1)
 334     {
 335         Info<< "Removing zone " << setName << " at index " << zoneID << endl;
 336         // Shuffle to last position
 337         labelList oldToNew(zones.size());
 338         label newI = 0;
 339         forAll(oldToNew, i)
 340         {
 341             if (i != zoneID)
 342             {
 343                 oldToNew[i] = newI++;
 344             }
 345         }
 346         oldToNew[zoneID] = newI;
 347         zones.reorder(oldToNew);
 348         // Remove last element
 349         zones.setSize(zones.size()-1);
 350         zones.clearAddressing();
 351         zones.write();
 352     }
 353 }
 354
 355
 356 // Physically remove a set
 357 void removeSet
 358 (
 359     const polyMesh& mesh,
 360     const word& setType,
 361     const word& setName
 362 )
 363 {
 364     // Remove the file
 365     IOobjectList objects
 366     (
 367         mesh,
 368         mesh.time().findInstance
 369         (
 370             polyMesh::meshSubDir/"sets",
 371             word::null,
 372             IOobject::READ_IF_PRESENT,
 373             mesh.facesInstance()
 374         ),
 375         polyMesh::meshSubDir/"sets"
 376     );
 377
 378     if (objects.found(setName))
 379     {
 380         // Remove file
 381         fileName object = objects[setName]->objectPath();
 382         Info<< "Removing file " << object << endl;
 383         rm(object);
 384     }
 385
 386     // See if zone
 387     if (setType == cellZoneSet::typeName)
 388     {
 389         removeZone
 390         (
 391             const_cast<cellZoneMesh&>(mesh.cellZones()),
 392             setName
 393         );
 394     }
 395     else if (setType == faceZoneSet::typeName)
 396     {
 397         removeZone
 398         (
 399             const_cast<faceZoneMesh&>(mesh.faceZones()),
 400             setName
 401         );
 402     }
 403     else if (setType == pointZoneSet::typeName)
 404     {
 405         removeZone
 406         (
 407             const_cast<pointZoneMesh&>(mesh.pointZones()),
 408             setName
 409         );
 410     }
 411 }
 412
 413
 414 // Read command and execute. Return true if ok, false otherwise.
 415 bool doCommand
 416 (
 417     const polyMesh& mesh,
 418     const word& setType,
 419     const word& setName,
 420     const word& actionName,
 421     const bool writeVTKFile,
 422     const bool writeCurrentTime,
 423     const bool noSync,
 424     Istream& is
 425 )
 426 {
 427     // Get some size estimate for set.
 428     const globalMeshData& parData = mesh.globalData();
 429
 430     label typSize =
 431         max
 432         (
 433             parData.nTotalCells(),
 434             max
 435             (
 436                 parData.nTotalFaces(),
 437                 parData.nTotalPoints()
 438             )
 439         )
 440       / (10*Pstream::nProcs());
 441
 442
 443     bool ok = true;
 444
 445     // Set to work on
 446     autoPtr<topoSet> currentSetPtr;
 447
 448     word sourceType;
 449
 450     try
 451     {
 452         topoSetSource::setAction action =
 453             topoSetSource::toAction(actionName);
 454
 455
 456         IOobject::readOption r;
 457
 458         if (action == topoSetSource::REMOVE)
 459         {
 460             removeSet(mesh, setType, setName);
 461         }
 462         else if
 463         (
 464             (action == topoSetSource::NEW)
 465          || (action == topoSetSource::CLEAR)
 466         )
 467         {
 468             r = IOobject::NO_READ;
 469             currentSetPtr = topoSet::New(setType, mesh, setName, typSize);
 470         }
 471         else
 472         {
 473             r = IOobject::MUST_READ;
 474             currentSetPtr = topoSet::New(setType, mesh, setName, r);
 475             topoSet& currentSet = currentSetPtr();
 476             // Presize it according to current mesh data.
 477             currentSet.resize(max(currentSet.size(), typSize));
 478         }
 479
 480         if (currentSetPtr.valid())
 481         {
 482             topoSet& currentSet = currentSetPtr();
 483
 484             Info<< "    Set:" << currentSet.name()
 485                 << "  Size:" << returnReduce(currentSet.size(), sumOp<label>())
 486                 << "  Action:" << actionName
 487                 << endl;
 488
 489             switch (action)
 490             {
 491                 case topoSetSource::CLEAR:
 492                 {
 493                     // Already handled above by not reading
 494                     break;
 495                 }
 496                 case topoSetSource::INVERT:
 497                 {
 498                     currentSet.invert(currentSet.maxSize(mesh));
 499                     break;
 500                 }
 501                 case topoSetSource::LIST:
 502                 {
 503                     currentSet.writeDebug(Pout, mesh, 100);
 504                     Pout<< endl;
 505                     break;
 506                 }
 507                 case topoSetSource::SUBSET:
 508                 {
 509                     if (is >> sourceType)
 510                     {
 511                         autoPtr<topoSetSource> setSource
 512                         (
 513                             topoSetSource::New
 514                             (
 515                                 sourceType,
 516                                 mesh,
 517                                 is
 518                             )
 519                         );
 520
 521                         // Backup current set.
 522                         autoPtr<topoSet> oldSet
 523                         (
 524                             topoSet::New
 525                             (
 526                                 setType,
 527                                 mesh,
 528                                 currentSet.name() + "_old2",
 529                                 currentSet
 530                             )
 531                         );
 532
 533                         currentSet.clear();
 534                         setSource().applyToSet(topoSetSource::NEW, currentSet);
 535
 536                         // Combine new value of currentSet with old one.
 537                         currentSet.subset(oldSet);
 538                     }
 539                     break;
 540                 }
 541                 default:
 542                 {
 543                     if (is >> sourceType)
 544                     {
 545                         autoPtr<topoSetSource> setSource
 546                         (
 547                             topoSetSource::New
 548                             (
 549                                 sourceType,
 550                                 mesh,
 551                                 is
 552                             )
 553                         );
 554
 555                         setSource().applyToSet(action, currentSet);
 556                     }
 557                 }
 558             }
 559
 560
 561             if (action != topoSetSource::LIST)
 562             {
 563                 // Set will have been modified.
 564
 565                 // Synchronize for coupled patches.
 566                 if (!noSync) currentSet.sync(mesh);
 567
 568                 // Write
 569                 if (writeVTKFile)
 570                 {
 571                     mkDir(mesh.time().path()/"VTK"/currentSet.name());
 572
 573                     fileName vtkName
 574                     (
 575                         "VTK"/currentSet.name()/currentSet.name()
 576                       + "_"
 577                       + name(mesh.time().timeIndex())
 578                       + ".vtk"
 579                     );
 580
 581                     Info<< "    Writing " << currentSet.name()
 582                         << " (size "
 583                         << returnReduce(currentSet.size(), sumOp<label>())
 584                         << ") to "
 585                         << currentSet.instance()/currentSet.local()
 586                            /currentSet.name()
 587                         << " and to vtk file " << vtkName << endl << endl;
 588
 589                     writeVTK(mesh, currentSet, vtkName);
 590                 }
 591                 else
 592                 {
 593                     Info<< "    Writing " << currentSet.name()
 594                         << " (size "
 595                         << returnReduce(currentSet.size(), sumOp<label>())
 596                         << ") to "
 597                         << currentSet.instance()/currentSet.local()
 598                            /currentSet.name() << endl << endl;
 599                 }
 600
 601                 if (writeCurrentTime)
 602                 {
 603                     currentSet.instance() = mesh.time().timeName();
 604                 }
 605                 currentSet.write();
 606             }
 607         }
 608     }
 609     catch (Foam::IOerror& fIOErr)
 610     {
 611         ok = false;
 612
 613         Pout<< fIOErr.message().c_str() << endl;
 614
 615         if (sourceType.size())
 616         {
 617             Pout<< topoSetSource::usage(sourceType).c_str();
 618         }
 619     }
 620     catch (Foam::error& fErr)
 621     {
 622         ok = false;
 623
 624         Pout<< fErr.message().c_str() << endl;
 625
 626         if (sourceType.size())
 627         {
 628             Pout<< topoSetSource::usage(sourceType).c_str();
 629         }
 630     }
 631
 632     return ok;
 633 }
 634
 635
 636 // Status returned from parsing the first token of the line
 637 enum commandStatus
 638 {
 639     QUIT,           // quit program
 640     INVALID,        // token is not a valid set manipulation command
 641     VALIDSETCMD,    // ,,    is a valid     ,,
 642     VALIDZONECMD    // ,,    is a valid     zone      ,,
 643 };
 644
 645
 646 void printMesh(const Time& runTime, const polyMesh& mesh)
 647 {
 648     Info<< "Time:" << runTime.timeName()
 649         << "  cells:" << mesh.globalData().nTotalCells()
 650         << "  faces:" << mesh.globalData().nTotalFaces()
 651         << "  points:" << mesh.globalData().nTotalPoints()
 652         << "  patches:" << mesh.boundaryMesh().size()
 653         << "  bb:" << mesh.bounds() << nl;
 654 }
 655
 656
 657 polyMesh::readUpdateState meshReadUpdate(polyMesh& mesh)
 658 {
 659     polyMesh::readUpdateState stat = mesh.readUpdate();
 660
 661     switch(stat)
 662     {
 663         case polyMesh::UNCHANGED:
 664         {
 665             Info<< "    mesh not changed." << endl;
 666             break;
 667         }
 668         case polyMesh::POINTS_MOVED:
 669         {
 670             Info<< "    points moved; topology unchanged." << endl;
 671             break;
 672         }
 673         case polyMesh::TOPO_CHANGE:
 674         {
 675             Info<< "    topology changed; patches unchanged." << nl
 676                 << "    ";
 677             printMesh(mesh.time(), mesh);
 678             break;
 679         }
 680         case polyMesh::TOPO_PATCH_CHANGE:
 681         {
 682             Info<< "    topology changed and patches changed." << nl
 683                 << "    ";
 684             printMesh(mesh.time(), mesh);
 685
 686             break;
 687         }
 688         default:
 689         {
 690             FatalErrorIn("meshReadUpdate(polyMesh&)")
 691                 << "Illegal mesh update state "
 692                 << stat  << abort(FatalError);
 693             break;
 694         }
 695     }
 696     return stat;
 697 }
 698
 699
 700 commandStatus parseType
 701 (
 702     Time& runTime,
 703     polyMesh& mesh,
 704     const word& setType,
 705     IStringStream& is
 706 )
 707 {
 708     if (setType.empty())
 709     {
 710         Info<< "Type 'help' for usage information" << endl;
 711
 712         return INVALID;
 713     }
 714     else if (setType == "help")
 715     {
 716         printHelp(Info);
 717
 718         return INVALID;
 719     }
 720     else if (setType == "list")
 721     {
 722         printAllSets(mesh, Info);
 723
 724         return INVALID;
 725     }
 726     else if (setType == "time")
 727     {
 728         scalar requestedTime = readScalar(is);
 729         instantList Times = runTime.times();
 730
 731         label nearestIndex = Time::findClosestTimeIndex(Times, requestedTime);
 732
 733         Info<< "Changing time from " << runTime.timeName()
 734             << " to " << Times[nearestIndex].name()
 735             << endl;
 736
 737         // Set time
 738         runTime.setTime(Times[nearestIndex], nearestIndex);
 739         // Optionally re-read mesh
 740         meshReadUpdate(mesh);
 741
 742         return INVALID;
 743     }
 744     else if (setType == "quit")
 745     {
 746         Info<< "Quitting ..." << endl;
 747
 748         return QUIT;
 749     }
 750     else if
 751     (
 752         setType == "cellSet"
 753      || setType == "faceSet"
 754      || setType == "pointSet"
 755     )
 756     {
 757         return VALIDSETCMD;
 758     }
 759     else if
 760     (
 761         setType == "cellZoneSet"
 762      || setType == "faceZoneSet"
 763      || setType == "pointZoneSet"
 764     )
 765     {
 766         return VALIDZONECMD;
 767     }
 768     else
 769     {
 770         SeriousErrorIn
 771         (
 772             "commandStatus parseType(Time&, polyMesh&, const word&"
 773             ", IStringStream&)"
 774         )   << "Illegal command " << setType << endl
 775             << "Should be one of 'help', 'list', 'time' or a set type :"
 776             << " 'cellSet', 'faceSet', 'pointSet', 'faceZoneSet'"
 777             << endl;
 778
 779         return INVALID;
 780     }
 781 }
 782
 783
 784 commandStatus parseAction(const word& actionName)
 785 {
 786     commandStatus stat = INVALID;
 787
 788     if (actionName.size())
 789     {
 790         try
 791         {
 792             (void)topoSetSource::toAction(actionName);
 793
 794             stat = VALIDSETCMD;
 795         }
 796         catch (Foam::IOerror& fIOErr)
 797         {
 798             stat = INVALID;
 799         }
 800         catch (Foam::error& fErr)
 801         {
 802             stat = INVALID;
 803         }
 804     }
 805     return stat;
 806 }
 807
 808
 809 // Main program:
 810
 811 int main(int argc, char *argv[])
 812 {
 813     timeSelector::addOptions(true, false);
 814 #   include "addRegionOption.H"
 815     argList::addBoolOption("noVTK", "do not write VTK files");
 816     argList::addBoolOption("loop", "execute batch commands for all timesteps");
 817     argList::addOption
 818     (
 819         "batch",
 820         "file",
 821         "process in batch mode, using input from specified file"
 822     );
 823     argList::addBoolOption
 824     (
 825         "noSync",
 826         "do not synchronise selection across coupled patches"
 827     );
 828
 829 #   include "setRootCase.H"
 830 #   include "createTime.H"
 831     instantList timeDirs = timeSelector::select0(runTime, args);
 832
 833     const bool writeVTK = !args.optionFound("noVTK");
 834     const bool loop = args.optionFound("loop");
 835     const bool batch = args.optionFound("batch");
 836     const bool noSync = args.optionFound("noSync");
 837
 838     if (loop && !batch)
 839     {
 840         FatalErrorIn(args.executable())
 841             << "Can only loop in batch mode."
 842             << exit(FatalError);
 843     }
 844
 845
 846 #   include "createNamedPolyMesh.H"
 847
 848     // Print some mesh info
 849     printMesh(runTime, mesh);
 850
 851     // Print current sets
 852     printAllSets(mesh, Info);
 853
 854     // Read history if interactive
 855 #   ifdef HAS_READLINE
 856     if (!batch && !read_history(historyFile))
 857     {
 858         Info<< "Successfully read history from " << historyFile << endl;
 859     }
 860 #   endif
 861
 862
 863     // Exit status
 864     int status = 0;
 865
 866
 867     forAll(timeDirs, timeI)
 868     {
 869         runTime.setTime(timeDirs[timeI], timeI);
 870         Info<< "Time = " << runTime.timeName() << endl;
 871
 872         // Handle geometry/topology changes
 873         meshReadUpdate(mesh);
 874
 875
 876         // Main command read & execute loop
 877         // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 878
 879         autoPtr<IFstream> fileStreamPtr(NULL);
 880
 881         if (batch)
 882         {
 883             const fileName batchFile = args["batch"];
 884
 885             Info<< "Reading commands from file " << batchFile << endl;
 886
 887             // we cannot handle .gz files
 888             if (!isFile(batchFile, false))
 889             {
 890                 FatalErrorIn(args.executable())
 891                     << "Cannot open file " << batchFile << exit(FatalError);
 892             }
 893
 894             fileStreamPtr.reset(new IFstream(batchFile));
 895         }
 896
 897         Info<< "Please type 'help', 'quit' or a set command after prompt."
 898             << endl;
 899
 900         // Whether to quit
 901         bool quit = false;
 902
 903         FatalError.throwExceptions();
 904         FatalIOError.throwExceptions();
 905
 906         do
 907         {
 908             string rawLine;
 909
 910             // Type: cellSet, faceSet, pointSet
 911             word setType;
 912             // Name of destination set.
 913             word setName;
 914             // Action (new, invert etc.)
 915             word actionName;
 916
 917             commandStatus stat = INVALID;
 918
 919             if (fileStreamPtr.valid())
 920             {
 921                 if (!fileStreamPtr().good())
 922                 {
 923                     Info<< "End of batch file" << endl;
 924                     // No error.
 925                     break;
 926                 }
 927
 928                 fileStreamPtr().getLine(rawLine);
 929
 930                 if (rawLine.size())
 931                 {
 932                     Info<< "Doing:" << rawLine << endl;
 933                 }
 934             }
 935             else
 936             {
 937 #               ifdef HAS_READLINE
 938                 {
 939                     char* linePtr = readline("readline>");
 940
 941                     if (linePtr)
 942                     {
 943                         rawLine = string(linePtr);
 944
 945                         if (*linePtr)
 946                         {
 947                             add_history(linePtr);
 948                             write_history(historyFile);
 949                         }
 950
 951                         free(linePtr);   // readline uses malloc, not new.
 952                     }
 953                     else
 954                     {
 955                         break;
 956                     }
 957                 }
 958 #               else
 959                 {
 960                     if (!std::cin.good())
 961                     {
 962                         Info<< "End of cin" << endl;
 963                         // No error.
 964                         break;
 965                     }
 966                     Info<< "Command>" << flush;
 967                     std::getline(std::cin, rawLine);
 968                 }
 969 #               endif
 970             }
 971
 972             // Strip off anything after #
 973             string::size_type i = rawLine.find_first_of("#");
 974             if (i != string::npos)
 975             {
 976                 rawLine = rawLine(0, i);
 977             }
 978
 979             if (rawLine.empty())
 980             {
 981                 continue;
 982             }
 983
 984             IStringStream is(rawLine + ' ');
 985
 986             // Type: cellSet, faceSet, pointSet, faceZoneSet
 987             is  >> setType;
 988
 989             stat = parseType(runTime, mesh, setType, is);
 990
 991             if (stat == VALIDSETCMD || stat == VALIDZONECMD)
 992             {
 993                 if (is >> setName)
 994                 {
 995                     if (is >> actionName)
 996                     {
 997                         stat = parseAction(actionName);
 998                     }
 999                 }
1000             }
1001
1002             if (stat == QUIT)
1003             {
1004                 // Make sure to quit
1005                 quit = true;
1006             }
1007             else if (stat == VALIDSETCMD || stat == VALIDZONECMD)
1008             {
1009                 bool ok = doCommand
1010                 (
1011                     mesh,
1012                     setType,
1013                     setName,
1014                     actionName,
1015                     writeVTK,
1016                     loop,   // if in looping mode dump sets to time directory
1017                     noSync,
1018                     is
1019                 );
1020
1021                 if (!ok)
1022                 {
1023                     // Exit with error.
1024                     quit = true;
1025                     status = 1;
1026                 }
1027             }
1028
1029         } while (!quit);
1030
1031         if (quit)
1032         {
1033             break;
1034         }
1035     }
1036
1037     Info<< "\nEnd\n" << endl;
1038
1039     return status;
1040 }
1041
1042
1043 // ************************************************************************* //