57 void Foam::polyTopoChange::renumberReverseMap
60 DynamicList<label>& elems
65 const label val = elems[elemI];
69 elems[elemI] = oldToNew[val];
73 const label mergedVal = -val-2;
74 elems[elemI] = -oldToNew[mergedVal]-2;
80 void Foam::polyTopoChange::renumber
88 for (
const label val : labels)
90 const label newVal = oldToNew[val];
98 labels.transfer(newSet);
103 void Foam::polyTopoChange::renumberCompact
111 for (
const label val : elems)
113 const label newVal = oldToNew[val];
117 elems[nElem++] = newVal;
120 elems.setSize(nElem);
124 void Foam::polyTopoChange::countMap
141 const label oldCelli = map[newCelli];
147 oldCelli < reverseMap.size()
148 && reverseMap[oldCelli] == newCelli
159 else if (oldCelli == -1)
171 forAll(reverseMap, oldCelli)
173 const label newCelli = reverseMap[oldCelli];
179 else if (newCelli == -1)
193 void Foam::polyTopoChange::writeMeshStats(
const polyMesh&
mesh, Ostream& os)
201 patchSizes[patchi] =
patches[patchi].size();
208 <<
" PatchSizes : " << patchSizes <<
nl
209 <<
" PatchStarts : " << patchStarts <<
nl
214 void Foam::polyTopoChange::getMergeSets
218 List<objectMap>& cellsFromCells
224 forAll(reverseCellMap, oldCelli)
226 const label newCelli = reverseCellMap[oldCelli];
230 label mergeCelli = -newCelli-2;
232 nMerged[mergeCelli]++;
237 labelList cellToMergeSet(cellMap.size(), -1);
243 if (nMerged[celli] > 1)
245 cellToMergeSet[celli] = nSets++;
255 cellsFromCells.setSize(nSets);
257 forAll(reverseCellMap, oldCelli)
259 const label newCelli = reverseCellMap[oldCelli];
263 const label mergeCelli = -newCelli-2;
267 const label setI = cellToMergeSet[mergeCelli];
269 objectMap& mergeSet = cellsFromCells[setI];
271 if (mergeSet.masterObjects().empty())
275 mergeSet.index() = mergeCelli;
276 mergeSet.masterObjects().setSize(nMerged[mergeCelli]);
279 mergeSet.masterObjects()[0] = cellMap[mergeCelli];
282 mergeSet.masterObjects()[1] = oldCelli;
284 nMerged[mergeCelli] = 2;
288 mergeSet.masterObjects()[nMerged[mergeCelli]++] = oldCelli;
295 bool Foam::polyTopoChange::hasValidPoints(
const face&
f)
const
297 for (
const label fp :
f)
299 if (fp < 0 || fp >= points_.size())
313 if (
f[fp] < 0 &&
f[fp] >= points_.size())
324 void Foam::polyTopoChange::checkFace
336 if (own == -1 && zoneI != -1)
340 else if (patchi == -1 || patchi >= nPatches_)
343 <<
"Face has no neighbour (so external) but does not have"
344 <<
" a valid patch" <<
nl
346 <<
" facei(-1 if added face):" << facei
347 <<
" own:" << own <<
" nei:" << nei
348 <<
" patchi:" << patchi <<
nl;
349 if (hasValidPoints(
f))
352 <<
"points (removed points marked with "
363 <<
"Cannot both have valid patchi and neighbour" <<
nl
365 <<
" facei(-1 if added face):" << facei
366 <<
" own:" << own <<
" nei:" << nei
367 <<
" patchi:" << patchi <<
nl;
368 if (hasValidPoints(
f))
371 <<
"points (removed points marked with "
380 <<
"Owner cell label should be less than neighbour cell label"
383 <<
" facei(-1 if added face):" << facei
384 <<
" own:" << own <<
" nei:" << nei
385 <<
" patchi:" << patchi <<
nl;
386 if (hasValidPoints(
f))
389 <<
"points (removed points marked with "
396 if (
f.size() < 3 ||
f.found(-1))
399 <<
"Illegal vertices in face"
402 <<
" facei(-1 if added face):" << facei
403 <<
" own:" << own <<
" nei:" << nei
404 <<
" patchi:" << patchi <<
nl;
405 if (hasValidPoints(
f))
408 <<
"points (removed points marked with "
413 if (facei >= 0 && facei < faces_.size() && faceRemoved(facei))
416 <<
"Face already marked for removal"
419 <<
" facei(-1 if added face):" << facei
420 <<
" own:" << own <<
" nei:" << nei
421 <<
" patchi:" << patchi <<
nl;
422 if (hasValidPoints(
f))
425 <<
"points (removed points marked with "
432 if (
f[fp] < points_.size() && pointRemoved(
f[fp]))
435 <<
"Face uses removed vertices"
438 <<
" facei(-1 if added face):" << facei
439 <<
" own:" << own <<
" nei:" << nei
440 <<
" patchi:" << patchi <<
nl;
441 if (hasValidPoints(
f))
444 <<
"points (removed points marked with "
453 void Foam::polyTopoChange::makeCells
455 const label nActiveFaces,
460 cellFaces.setSize(2*nActiveFaces);
461 cellFaceOffsets.setSize(cellMap_.size() + 1);
468 for (label facei = 0; facei < nActiveFaces; facei++)
470 if (faceOwner_[facei] < 0)
473 if (facei < faces_.size())
475 const face&
f = faces_[facei];
479 if (
f[fp] < points_.size())
481 newPoints[fp] = points_[
f[fp]];
488 <<
"Face " << facei <<
" is active but its owner has"
489 <<
" been deleted. This is usually due to deleting cells"
490 <<
" without modifying exposed faces to be boundary faces."
493 nNbrs[faceOwner_[facei]]++;
495 for (label facei = 0; facei < nActiveFaces; facei++)
497 if (faceNeighbour_[facei] >= 0)
499 nNbrs[faceNeighbour_[facei]]++;
505 cellFaceOffsets[0] = 0;
508 cellFaceOffsets[celli+1] = cellFaceOffsets[celli] + nNbrs[celli];
516 for (label facei = 0; facei < nActiveFaces; facei++)
518 label celli = faceOwner_[facei];
520 cellFaces[cellFaceOffsets[celli] + nNbrs[celli]++] = facei;
523 for (label facei = 0; facei < nActiveFaces; facei++)
525 label celli = faceNeighbour_[facei];
529 cellFaces[cellFaceOffsets[celli] + nNbrs[celli]++] = facei;
534 cellFaces.setSize(cellFaceOffsets[cellMap_.size()]);
540 void Foam::polyTopoChange::makeCellCells
542 const label nActiveFaces,
543 CompactListList<label>& cellCells
551 for (label facei = 0; facei < nActiveFaces; facei++)
553 if (faceNeighbour_[facei] >= 0)
555 nNbrs[faceOwner_[facei]]++;
556 nNbrs[faceNeighbour_[facei]]++;
569 for (label facei = 0; facei < nActiveFaces; facei++)
571 label nei = faceNeighbour_[facei];
575 label own = faceOwner_[facei];
576 cellCells.m()[cellCells.index(own, nNbrs[own]++)] = nei;
577 cellCells.m()[cellCells.index(nei, nNbrs[nei]++)] = own;
585 Foam::label Foam::polyTopoChange::getCellOrder
587 const CompactListList<label, labelList>& cellCellAddressing,
591 labelList newOrder(cellCellAddressing.size());
594 SLList<label> nextCell;
597 bitSet visited(cellCellAddressing.size());
599 label cellInOrder = 0;
604 DynamicList<label> nbrs;
606 DynamicList<label> weights;
616 label currentCell = -1;
622 if (!cellRemoved(celli) && !visited.test(celli))
624 if (cellCellAddressing[celli].size() < minWeight)
626 minWeight = cellCellAddressing[celli].size();
633 if (currentCell == -1)
643 nextCell.append(currentCell);
650 while (nextCell.size())
652 currentCell = nextCell.removeHead();
654 if (visited.set(currentCell))
659 newOrder[cellInOrder] = currentCell;
663 const labelUList neighbours = cellCellAddressing[currentCell];
671 for (
const label nbr : neighbours)
673 if (!cellRemoved(nbr) && !visited.test(nbr))
677 weights.append(cellCellAddressing[nbr].size());
683 for (
const label nbri : order)
685 nextCell.append(nbrs[nbri]);
692 newOrder.setSize(cellInOrder);
695 oldToNew =
invert(cellCellAddressing.size(), newOrder);
704 void Foam::polyTopoChange::getFaceOrder
706 const label nActiveFaces,
715 oldToNew.setSize(faceOwner_.size());
726 label startOfCell = cellFaceOffsets[celli];
727 label nFaces = cellFaceOffsets[celli+1] - startOfCell;
732 for (label i = 0; i < nFaces; i++)
734 label facei = cellFaces[startOfCell + i];
736 label nbrCelli = faceNeighbour_[facei];
738 if (facei >= nActiveFaces)
743 else if (nbrCelli != -1)
746 if (nbrCelli == celli)
748 nbrCelli = faceOwner_[facei];
751 if (celli < nbrCelli)
771 for (
const label index : order)
773 if (nbr[index] != -1)
775 oldToNew[cellFaces[startOfCell + index]] = newFacei++;
782 patchStarts.setSize(nPatches_);
784 patchSizes.setSize(nPatches_);
789 patchStarts[0] = newFacei;
791 for (label facei = 0; facei < nActiveFaces; facei++)
793 if (region_[facei] >= 0)
795 patchSizes[region_[facei]]++;
799 label facei = patchStarts[0];
801 forAll(patchStarts, patchi)
803 patchStarts[patchi] = facei;
804 facei += patchSizes[patchi];
816 for (label facei = 0; facei < nActiveFaces; facei++)
818 if (region_[facei] >= 0)
820 oldToNew[facei] = workPatchStarts[region_[facei]]++;
825 for (label facei = nActiveFaces; facei < oldToNew.size(); facei++)
827 oldToNew[facei] = facei;
833 if (oldToNew[facei] == -1)
836 <<
"Did not determine new position"
837 <<
" for face " << facei
838 <<
" owner " << faceOwner_[facei]
839 <<
" neighbour " << faceNeighbour_[facei]
840 <<
" region " << region_[facei] <<
endl
841 <<
"This is usually caused by not specifying a patch for"
842 <<
" a boundary face." <<
nl
843 <<
"Switch on the polyTopoChange::debug flag to catch"
844 <<
" this error earlier." <<
nl;
845 if (hasValidPoints(faces_[facei]))
848 <<
"points (removed points marked with "
849 <<
vector::max <<
") " << facePoints(faces_[facei]);
858 void Foam::polyTopoChange::reorderCompactFaces
865 faces_.setCapacity(newSize);
868 region_.setCapacity(newSize);
871 faceOwner_.setCapacity(newSize);
873 reorder(oldToNew, faceNeighbour_);
874 faceNeighbour_.setCapacity(newSize);
878 faceMap_.setCapacity(newSize);
880 renumberReverseMap(oldToNew, reverseFaceMap_);
882 renumberKey(oldToNew, faceFromPoint_);
883 renumberKey(oldToNew, faceFromEdge_);
885 flipFaceFlux_.setCapacity(newSize);
886 renumberKey(oldToNew, faceZone_);
888 faceZoneFlip_.setCapacity(newSize);
896 void Foam::polyTopoChange::compact
898 const bool orderCells,
899 const bool orderPoints,
900 label& nInternalPoints,
907 reversePointMap_.shrink();
912 faceNeighbour_.shrink();
914 reverseFaceMap_.shrink();
917 reverseCellMap_.shrink();
922 label nActivePoints = 0;
924 labelList localPointMap(points_.size(), -1);
929 nInternalPoints = -1;
933 if (!pointRemoved(pointi) && !retiredPoints_.found(pointi))
944 if (!pointRemoved(pointi) && !retiredPoints_.found(pointi))
956 && faceOwner_[facei] >= 0
957 && faceNeighbour_[facei] < 0
961 const face&
f = faces_[facei];
965 label pointi =
f[fp];
967 if (localPointMap[pointi] == -1)
972 || retiredPoints_.found(pointi)
976 <<
"Removed or retired point " << pointi
978 <<
" at position " << facei <<
endl
979 <<
"Probably face has not been adapted for"
989 nInternalPoints = nActivePoints - nBoundaryPoints;
992 forAll(localPointMap, pointi)
994 if (localPointMap[pointi] != -1)
996 localPointMap[pointi] += nInternalPoints;
1003 forAll(faceOwner_, facei)
1008 && faceOwner_[facei] >= 0
1009 && faceNeighbour_[facei] >= 0
1013 const face&
f = faces_[facei];
1017 label pointi =
f[fp];
1019 if (localPointMap[pointi] == -1)
1023 pointRemoved(pointi)
1024 || retiredPoints_.found(pointi)
1028 <<
"Removed or retired point " << pointi
1030 <<
" at position " << facei <<
endl
1031 <<
"Probably face has not been adapted for"
1048 for (
const label pointi : retiredPoints_)
1056 Pout<<
"Points : active:" << nActivePoints
1060 reorder(localPointMap, points_);
1064 reorder(localPointMap, pointMap_);
1066 renumberReverseMap(localPointMap, reversePointMap_);
1068 renumberKey(localPointMap, pointZone_);
1069 renumber(localPointMap, retiredPoints_);
1074 face&
f = faces_[facei];
1077 renumberCompact(localPointMap,
f);
1079 if (!faceRemoved(facei) &&
f.size() < 3)
1082 <<
"Created illegal face " <<
f
1084 <<
" at position:" << facei
1085 <<
" when filtering removed points"
1094 labelList localFaceMap(faces_.size(), -1);
1099 if (!faceRemoved(facei) && faceOwner_[facei] >= 0)
1101 localFaceMap[facei] = newFacei++;
1104 nActiveFaces_ = newFacei;
1108 if (!faceRemoved(facei) && faceOwner_[facei] < 0)
1111 localFaceMap[facei] = newFacei++;
1117 Pout<<
"Faces : active:" << nActiveFaces_
1118 <<
" removed:" << faces_.size()-newFacei <<
endl;
1122 reorderCompactFaces(newFacei, localFaceMap);
1133 CompactListList<label> cellCells;
1134 makeCellCells(nActiveFaces_, cellCells);
1137 newCelli = getCellOrder(cellCells, localCellMap);
1142 localCellMap.setSize(cellMap_.size());
1148 if (!cellRemoved(celli))
1150 localCellMap[celli] = newCelli++;
1157 Pout<<
"Cells : active:" << newCelli
1158 <<
" removed:" << cellMap_.size()-newCelli <<
endl;
1162 if (orderCells || (newCelli != cellMap_.size()))
1164 reorder(localCellMap, cellMap_);
1165 cellMap_.setCapacity(newCelli);
1166 renumberReverseMap(localCellMap, reverseCellMap_);
1168 reorder(localCellMap, cellZone_);
1169 cellZone_.setCapacity(newCelli);
1171 renumberKey(localCellMap, cellFromPoint_);
1172 renumberKey(localCellMap, cellFromEdge_);
1173 renumberKey(localCellMap, cellFromFace_);
1177 forAll(faceOwner_, facei)
1179 label own = faceOwner_[facei];
1180 label nei = faceNeighbour_[facei];
1185 faceOwner_[facei] = localCellMap[own];
1190 faceNeighbour_[facei] = localCellMap[nei];
1195 faceNeighbour_[facei] >= 0
1196 && faceNeighbour_[facei] < faceOwner_[facei]
1199 faces_[facei].flip();
1200 Swap(faceOwner_[facei], faceNeighbour_[facei]);
1201 flipFaceFlux_.flip(facei);
1202 faceZoneFlip_.flip(facei);
1209 faceNeighbour_[facei] = localCellMap[nei];
1220 makeCells(nActiveFaces_, cellFaces, cellFaceOffsets);
1236 reorderCompactFaces(localFaceMap.size(), localFaceMap);
1250 const primitiveMesh&
mesh,
1252 const bool internalFacesOnly
1259 label facei = faceLabels[i];
1273 collectedFaces = faceLabels;
1277 collectedFaces.
setSize(nFaces);
1283 label facei = faceLabels[i];
1287 collectedFaces[nFaces++] = facei;
1292 return collectedFaces;
1298 void Foam::polyTopoChange::calcPatchPointMap
1300 const UList<Map<label>>& oldPatchMeshPointMaps,
1305 patchPointMap.setSize(
boundary.size());
1311 const Map<label>& oldMeshPointMap = oldPatchMeshPointMaps[patchi];
1313 labelList& curPatchPointRnb = patchPointMap[patchi];
1315 curPatchPointRnb.
setSize(meshPoints.size());
1319 if (meshPoints[i] < pointMap_.size())
1322 curPatchPointRnb[i] = oldMeshPointMap.lookup
1324 pointMap_[meshPoints[i]],
1330 curPatchPointRnb[i] = -1;
1337 void Foam::polyTopoChange::calcFaceInflationMaps
1339 const polyMesh&
mesh,
1340 List<objectMap>& facesFromPoints,
1341 List<objectMap>& facesFromEdges,
1342 List<objectMap>& facesFromFaces
1348 facesFromPoints.setSize(faceFromPoint_.size());
1350 if (faceFromPoint_.size())
1352 label nFacesFromPoints = 0;
1357 const label facei = iter.key();
1358 const label pointi = iter.val();
1361 const bool internal = (region_[facei] == -1);
1363 facesFromPoints[nFacesFromPoints++] = objectMap
1380 facesFromEdges.setSize(faceFromEdge_.size());
1382 if (faceFromEdge_.size())
1384 label nFacesFromEdges = 0;
1389 const label facei = iter.key();
1390 const label edgei = iter.val();
1393 const bool internal = (region_[facei] == -1);
1395 facesFromEdges[nFacesFromEdges++] = objectMap
1421 void Foam::polyTopoChange::calcCellInflationMaps
1423 const polyMesh&
mesh,
1424 List<objectMap>& cellsFromPoints,
1425 List<objectMap>& cellsFromEdges,
1426 List<objectMap>& cellsFromFaces,
1427 List<objectMap>& cellsFromCells
1430 cellsFromPoints.setSize(cellFromPoint_.size());
1432 if (cellFromPoint_.size())
1434 label nCellsFromPoints = 0;
1439 const label celli = iter.key();
1440 const label pointi = iter.val();
1442 cellsFromPoints[nCellsFromPoints++] = objectMap
1451 cellsFromEdges.setSize(cellFromEdge_.size());
1453 if (cellFromEdge_.size())
1455 label nCellsFromEdges = 0;
1460 const label celli = iter.key();
1461 const label edgei = iter.val();
1463 cellsFromEdges[nCellsFromEdges++] = objectMap
1472 cellsFromFaces.setSize(cellFromFace_.size());
1474 if (cellFromFace_.size())
1476 label nCellsFromFaces = 0;
1483 const label celli = iter.key();
1484 const label oldFacei = iter.val();
1490 cellsFromFaces[nCellsFromFaces++] = objectMap
1498 cellsFromFaces[nCellsFromFaces++] = objectMap
1520 void Foam::polyTopoChange::resetZones
1522 const polyMesh&
mesh,
1542 const label pointi = iter.key();
1543 const label zonei = iter.val();
1545 if (zonei < 0 || zonei >= pointZones.size())
1548 <<
"Illegal zoneID " << zonei <<
" for point "
1549 << pointi <<
" coord " <<
mesh.
points()[pointi]
1558 forAll(addressing, zonei)
1560 addressing[zonei].setSize(
nPoints[zonei]);
1566 const label pointi = iter.key();
1567 const label zonei = iter.val();
1569 addressing[zonei][
nPoints[zonei]++] = pointi;
1572 forAll(addressing, zonei)
1579 forAll(addressing, zonei)
1581 const pointZone& oldZone = pointZones[zonei];
1582 const labelList& newZoneAddr = addressing[zonei];
1584 labelList& curPzRnb = pointZoneMap[zonei];
1585 curPzRnb.
setSize(newZoneAddr.size());
1589 if (newZoneAddr[i] < pointMap_.size())
1591 curPzRnb[i] = oldZone.whichPoint(pointMap_[newZoneAddr[i]]);
1601 newMesh.pointZones().clearAddressing();
1602 forAll(newMesh.pointZones(), zonei)
1606 Pout<<
"pointZone:" << zonei
1607 <<
" name:" << newMesh.pointZones()[zonei].
name()
1608 <<
" size:" << addressing[zonei].size()
1612 newMesh.pointZones()[zonei] = addressing[zonei];
1628 const label facei = iter.key();
1629 const label zonei = iter.val();
1631 if (zonei < 0 || zonei >= faceZones.size())
1634 <<
"Illegal zoneID " << zonei <<
" for face "
1644 forAll(addressing, zonei)
1646 addressing[zonei].setSize(nFaces[zonei]);
1647 flipMode[zonei].setSize(nFaces[zonei]);
1653 const label facei = iter.key();
1654 const label zonei = iter.val();
1656 const label index = nFaces[zonei]++;
1658 addressing[zonei][index] = facei;
1659 flipMode[zonei][index] = faceZoneFlip_.test(facei);
1663 forAll(addressing, zonei)
1667 labelList newAddressing(addressing[zonei].size());
1670 newAddressing[i] = addressing[zonei][newToOld[i]];
1672 addressing[zonei].transfer(newAddressing);
1675 boolList newFlipMode(flipMode[zonei].size());
1678 newFlipMode[i] = flipMode[zonei][newToOld[i]];
1680 flipMode[zonei].transfer(newFlipMode);
1686 forAll(addressing, zonei)
1688 const faceZone& oldZone = faceZones[zonei];
1689 const labelList& newZoneAddr = addressing[zonei];
1691 labelList& curFzFaceRnb = faceZoneFaceMap[zonei];
1693 curFzFaceRnb.
setSize(newZoneAddr.size());
1697 if (newZoneAddr[i] < faceMap_.size())
1700 oldZone.whichFace(faceMap_[newZoneAddr[i]]);
1704 curFzFaceRnb[i] = -1;
1711 newMesh.faceZones().clearAddressing();
1712 forAll(newMesh.faceZones(), zonei)
1716 Pout<<
"faceZone:" << zonei
1717 <<
" name:" << newMesh.faceZones()[zonei].
name()
1718 <<
" size:" << addressing[zonei].size()
1722 newMesh.faceZones()[zonei].resetAddressing
1742 const label zonei = cellZone_[celli];
1744 if (zonei >= cellZones.size())
1747 <<
"Illegal zoneID " << zonei <<
" for cell "
1758 forAll(addressing, zonei)
1760 addressing[zonei].setSize(nCells[zonei]);
1766 const label zonei = cellZone_[celli];
1770 addressing[zonei][nCells[zonei]++] = celli;
1774 forAll(addressing, zonei)
1781 forAll(addressing, zonei)
1783 const cellZone& oldZone = cellZones[zonei];
1784 const labelList& newZoneAddr = addressing[zonei];
1786 labelList& curCellRnb = cellZoneMap[zonei];
1788 curCellRnb.
setSize(newZoneAddr.size());
1792 if (newZoneAddr[i] < cellMap_.size())
1795 oldZone.whichCell(cellMap_[newZoneAddr[i]]);
1805 newMesh.cellZones().clearAddressing();
1806 forAll(newMesh.cellZones(), zonei)
1810 Pout<<
"cellZone:" << zonei
1811 <<
" name:" << newMesh.cellZones()[zonei].
name()
1812 <<
" size:" << addressing[zonei].size()
1816 newMesh.cellZones()[zonei] = addressing[zonei];
1822 void Foam::polyTopoChange::calcFaceZonePointMap
1824 const polyMesh&
mesh,
1825 const UList<Map<label>>& oldFaceZoneMeshPointMaps,
1831 faceZonePointMap.setSize(faceZones.size());
1835 const faceZone& newZone = faceZones[zonei];
1837 const labelList& newZoneMeshPoints = newZone().meshPoints();
1839 const Map<label>& oldZoneMeshPointMap = oldFaceZoneMeshPointMaps[zonei];
1841 labelList& curFzPointRnb = faceZonePointMap[zonei];
1843 curFzPointRnb.
setSize(newZoneMeshPoints.size());
1845 forAll(newZoneMeshPoints, pointi)
1847 if (newZoneMeshPoints[pointi] < pointMap_.size())
1849 curFzPointRnb[pointi] =
1850 oldZoneMeshPointMap.lookup
1852 pointMap_[newZoneMeshPoints[pointi]],
1858 curFzPointRnb[pointi] = -1;
1865 void Foam::polyTopoChange::reorderCoupledFaces
1867 const bool syncParallel,
1886 if (syncParallel || !isA<processorPolyPatch>(
boundary[patchi]))
1907 pBufs.finishedSends();
1912 bool anyChanged =
false;
1916 if (syncParallel || !isA<processorPolyPatch>(
boundary[patchi]))
1918 labelList patchFaceMap(patchSizes[patchi], -1);
1921 const bool changed =
boundary[patchi].order
1941 const label start = patchStarts[patchi];
1943 forAll(patchFaceMap, patchFacei)
1945 oldToNew[patchFacei + start] =
1946 start + patchFaceMap[patchFacei];
1949 forAll(patchFaceRotation, patchFacei)
1951 rotation[patchFacei + start] =
1952 patchFaceRotation[patchFacei];
1962 reduce(anyChanged, orOp<bool>());
1968 reorderCompactFaces(oldToNew.size(), oldToNew);
1973 if (rotation[facei] != 0)
1975 inplaceRotateList<List, label>(faces_[facei], rotation[facei]);
1982 void Foam::polyTopoChange::compactAndReorder
1984 const polyMesh&
mesh,
1985 const bool syncParallel,
1986 const bool orderCells,
1987 const bool orderPoints,
1989 label& nInternalPoints,
1993 List<objectMap>& pointsFromPoints,
1994 List<objectMap>& facesFromPoints,
1995 List<objectMap>& facesFromEdges,
1996 List<objectMap>& facesFromFaces,
1997 List<objectMap>& cellsFromPoints,
1998 List<objectMap>& cellsFromEdges,
1999 List<objectMap>& cellsFromFaces,
2000 List<objectMap>& cellsFromCells,
2001 List<Map<label>>& oldPatchMeshPointMaps,
2004 List<Map<label>>& oldFaceZoneMeshPointMaps
2010 <<
"polyTopoChange was constructed with a mesh with "
2011 << nPatches_ <<
" patches." <<
endl
2012 <<
"The mesh now provided has a different number of patches "
2014 <<
" which is illegal" <<
endl
2020 compact(orderCells, orderPoints, nInternalPoints, patchSizes, patchStarts);
2024 newPoints.transfer(points_);
2056 calcFaceInflationMaps
2064 calcCellInflationMaps
2075 faceFromPoint_.clearStorage();
2076 faceFromEdge_.clearStorage();
2078 cellFromPoint_.clearStorage();
2079 cellFromEdge_.clearStorage();
2080 cellFromFace_.clearStorage();
2088 oldPatchNMeshPoints.setSize(
boundary.size());
2089 oldPatchStarts.setSize(
boundary.size());
2094 oldPatchMeshPointMaps[patchi] =
boundary[patchi].meshPointMap();
2095 oldPatchNMeshPoints[patchi] =
boundary[patchi].meshPoints().size();
2096 oldPatchStarts[patchi] =
boundary[patchi].start();
2108 oldFaceZoneMeshPointMaps[zonei] = oldZone().meshPointMap();
2122 reversePointMap_(0),
2157 reversePointMap_(0),
2194 points_.clearStorage();
2195 pointMap_.clearStorage();
2196 reversePointMap_.clearStorage();
2197 pointZone_.clearStorage();
2198 retiredPoints_.clearStorage();
2200 faces_.clearStorage();
2201 region_.clearStorage();
2202 faceOwner_.clearStorage();
2203 faceNeighbour_.clearStorage();
2204 faceMap_.clearStorage();
2205 reverseFaceMap_.clearStorage();
2206 faceFromPoint_.clearStorage();
2207 faceFromEdge_.clearStorage();
2208 flipFaceFlux_.clearStorage();
2209 faceZone_.clearStorage();
2210 faceZoneFlip_.clearStorage();
2213 cellMap_.clearStorage();
2214 reverseCellMap_.clearStorage();
2215 cellZone_.clearStorage();
2216 cellFromPoint_.clearStorage();
2217 cellFromEdge_.clearStorage();
2218 cellFromFace_.clearStorage();
2231 label maxRegion = nPatches_ - 1;
2232 for (
const label regioni : patchMap)
2234 maxRegion =
max(maxRegion, regioni);
2236 nPatches_ = maxRegion + 1;
2245 points_.setCapacity(points_.size() +
points.size());
2246 pointMap_.setCapacity(pointMap_.size() +
points.size());
2247 reversePointMap_.setCapacity(reversePointMap_.size() +
points.size());
2248 pointZone_.resize(pointZone_.size() +
points.size()/100);
2253 forAll(pointZones, zonei)
2259 newZoneID[pointi] = pointZoneMap[zonei];
2264 for (label pointi = 0; pointi <
mesh.
nPoints(); pointi++)
2286 cellMap_.setCapacity(cellMap_.size() + nAllCells);
2287 reverseCellMap_.setCapacity(reverseCellMap_.size() + nAllCells);
2288 cellFromPoint_.resize(cellFromPoint_.size() + nAllCells/100);
2289 cellFromEdge_.resize(cellFromEdge_.size() + nAllCells/100);
2290 cellFromFace_.resize(cellFromFace_.size() + nAllCells/100);
2291 cellZone_.setCapacity(cellZone_.size() + nAllCells);
2299 const labelList& cellLabels = cellZones[zonei];
2301 for (
const label celli : cellLabels)
2303 if (newZoneID[celli] != -1)
2308 <<
" is in two zones:"
2309 << cellZones[newZoneID[celli]].name()
2310 <<
" and " << cellZones[zonei].name() <<
endl
2311 <<
" This is not supported."
2312 <<
" Continuing with first zone only." <<
endl;
2316 newZoneID[celli] = cellZoneMap[zonei];
2322 for (label celli = 0; celli < nAllCells; celli++)
2325 addCell(-1, -1, -1, celli, newZoneID[celli]);
2340 faces_.setCapacity(faces_.size() + nAllFaces);
2341 region_.setCapacity(region_.size() + nAllFaces);
2342 faceOwner_.setCapacity(faceOwner_.size() + nAllFaces);
2343 faceNeighbour_.setCapacity(faceNeighbour_.size() + nAllFaces);
2344 faceMap_.setCapacity(faceMap_.size() + nAllFaces);
2345 reverseFaceMap_.setCapacity(reverseFaceMap_.size() + nAllFaces);
2346 faceFromPoint_.
resize(faceFromPoint_.size() + nAllFaces/100);
2347 faceFromEdge_.resize(faceFromEdge_.size() + nAllFaces/100);
2348 flipFaceFlux_.setCapacity(faces_.size() + nAllFaces);
2349 faceZone_.resize(faceZone_.size() + nAllFaces/100);
2350 faceZoneFlip_.setCapacity(faces_.size() + nAllFaces);
2355 boolList zoneFlip(nAllFaces,
false);
2359 const labelList& faceLabels = faceZones[zonei];
2360 const boolList& flipMap = faceZones[zonei].flipMap();
2362 forAll(faceLabels, facei)
2364 newZoneID[faceLabels[facei]] = faceZoneMap[zonei];
2365 zoneFlip[faceLabels[facei]] = flipMap[facei];
2378 faceNeighbour[facei],
2394 if (pp.
start() != faces_.size())
2398 <<
"Patch " << pp.
name() <<
" starts at " << pp.
start()
2400 <<
"Current face counter at " << faces_.size() <<
endl
2401 <<
"Are patches in incremental order?"
2406 const label facei = pp.
start() + patchFacei;
2435 pointMap_.setCapacity(
nPoints);
2436 reversePointMap_.setCapacity(
nPoints);
2437 pointZone_.resize(pointZone_.size() +
nPoints/100);
2439 faces_.setCapacity(nFaces);
2440 region_.setCapacity(nFaces);
2441 faceOwner_.setCapacity(nFaces);
2442 faceNeighbour_.setCapacity(nFaces);
2443 faceMap_.setCapacity(nFaces);
2444 reverseFaceMap_.setCapacity(nFaces);
2445 faceFromPoint_.resize(faceFromPoint_.size() + nFaces/100);
2446 faceFromEdge_.resize(faceFromEdge_.size() + nFaces/100);
2447 flipFaceFlux_.setCapacity(nFaces);
2448 faceZone_.resize(faceZone_.size() + nFaces/100);
2449 faceZoneFlip_.setCapacity(nFaces);
2451 cellMap_.setCapacity(nCells);
2452 reverseCellMap_.setCapacity(nCells);
2453 cellFromPoint_.resize(cellFromPoint_.size() + nCells/100);
2454 cellFromEdge_.resize(cellFromEdge_.size() + nCells/100);
2455 cellFromFace_.resize(cellFromFace_.size() + nCells/100);
2456 cellZone_.setCapacity(nCells);
2462 if (isType<polyAddPoint>(action))
2464 const polyAddPoint& pap = refCast<const polyAddPoint>(action);
2474 else if (isType<polyModifyPoint>(action))
2488 else if (isType<polyRemovePoint>(action))
2496 else if (isType<polyAddFace>(action))
2498 const polyAddFace& paf = refCast<const polyAddFace>(action);
2514 else if (isType<polyModifyFace>(action))
2516 const polyModifyFace& pmf = refCast<const polyModifyFace>(action);
2532 else if (isType<polyRemoveFace>(action))
2534 const polyRemoveFace& prf = refCast<const polyRemoveFace>(action);
2540 else if (isType<polyAddCell>(action))
2542 const polyAddCell& pac = refCast<const polyAddCell>(action);
2553 else if (isType<polyModifyCell>(action))
2555 const polyModifyCell& pmc = refCast<const polyModifyCell>(action);
2559 modifyCell(pmc.
cellID(), -1);
2568 else if (isType<polyRemoveCell>(action))
2570 const polyRemoveCell& prc = refCast<const polyRemoveCell>(action);
2579 <<
"Unknown type of topoChange: " << action.type()
2591 const label masterPointID,
2596 const label pointi = points_.size();
2599 pointMap_.append(masterPointID);
2600 reversePointMap_.append(pointi);
2604 pointZone_.insert(pointi,
zoneID);
2609 retiredPoints_.insert(pointi);
2624 if (pointi < 0 || pointi >= points_.size())
2627 <<
"illegal point label " << pointi <<
endl
2628 <<
"Valid point labels are 0 .. " << points_.size()-1
2631 if (pointRemoved(pointi) || pointMap_[pointi] == -1)
2634 <<
"point " << pointi <<
" already marked for removal"
2637 points_[pointi] = pt;
2641 pointZone_.set(pointi,
zoneID);
2645 pointZone_.erase(pointi);
2650 retiredPoints_.erase(pointi);
2654 retiredPoints_.insert(pointi);
2661 if (newPoints.size() != points_.size())
2664 <<
"illegal pointField size." <<
endl
2665 <<
"Size:" << newPoints.size() <<
endl
2666 <<
"Points in mesh:" << points_.size()
2672 points_[pointi] = newPoints[pointi];
2680 const label mergePointi
2683 if (pointi < 0 || pointi >= points_.size())
2686 <<
"illegal point label " << pointi <<
endl
2687 <<
"Valid point labels are 0 .. " << points_.size()-1
2694 && (pointRemoved(pointi) || pointMap_[pointi] == -1)
2698 <<
"point " << pointi <<
" already marked for removal" <<
nl
2699 <<
"Point:" << points_[pointi] <<
" pointMap:" << pointMap_[pointi]
2703 if (pointi == mergePointi)
2706 <<
"Cannot remove/merge point " << pointi <<
" onto itself."
2711 pointMap_[pointi] = -1;
2712 if (mergePointi >= 0)
2714 reversePointMap_[pointi] = -mergePointi-2;
2718 reversePointMap_[pointi] = -1;
2720 pointZone_.erase(pointi);
2721 retiredPoints_.erase(pointi);
2730 const label masterPointID,
2731 const label masterEdgeID,
2732 const label masterFaceID,
2733 const bool flipFaceFlux,
2745 label facei = faces_.size();
2749 faceOwner_.append(own);
2750 faceNeighbour_.append(nei);
2752 if (masterPointID >= 0)
2754 faceMap_.append(-1);
2755 faceFromPoint_.insert(facei, masterPointID);
2757 else if (masterEdgeID >= 0)
2759 faceMap_.append(-1);
2760 faceFromEdge_.insert(facei, masterEdgeID);
2762 else if (masterFaceID >= 0)
2764 faceMap_.append(masterFaceID);
2773 faceMap_.append(-1);
2775 reverseFaceMap_.append(facei);
2777 flipFaceFlux_.set(facei, flipFaceFlux);
2781 faceZone_.insert(facei,
zoneID);
2783 faceZoneFlip_.set(facei, zoneFlip);
2795 const bool flipFaceFlux,
2808 faceOwner_[facei] = own;
2809 faceNeighbour_[facei] = nei;
2812 flipFaceFlux_.set(facei, flipFaceFlux);
2813 faceZoneFlip_.set(facei, zoneFlip);
2817 faceZone_.set(facei,
zoneID);
2821 faceZone_.erase(facei);
2829 const label mergeFacei
2832 if (facei < 0 || facei >= faces_.size())
2835 <<
"illegal face label " << facei <<
endl
2836 <<
"Valid face labels are 0 .. " << faces_.size()-1
2843 && (faceRemoved(facei) || faceMap_[facei] == -1)
2848 <<
" already marked for removal"
2852 faces_[facei].setSize(0);
2853 region_[facei] = -1;
2854 faceOwner_[facei] = -1;
2855 faceNeighbour_[facei] = -1;
2856 faceMap_[facei] = -1;
2857 if (mergeFacei >= 0)
2859 reverseFaceMap_[facei] = -mergeFacei-2;
2863 reverseFaceMap_[facei] = -1;
2865 faceFromEdge_.erase(facei);
2866 faceFromPoint_.erase(facei);
2867 flipFaceFlux_.unset(facei);
2868 faceZoneFlip_.unset(facei);
2869 faceZone_.erase(facei);
2875 const label masterPointID,
2876 const label masterEdgeID,
2877 const label masterFaceID,
2878 const label masterCellID,
2882 label celli = cellMap_.size();
2884 if (masterPointID >= 0)
2886 cellMap_.append(-1);
2887 cellFromPoint_.insert(celli, masterPointID);
2889 else if (masterEdgeID >= 0)
2891 cellMap_.append(-1);
2892 cellFromEdge_.insert(celli, masterEdgeID);
2894 else if (masterFaceID >= 0)
2896 cellMap_.append(-1);
2897 cellFromFace_.insert(celli, masterFaceID);
2901 cellMap_.append(masterCellID);
2903 reverseCellMap_.append(celli);
2904 cellZone_.append(
zoneID);
2916 cellZone_[celli] =
zoneID;
2923 const label mergeCelli
2926 if (celli < 0 || celli >= cellMap_.size())
2929 <<
"illegal cell label " << celli <<
endl
2930 <<
"Valid cell labels are 0 .. " << cellMap_.size()-1
2934 if (strict_ && cellMap_[celli] == -2)
2938 <<
" already marked for removal"
2942 cellMap_[celli] = -2;
2943 if (mergeCelli >= 0)
2945 reverseCellMap_[celli] = -mergeCelli-2;
2949 reverseCellMap_[celli] = -1;
2951 cellFromPoint_.erase(celli);
2952 cellFromEdge_.erase(celli);
2953 cellFromFace_.erase(celli);
2954 cellZone_[celli] = -1;
2962 const bool syncParallel,
2963 const bool orderCells,
2964 const bool orderPoints
2969 Pout<<
"polyTopoChange::changeMesh"
2970 <<
"(polyMesh&, const bool, const bool, const bool, const bool)"
2976 Pout<<
"Old mesh:" <<
nl;
2983 label nInternalPoints;
3022 oldPatchMeshPointMaps,
3023 oldPatchNMeshPoints,
3025 oldFaceZoneMeshPointMaps
3043 pointField renumberedMeshPoints(newPoints.size());
3047 const label oldPointi = pointMap_[
newPointi];
3092 retiredPoints_.clearStorage();
3093 region_.clearStorage();
3100 label nAdd, nInflate, nMerge, nRemove;
3101 countMap(pointMap_, reversePointMap_, nAdd, nInflate, nMerge, nRemove);
3103 <<
" added(from point):" << nAdd
3104 <<
" added(from nothing):" << nInflate
3105 <<
" merged(into other point):" << nMerge
3106 <<
" removed:" << nRemove
3109 countMap(faceMap_, reverseFaceMap_, nAdd, nInflate, nMerge, nRemove);
3111 <<
" added(from face):" << nAdd
3112 <<
" added(inflated):" << nInflate
3113 <<
" merged(into other face):" << nMerge
3114 <<
" removed:" << nRemove
3117 countMap(cellMap_, reverseCellMap_, nAdd, nInflate, nMerge, nRemove);
3119 <<
" added(from cell):" << nAdd
3120 <<
" added(inflated):" << nInflate
3121 <<
" merged(into other cell):" << nMerge
3122 <<
" removed:" << nRemove
3129 Pout<<
"New mesh:" <<
nl;
3144 resetZones(
mesh,
mesh, pointZoneMap, faceZoneFaceMap, cellZoneMap);
3148 pointZone_.clearStorage();
3149 faceZone_.clearStorage();
3150 faceZoneFlip_.clearStorage();
3151 cellZone_.clearStorage();
3161 oldPatchMeshPointMaps,
3168 calcFaceZonePointMap(
mesh, oldFaceZoneMeshPointMaps, faceZonePointMap);
3209 oldPatchNMeshPoints,