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" 302 void Foam::holeToFace::writeFaces
305 const bitSet& faceFld
310 Pout<<
"Writing " << faceFld.count() <<
" faces to " << str.
name() <<
endl;
312 for (
const label facei : faceFld)
319 void Foam::holeToFace::calculateDistance
322 const bitSet& isBlockedCell,
328 if (isBlockedCell.size() != mesh_.nCells())
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);
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();
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);
729 bitSet isFrontFace(frontFaces(isSurfaceFace, locationFaces, isHoleCell));
733 label surfaceDist =
gMax(surfaceCellDist);
739 while (surfaceDist >= 0)
749 for (
const label facei : isFrontFace)
751 const label own = faceOwner[facei];
754 const label ownDist = surfaceCellDist[own];
755 if (ownDist >= surfaceDist)
758 isHoleCell.unset(own);
761 const cell& cFaces = mesh_.cells()[own];
763 const unsigned int mask = locationFaces[facei];
764 for (
const label fi : cFaces)
767 locationFaces[fi] |= mask;
771 else if (mesh_.isInternalFace(facei))
773 const label nei = faceNeighbour[facei];
774 const label neiDist = surfaceCellDist[nei];
776 if (isHoleCell[nei] && neiDist >= surfaceDist)
779 isHoleCell[nei] =
false;
782 const cell& cFaces = mesh_.cells()[nei];
784 const unsigned int mask = locationFaces[facei];
785 for (
const label fi : cFaces)
788 locationFaces[fi] |= mask;
794 reduce(nChanged, sumOp<label>());
799 Pout<<
"holeToFace::findClosure :" 800 <<
" surfaceDist:" << surfaceDist
801 <<
" front:" << isFrontFace.count()
802 <<
" nChangedCells:" << nChanged
822 bitOrEqOp<unsigned int>()
827 isFrontFace = frontFaces(isSurfaceFace, locationFaces, isHoleCell);
853 bitSet isCommonFace(mesh_.nFaces());
854 forAll(locationFaces, facei)
856 unsigned int masks = locationFaces[facei];
863 isCommonFace.set(facei);
869 for (
const label facei : isCommonFace)
871 if (surfaceFaceDist[facei] == 0)
873 isCommonFace.unset(facei);
879 Pout<<
"holeToFace::findClosure :" 880 <<
" closure faces:" << isCommonFace.count() <<
endl;
889 const bitSet& isSurfaceFace,
890 const bitSet& isSetFace
895 const labelList& faceOwner = mesh_.faceOwner();
896 const labelList& faceNeighbour = mesh_.faceNeighbour();
898 bitSet isSetCell(mesh_.nCells());
899 for (
const label facei : isSetFace)
901 isSetCell.
set(faceOwner[facei]);
902 if (mesh_.isInternalFace(facei))
904 isSetCell.set(faceNeighbour[facei]);
910 bitSet erodedSet(isSetFace);
911 for (
const label celli : isSetCell)
913 const cell& cFaces = mesh_.cells()[celli];
915 label nBlockedFaces = 0;
916 label nSurfaceFaces = 0;
917 for (
const label facei : cFaces)
919 if (erodedSet[facei])
923 else if (isSurfaceFace[facei])
929 if ((nSurfaceFaces + nBlockedFaces) == cFaces.size())
932 for (
const label facei : cFaces)
934 erodedSet.unset(facei);
942 andEqOp<unsigned int>()
946 for (
const label celli : isSetCell)
948 const cell& cFaces = mesh_.cells()[celli];
950 label nBlockedFaces = 0;
951 for (
const label facei : cFaces)
953 if (erodedSet[facei])
958 if (nBlockedFaces >= cFaces.size()-2)
960 for (
const label facei : cFaces)
962 erodedSet.flip(facei);
970 andEqOp<unsigned int>()
975 Pout<<
"holeToFace::erodeSet :" 976 <<
" starting set:" << isSetFace.count()
977 <<
" eroded set:" << erodedSet.count() <<
endl;
988 const bitSet& isSurfaceFace,
994 forAll(zonePoints_, zonei)
996 const pointField& zoneLocations = zonePoints_[zonei];
997 labelList& zoneCells = locationCells[zonei];
1001 const label celli = mesh_.findCell(zoneLocations[i]);
1002 zoneCells[i] = celli;
1008 const cell& cFaces = mesh_.cells()[celli];
1009 bool hasUnblocked =
false;
1010 for (
const label facei : cFaces)
1012 if (!isSurfaceFace[facei])
1014 hasUnblocked =
true;
1022 <<
"problem : location:" << zoneLocations[i]
1023 <<
" in zone:" << zonei
1024 <<
" is found in cell at:" << celli
1025 << mesh_.cellCentres()[celli]
1026 <<
" which is completely surrounded by blocked faces" 1033 bitSet isClosingFace
1045 isClosingFace = erodeSet(isSurfaceFace, isClosingFace);
1050 writeFaces(
"isClosingFace.obj", isClosingFace);
1054 for (
const label facei : isClosingFace)
1056 addOrDelete(
set, facei,
add);
1063 List<pointField> allPts(
pts.
size());
1085 zonePoints_(zonePoints),
1086 blockedFaceNames_(blockedFaceNames),
1087 blockedCellNames_(blockedCellNames),
1100 blockedFaceNames_(),
1101 blockedCellNames_(),
1102 erode_(
dict.getOrDefault<bool>(
"erode", false))
1106 if (!
dict.readIfPresent(
"faceSets", blockedFaceNames_))
1108 if (
dict.readEntry(
"faceSet", setName))
1110 blockedFaceNames_.
resize(1);
1111 blockedFaceNames_.
front() = std::move(setName);
1114 if (!
dict.readIfPresent(
"cellSets", blockedCellNames_))
1116 if (
dict.readEntry(
"cellSet", setName))
1118 blockedCellNames_.
resize(1);
1119 blockedCellNames_.
front() = std::move(setName);
1133 blockedFaceNames_(
one{}, word(checkIs(is))),
1134 blockedCellNames_(),
1149 for (
const word& setName : blockedFaceNames_)
1156 bitSet isCandidateCell(mesh_.nCells());
1157 if (blockedFaceNames_.size())
1159 for (
const word& setName : blockedCellNames_)
1162 isCandidateCell.setMany(loadedSet.
begin(), loadedSet.
end());
1167 isCandidateCell =
true;
1174 Info<<
" Adding all faces to disconnect regions: " 1184 Info<<
" Removing all faces to disconnect regions: " 1208 <<
" globalBlockedFaces:" << globalBlockedFaces.
localSize()
1215 const holeToFace faceSetSource
1223 faceSet closureFaceSet(
mesh,
"calcClosure", 256);
1226 const bitSet isActiveCell(
mesh.
nCells(),
true);
1228 faceSetSource.combine
1236 closureFaces = closureFaceSet.sortedToc();
1242 closureToBlocked.
clear();
1251 IndirectList<face>(
mesh.
faces(), closureFaces),
1260 EdgeMap<label> edgeMap(
pp.
nEdges());
1263 const label globalBlockedi = globalBlockedFaces.
toGlobal(i);
1264 const label facei = blockedFaces[i];
1269 edgeMap.insert(edge(
f[fp],
f[nextFp]), globalBlockedi);
1277 DynamicList<label> initialEdges(2*nBndEdges);
1278 DynamicList<edgeTopoDistanceData<label, indirectPrimitivePatch>>
1279 initialEdgesInfo(2*nBndEdges);
1282 const edge&
e = edges[edgei];
1283 const edge meshE = edge(
mp[
e[0]],
mp[
e[1]]);
1285 auto iter = edgeMap.cfind(meshE);
1291 initialEdges.append(edgei);
1292 initialEdgesInfo.append
1294 edgeTopoDistanceData<label, indirectPrimitivePatch>
1304 List<edgeTopoDistanceData<label, indirectPrimitivePatch>> allEdgeInfo
1308 List<edgeTopoDistanceData<label, indirectPrimitivePatch>>
allFaceInfo 1317 edgeTopoDistanceData<label, indirectPrimitivePatch>
1332 closureToBlocked = -1;
1337 closureToBlocked[facei] =
allFaceInfo[facei].data();
1344 UIndirectList<label>(syncFld,
pp.addressing()) = closureToBlocked;
1346 closureToBlocked = UIndirectList<label>(syncFld,
pp.addressing());
1349 List<Map<label>> compactMap;
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
const polyBoundaryMesh & pbm
void size(const label n)
Older name for setAddressableSize.
label setMany(InputIter first, InputIter last)
Set the locations listed by the iterator range, auto-vivify entries if needed.
errorManipArg< error, int > exit(error &err, const int errNo=1)
void resize(const label len)
Adjust allocated size of list.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
fileName timePath() const
Return current time path = path/timeName.
void append(const T &val)
Append an element at the end of the list.
List< edge > edgeList
List of edge.
const bitSet isBlockedFace(intersectedFaces())
Create a new set and ADD elements to it.
Add elements to current set.
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type set(const label i, bool val=true)
A bitSet::set() method for a list of bool.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
constexpr char nl
The newline '\n' character (0x0a)
bool empty() const noexcept
True if List is empty (ie, size() is zero)
T & front()
Access first element of the list, position [0].
Ostream & endl(Ostream &os)
Add newline and flush stream.
label nInternalEdges() const
Number of internal edges.
void resize_nocopy(const label len)
Adjust allocated size of list without necessarily.
label nFaces() const noexcept
Number of mesh faces.
PrimitivePatch< IndirectList< face >, const pointField & > indirectPrimitivePatch
A PrimitivePatch with an IndirectList for the faces, const reference for the point field...
virtual const fileName & name() const override
Get the name of the output serial stream. (eg, the name of the Fstream file name) ...
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Perform reduction on a copy, using specified binary operation.
List< labelList > labelListList
List of labelList.
Macros for easy insertion into run-time selection tables.
virtual void applyToSet(const topoSetSource::setAction action, topoSet &) const
Apply specified action to the topoSet.
UList< label > labelUList
A UList of labels.
const labelList & meshPoints() const
Return labelList of mesh points in patch.
addNamedToRunTimeSelectionTable(topoSetCellSource, badQualityToCell, word, badQuality)
virtual const pointField & points() const
Return raw points.
#define forAll(list, i)
Loop across all elements in list.
label fcIndex(const label i) const noexcept
The forward circular index. The next index in the list which returns to the first at the end of the l...
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
AccessType combine(const UList< T > &lists, AccessOp aop=accessOp< T >())
Combines sub-lists into a single list.
bool returnReduceAnd(const bool value, const label comm=UPstream::worldComm)
Perform logical (and) MPI Allreduce on a copy. Uses UPstream::reduceAnd.
vectorField pointField
pointField is a vectorField.
const dimensionedScalar e
Elementary charge.
void setSize(const label n)
Alias for resize()
bool mkDir(const fileName &pathName, mode_t mode=0777)
Make a directory and return an error if it could not be created.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
void clear()
Clear the list, i.e. set size to zero.
A class for handling words, derived from Foam::string.
const Time & time() const noexcept
Return time registry.
void combine(topoSet &set, const bitSet &isBlockedFace, const bitSet &isActiveCell, const bool add) const
Optional direct use to generate a faceSet.
The topoSetFaceSource is a intermediate class for handling topoSet sources for selecting faces...
setAction
Enumeration defining various actions.
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
label localSize(const label proci) const
Size of proci data.
virtual const faceList & faces() const
Return raw faces.
const polyMesh & mesh_
Reference to the mesh.
label nEdges() const
Number of edges in patch.
void add(FieldField< Field1, typename typeOfSum< Type1, Type2 >::type > &f, const FieldField< Field1, Type1 > &f1, const FieldField< Field2, Type2 > &f2)
int debug
Static debugging option.
Type gMax(const FieldField< Field, Type > &f)
defineTypeNameAndDebug(combustionModel, 0)
General set of labels of mesh quantity (points, cells, faces).
List< word > wordList
List of word.
Subtract elements from current set.
label toGlobal(const label proci, const label i) const
From local to global on proci.
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
Class with constructor to add usage string to table.
label nCells() const noexcept
Number of mesh cells.
A cell is defined as a list of faces with extra functionality.
A collection of cell labels.
unsigned int bit_count(UIntType x)
Count arbitrary number of bits (of an integral type)
messageStream Info
Information stream (stdout output on master, null elsewhere)
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.
Mesh consisting of general polyhedral cells.
List< label > labelList
A List of labels.
static autoPtr< T > New(Args &&... args)
Construct autoPtr with forwarding arguments.
holeToFace(const polyMesh &mesh, const List< pointField > &zonePoints, const wordList &blockedFaceNames, const wordList &blockedCellNames, const bool erode)
Construct from components.
List< wallPoints > allFaceInfo(mesh_.nFaces())
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
label size() const noexcept
Number of entries.
Do not request registration (bool: false)
void reduce(T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce inplace (cf. MPI Allreduce) using linear/tree communication schedule.
const dimensionedScalar mp
Proton mass.
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
addToRunTimeSelectionTable(functionObject, pointHistory, dictionary)
A class representing the concept of 1 (one) that can be used to avoid manipulating objects known to b...
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
string expand(const std::string &s, const HashTable< string > &mapping, const char sigil='$')
Expand occurrences of variables according to the mapping and return the expanded string.