Files
openfoam/applications/utilities/mesh/manipulation/setSet/setSet.C
mattijs 8b7f6666ca ENH: sets are now searched for.
Some tools now output sets at the pointsInstance (if relating to geometry),
some at the facesInstance (if relating to topology).
So when loading them we now search for them from current time down
to facesInstance. They cannot be beyond facesInstance.
Adapted Time::findInstance to take optional limiting instance. Adapted
topoSet to use this. Adapted all uses of IOobjectList on sets to do the
same.
2010-01-29 10:25:13 +00:00

1014 lines
26 KiB
C

/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Description
Manipulate a cell/face/point/ set or zone interactively.
\*---------------------------------------------------------------------------*/
#include "argList.H"
#include "Time.H"
#include "polyMesh.H"
#include "globalMeshData.H"
#include "IStringStream.H"
#include "cellSet.H"
#include "faceSet.H"
#include "pointSet.H"
#include "topoSetSource.H"
#include "OFstream.H"
#include "IFstream.H"
#include "demandDrivenData.H"
#include "writePatch.H"
#include "writePointSet.H"
#include "IOobjectList.H"
#include "cellZoneSet.H"
#include "faceZoneSet.H"
#include "pointZoneSet.H"
#include <stdio.h>
#if READLINE != 0
# include <readline/readline.h>
# include <readline/history.h>
#endif
using namespace Foam;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#if READLINE != 0
static const char* historyFile = ".setSet";
#endif
Istream& selectStream(Istream* is0Ptr, Istream* is1Ptr)
{
if (is0Ptr)
{
return *is0Ptr;
}
else if (is1Ptr)
{
return *is1Ptr;
}
else
{
FatalErrorIn("selectStream(Istream*, Istream*)")
<< "No valid stream opened" << abort(FatalError);
return *is0Ptr;
}
}
// Copy set
void backup
(
const word& setType,
const polyMesh& mesh,
const word& fromName,
const topoSet& fromSet,
const word& toName
)
{
if (fromSet.size())
{
Info<< " Backing up " << fromName << " into " << toName << endl;
topoSet::New(setType, mesh, toName, fromSet)().write();
}
}
// Read and copy set
void backup
(
const word& setType,
const polyMesh& mesh,
const word& fromName,
const word& toName
)
{
autoPtr<topoSet> fromSet = topoSet::New
(
setType,
mesh,
fromName,
IOobject::READ_IF_PRESENT
);
backup(setType, mesh, fromName, fromSet(), toName);
}
// Write set to VTK readable files
void writeVTK
(
const polyMesh& mesh,
const topoSet& currentSet,
const fileName& vtkName
)
{
if (isA<faceSet>(currentSet))
{
// Faces of set with OpenFOAM faceID as value
faceList setFaces(currentSet.size());
labelList faceValues(currentSet.size());
label setFaceI = 0;
forAllConstIter(topoSet, currentSet, iter)
{
setFaces[setFaceI] = mesh.faces()[iter.key()];
faceValues[setFaceI] = iter.key();
setFaceI++;
}
primitiveFacePatch fp(setFaces, mesh.points());
writePatch
(
true,
currentSet.name(),
fp,
"faceID",
faceValues,
mesh.time().path()/vtkName
);
}
else if (isA<cellSet>(currentSet))
{
// External faces of cellset with OpenFOAM cellID as value
Map<label> cellFaces(currentSet.size());
forAllConstIter(cellSet, currentSet, iter)
{
label cellI = iter.key();
const cell& cFaces = mesh.cells()[cellI];
forAll(cFaces, i)
{
label faceI = cFaces[i];
if (mesh.isInternalFace(faceI))
{
label otherCellI = mesh.faceOwner()[faceI];
if (otherCellI == cellI)
{
otherCellI = mesh.faceNeighbour()[faceI];
}
if (!currentSet.found(otherCellI))
{
cellFaces.insert(faceI, cellI);
}
}
else
{
cellFaces.insert(faceI, cellI);
}
}
}
faceList setFaces(cellFaces.size());
labelList faceValues(cellFaces.size());
label setFaceI = 0;
forAllConstIter(Map<label>, cellFaces, iter)
{
setFaces[setFaceI] = mesh.faces()[iter.key()];
faceValues[setFaceI] = iter(); // Cell ID
setFaceI++;
}
primitiveFacePatch fp(setFaces, mesh.points());
writePatch
(
true,
currentSet.name(),
fp,
"cellID",
faceValues,
mesh.time().path()/vtkName
);
}
else if (isA<pointSet>(currentSet))
{
writePointSet
(
true,
mesh,
currentSet,
mesh.time().path()/vtkName
);
}
else
{
WarningIn
(
"void writeVTK"
"(const polyMesh& mesh, const topoSet& currentSet,"
"const fileName& vtkName)"
) << "Don't know how to handle set of type " << currentSet.type()
<< endl;
}
}
void printHelp(Ostream& os)
{
os << "Please type 'help', 'list', 'quit', 'time ddd'"
<< " or a set command after prompt." << endl
<< "'list' will show all current cell/face/point sets." << endl
<< "'time ddd' will change the current time." << endl
<< endl
<< "A set command should be of the following form" << endl
<< endl
<< " cellSet|faceSet|pointSet <setName> <action> <source>"
<< endl
<< endl
<< "The <action> is one of" << endl
<< " list - prints the contents of the set" << endl
<< " clear - clears the set" << endl
<< " invert - inverts the set" << endl
<< " remove - remove the set" << endl
<< " new <source> - sets to set to the source set" << endl
<< " add <source> - adds all elements from the source set" << endl
<< " delete <source> - deletes ,," << endl
<< " subset <source> - combines current set with the source set"
<< endl
<< endl
<< "The sources come in various forms. Type a wrong source"
<< " to see all the types available." << endl
<< endl
<< "Example: pick up all cells connected by point or face to patch"
<< " movingWall" << endl
<< endl
<< "Pick up all faces of patch:" << endl
<< " faceSet f0 new patchToFace movingWall" << endl
<< "Add faces 0,1,2:" << endl
<< " faceSet f0 add labelToFace (0 1 2)" << endl
<< "Pick up all points used by faces in faceSet f0:" << endl
<< " pointSet p0 new faceToPoint f0 all" << endl
<< "Pick up cell which has any face in f0:" << endl
<< " cellSet c0 new faceToCell f0 any" << endl
<< "Add cells which have any point in p0:" << endl
<< " cellSet c0 add pointToCell p0 any" << endl
<< "List set:" << endl
<< " cellSet c0 list" << endl
<< endl
<< "Zones can be set using zoneSets from corresponding sets:" << endl
<< " cellZoneSet c0Zone new setToCellZone c0" << endl
<< " faceZoneSet f0Zone new setToFaceZone f0" << endl
<< endl
<< "or if orientation is important:" << endl
<< " faceZoneSet f0Zone new setsToFaceZone f0 c0" << endl
<< endl
<< "ZoneSets can be manipulated using the general actions:" << endl
<< " list - prints the contents of the set" << endl
<< " clear - clears the set" << endl
<< " invert - inverts the set (undefined orientation)"
<< endl
<< " remove - remove the set" << endl
<< endl;
}
void printAllSets(const polyMesh& mesh, Ostream& os)
{
IOobjectList objects
(
mesh,
mesh.time().findInstance
(
polyMesh::meshSubDir/"sets",
word::null,
IOobject::READ_IF_PRESENT,
mesh.facesInstance()
),
polyMesh::meshSubDir/"sets"
);
IOobjectList cellSets(objects.lookupClass(cellSet::typeName));
if (cellSets.size())
{
os << "cellSets:" << endl;
forAllConstIter(IOobjectList, cellSets, iter)
{
cellSet set(*iter());
os << '\t' << set.name() << "\tsize:" << set.size() << endl;
}
}
IOobjectList faceSets(objects.lookupClass(faceSet::typeName));
if (faceSets.size())
{
os << "faceSets:" << endl;
forAllConstIter(IOobjectList, faceSets, iter)
{
faceSet set(*iter());
os << '\t' << set.name() << "\tsize:" << set.size() << endl;
}
}
IOobjectList pointSets(objects.lookupClass(pointSet::typeName));
if (pointSets.size())
{
os << "pointSets:" << endl;
forAllConstIter(IOobjectList, pointSets, iter)
{
pointSet set(*iter());
os << '\t' << set.name() << "\tsize:" << set.size() << endl;
}
}
const cellZoneMesh& cellZones = mesh.cellZones();
if (cellZones.size())
{
os << "cellZones:" << endl;
forAll(cellZones, i)
{
const cellZone& zone = cellZones[i];
os << '\t' << zone.name() << "\tsize:" << zone.size() << endl;
}
}
const faceZoneMesh& faceZones = mesh.faceZones();
if (faceZones.size())
{
os << "faceZones:" << endl;
forAll(faceZones, i)
{
const faceZone& zone = faceZones[i];
os << '\t' << zone.name() << "\tsize:" << zone.size() << endl;
}
}
const pointZoneMesh& pointZones = mesh.pointZones();
if (pointZones.size())
{
os << "pointZones:" << endl;
forAll(pointZones, i)
{
const pointZone& zone = pointZones[i];
os << '\t' << zone.name() << "\tsize:" << zone.size() << endl;
}
}
os << endl;
}
template<class ZoneType>
void removeZone
(
ZoneMesh<ZoneType, polyMesh>& zones,
const word& setName
)
{
label zoneID = zones.findZoneID(setName);
if (zoneID != -1)
{
Info<< "Removing zone " << setName << " at index " << zoneID << endl;
// Shuffle to last position
labelList oldToNew(zones.size());
label newI = 0;
forAll(oldToNew, i)
{
if (i != zoneID)
{
oldToNew[i] = newI++;
}
}
oldToNew[zoneID] = newI;
zones.reorder(oldToNew);
// Remove last element
zones.setSize(zones.size()-1);
zones.clearAddressing();
zones.write();
}
}
// Physically remove a set
void removeSet
(
const polyMesh& mesh,
const word& setType,
const word& setName
)
{
// Remove the file
IOobjectList objects
(
mesh,
mesh.time().findInstance
(
polyMesh::meshSubDir/"sets",
word::null,
IOobject::READ_IF_PRESENT,
mesh.facesInstance()
),
polyMesh::meshSubDir/"sets"
);
if (objects.found(setName))
{
// Remove file
fileName object = objects[setName]->objectPath();
Info<< "Removing file " << object << endl;
rm(object);
}
// See if zone
if (setType == cellZoneSet::typeName)
{
removeZone
(
const_cast<cellZoneMesh&>(mesh.cellZones()),
setName
);
}
else if (setType == faceZoneSet::typeName)
{
removeZone
(
const_cast<faceZoneMesh&>(mesh.faceZones()),
setName
);
}
else if (setType == pointZoneSet::typeName)
{
removeZone
(
const_cast<pointZoneMesh&>(mesh.pointZones()),
setName
);
}
}
// Read command and execute. Return true if ok, false otherwise.
bool doCommand
(
const polyMesh& mesh,
const word& setType,
const word& setName,
const word& actionName,
const bool writeVTKFile,
Istream& is
)
{
// Get some size estimate for set.
const globalMeshData& parData = mesh.globalData();
label typSize =
max
(
parData.nTotalCells(),
max
(
parData.nTotalFaces(),
parData.nTotalPoints()
)
)
/ (10*Pstream::nProcs());
bool ok = true;
// Set to work on
autoPtr<topoSet> currentSetPtr;
word sourceType;
try
{
topoSetSource::setAction action =
topoSetSource::toAction(actionName);
IOobject::readOption r;
if (action == topoSetSource::REMOVE)
{
removeSet(mesh, setType, setName);
}
else if
(
(action == topoSetSource::NEW)
|| (action == topoSetSource::CLEAR)
)
{
r = IOobject::NO_READ;
currentSetPtr = topoSet::New(setType, mesh, setName, typSize);
}
else
{
r = IOobject::MUST_READ;
currentSetPtr = topoSet::New(setType, mesh, setName, r);
topoSet& currentSet = currentSetPtr();
// Presize it according to current mesh data.
currentSet.resize(max(currentSet.size(), typSize));
}
if (currentSetPtr.valid())
{
topoSet& currentSet = currentSetPtr();
Info<< " Set:" << currentSet.name()
<< " Size:" << currentSet.size()
<< " Action:" << actionName
<< endl;
switch (action)
{
case topoSetSource::CLEAR:
{
// Already handled above by not reading
break;
}
case topoSetSource::INVERT:
{
currentSet.invert(currentSet.maxSize(mesh));
break;
}
case topoSetSource::LIST:
{
currentSet.writeDebug(Pout, mesh, 100);
Pout<< endl;
break;
}
case topoSetSource::SUBSET:
{
if (is >> sourceType)
{
autoPtr<topoSetSource> setSource
(
topoSetSource::New
(
sourceType,
mesh,
is
)
);
// Backup current set.
autoPtr<topoSet> oldSet
(
topoSet::New
(
setType,
mesh,
currentSet.name() + "_old2",
currentSet
)
);
currentSet.clear();
setSource().applyToSet(topoSetSource::NEW, currentSet);
// Combine new value of currentSet with old one.
currentSet.subset(oldSet);
}
break;
}
default:
{
if (is >> sourceType)
{
autoPtr<topoSetSource> setSource
(
topoSetSource::New
(
sourceType,
mesh,
is
)
);
setSource().applyToSet(action, currentSet);
}
}
}
if (action != topoSetSource::LIST)
{
// Set will have been modified.
// Synchronize for coupled patches.
currentSet.sync(mesh);
// Write
if (writeVTKFile)
{
mkDir(mesh.time().path()/"VTK"/currentSet.name());
fileName vtkName
(
"VTK"/currentSet.name()/currentSet.name()
+ "_"
+ name(mesh.time().timeIndex())
+ ".vtk"
);
Info<< " Writing " << currentSet.name()
<< " (size " << currentSet.size() << ") to "
<< currentSet.instance()/currentSet.local()
/currentSet.name()
<< " and to vtk file " << vtkName << endl << endl;
currentSet.write();
writeVTK(mesh, currentSet, vtkName);
}
else
{
Info<< " Writing " << currentSet.name()
<< " (size " << currentSet.size() << ") to "
<< currentSet.instance()/currentSet.local()
/currentSet.name() << endl << endl;
currentSet.write();
}
}
}
}
catch (Foam::IOerror& fIOErr)
{
ok = false;
Pout<< fIOErr.message().c_str() << endl;
if (sourceType.size())
{
Pout<< topoSetSource::usage(sourceType).c_str();
}
}
catch (Foam::error& fErr)
{
ok = false;
Pout<< fErr.message().c_str() << endl;
if (sourceType.size())
{
Pout<< topoSetSource::usage(sourceType).c_str();
}
}
return ok;
}
// Status returned from parsing the first token of the line
enum commandStatus
{
QUIT, // quit program
INVALID, // token is not a valid set manipulation command
VALIDSETCMD, // ,, is a valid ,,
VALIDZONECMD // ,, is a valid zone ,,
};
void printMesh(const Time& runTime, const polyMesh& mesh)
{
Info<< "Time:" << runTime.timeName()
<< " cells:" << mesh.nCells()
<< " faces:" << mesh.nFaces()
<< " points:" << mesh.nPoints()
<< " patches:" << mesh.boundaryMesh().size()
<< " bb:" << mesh.bounds() << nl;
}
commandStatus parseType
(
Time& runTime,
polyMesh& mesh,
const word& setType,
IStringStream& is
)
{
if (setType.empty())
{
Info<< "Type 'help' for usage information" << endl;
return INVALID;
}
else if (setType == "help")
{
printHelp(Info);
return INVALID;
}
else if (setType == "list")
{
printAllSets(mesh, Info);
return INVALID;
}
else if (setType == "time")
{
scalar requestedTime = readScalar(is);
instantList Times = runTime.times();
label nearestIndex = Time::findClosestTimeIndex(Times, requestedTime);
Info<< "Changing time from " << runTime.timeName()
<< " to " << Times[nearestIndex].name()
<< endl;
runTime.setTime(Times[nearestIndex], nearestIndex);
polyMesh::readUpdateState stat = mesh.readUpdate();
switch(stat)
{
case polyMesh::UNCHANGED:
{
Info<< " mesh not changed." << endl;
break;
}
case polyMesh::POINTS_MOVED:
{
Info<< " points moved; topology unchanged." << endl;
break;
}
case polyMesh::TOPO_CHANGE:
{
Info<< " topology changed; patches unchanged." << nl
<< " ";
printMesh(runTime, mesh);
break;
}
case polyMesh::TOPO_PATCH_CHANGE:
{
Info<< " topology changed and patches changed." << nl
<< " ";
printMesh(runTime, mesh);
break;
}
default:
{
FatalErrorIn("parseType")
<< "Illegal mesh update state "
<< stat << abort(FatalError);
break;
}
}
return INVALID;
}
else if (setType == "quit")
{
Info<< "Quitting ..." << endl;
return QUIT;
}
else if
(
setType == "cellSet"
|| setType == "faceSet"
|| setType == "pointSet"
)
{
return VALIDSETCMD;
}
else if
(
setType == "cellZoneSet"
|| setType == "faceZoneSet"
|| setType == "pointZoneSet"
)
{
return VALIDZONECMD;
}
else
{
SeriousErrorIn
(
"commandStatus parseType(Time&, polyMesh&, const word&"
", IStringStream&)"
) << "Illegal command " << setType << endl
<< "Should be one of 'help', 'list', 'time' or a set type :"
<< " 'cellSet', 'faceSet', 'pointSet', 'faceZoneSet'"
<< endl;
return INVALID;
}
}
commandStatus parseAction(const word& actionName)
{
commandStatus stat = INVALID;
if (actionName.size())
{
try
{
(void)topoSetSource::toAction(actionName);
stat = VALIDSETCMD;
}
catch (Foam::IOerror& fIOErr)
{
stat = INVALID;
}
catch (Foam::error& fErr)
{
stat = INVALID;
}
}
return stat;
}
// Main program:
int main(int argc, char *argv[])
{
# include "addRegionOption.H"
# include "addTimeOptions.H"
argList::addBoolOption("noVTK");
argList::addOption("batch", "file");
# include "setRootCase.H"
# include "createTime.H"
bool writeVTK = !args.optionFound("noVTK");
// Get times list
instantList Times = runTime.times();
# include "checkTimeOptions.H"
runTime.setTime(Times[startTime], startTime);
# include "createNamedPolyMesh.H"
// Print some mesh info
printMesh(runTime, mesh);
// Print current sets
printAllSets(mesh, Info);
std::ifstream* fileStreamPtr(NULL);
if (args.optionFound("batch"))
{
fileName batchFile(args.option("batch"));
Info<< "Reading commands from file " << batchFile << endl;
// we cannot handle .gz files
if (!isFile(batchFile, false))
{
FatalErrorIn(args.executable())
<< "Cannot open file " << batchFile << exit(FatalError);
}
fileStreamPtr = new std::ifstream(batchFile.c_str());
}
#if READLINE != 0
else if (!read_history(historyFile))
{
Info<< "Successfully read history from " << historyFile << endl;
}
#endif
Info<< "Please type 'help', 'quit' or a set command after prompt." << endl;
bool ok = true;
FatalError.throwExceptions();
FatalIOError.throwExceptions();
do
{
string rawLine;
// Type: cellSet, faceSet, pointSet
word setType;
// Name of destination set.
word setName;
// Action (new, invert etc.)
word actionName;
commandStatus stat = INVALID;
if (fileStreamPtr)
{
if (!fileStreamPtr->good())
{
Info<< "End of batch file" << endl;
break;
}
std::getline(*fileStreamPtr, rawLine);
if (rawLine.size())
{
Info<< "Doing:" << rawLine << endl;
}
}
else
{
# if READLINE != 0
{
char* linePtr = readline("readline>");
rawLine = string(linePtr);
if (*linePtr)
{
add_history(linePtr);
write_history(historyFile);
}
free(linePtr); // readline uses malloc, not new.
}
# else
{
Info<< "Command>" << flush;
std::getline(std::cin, rawLine);
}
# endif
}
if (rawLine.empty() || rawLine[0] == '#')
{
continue;
}
IStringStream is(rawLine + ' ');
// Type: cellSet, faceSet, pointSet, faceZoneSet
is >> setType;
stat = parseType(runTime, mesh, setType, is);
if (stat == VALIDSETCMD || stat == VALIDZONECMD)
{
if (is >> setName)
{
if (is >> actionName)
{
stat = parseAction(actionName);
}
}
}
ok = true;
if (stat == QUIT)
{
break;
}
else if (stat == VALIDSETCMD || stat == VALIDZONECMD)
{
ok = doCommand(mesh, setType, setName, actionName, writeVTK, is);
}
} while (ok);
if (fileStreamPtr)
{
delete fileStreamPtr;
}
Info<< "\nEnd\n" << endl;
return 0;
}
// ************************************************************************* //