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];
145 if (reverseMap[oldCelli] == newCelli)
155 else if (oldCelli == -1)
167 forAll(reverseMap, oldCelli)
169 const label newCelli = reverseMap[oldCelli];
175 else if (newCelli == -1)
189 void Foam::polyTopoChange::writeMeshStats(
const polyMesh&
mesh, Ostream& os)
197 patchSizes[patchi] =
patches[patchi].size();
204 <<
" PatchSizes : " << patchSizes <<
nl
205 <<
" PatchStarts : " << patchStarts <<
nl
210 void Foam::polyTopoChange::getMergeSets
214 List<objectMap>& cellsFromCells
220 forAll(reverseCellMap, oldCelli)
222 const label newCelli = reverseCellMap[oldCelli];
226 label mergeCelli = -newCelli-2;
228 nMerged[mergeCelli]++;
233 labelList cellToMergeSet(cellMap.size(), -1);
239 if (nMerged[celli] > 1)
241 cellToMergeSet[celli] = nSets++;
251 cellsFromCells.setSize(nSets);
253 forAll(reverseCellMap, oldCelli)
255 const label newCelli = reverseCellMap[oldCelli];
259 const label mergeCelli = -newCelli-2;
263 const label setI = cellToMergeSet[mergeCelli];
265 objectMap& mergeSet = cellsFromCells[setI];
267 if (mergeSet.masterObjects().empty())
271 mergeSet.index() = mergeCelli;
272 mergeSet.masterObjects().setSize(nMerged[mergeCelli]);
275 mergeSet.masterObjects()[0] = cellMap[mergeCelli];
278 mergeSet.masterObjects()[1] = oldCelli;
280 nMerged[mergeCelli] = 2;
284 mergeSet.masterObjects()[nMerged[mergeCelli]++] = oldCelli;
291 bool Foam::polyTopoChange::hasValidPoints(
const face&
f)
const
293 for (
const label fp :
f)
295 if (fp < 0 || fp >= points_.size())
309 if (
f[fp] < 0 &&
f[fp] >= points_.size())
320 void Foam::polyTopoChange::checkFace
332 if (own == -1 && zoneI != -1)
336 else if (patchi == -1 || patchi >= nPatches_)
339 <<
"Face has no neighbour (so external) but does not have"
340 <<
" a valid patch" <<
nl
342 <<
" facei(-1 if added face):" << facei
343 <<
" own:" << own <<
" nei:" << nei
344 <<
" patchi:" << patchi <<
nl;
345 if (hasValidPoints(
f))
348 <<
"points (removed points marked with "
359 <<
"Cannot both have valid patchi and neighbour" <<
nl
361 <<
" facei(-1 if added face):" << facei
362 <<
" own:" << own <<
" nei:" << nei
363 <<
" patchi:" << patchi <<
nl;
364 if (hasValidPoints(
f))
367 <<
"points (removed points marked with "
376 <<
"Owner cell label should be less than neighbour cell label"
379 <<
" facei(-1 if added face):" << facei
380 <<
" own:" << own <<
" nei:" << nei
381 <<
" patchi:" << patchi <<
nl;
382 if (hasValidPoints(
f))
385 <<
"points (removed points marked with "
392 if (
f.size() < 3 ||
f.found(-1))
395 <<
"Illegal vertices in face"
398 <<
" facei(-1 if added face):" << facei
399 <<
" own:" << own <<
" nei:" << nei
400 <<
" patchi:" << patchi <<
nl;
401 if (hasValidPoints(
f))
404 <<
"points (removed points marked with "
409 if (facei >= 0 && facei < faces_.size() && faceRemoved(facei))
412 <<
"Face already marked for removal"
415 <<
" facei(-1 if added face):" << facei
416 <<
" own:" << own <<
" nei:" << nei
417 <<
" patchi:" << patchi <<
nl;
418 if (hasValidPoints(
f))
421 <<
"points (removed points marked with "
428 if (
f[fp] < points_.size() && pointRemoved(
f[fp]))
431 <<
"Face uses removed vertices"
434 <<
" facei(-1 if added face):" << facei
435 <<
" own:" << own <<
" nei:" << nei
436 <<
" patchi:" << patchi <<
nl;
437 if (hasValidPoints(
f))
440 <<
"points (removed points marked with "
449 void Foam::polyTopoChange::makeCells
451 const label nActiveFaces,
456 cellFaces.setSize(2*nActiveFaces);
457 cellFaceOffsets.setSize(cellMap_.size() + 1);
464 for (label facei = 0; facei < nActiveFaces; facei++)
466 if (faceOwner_[facei] < 0)
469 if (facei < faces_.size())
471 const face&
f = faces_[facei];
475 if (
f[fp] < points_.size())
477 newPoints[fp] = points_[
f[fp]];
484 <<
"Face " << facei <<
" is active but its owner has"
485 <<
" been deleted. This is usually due to deleting cells"
486 <<
" without modifying exposed faces to be boundary faces."
489 nNbrs[faceOwner_[facei]]++;
491 for (label facei = 0; facei < nActiveFaces; facei++)
493 if (faceNeighbour_[facei] >= 0)
495 nNbrs[faceNeighbour_[facei]]++;
501 cellFaceOffsets[0] = 0;
504 cellFaceOffsets[celli+1] = cellFaceOffsets[celli] + nNbrs[celli];
512 for (label facei = 0; facei < nActiveFaces; facei++)
514 label celli = faceOwner_[facei];
516 cellFaces[cellFaceOffsets[celli] + nNbrs[celli]++] = facei;
519 for (label facei = 0; facei < nActiveFaces; facei++)
521 label celli = faceNeighbour_[facei];
525 cellFaces[cellFaceOffsets[celli] + nNbrs[celli]++] = facei;
530 cellFaces.setSize(cellFaceOffsets[cellMap_.size()]);
536 void Foam::polyTopoChange::makeCellCells
538 const label nActiveFaces,
539 CompactListList<label>& cellCells
547 for (label facei = 0; facei < nActiveFaces; facei++)
549 if (faceNeighbour_[facei] >= 0)
551 nNbrs[faceOwner_[facei]]++;
552 nNbrs[faceNeighbour_[facei]]++;
565 for (label facei = 0; facei < nActiveFaces; facei++)
567 label nei = faceNeighbour_[facei];
571 label own = faceOwner_[facei];
572 cellCells.m()[cellCells.index(own, nNbrs[own]++)] = nei;
573 cellCells.m()[cellCells.index(nei, nNbrs[nei]++)] = own;
581 Foam::label Foam::polyTopoChange::getCellOrder
583 const CompactListList<label, labelList>& cellCellAddressing,
587 labelList newOrder(cellCellAddressing.size());
590 SLList<label> nextCell;
593 bitSet visited(cellCellAddressing.size());
595 label cellInOrder = 0;
600 DynamicList<label> nbrs;
602 DynamicList<label> weights;
612 label currentCell = -1;
618 if (!cellRemoved(celli) && !visited.test(celli))
620 if (cellCellAddressing[celli].size() < minWeight)
622 minWeight = cellCellAddressing[celli].size();
629 if (currentCell == -1)
639 nextCell.append(currentCell);
646 while (nextCell.size())
648 currentCell = nextCell.removeHead();
650 if (visited.set(currentCell))
655 newOrder[cellInOrder] = currentCell;
659 const labelUList neighbours = cellCellAddressing[currentCell];
667 for (
const label nbr : neighbours)
669 if (!cellRemoved(nbr) && !visited.test(nbr))
673 weights.append(cellCellAddressing[nbr].size());
679 for (
const label nbri : order)
681 nextCell.append(nbrs[nbri]);
688 newOrder.setSize(cellInOrder);
691 oldToNew =
invert(cellCellAddressing.size(), newOrder);
700 void Foam::polyTopoChange::getFaceOrder
702 const label nActiveFaces,
711 oldToNew.setSize(faceOwner_.size());
722 label startOfCell = cellFaceOffsets[celli];
723 label nFaces = cellFaceOffsets[celli+1] - startOfCell;
728 for (label i = 0; i < nFaces; i++)
730 label facei = cellFaces[startOfCell + i];
732 label nbrCelli = faceNeighbour_[facei];
734 if (facei >= nActiveFaces)
739 else if (nbrCelli != -1)
742 if (nbrCelli == celli)
744 nbrCelli = faceOwner_[facei];
747 if (celli < nbrCelli)
767 for (
const label index : order)
769 if (nbr[index] != -1)
771 oldToNew[cellFaces[startOfCell + index]] = newFacei++;
778 patchStarts.setSize(nPatches_);
780 patchSizes.setSize(nPatches_);
785 patchStarts[0] = newFacei;
787 for (label facei = 0; facei < nActiveFaces; facei++)
789 if (region_[facei] >= 0)
791 patchSizes[region_[facei]]++;
795 label facei = patchStarts[0];
797 forAll(patchStarts, patchi)
799 patchStarts[patchi] = facei;
800 facei += patchSizes[patchi];
812 for (label facei = 0; facei < nActiveFaces; facei++)
814 if (region_[facei] >= 0)
816 oldToNew[facei] = workPatchStarts[region_[facei]]++;
821 for (label facei = nActiveFaces; facei < oldToNew.size(); facei++)
823 oldToNew[facei] = facei;
829 if (oldToNew[facei] == -1)
832 <<
"Did not determine new position"
833 <<
" for face " << facei
834 <<
" owner " << faceOwner_[facei]
835 <<
" neighbour " << faceNeighbour_[facei]
836 <<
" region " << region_[facei] <<
endl
837 <<
"This is usually caused by not specifying a patch for"
838 <<
" a boundary face." <<
nl
839 <<
"Switch on the polyTopoChange::debug flag to catch"
840 <<
" this error earlier." <<
nl;
841 if (hasValidPoints(faces_[facei]))
844 <<
"points (removed points marked with "
845 <<
vector::max <<
") " << facePoints(faces_[facei]);
854 void Foam::polyTopoChange::reorderCompactFaces
861 faces_.setCapacity(newSize);
864 region_.setCapacity(newSize);
867 faceOwner_.setCapacity(newSize);
869 reorder(oldToNew, faceNeighbour_);
870 faceNeighbour_.setCapacity(newSize);
874 faceMap_.setCapacity(newSize);
876 renumberReverseMap(oldToNew, reverseFaceMap_);
878 renumberKey(oldToNew, faceFromPoint_);
879 renumberKey(oldToNew, faceFromEdge_);
881 flipFaceFlux_.setCapacity(newSize);
882 renumberKey(oldToNew, faceZone_);
884 faceZoneFlip_.setCapacity(newSize);
892 void Foam::polyTopoChange::compact
894 const bool orderCells,
895 const bool orderPoints,
896 label& nInternalPoints,
903 reversePointMap_.shrink();
908 faceNeighbour_.shrink();
910 reverseFaceMap_.shrink();
913 reverseCellMap_.shrink();
918 label nActivePoints = 0;
920 labelList localPointMap(points_.size(), -1);
925 nInternalPoints = -1;
929 if (!pointRemoved(pointi) && !retiredPoints_.found(pointi))
940 if (!pointRemoved(pointi) && !retiredPoints_.found(pointi))
952 && faceOwner_[facei] >= 0
953 && faceNeighbour_[facei] < 0
957 const face&
f = faces_[facei];
961 label pointi =
f[fp];
963 if (localPointMap[pointi] == -1)
968 || retiredPoints_.found(pointi)
972 <<
"Removed or retired point " << pointi
974 <<
" at position " << facei <<
endl
975 <<
"Probably face has not been adapted for"
985 nInternalPoints = nActivePoints - nBoundaryPoints;
988 forAll(localPointMap, pointi)
990 if (localPointMap[pointi] != -1)
992 localPointMap[pointi] += nInternalPoints;
1004 && faceOwner_[facei] >= 0
1005 && faceNeighbour_[facei] >= 0
1009 const face&
f = faces_[facei];
1013 label pointi =
f[fp];
1015 if (localPointMap[pointi] == -1)
1019 pointRemoved(pointi)
1020 || retiredPoints_.found(pointi)
1024 <<
"Removed or retired point " << pointi
1026 <<
" at position " << facei <<
endl
1027 <<
"Probably face has not been adapted for"
1044 for (
const label pointi : retiredPoints_)
1052 Pout<<
"Points : active:" << nActivePoints
1056 reorder(localPointMap, points_);
1060 reorder(localPointMap, pointMap_);
1062 renumberReverseMap(localPointMap, reversePointMap_);
1064 renumberKey(localPointMap, pointZone_);
1065 renumber(localPointMap, retiredPoints_);
1070 face&
f = faces_[facei];
1073 renumberCompact(localPointMap,
f);
1075 if (!faceRemoved(facei) &&
f.size() < 3)
1078 <<
"Created illegal face " <<
f
1080 <<
" at position:" << facei
1081 <<
" when filtering removed points"
1090 labelList localFaceMap(faces_.size(), -1);
1095 if (!faceRemoved(facei) && faceOwner_[facei] >= 0)
1097 localFaceMap[facei] = newFacei++;
1100 nActiveFaces_ = newFacei;
1104 if (!faceRemoved(facei) && faceOwner_[facei] < 0)
1107 localFaceMap[facei] = newFacei++;
1113 Pout<<
"Faces : active:" << nActiveFaces_
1114 <<
" removed:" << faces_.size()-newFacei <<
endl;
1118 reorderCompactFaces(newFacei, localFaceMap);
1129 CompactListList<label> cellCells;
1130 makeCellCells(nActiveFaces_, cellCells);
1133 newCelli = getCellOrder(cellCells, localCellMap);
1138 localCellMap.setSize(cellMap_.size());
1144 if (!cellRemoved(celli))
1146 localCellMap[celli] = newCelli++;
1153 Pout<<
"Cells : active:" << newCelli
1154 <<
" removed:" << cellMap_.size()-newCelli <<
endl;
1158 if (orderCells || (newCelli != cellMap_.size()))
1160 reorder(localCellMap, cellMap_);
1161 cellMap_.setCapacity(newCelli);
1162 renumberReverseMap(localCellMap, reverseCellMap_);
1164 reorder(localCellMap, cellZone_);
1165 cellZone_.setCapacity(newCelli);
1167 renumberKey(localCellMap, cellFromPoint_);
1168 renumberKey(localCellMap, cellFromEdge_);
1169 renumberKey(localCellMap, cellFromFace_);
1173 forAll(faceOwner_, facei)
1175 label own = faceOwner_[facei];
1176 label nei = faceNeighbour_[facei];
1181 faceOwner_[facei] = localCellMap[own];
1186 faceNeighbour_[facei] = localCellMap[nei];
1191 faceNeighbour_[facei] >= 0
1192 && faceNeighbour_[facei] < faceOwner_[facei]
1195 faces_[facei].flip();
1196 Swap(faceOwner_[facei], faceNeighbour_[facei]);
1197 flipFaceFlux_.flip(facei);
1198 faceZoneFlip_.flip(facei);
1205 faceNeighbour_[facei] = localCellMap[nei];
1216 makeCells(nActiveFaces_, cellFaces, cellFaceOffsets);
1232 reorderCompactFaces(localFaceMap.size(), localFaceMap);
1246 const primitiveMesh&
mesh,
1248 const bool internalFacesOnly
1255 label facei = faceLabels[i];
1269 collectedFaces = faceLabels;
1273 collectedFaces.
setSize(nFaces);
1279 label facei = faceLabels[i];
1283 collectedFaces[nFaces++] = facei;
1288 return collectedFaces;
1294 void Foam::polyTopoChange::calcPatchPointMap
1296 const UList<Map<label>>& oldPatchMeshPointMaps,
1301 patchPointMap.setSize(
boundary.size());
1307 const Map<label>& oldMeshPointMap = oldPatchMeshPointMaps[patchi];
1309 labelList& curPatchPointRnb = patchPointMap[patchi];
1311 curPatchPointRnb.
setSize(meshPoints.size());
1315 if (meshPoints[i] < pointMap_.size())
1318 curPatchPointRnb[i] = oldMeshPointMap.lookup
1320 pointMap_[meshPoints[i]],
1326 curPatchPointRnb[i] = -1;
1333 void Foam::polyTopoChange::calcFaceInflationMaps
1335 const polyMesh&
mesh,
1336 List<objectMap>& facesFromPoints,
1337 List<objectMap>& facesFromEdges,
1338 List<objectMap>& facesFromFaces
1344 facesFromPoints.setSize(faceFromPoint_.size());
1346 if (faceFromPoint_.size())
1348 label nFacesFromPoints = 0;
1353 const label facei = iter.key();
1354 const label pointi = iter.val();
1357 const bool internal = (region_[facei] == -1);
1359 facesFromPoints[nFacesFromPoints++] = objectMap
1376 facesFromEdges.setSize(faceFromEdge_.size());
1378 if (faceFromEdge_.size())
1380 label nFacesFromEdges = 0;
1385 const label facei = iter.key();
1386 const label edgei = iter.val();
1389 const bool internal = (region_[facei] == -1);
1391 facesFromEdges[nFacesFromEdges++] = objectMap
1417 void Foam::polyTopoChange::calcCellInflationMaps
1419 const polyMesh&
mesh,
1420 List<objectMap>& cellsFromPoints,
1421 List<objectMap>& cellsFromEdges,
1422 List<objectMap>& cellsFromFaces,
1423 List<objectMap>& cellsFromCells
1426 cellsFromPoints.setSize(cellFromPoint_.size());
1428 if (cellFromPoint_.size())
1430 label nCellsFromPoints = 0;
1435 const label celli = iter.key();
1436 const label pointi = iter.val();
1438 cellsFromPoints[nCellsFromPoints++] = objectMap
1447 cellsFromEdges.setSize(cellFromEdge_.size());
1449 if (cellFromEdge_.size())
1451 label nCellsFromEdges = 0;
1456 const label celli = iter.key();
1457 const label edgei = iter.val();
1459 cellsFromEdges[nCellsFromEdges++] = objectMap
1468 cellsFromFaces.setSize(cellFromFace_.size());
1470 if (cellFromFace_.size())
1472 label nCellsFromFaces = 0;
1479 const label celli = iter.key();
1480 const label oldFacei = iter.val();
1486 cellsFromFaces[nCellsFromFaces++] = objectMap
1494 cellsFromFaces[nCellsFromFaces++] = objectMap
1516 void Foam::polyTopoChange::resetZones
1518 const polyMesh&
mesh,
1538 const label pointi = iter.key();
1539 const label zonei = iter.val();
1541 if (zonei < 0 || zonei >= pointZones.size())
1544 <<
"Illegal zoneID " << zonei <<
" for point "
1545 << pointi <<
" coord " <<
mesh.
points()[pointi]
1554 forAll(addressing, zonei)
1556 addressing[zonei].setSize(
nPoints[zonei]);
1562 const label pointi = iter.key();
1563 const label zonei = iter.val();
1565 addressing[zonei][
nPoints[zonei]++] = pointi;
1568 forAll(addressing, zonei)
1575 forAll(addressing, zonei)
1577 const pointZone& oldZone = pointZones[zonei];
1578 const labelList& newZoneAddr = addressing[zonei];
1580 labelList& curPzRnb = pointZoneMap[zonei];
1581 curPzRnb.
setSize(newZoneAddr.size());
1585 if (newZoneAddr[i] < pointMap_.size())
1587 curPzRnb[i] = oldZone.whichPoint(pointMap_[newZoneAddr[i]]);
1597 newMesh.pointZones().clearAddressing();
1598 forAll(newMesh.pointZones(), zonei)
1602 Pout<<
"pointZone:" << zonei
1603 <<
" name:" << newMesh.pointZones()[zonei].
name()
1604 <<
" size:" << addressing[zonei].size()
1608 newMesh.pointZones()[zonei] = addressing[zonei];
1624 const label facei = iter.key();
1625 const label zonei = iter.val();
1627 if (zonei < 0 || zonei >= faceZones.size())
1630 <<
"Illegal zoneID " << zonei <<
" for face "
1640 forAll(addressing, zonei)
1642 addressing[zonei].setSize(nFaces[zonei]);
1643 flipMode[zonei].setSize(nFaces[zonei]);
1649 const label facei = iter.key();
1650 const label zonei = iter.val();
1652 const label index = nFaces[zonei]++;
1654 addressing[zonei][index] = facei;
1655 flipMode[zonei][index] = faceZoneFlip_.test(facei);
1659 forAll(addressing, zonei)
1663 labelList newAddressing(addressing[zonei].size());
1666 newAddressing[i] = addressing[zonei][newToOld[i]];
1668 addressing[zonei].transfer(newAddressing);
1671 boolList newFlipMode(flipMode[zonei].size());
1674 newFlipMode[i] = flipMode[zonei][newToOld[i]];
1676 flipMode[zonei].transfer(newFlipMode);
1682 forAll(addressing, zonei)
1684 const faceZone& oldZone = faceZones[zonei];
1685 const labelList& newZoneAddr = addressing[zonei];
1687 labelList& curFzFaceRnb = faceZoneFaceMap[zonei];
1689 curFzFaceRnb.
setSize(newZoneAddr.size());
1693 if (newZoneAddr[i] < faceMap_.size())
1696 oldZone.whichFace(faceMap_[newZoneAddr[i]]);
1700 curFzFaceRnb[i] = -1;
1707 newMesh.faceZones().clearAddressing();
1708 forAll(newMesh.faceZones(), zonei)
1712 Pout<<
"faceZone:" << zonei
1713 <<
" name:" << newMesh.faceZones()[zonei].
name()
1714 <<
" size:" << addressing[zonei].size()
1718 newMesh.faceZones()[zonei].resetAddressing
1738 const label zonei = cellZone_[celli];
1740 if (zonei >= cellZones.size())
1743 <<
"Illegal zoneID " << zonei <<
" for cell "
1754 forAll(addressing, zonei)
1756 addressing[zonei].setSize(nCells[zonei]);
1762 const label zonei = cellZone_[celli];
1766 addressing[zonei][nCells[zonei]++] = celli;
1770 forAll(addressing, zonei)
1777 forAll(addressing, zonei)
1779 const cellZone& oldZone = cellZones[zonei];
1780 const labelList& newZoneAddr = addressing[zonei];
1782 labelList& curCellRnb = cellZoneMap[zonei];
1784 curCellRnb.
setSize(newZoneAddr.size());
1788 if (newZoneAddr[i] < cellMap_.size())
1791 oldZone.whichCell(cellMap_[newZoneAddr[i]]);
1801 newMesh.cellZones().clearAddressing();
1802 forAll(newMesh.cellZones(), zonei)
1806 Pout<<
"cellZone:" << zonei
1807 <<
" name:" << newMesh.cellZones()[zonei].
name()
1808 <<
" size:" << addressing[zonei].size()
1812 newMesh.cellZones()[zonei] = addressing[zonei];
1818 void Foam::polyTopoChange::calcFaceZonePointMap
1820 const polyMesh&
mesh,
1821 const UList<Map<label>>& oldFaceZoneMeshPointMaps,
1827 faceZonePointMap.setSize(faceZones.size());
1831 const faceZone& newZone = faceZones[zonei];
1833 const labelList& newZoneMeshPoints = newZone().meshPoints();
1835 const Map<label>& oldZoneMeshPointMap = oldFaceZoneMeshPointMaps[zonei];
1837 labelList& curFzPointRnb = faceZonePointMap[zonei];
1839 curFzPointRnb.
setSize(newZoneMeshPoints.size());
1841 forAll(newZoneMeshPoints, pointi)
1843 if (newZoneMeshPoints[pointi] < pointMap_.size())
1845 curFzPointRnb[pointi] =
1846 oldZoneMeshPointMap.lookup
1848 pointMap_[newZoneMeshPoints[pointi]],
1854 curFzPointRnb[pointi] = -1;
1861 void Foam::polyTopoChange::reorderCoupledFaces
1863 const bool syncParallel,
1882 if (syncParallel || !isA<processorPolyPatch>(
boundary[patchi]))
1903 pBufs.finishedSends();
1908 bool anyChanged =
false;
1912 if (syncParallel || !isA<processorPolyPatch>(
boundary[patchi]))
1914 labelList patchFaceMap(patchSizes[patchi], -1);
1917 const bool changed =
boundary[patchi].order
1937 const label start = patchStarts[patchi];
1939 forAll(patchFaceMap, patchFacei)
1941 oldToNew[patchFacei + start] =
1942 start + patchFaceMap[patchFacei];
1945 forAll(patchFaceRotation, patchFacei)
1947 rotation[patchFacei + start] =
1948 patchFaceRotation[patchFacei];
1958 reduce(anyChanged, orOp<bool>());
1964 reorderCompactFaces(oldToNew.size(), oldToNew);
1969 if (rotation[facei] != 0)
1971 inplaceRotateList<List, label>(faces_[facei], rotation[facei]);
1978 void Foam::polyTopoChange::compactAndReorder
1980 const polyMesh&
mesh,
1981 const bool syncParallel,
1982 const bool orderCells,
1983 const bool orderPoints,
1985 label& nInternalPoints,
1989 List<objectMap>& pointsFromPoints,
1990 List<objectMap>& facesFromPoints,
1991 List<objectMap>& facesFromEdges,
1992 List<objectMap>& facesFromFaces,
1993 List<objectMap>& cellsFromPoints,
1994 List<objectMap>& cellsFromEdges,
1995 List<objectMap>& cellsFromFaces,
1996 List<objectMap>& cellsFromCells,
1997 List<Map<label>>& oldPatchMeshPointMaps,
2000 List<Map<label>>& oldFaceZoneMeshPointMaps
2006 <<
"polyTopoChange was constructed with a mesh with "
2007 << nPatches_ <<
" patches." <<
endl
2008 <<
"The mesh now provided has a different number of patches "
2010 <<
" which is illegal" <<
endl
2016 compact(orderCells, orderPoints, nInternalPoints, patchSizes, patchStarts);
2020 newPoints.transfer(points_);
2052 calcFaceInflationMaps
2060 calcCellInflationMaps
2071 faceFromPoint_.clearStorage();
2072 faceFromEdge_.clearStorage();
2074 cellFromPoint_.clearStorage();
2075 cellFromEdge_.clearStorage();
2076 cellFromFace_.clearStorage();
2084 oldPatchNMeshPoints.setSize(
boundary.size());
2085 oldPatchStarts.setSize(
boundary.size());
2090 oldPatchMeshPointMaps[patchi] =
boundary[patchi].meshPointMap();
2091 oldPatchNMeshPoints[patchi] =
boundary[patchi].meshPoints().size();
2092 oldPatchStarts[patchi] =
boundary[patchi].start();
2104 oldFaceZoneMeshPointMaps[zonei] = oldZone().meshPointMap();
2118 reversePointMap_(0),
2153 reversePointMap_(0),
2190 points_.clearStorage();
2191 pointMap_.clearStorage();
2192 reversePointMap_.clearStorage();
2193 pointZone_.clearStorage();
2194 retiredPoints_.clearStorage();
2196 faces_.clearStorage();
2197 region_.clearStorage();
2198 faceOwner_.clearStorage();
2199 faceNeighbour_.clearStorage();
2200 faceMap_.clearStorage();
2201 reverseFaceMap_.clearStorage();
2202 faceFromPoint_.clearStorage();
2203 faceFromEdge_.clearStorage();
2204 flipFaceFlux_.clearStorage();
2205 faceZone_.clearStorage();
2206 faceZoneFlip_.clearStorage();
2209 cellMap_.clearStorage();
2210 reverseCellMap_.clearStorage();
2211 cellZone_.clearStorage();
2212 cellFromPoint_.clearStorage();
2213 cellFromEdge_.clearStorage();
2214 cellFromFace_.clearStorage();
2227 label maxRegion = nPatches_ - 1;
2228 for (
const label regioni : patchMap)
2230 maxRegion =
max(maxRegion, regioni);
2232 nPatches_ = maxRegion + 1;
2241 points_.setCapacity(points_.size() +
points.size());
2242 pointMap_.setCapacity(pointMap_.size() +
points.size());
2243 reversePointMap_.setCapacity(reversePointMap_.size() +
points.size());
2244 pointZone_.resize(pointZone_.size() +
points.size()/100);
2249 forAll(pointZones, zonei)
2255 newZoneID[pointi] = pointZoneMap[zonei];
2260 for (label pointi = 0; pointi <
mesh.
nPoints(); pointi++)
2282 cellMap_.setCapacity(cellMap_.size() + nAllCells);
2283 reverseCellMap_.setCapacity(reverseCellMap_.size() + nAllCells);
2284 cellFromPoint_.resize(cellFromPoint_.size() + nAllCells/100);
2285 cellFromEdge_.resize(cellFromEdge_.size() + nAllCells/100);
2286 cellFromFace_.resize(cellFromFace_.size() + nAllCells/100);
2287 cellZone_.setCapacity(cellZone_.size() + nAllCells);
2295 const labelList& cellLabels = cellZones[zonei];
2297 for (
const label celli : cellLabels)
2299 if (newZoneID[celli] != -1)
2304 <<
" is in two zones:"
2305 << cellZones[newZoneID[celli]].name()
2306 <<
" and " << cellZones[zonei].name() <<
endl
2307 <<
" This is not supported."
2308 <<
" Continuing with first zone only." <<
endl;
2312 newZoneID[celli] = cellZoneMap[zonei];
2318 for (label celli = 0; celli < nAllCells; celli++)
2321 addCell(-1, -1, -1, celli, newZoneID[celli]);
2336 faces_.setCapacity(faces_.size() + nAllFaces);
2337 region_.setCapacity(region_.size() + nAllFaces);
2338 faceOwner_.setCapacity(faceOwner_.size() + nAllFaces);
2339 faceNeighbour_.setCapacity(faceNeighbour_.size() + nAllFaces);
2340 faceMap_.setCapacity(faceMap_.size() + nAllFaces);
2341 reverseFaceMap_.setCapacity(reverseFaceMap_.size() + nAllFaces);
2342 faceFromPoint_.
resize(faceFromPoint_.size() + nAllFaces/100);
2343 faceFromEdge_.resize(faceFromEdge_.size() + nAllFaces/100);
2344 flipFaceFlux_.setCapacity(faces_.size() + nAllFaces);
2345 faceZone_.resize(faceZone_.size() + nAllFaces/100);
2346 faceZoneFlip_.setCapacity(faces_.size() + nAllFaces);
2351 boolList zoneFlip(nAllFaces,
false);
2355 const labelList& faceLabels = faceZones[zonei];
2356 const boolList& flipMap = faceZones[zonei].flipMap();
2358 forAll(faceLabels, facei)
2360 newZoneID[faceLabels[facei]] = faceZoneMap[zonei];
2361 zoneFlip[faceLabels[facei]] = flipMap[facei];
2374 faceNeighbour[facei],
2390 if (pp.
start() != faces_.size())
2394 <<
"Patch " << pp.
name() <<
" starts at " << pp.
start()
2396 <<
"Current face counter at " << faces_.size() <<
endl
2397 <<
"Are patches in incremental order?"
2402 const label facei = pp.
start() + patchFacei;
2431 pointMap_.setCapacity(
nPoints);
2432 reversePointMap_.setCapacity(
nPoints);
2433 pointZone_.resize(pointZone_.size() +
nPoints/100);
2435 faces_.setCapacity(nFaces);
2436 region_.setCapacity(nFaces);
2437 faceOwner_.setCapacity(nFaces);
2438 faceNeighbour_.setCapacity(nFaces);
2439 faceMap_.setCapacity(nFaces);
2440 reverseFaceMap_.setCapacity(nFaces);
2441 faceFromPoint_.resize(faceFromPoint_.size() + nFaces/100);
2442 faceFromEdge_.resize(faceFromEdge_.size() + nFaces/100);
2443 flipFaceFlux_.setCapacity(nFaces);
2444 faceZone_.resize(faceZone_.size() + nFaces/100);
2445 faceZoneFlip_.setCapacity(nFaces);
2447 cellMap_.setCapacity(nCells);
2448 reverseCellMap_.setCapacity(nCells);
2449 cellFromPoint_.resize(cellFromPoint_.size() + nCells/100);
2450 cellFromEdge_.resize(cellFromEdge_.size() + nCells/100);
2451 cellFromFace_.resize(cellFromFace_.size() + nCells/100);
2452 cellZone_.setCapacity(nCells);
2458 if (isType<polyAddPoint>(action))
2460 const polyAddPoint& pap = refCast<const polyAddPoint>(action);
2470 else if (isType<polyModifyPoint>(action))
2484 else if (isType<polyRemovePoint>(action))
2492 else if (isType<polyAddFace>(action))
2494 const polyAddFace& paf = refCast<const polyAddFace>(action);
2510 else if (isType<polyModifyFace>(action))
2512 const polyModifyFace& pmf = refCast<const polyModifyFace>(action);
2528 else if (isType<polyRemoveFace>(action))
2530 const polyRemoveFace& prf = refCast<const polyRemoveFace>(action);
2536 else if (isType<polyAddCell>(action))
2538 const polyAddCell& pac = refCast<const polyAddCell>(action);
2549 else if (isType<polyModifyCell>(action))
2551 const polyModifyCell& pmc = refCast<const polyModifyCell>(action);
2555 modifyCell(pmc.
cellID(), -1);
2564 else if (isType<polyRemoveCell>(action))
2566 const polyRemoveCell& prc = refCast<const polyRemoveCell>(action);
2575 <<
"Unknown type of topoChange: " << action.type()
2587 const label masterPointID,
2592 const label pointi = points_.size();
2595 pointMap_.append(masterPointID);
2596 reversePointMap_.append(pointi);
2600 pointZone_.insert(pointi,
zoneID);
2605 retiredPoints_.insert(pointi);
2620 if (pointi < 0 || pointi >= points_.size())
2623 <<
"illegal point label " << pointi <<
endl
2624 <<
"Valid point labels are 0 .. " << points_.size()-1
2627 if (pointRemoved(pointi) || pointMap_[pointi] == -1)
2630 <<
"point " << pointi <<
" already marked for removal"
2633 points_[pointi] = pt;
2637 pointZone_.set(pointi,
zoneID);
2641 pointZone_.erase(pointi);
2646 retiredPoints_.erase(pointi);
2650 retiredPoints_.insert(pointi);
2657 if (newPoints.size() != points_.size())
2660 <<
"illegal pointField size." <<
endl
2661 <<
"Size:" << newPoints.size() <<
endl
2662 <<
"Points in mesh:" << points_.size()
2668 points_[pointi] = newPoints[pointi];
2676 const label mergePointi
2679 if (pointi < 0 || pointi >= points_.size())
2682 <<
"illegal point label " << pointi <<
endl
2683 <<
"Valid point labels are 0 .. " << points_.size()-1
2690 && (pointRemoved(pointi) || pointMap_[pointi] == -1)
2694 <<
"point " << pointi <<
" already marked for removal" <<
nl
2695 <<
"Point:" << points_[pointi] <<
" pointMap:" << pointMap_[pointi]
2699 if (pointi == mergePointi)
2702 <<
"Cannot remove/merge point " << pointi <<
" onto itself."
2707 pointMap_[pointi] = -1;
2708 if (mergePointi >= 0)
2710 reversePointMap_[pointi] = -mergePointi-2;
2714 reversePointMap_[pointi] = -1;
2716 pointZone_.erase(pointi);
2717 retiredPoints_.erase(pointi);
2726 const label masterPointID,
2727 const label masterEdgeID,
2728 const label masterFaceID,
2729 const bool flipFaceFlux,
2741 label facei = faces_.size();
2745 faceOwner_.append(own);
2746 faceNeighbour_.append(nei);
2748 if (masterPointID >= 0)
2750 faceMap_.append(-1);
2751 faceFromPoint_.insert(facei, masterPointID);
2753 else if (masterEdgeID >= 0)
2755 faceMap_.append(-1);
2756 faceFromEdge_.insert(facei, masterEdgeID);
2758 else if (masterFaceID >= 0)
2760 faceMap_.append(masterFaceID);
2769 faceMap_.append(-1);
2771 reverseFaceMap_.append(facei);
2773 flipFaceFlux_.set(facei, flipFaceFlux);
2777 faceZone_.insert(facei,
zoneID);
2779 faceZoneFlip_.set(facei, zoneFlip);
2791 const bool flipFaceFlux,
2804 faceOwner_[facei] = own;
2805 faceNeighbour_[facei] = nei;
2808 flipFaceFlux_.set(facei, flipFaceFlux);
2809 faceZoneFlip_.set(facei, zoneFlip);
2813 faceZone_.set(facei,
zoneID);
2817 faceZone_.erase(facei);
2825 const label mergeFacei
2828 if (facei < 0 || facei >= faces_.size())
2831 <<
"illegal face label " << facei <<
endl
2832 <<
"Valid face labels are 0 .. " << faces_.size()-1
2839 && (faceRemoved(facei) || faceMap_[facei] == -1)
2844 <<
" already marked for removal"
2848 faces_[facei].setSize(0);
2849 region_[facei] = -1;
2850 faceOwner_[facei] = -1;
2851 faceNeighbour_[facei] = -1;
2852 faceMap_[facei] = -1;
2853 if (mergeFacei >= 0)
2855 reverseFaceMap_[facei] = -mergeFacei-2;
2859 reverseFaceMap_[facei] = -1;
2861 faceFromEdge_.erase(facei);
2862 faceFromPoint_.erase(facei);
2863 flipFaceFlux_.unset(facei);
2864 faceZoneFlip_.unset(facei);
2865 faceZone_.erase(facei);
2871 const label masterPointID,
2872 const label masterEdgeID,
2873 const label masterFaceID,
2874 const label masterCellID,
2878 label celli = cellMap_.size();
2880 if (masterPointID >= 0)
2882 cellMap_.append(-1);
2883 cellFromPoint_.insert(celli, masterPointID);
2885 else if (masterEdgeID >= 0)
2887 cellMap_.append(-1);
2888 cellFromEdge_.insert(celli, masterEdgeID);
2890 else if (masterFaceID >= 0)
2892 cellMap_.append(-1);
2893 cellFromFace_.insert(celli, masterFaceID);
2897 cellMap_.append(masterCellID);
2899 reverseCellMap_.append(celli);
2900 cellZone_.append(
zoneID);
2912 cellZone_[celli] =
zoneID;
2919 const label mergeCelli
2922 if (celli < 0 || celli >= cellMap_.size())
2925 <<
"illegal cell label " << celli <<
endl
2926 <<
"Valid cell labels are 0 .. " << cellMap_.size()-1
2930 if (strict_ && cellMap_[celli] == -2)
2934 <<
" already marked for removal"
2938 cellMap_[celli] = -2;
2939 if (mergeCelli >= 0)
2941 reverseCellMap_[celli] = -mergeCelli-2;
2945 reverseCellMap_[celli] = -1;
2947 cellFromPoint_.erase(celli);
2948 cellFromEdge_.erase(celli);
2949 cellFromFace_.erase(celli);
2950 cellZone_[celli] = -1;
2958 const bool syncParallel,
2959 const bool orderCells,
2960 const bool orderPoints
2965 Pout<<
"polyTopoChange::changeMesh"
2966 <<
"(polyMesh&, const bool, const bool, const bool, const bool)"
2972 Pout<<
"Old mesh:" <<
nl;
2979 label nInternalPoints;
3018 oldPatchMeshPointMaps,
3019 oldPatchNMeshPoints,
3021 oldFaceZoneMeshPointMaps
3039 pointField renumberedMeshPoints(newPoints.size());
3043 const label oldPointi = pointMap_[
newPointi];
3088 retiredPoints_.clearStorage();
3089 region_.clearStorage();
3096 label nAdd, nInflate, nMerge, nRemove;
3097 countMap(pointMap_, reversePointMap_, nAdd, nInflate, nMerge, nRemove);
3099 <<
" added(from point):" << nAdd
3100 <<
" added(from nothing):" << nInflate
3101 <<
" merged(into other point):" << nMerge
3102 <<
" removed:" << nRemove
3105 countMap(faceMap_, reverseFaceMap_, nAdd, nInflate, nMerge, nRemove);
3107 <<
" added(from face):" << nAdd
3108 <<
" added(inflated):" << nInflate
3109 <<
" merged(into other face):" << nMerge
3110 <<
" removed:" << nRemove
3113 countMap(cellMap_, reverseCellMap_, nAdd, nInflate, nMerge, nRemove);
3115 <<
" added(from cell):" << nAdd
3116 <<
" added(inflated):" << nInflate
3117 <<
" merged(into other cell):" << nMerge
3118 <<
" removed:" << nRemove
3125 Pout<<
"New mesh:" <<
nl;
3140 resetZones(
mesh,
mesh, pointZoneMap, faceZoneFaceMap, cellZoneMap);
3144 pointZone_.clearStorage();
3145 faceZone_.clearStorage();
3146 faceZoneFlip_.clearStorage();
3147 cellZone_.clearStorage();
3157 oldPatchMeshPointMaps,
3164 calcFaceZonePointMap(
mesh, oldFaceZoneMeshPointMaps, faceZonePointMap);
3205 oldPatchNMeshPoints,