53 for (
auto idx : locations)
73 if (selectedElements.
empty() || subsetMap.
empty())
80 const bitSet selected(nElems, selectedElements);
86 if (selected[subsetMap[i]])
98 if (selected[subsetMap[i]])
100 subsettedElements[
n] = i;
105 return subsettedElements;
118 <<
"Mesh is not subsetted!" <<
nl 128 void Foam::fvMeshSubset::calcFaceFlipMap()
const 131 const labelList& subToBaseCell = cellMap();
133 faceFlipMapPtr_.reset(
new labelList(subToBaseFace.size()));
134 auto& faceFlipMap = *faceFlipMapPtr_;
138 const label subInt = subMesh().nInternalFaces();
140 const labelList& subOwn = subMesh().faceOwner();
141 const labelList& own = baseMesh_.faceOwner();
143 for (label subFaceI = 0; subFaceI < subInt; ++subFaceI)
145 faceFlipMap[subFaceI] = subToBaseFace[subFaceI]+1;
147 for (label subFaceI = subInt; subFaceI < subOwn.size(); ++subFaceI)
149 const label faceI = subToBaseFace[subFaceI];
150 if (subToBaseCell[subOwn[subFaceI]] == own[faceI])
152 faceFlipMap[subFaceI] = faceI+1;
156 faceFlipMap[subFaceI] = -faceI-1;
162 void Foam::fvMeshSubset::doCoupledPatches
171 const polyBoundaryMesh& oldPatches = baseMesh().
boundaryMesh();
173 label nUncoupled = 0;
180 if (!nCellsUsingFace.empty())
182 for (
const polyPatch&
pp : oldPatches)
184 const auto* procPatch = isA<processorPolyPatch>(
pp);
188 const label nbrProci = procPatch->neighbProcNo();
190 UOPstream toNbr(nbrProci, pBufs);
192 SubList<label>(nCellsUsingFace,
pp.size(),
pp.start());
197 pBufs.finishedSends();
200 for (
const polyPatch&
pp : oldPatches)
202 const auto* procPatch = isA<processorPolyPatch>(
pp);
206 const label nbrProci = procPatch->neighbProcNo();
208 if (!pBufs.recvDataCount(nbrProci))
216 UIPstream fromNbr(nbrProci, pBufs);
217 fromNbr >> nbrCellsUsingFace;
220 if (!nCellsUsingFace.empty() && !nbrCellsUsingFace.empty())
228 nCellsUsingFace[
pp.start()+i] == 1
229 && nbrCellsUsingFace[i] == 0
234 nCellsUsingFace[
pp.start()+i] = 3;
244 for (
const polyPatch&
pp : oldPatches)
246 const cyclicPolyPatch* cpp = isA<cyclicPolyPatch>(
pp);
248 if (cpp && !nCellsUsingFace.empty())
250 const auto& cycPatch = *cpp;
254 label thisFacei = cycPatch.start() + i;
255 label otherFacei = cycPatch.transformGlobalFace(thisFacei);
259 nCellsUsingFace[thisFacei] == 1
260 && nCellsUsingFace[otherFacei] == 0
263 nCellsUsingFace[thisFacei] = 3;
272 reduce(nUncoupled, sumOp<label>());
277 Info<<
"Uncoupled " << nUncoupled <<
" faces on coupled patches. " 278 <<
"(processorPolyPatch, cyclicPolyPatch)" <<
endl;
283 void Foam::fvMeshSubset::subsetZones()
291 auto& newSubMesh = subMeshPtr_();
297 List<pointZone*>
pZones(pointZones.size());
301 const pointZone& pz = pointZones[zonei];
303 pZones[zonei] =
new pointZone
308 newSubMesh.pointZones()
319 List<faceZone*> fZones(faceZones.size());
323 const faceZone& fz = faceZones[zonei];
332 zone[fz[j]] = (fz.flipMap()[j] ? 1 : -1);
349 const label meshFacei =
faceMap()[subFacei];
350 if (zone[meshFacei] != 0)
352 subAddressing[nSub] = subFacei;
353 const label subOwner = subMesh().faceOwner()[subFacei];
354 const label baseOwner = baseMesh().
faceOwner()[meshFacei];
356 const bool sameOwner = (cellMap()[subOwner] == baseOwner);
357 const bool flip = (zone[meshFacei] == 1);
358 subFlipStatus[nSub] = (sameOwner == flip);
364 fZones[zonei] =
new faceZone
370 newSubMesh.faceZones()
379 List<cellZone*> cZones(cellZones.size());
383 const cellZone& cz = cellZones[zonei];
385 cZones[zonei] =
new cellZone
390 newSubMesh.cellZones()
395 newSubMesh.addZones(
pZones, fZones, cZones);
404 subMeshPtr_(nullptr),
405 faceFlipMapPtr_(nullptr),
424 const bitSet& selectedCells,
431 reset(selectedCells,
patchID, syncPar);
445 reset(selectedCells,
patchID, syncPar);
459 reset(selectedCells,
patchID, syncPar);
474 reset(regioni, regions,
patchID, syncPar);
482 subMeshPtr_.reset(
nullptr);
483 faceFlipMapPtr_.reset(
nullptr);
507 subMeshPtr_.reset(std::move(subMeshPtr));
508 faceFlipMapPtr_.reset(
nullptr);
510 pointMap_ = std::move(pointMap);
512 cellMap_ = std::move(cellMap);
513 patchMap_ = std::move(patchMap);
535 baseMesh_.time().timeName(),
545 auto& newSubMesh = subMeshPtr_();
550 const polyBoundaryMesh& oldBoundary = baseMesh_.boundaryMesh();
551 const polyBoundaryMesh& newBoundary = newSubMesh.boundaryMesh();
555 patchMap_ =
identity(newPatches.size());
557 forAll(newPatches, patchi)
562 oldBoundary[patchi].clone
572 newSubMesh.addFvPatches(newPatches);
583 const bitSet& selectedCells,
603 if (wantedPatchID == -1)
607 wantedPatchID = oldPatches.
findPatchID(exposedPatchName);
609 else if (wantedPatchID < 0 || wantedPatchID >= oldPatches.
size())
612 <<
"Non-existing patch index " << wantedPatchID <<
endl 613 <<
"Should be between 0 and " << oldPatches.
size()-1
623 label len = cellMap_.
size();
627 i >= 0 && (cellMap_[i] >= oldCells.
size());
633 cellMap_.resize(len);
649 label nFacesInSet = 0;
651 forAll(oldFaces, oldFacei)
653 bool faceUsed =
false;
655 if (selectedCells.
test(oldOwner[oldFacei]))
657 ++nCellsUsingFace[oldFacei];
664 && selectedCells.
test(oldNeighbour[oldFacei])
667 ++nCellsUsingFace[oldFacei];
676 faceMap_.resize(nFacesInSet);
679 doCoupledPatches(syncPar, nCellsUsingFace);
683 label oldInternalPatchID = 0;
686 label nextPatchID = oldPatches.
size();
692 label nbSize = oldPatches.
size();
694 if (wantedPatchID == -1)
700 forAll(oldPatches, patchi)
702 if (isA<processorPolyPatch>(oldPatches[patchi]))
704 nextPatchID = patchi;
707 ++oldInternalPatchID;
713 for (label oldPatchi = 0; oldPatchi < nextPatchID; oldPatchi++)
715 globalPatchMap[oldPatchi] = oldPatchi;
719 label oldPatchi = nextPatchID;
720 oldPatchi < oldPatches.
size();
724 globalPatchMap[oldPatchi] = oldPatchi+1;
729 oldInternalPatchID = wantedPatchID;
730 nextPatchID = wantedPatchID+1;
746 for (label oldFacei = 0; oldFacei < oldNInternalFaces; ++oldFacei)
748 if (nCellsUsingFace[oldFacei] == 2)
750 globalFaceMap[oldFacei] = facei;
751 faceMap_[facei++] = oldFacei;
754 markUsed(oldFaces[oldFacei], globalPointMap);
759 const label nInternalFaces = facei;
765 oldPatchi < oldPatches.
size()
766 && oldPatchi < nextPatchID;
770 const polyPatch& oldPatch = oldPatches[oldPatchi];
772 label oldFacei = oldPatch.
start();
776 if (nCellsUsingFace[oldFacei] == 1)
781 globalFaceMap[oldFacei] = facei;
782 faceMap_[facei++] = oldFacei;
785 markUsed(oldFaces[oldFacei], globalPointMap);
788 ++boundaryPatchSizes[globalPatchMap[oldPatchi]];
795 for (label oldFacei = 0; oldFacei < oldNInternalFaces; ++oldFacei)
797 if (nCellsUsingFace[oldFacei] == 1)
799 globalFaceMap[oldFacei] = facei;
800 faceMap_[facei++] = oldFacei;
803 markUsed(oldFaces[oldFacei], globalPointMap);
806 ++boundaryPatchSizes[oldInternalPatchID];
813 label oldFacei = oldNInternalFaces;
814 oldFacei < oldFaces.
size();
818 if (nCellsUsingFace[oldFacei] == 3)
820 globalFaceMap[oldFacei] = facei;
821 faceMap_[facei++] = oldFacei;
824 markUsed(oldFaces[oldFacei], globalPointMap);
827 ++boundaryPatchSizes[oldInternalPatchID];
834 label oldPatchi = nextPatchID;
835 oldPatchi < oldPatches.
size();
839 const polyPatch& oldPatch = oldPatches[oldPatchi];
841 label oldFacei = oldPatch.
start();
845 if (nCellsUsingFace[oldFacei] == 1)
850 globalFaceMap[oldFacei] = facei;
851 faceMap_[facei++] = oldFacei;
854 markUsed(oldFaces[oldFacei], globalPointMap);
857 ++boundaryPatchSizes[globalPatchMap[oldPatchi]];
863 if (facei != nFacesInSet)
871 label nPointsInSet = 0;
873 forAll(globalPointMap, pointi)
875 if (globalPointMap[pointi] != -1)
880 pointMap_.setSize(nPointsInSet);
884 forAll(globalPointMap, pointi)
886 if (globalPointMap[pointi] != -1)
888 pointMap_[nPointsInSet] = pointi;
889 globalPointMap[pointi] = nPointsInSet;
915 auto iter = newFaces.
begin();
916 const auto& renumbering = globalPointMap;
919 for (label facei = 0; facei < nInternalFaces; ++facei)
921 face& newItem = *iter;
924 const face& oldItem = oldFaces[faceMap_[facei]];
927 newItem.
resize(oldItem.size());
931 newItem[i] = renumbering[oldItem[i]];
936 for (label facei = nInternalFaces; facei < faceMap_.size(); ++facei)
938 const label oldFacei = faceMap_[facei];
940 face& newItem = *iter;
947 && selectedCells.
test(oldNeighbour[oldFacei])
948 && !selectedCells.
test(oldOwner[oldFacei])
951 ? oldFaces[oldFacei].reverseFace()
956 newItem.resize(oldItem.size());
960 newItem[i] = renumbering[oldItem[i]];
972 auto iter = newCells.
begin();
973 const auto& renumbering = globalFaceMap;
975 for (
const label oldCelli : cellMap_)
977 cell& newItem = *iter;
980 const labelList& oldItem = oldCells[oldCelli];
983 newItem.
resize(oldItem.size());
987 newItem[i] = renumbering[oldItem[i]];
1005 baseMesh_.time().timeName(),
1011 std::move(newPoints),
1012 std::move(newFaces),
1013 std::move(newCells),
1020 patchMap_.resize(nbSize);
1021 label nNewPatches = 0;
1022 label patchStart = nInternalFaces;
1030 labelList globalPatchSizes(boundaryPatchSizes);
1031 globalPatchSizes.setSize(nextPatchID);
1049 bool samePatches =
true;
1055 samePatches =
false;
1073 label oldPatchi = 0;
1074 oldPatchi < oldPatches.
size()
1075 && oldPatchi < nextPatchID;
1079 const label newSize = boundaryPatchSizes[globalPatchMap[oldPatchi]];
1081 if (oldInternalPatchID != oldPatchi)
1086 for (label patchFacei = 0; patchFacei < newSize; patchFacei++)
1088 const label facei = patchStart+patchFacei;
1089 const label oldFacei = faceMap_[facei];
1090 map[patchFacei] = oldPatches[oldPatchi].whichFace(oldFacei);
1096 oldPatches[oldPatchi].clone
1098 subMeshPtr_().boundaryMesh(),
1111 oldPatches[oldPatchi].clone
1113 subMeshPtr_().boundaryMesh(),
1121 patchStart += newSize;
1122 patchMap_[nNewPatches] = oldPatchi;
1128 if (wantedPatchID == -1)
1130 label oldInternalSize = boundaryPatchSizes[oldInternalPatchID];
1134 reduce(oldInternalSize, sumOp<label>());
1138 if (oldInternalSize > 0)
1146 boundaryPatchSizes[oldInternalPatchID],
1149 subMeshPtr_().boundaryMesh(),
1150 emptyPolyPatch::typeName
1159 patchStart += boundaryPatchSizes[oldInternalPatchID];
1160 patchMap_[nNewPatches] = -1;
1169 label oldPatchi = nextPatchID;
1170 oldPatchi < oldPatches.
size();
1174 const label newSize = boundaryPatchSizes[globalPatchMap[oldPatchi]];
1176 if (oldInternalPatchID != oldPatchi)
1181 for (label patchFacei = 0; patchFacei < newSize; patchFacei++)
1183 const label facei = patchStart+patchFacei;
1184 const label oldFacei = faceMap_[facei];
1185 map[patchFacei] = oldPatches[oldPatchi].whichFace(oldFacei);
1191 oldPatches[oldPatchi].clone
1193 subMeshPtr_().boundaryMesh(),
1206 oldPatches[oldPatchi].clone
1208 subMeshPtr_().boundaryMesh(),
1219 patchStart += newSize;
1220 patchMap_[nNewPatches] = oldPatchi;
1226 newBoundary.resize(nNewPatches);
1227 patchMap_.resize(nNewPatches);
1231 subMeshPtr_().addFvPatches(newBoundary, syncPar);
1272 const label regioni,
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
void size(const label n)
Older name for setAddressableSize.
ZoneMesh< faceZone, polyMesh > faceZoneMesh
A ZoneMesh with the type faceZone.
List< cell > cellList
List of cell.
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...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
virtual const labelList & faceNeighbour() const
Return face neighbour.
constexpr char nl
The newline '\n' character (0x0a)
label start() const noexcept
The start label of boundary faces in the polyMesh face list.
bool empty() const noexcept
True if List is empty (ie, size() is zero)
labelList sortedToc() const
The indices of the on bits as a sorted labelList.
Ostream & endl(Ostream &os)
Add newline and flush stream.
static bool & parRun() noexcept
Test if this a parallel run.
const cellList & cells() const
Ignore writing from objectRegistry::writeObject()
static int myProcNo(const label communicator=worldComm)
Rank of this process in the communicator (starting from masterNo()). Can be negative if the process i...
bool isInternalFace(const label faceIndex) const noexcept
Return true if given face label is internal to the mesh.
UList< label > labelUList
A UList of labels.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
virtual const pointField & points() const
Return raw points.
#define forAll(list, i)
Loop across all elements in list.
static void allGatherList(List< T > &values, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Gather data, but keep individual values separate. Uses linear/tree communication. ...
List< face > faceList
List of faces.
PtrList< polyPatch > polyPatchList
Store lists of polyPatch as a PtrList.
static label nProcs(const label communicator=worldComm)
Number of ranks in parallel run (for given communicator). It is 1 for serial run. ...
vectorField pointField
pointField is a vectorField.
void setSize(const label n)
Alias for resize()
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i), works like std::iota() but returning a...
bool checkHasSubMesh() const
FatalError if subset has not been performed.
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
A class for handling words, derived from Foam::string.
wordList names() const
Return a list of patch names.
label size() const noexcept
The number of entries in the list.
wordList patchNames(nPatches)
virtual const labelList & faceOwner() const
Return face owner.
void reset()
Reset subMesh and all maps. Same as clear()
label nInternalFaces() const noexcept
Number of internal faces.
ZoneMesh< pointZone, polyMesh > pointZoneMesh
A ZoneMesh with the type pointZone.
virtual const faceList & faces() const
Return raw faces.
errorManip< error > abort(error &err)
const T * set(const label i) const
Return const pointer to element (can be nullptr), or nullptr for out-of-range access (ie...
iterator begin() noexcept
Return an iterator to begin traversing the UList.
bool test(const label pos) const
Test for True value at specified position, never auto-vivify entries.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO...
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
Holds a reference to the original mesh (the baseMesh) and optionally to a subset of that mesh (the su...
const pointZoneMesh & pointZones() const noexcept
Return point zone mesh.
static word exposedPatchName
Name for exposed internal faces (default: oldInternalFaces)
fvMeshSubset(const fvMeshSubset &)=delete
No copy construct.
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
Mesh data needed to do the Finite Volume discretisation.
static labelList subsetSubset(const label nElems, const labelUList &selectedElements, const labelUList &subsetMap)
const cellZoneMesh & cellZones() const noexcept
Return cell zone mesh.
void reduce(const List< UPstream::commsStruct > &comms, T &value, const BinaryOp &bop, const int tag, const label comm)
Reduce inplace (cf. MPI Allreduce) using specified communication schedule.
"nonBlocking" : (MPI_Isend, MPI_Irecv)
messageStream Info
Information stream (stdout output on master, null elsewhere)
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
List< label > labelList
A List of labels.
static autoPtr< T > New(Args &&... args)
Construct autoPtr with forwarding arguments.
bitSet create(const label n, const labelHashSet &locations, const bool on=true)
Create a bitSet with length n with the specified on locations.
List< bool > boolList
A List of bools.
ZoneMesh< cellZone, polyMesh > cellZoneMesh
A ZoneMesh with the type cellZone.
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
IOporosityModelList pZones(mesh)
void clear()
Reset subMesh and all maps.
static void listCombineReduce(List< T > &values, const CombineOp &cop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
After completion all processors have the same data.
static constexpr const zero Zero
Global zero (0)