74 "\n Usage: holeToFace <faceSet> ((x0 y0 z0) (x1 y1 z1))\n\n"
75 " Select faces disconnecting the individual regions"
76 " (each indicated by a locations).\n"
302void Foam::holeToFace::writeFaces
305 const bitSet& faceFld
310 Pout<<
"Writing " << faceFld.count() <<
" faces to " << str.
name() <<
endl;
312 for (
const label facei : faceFld)
319void Foam::holeToFace::calculateDistance
322 const bitSet& isBlockedCell,
323 const bitSet& isBlockedFace,
328 if (isBlockedCell.size() != mesh_.nCells())
332 if (isBlockedFace.size() != mesh_.nFaces())
340 List<topoDistanceData<label>> cellData(mesh_.nCells());
341 List<topoDistanceData<label>> faceData(mesh_.nFaces());
344 List<topoDistanceData<label>> seedData
347 topoDistanceData<label>(0, 123)
357 for (
const label celli : isBlockedCell)
359 cellData[celli] = topoDistanceData<label>(0, 0);
363 for (
const label facei : isBlockedFace)
365 faceData[facei] = topoDistanceData<label>(0, 0);
369 FaceCellWave<topoDistanceData<label>> deltaCalc
376 mesh_.globalData().nTotalCells()+1
383 if (!isBlockedCell[celli])
385 if (!cellData[celli].valid(deltaCalc.data()))
397 cellDist[celli] = cellData[celli].distance();
404 if (!isBlockedFace[facei])
406 if (!faceData[facei].valid(deltaCalc.data()))
418 faceDist[facei] = faceData[facei].distance();
427 const bitSet& isSurfaceFace,
428 const List<unsigned int>& locationFaces,
429 const bitSet& isHoleCell
432 const labelList& faceOwner = mesh_.faceOwner();
433 const labelList& faceNeighbour = mesh_.faceNeighbour();
434 const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
436 bitSet isFrontFace(mesh_.nFaces());
437 for (label facei = 0; facei < mesh_.nInternalFaces(); facei++)
439 if (!isSurfaceFace[facei])
441 const label ownHole = isHoleCell[faceOwner[facei]];
442 const label neiHole = isHoleCell[faceNeighbour[facei]];
444 if (ownHole != neiHole)
446 unsigned int masks = locationFaces[facei];
450 <<
" at:" << mesh_.faceCentres()[facei]
459 isFrontFace.set(facei);
466 bitSet isHoleNeiCell(mesh_.nBoundaryFaces());
468 for (
const polyPatch& pp : pbm)
470 label bFacei = pp.start()-mesh_.nInternalFaces();
473 for (
const label celli : faceCells)
475 isHoleNeiCell[bFacei] = isHoleCell[celli];
483 for (label facei = mesh_.nInternalFaces(); facei < mesh_.nFaces(); facei++)
485 if (!isSurfaceFace[facei])
487 const label ownHole = isHoleCell[faceOwner[facei]];
488 const label neiHole = isHoleNeiCell[facei-mesh_.nInternalFaces()];
490 if (ownHole != neiHole)
492 unsigned int masks = locationFaces[facei];
496 <<
" at:" << mesh_.faceCentres()[facei]
505 isFrontFace.set(facei);
516 const bitSet& isSurfaceFace,
517 const bitSet& isCandidateHoleCell,
521 const labelList& faceOwner = mesh_.faceOwner();
522 const labelList& faceNeighbour = mesh_.faceNeighbour();
524 if (zonePoints_.size() < 2)
530 if (zonePoints_.size() > 31)
539 bitSet isHoleCell(isCandidateHoleCell);
540 for (
const labelList& zoneCells : locationCells)
542 for (
const label celli : zoneCells)
546 isHoleCell.unset(celli);
551 bitSet notHoleCell(isHoleCell);
556 Pout<<
"holeToFace::findClosure :"
558 <<
"holeToFace::findClosure :"
559 <<
" initial blocked faces:" << isSurfaceFace.count()
560 <<
" candidate closure cells:" << isHoleCell.count()
565 labelList surfaceCellDist(mesh_.nCells(), -1);
566 labelList surfaceNeiCellDist(mesh_.nBoundaryFaces(), -1);
567 labelList surfaceFaceDist(mesh_.nFaces(), -1);
573 bitSet(mesh_.nFaces()),
590 Pout<<
"holeToFace::findClosure :"
591 <<
" calculated topological distance to initial blocked faces."
592 <<
" max distance:" <<
gMax(surfaceCellDist)
603 List<unsigned int> locationFaces(mesh_.nFaces(), 0u);
604 forAll(locationCells, zonei)
611 const labelList& zoneCells = locationCells[zonei];
612 for (
const label celli : zoneCells)
617 const cell& cFaces = mesh_.cells()[celli];
619 UIndirectList<label>(faceDist, cFaces) = 0;
625 bitSet isSeedFace(mesh_.nFaces(), seedFaces);
630 orEqOp<unsigned int>()
632 seedFaces = isSeedFace.toc();
649 const unsigned int mask = (1u << zonei);
652 if (faceDist[facei] >= 0)
654 locationFaces[facei] |= mask;
661 writeFaces(
"isSurfaceFace.obj", isSurfaceFace);
730 bitSet isFrontFace(frontFaces(isSurfaceFace, locationFaces, isHoleCell));
734 label surfaceDist =
gMax(surfaceCellDist);
740 while (surfaceDist >= 0)
750 for (
const label facei : isFrontFace)
752 const label own = faceOwner[facei];
755 const label ownDist = surfaceCellDist[own];
756 if (ownDist >= surfaceDist)
759 isHoleCell.unset(own);
762 const cell& cFaces = mesh_.cells()[own];
764 const unsigned int mask = locationFaces[facei];
765 for (
const label fi : cFaces)
768 locationFaces[fi] |= mask;
772 else if (mesh_.isInternalFace(facei))
774 const label nei = faceNeighbour[facei];
775 const label neiDist = surfaceCellDist[nei];
777 if (isHoleCell[nei] && neiDist >= surfaceDist)
780 isHoleCell[nei] =
false;
783 const cell& cFaces = mesh_.cells()[nei];
785 const unsigned int mask = locationFaces[facei];
786 for (
const label fi : cFaces)
789 locationFaces[fi] |= mask;
795 reduce(nChanged, sumOp<label>());
800 Pout<<
"holeToFace::findClosure :"
801 <<
" surfaceDist:" << surfaceDist
802 <<
" front:" << isFrontFace.count()
803 <<
" nChangedCells:" << nChanged
823 bitOrEqOp<unsigned int>()
828 isFrontFace = frontFaces(isSurfaceFace, locationFaces, isHoleCell);
854 bitSet isCommonFace(mesh_.nFaces());
855 forAll(locationFaces, facei)
857 unsigned int masks = locationFaces[facei];
864 isCommonFace.set(facei);
870 for (
const label facei : isCommonFace)
872 if (surfaceFaceDist[facei] == 0)
874 isCommonFace.unset(facei);
880 Pout<<
"holeToFace::findClosure :"
881 <<
" closure faces:" << isCommonFace.count() <<
endl;
890 const bitSet& isSurfaceFace,
891 const bitSet& isSetFace
896 const labelList& faceOwner = mesh_.faceOwner();
897 const labelList& faceNeighbour = mesh_.faceNeighbour();
899 bitSet isSetCell(mesh_.nCells());
900 for (
const label facei : isSetFace)
902 isSetCell.set(faceOwner[facei]);
903 if (mesh_.isInternalFace(facei))
905 isSetCell.set(faceNeighbour[facei]);
911 bitSet erodedSet(isSetFace);
912 for (
const label celli : isSetCell)
914 const cell& cFaces = mesh_.cells()[celli];
916 label nBlockedFaces = 0;
917 label nSurfaceFaces = 0;
918 for (
const label facei : cFaces)
920 if (erodedSet[facei])
924 else if (isSurfaceFace[facei])
930 if ((nSurfaceFaces + nBlockedFaces) == cFaces.size())
933 for (
const label facei : cFaces)
935 erodedSet.unset(facei);
943 andEqOp<unsigned int>()
947 for (
const label celli : isSetCell)
949 const cell& cFaces = mesh_.cells()[celli];
951 label nBlockedFaces = 0;
952 for (
const label facei : cFaces)
954 if (erodedSet[facei])
959 if (nBlockedFaces >= cFaces.size()-2)
961 for (
const label facei : cFaces)
963 erodedSet.flip(facei);
971 andEqOp<unsigned int>()
976 Pout<<
"holeToFace::erodeSet :"
977 <<
" starting set:" << isSetFace.count()
978 <<
" eroded set:" << erodedSet.count() <<
endl;
989 const bitSet& isSurfaceFace,
995 forAll(zonePoints_, zonei)
997 const pointField& zoneLocations = zonePoints_[zonei];
998 labelList& zoneCells = locationCells[zonei];
1002 const label celli = mesh_.findCell(zoneLocations[i]);
1003 zoneCells[i] = celli;
1009 const cell& cFaces = mesh_.cells()[celli];
1010 bool hasUnblocked =
false;
1011 for (
const label facei : cFaces)
1013 if (!isSurfaceFace[facei])
1015 hasUnblocked =
true;
1023 <<
"problem : location:" << zoneLocations[i]
1024 <<
" in zone:" << zonei
1025 <<
" is found in cell at:" << celli
1026 << mesh_.cellCentres()[celli]
1027 <<
" which is completely surrounded by blocked faces"
1046 isClosingFace = erodeSet(isSurfaceFace, isClosingFace);
1051 writeFaces(
"isClosingFace.obj", isClosingFace);
1055 for (
const label facei : isClosingFace)
1057 addOrDelete(set, facei,
add);
1086 zonePoints_(zonePoints),
1087 blockedFaceNames_(blockedFaceNames),
1088 blockedCellNames_(blockedCellNames),
1101 blockedFaceNames_(),
1102 blockedCellNames_(),
1103 erode_(
dict.getOrDefault<
bool>(
"erode", false))
1108 blockedFaceNames_.
resize(1);
1111 blockedFaceNames_.
clear();
1116 blockedCellNames_.
resize(1);
1119 blockedCellNames_.
clear();
1133 blockedFaceNames_(
one(),
word(checkIs(is))),
1134 blockedCellNames_(),
1148 bitSet isBlockedFace(mesh_.nFaces());
1149 for (
const word& setName : blockedFaceNames_)
1151 const faceSet loadedSet(mesh_, setName);
1152 isBlockedFace.
set(loadedSet.
toc());
1156 bitSet isCandidateCell(mesh_.nCells());
1157 if (blockedFaceNames_.size())
1159 for (
const word& setName : blockedCellNames_)
1161 const cellSet loadedSet(mesh_, setName);
1162 isCandidateCell.
set(loadedSet.
toc());
1167 isCandidateCell =
true;
1174 Info<<
" Adding all faces to disconnect regions "
1178 combine(set, isBlockedFace, isCandidateCell,
true);
1184 Info<<
" Removing all faces to disconnect regions "
1188 combine(set, isBlockedFace, isCandidateCell,
false);
1208 <<
" globalBlockedFaces:" << globalBlockedFaces.
localSize()
1236 closureFaces = closureFaceSet.
sortedToc();
1242 closureToBlocked.
clear();
1263 const label globalBlockedi = globalBlockedFaces.
toGlobal(i);
1264 const label facei = blockedFaces[i];
1269 edgeMap.
insert(
edge(
f[fp],
f[nextFp]), globalBlockedi);
1282 const edge&
e = edges[edgei];
1285 auto iter = edgeMap.
cfind(meshE);
1291 initialEdges.
append(edgei);
1331 closureToBlocked.
setSize(pp.size());
1332 closureToBlocked = -1;
1333 forAll(allFaceInfo, facei)
1335 if (allFaceInfo[facei].valid(calc.data()))
1337 closureToBlocked[facei] = allFaceInfo[facei].
data();
Macros for easy insertion into run-time selection tables.
#define addNamedToRunTimeSelectionTable(baseType, thisType, argNames, lookupName)
Add to construction table with 'lookupName' as the key.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
void append(const T &val)
Copy append an element to the end of this list.
Map from edge (expressed as its endpoints) to value. For easier forward declaration it is currently i...
List< Key > sortedToc() const
The table of contents (the keys) in sorted order.
List< Key > toc() const
The table of contents (the keys) in unsorted order.
const_iterator cfind(const Key &key) const
Find and return an const_iterator set at the hashed entry.
bool insert(const Key &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
A List with indirect addressing.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
void setSize(const label n)
Alias for resize()
void append(const T &val)
Append an element at the end of the list.
void resize(const label len)
Adjust allocated size of list.
void clear()
Clear the list, i.e. set size to zero.
virtual const fileName & name() const
Get the name of the stream.
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
Wave propagation of information along patch. Every iteration information goes through one layer of fa...
A list of faces which address into the list of points.
label nEdges() const
Number of edges in patch.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
label nInternalEdges() const
Number of internal edges.
const labelList & meshPoints() const
Return labelList of mesh points in patch.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
fileName timePath() const
Return current time path.
A List with indirect addressing. Like IndirectList but does not store addressing.
T & first()
Return the first element of the list.
T * data() noexcept
Return pointer to the underlying array serving as data storage.
void size(const label n)
Older name for setAddressableSize.
label fcIndex(const label i) const noexcept
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
void set(const bitSet &bitset)
Set specified bits from another bitset.
A collection of cell labels.
A cell is defined as a list of faces with extra functionality.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
For use with PatchEdgeFaceWave. Determines topological distance to starting edges....
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
A face is a list of labels corresponding to mesh vertices.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
label localSize() const
My local size.
label toGlobal(const label i) const
From local to global index.
A topoSetFaceSource to select a set of faces that closes a hole i.e. disconnects zones (specified by ...
virtual void applyToSet(const topoSetSource::setAction action, topoSet &) const
Apply specified action to the topoSet.
static autoPtr< mapDistribute > calcClosure(const polyMesh &mesh, const List< pointField > &zonePoints, const labelList &blockedFaces, const globalIndex &globalBlockedFaces, const bool erode, labelList &closureFaces, labelList &closureToBlocked)
Optional direct use to generate the set of faces and the method to.
void combine(topoSet &set, const bitSet &isBlockedFace, const bitSet &isActiveCell, const bool add) const
Optional direct use to generate a faceSet.
const Time & time() const noexcept
Return time registry.
A class representing the concept of 1 (one) that can be used to avoid manipulating objects known to b...
Mesh consisting of general polyhedral cells.
virtual const faceList & faces() const
Return raw faces.
virtual const pointField & points() const
Return raw points.
label nCells() const noexcept
Number of mesh cells.
label nFaces() const noexcept
Number of mesh faces.
The topoSetFaceSource is a intermediate class for handling topoSet sources for selecting faces.
Class with constructor to add usage string to table.
Base class of a source for a topoSet.
setAction
Enumeration defining various actions.
@ SUBTRACT
Subtract elements from current set.
@ ADD
Add elements to current set.
@ NEW
Create a new set and ADD elements to it.
const polyMesh & mesh_
Reference to the mesh.
General set of labels of mesh quantity (points, cells, faces).
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const labelList nBndEdges(UPstream::listGatherValues< label >(aMesh.nBoundaryEdges() - nLocalProcEdges))
unsigned int bit_count(UIntType x)
Count arbitrary number of bits (of an integral type)
List< word > wordList
A List of words.
List< label > labelList
A List of labels.
bool mkDir(const fileName &pathName, mode_t mode=0777)
Make a directory and return an error if it could not be created.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
List< labelList > labelListList
A List of labelList.
void add(FieldField< Field1, typename typeOfSum< Type1, Type2 >::type > &f, const FieldField< Field1, Type1 > &f1, const FieldField< Field2, Type2 > &f2)
void reduce(const List< UPstream::commsStruct > &comms, T &value, const BinaryOp &bop, const int tag, const label comm)
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
PrimitivePatch< IndirectList< face >, const pointField & > indirectPrimitivePatch
A PrimitivePatch with an IndirectList for the faces, const reference for the point field.
errorManipArg< error, int > exit(error &err, const int errNo=1)
UList< label > labelUList
A UList of labels.
Type gMax(const FieldField< Field, Type > &f)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.