89 return less(values_[a], values_[
b]);
126 const mapPolyMesh& map,
127 const labelList& oldCellsToRefine
130 const polyMesh&
mesh = map.mesh();
135 label nMasterChanged = 0;
153 bitSet oldRefineCell(map.nOldCells(), oldCellsToRefine);
156 bitSet refinedInternalFace(nInternalFaces);
160 for (label faceI = 0; faceI < nInternalFaces; faceI++)
162 label oldOwn = map.cellMap()[faceOwner[faceI]];
163 label oldNei = map.cellMap()[faceNeighbour[faceI]];
167 oldOwn >= 0 && !oldRefineCell.test(oldOwn)
168 && oldNei >= 0 && !oldRefineCell.test(oldNei)
175 refinedInternalFace.set(faceI);
188 label faceI = pp.
start();
192 label oldOwn = map.cellMap()[faceOwner[faceI]];
194 if (oldOwn >= 0 && !oldRefineCell.test(oldOwn))
200 refinedBoundaryFace[faceI-nInternalFaces] =
true;
220 forAll(refinedInternalFace, faceI)
222 if (refinedInternalFace.test(faceI))
224 const cell& ownFaces =
cells[faceOwner[faceI]];
227 changedFace[ownFaces[ownI]] =
true;
229 const cell& neiFaces =
cells[faceNeighbour[faceI]];
232 changedFace[neiFaces[neiI]] =
true;
237 forAll(refinedBoundaryFace, i)
239 if (refinedBoundaryFace[i])
241 const cell& ownFaces =
cells[faceOwner[i+nInternalFaces]];
244 changedFace[ownFaces[ownI]] =
true;
265 forAll(changedFace, faceI)
267 if (changedFace[faceI] && isMasterFace.test(faceI))
277 Pout<<
"getChangedFaces : Detected "
278 <<
" local:" << changedFaces.size()
279 <<
" global:" <<
returnReduce(nMasterChanged, sumOp<label>())
283 faceSet changedFacesSet(
mesh,
"changedFaces", changedFaces);
284 Pout<<
"getChangedFaces : Writing " << changedFaces.size()
285 <<
" changed faces to faceSet " << changedFacesSet.
name()
287 changedFacesSet.
write();
297bool Foam::meshRefinement::markForRefine
299 const label markValue,
300 const label nAllowRefine,
308 cellValue = markValue;
312 return nRefine <= nAllowRefine;
316void Foam::meshRefinement::markFeatureCellLevel
338 Cloud<trackedParticle> startPointCloud
342 IDLList<trackedParticle>()
350 for (
const point& keepPoint : keepPoints)
352 const label celli = mesh_.cellTree().findInside(keepPoint);
360 const edgeMesh& featureMesh = features_[feati];
361 const label featureLevel = features_.levels()[feati][0];
366 label nRegions = featureMesh.regions(edgeRegion);
369 bitSet regionVisited(nRegions);
375 forAll(pointEdges, pointi)
377 if (pointEdges[pointi].size() != 2)
381 Pout<<
"Adding particle from point:" << pointi
382 <<
" coord:" << featureMesh.points()[pointi]
383 <<
" since number of emanating edges:"
384 << pointEdges[pointi].size()
389 startPointCloud.addParticle
396 featureMesh.points()[pointi],
405 if (pointEdges[pointi].size() > 0)
407 label e0 = pointEdges[pointi][0];
408 label regioni = edgeRegion[e0];
409 regionVisited.set(regioni);
417 forAll(featureMesh.edges(), edgei)
419 if (regionVisited.set(edgeRegion[edgei]))
421 const edge&
e = featureMesh.edges()[edgei];
422 label pointi =
e.start();
425 Pout<<
"Adding particle from point:" << pointi
426 <<
" coord:" << featureMesh.points()[pointi]
427 <<
" on circular region:" << edgeRegion[edgei]
432 startPointCloud.addParticle
439 featureMesh.points()[pointi],
454 maxFeatureLevel =
labelList(mesh_.nCells(), -1);
457 List<bitSet> featureEdgeVisited(features_.size());
461 featureEdgeVisited[featI].setSize(features_[featI].edges().size());
462 featureEdgeVisited[featI] =
false;
466 trackedParticle::trackingData td
477 scalar maxTrackLen = 2.0*mesh_.bounds().mag();
481 Pout<<
"Tracking " << startPointCloud.size()
482 <<
" particles over distance " << maxTrackLen
483 <<
" to find the starting cell" <<
endl;
485 startPointCloud.move(startPointCloud, td, maxTrackLen);
489 maxFeatureLevel = -1;
492 featureEdgeVisited[featI] =
false;
496 Cloud<trackedParticle> cloud
500 IDLList<trackedParticle>()
505 Pout<<
"Constructing cloud for cell marking" <<
endl;
508 for (
const trackedParticle& startTp : startPointCloud)
510 const label featI = startTp.i();
511 const label pointI = startTp.j();
513 const edgeMesh& featureMesh = features_[featI];
514 const labelList& pEdges = featureMesh.pointEdges()[pointI];
519 label edgeI = pEdges[pEdgeI];
521 if (featureEdgeVisited[featI].
set(edgeI))
526 const edge&
e = featureMesh.edges()[edgeI];
527 label otherPointi =
e.otherVertex(pointI);
529 trackedParticle* tp(
new trackedParticle(startTp));
530 tp->start() = tp->position();
531 tp->end() = featureMesh.points()[otherPointi];
532 tp->j() = otherPointi;
537 Pout<<
"Adding particle for point:" << pointI
538 <<
" coord:" << tp->position()
539 <<
" feature:" << featI
540 <<
" to track to:" << tp->end()
544 cloud.addParticle(tp);
549 startPointCloud.clear();
557 Pout<<
"Tracking " << cloud.size()
558 <<
" particles over distance " << maxTrackLen
559 <<
" to mark cells" <<
endl;
561 cloud.move(cloud, td, maxTrackLen);
564 for (trackedParticle& tp : cloud)
566 const label featI = tp.i();
567 const label pointI = tp.j();
569 const edgeMesh& featureMesh = features_[featI];
570 const labelList& pEdges = featureMesh.pointEdges()[pointI];
575 bool keepParticle =
false;
579 label edgeI = pEdges[i];
581 if (featureEdgeVisited[featI].
set(edgeI))
586 const edge&
e = featureMesh.edges()[edgeI];
587 label otherPointi =
e.otherVertex(pointI);
589 tp.start() = tp.position();
590 tp.
end() = featureMesh.points()[otherPointi];
591 tp.j() = otherPointi;
602 cloud.deleteParticle(tp);
609 Pout<<
"Remaining particles " << cloud.size() <<
endl;
637Foam::label Foam::meshRefinement::markFeatureRefinement
640 const label nAllowRefine,
648 markFeatureCellLevel(keepPoints, maxFeatureLevel);
653 const labelList& cellLevel = meshCutter_.cellLevel();
655 label oldNRefine = nRefine;
657 forAll(maxFeatureLevel, cellI)
659 if (maxFeatureLevel[cellI] > cellLevel[cellI])
685 Info<<
"Reached refinement limit." <<
endl;
688 return returnReduce(nRefine-oldNRefine, sumOp<label>());
693Foam::label Foam::meshRefinement::markInternalDistanceToFeatureRefinement
695 const label nAllowRefine,
701 const labelList& cellLevel = meshCutter_.cellLevel();
702 const pointField& cellCentres = mesh_.cellCentres();
705 if (features_.maxDistance() <= 0.0)
710 label oldNRefine = nRefine;
714 labelList testLevels(cellLevel.size()-nRefine);
719 if (refineCell[cellI] == -1)
721 testCc[testI] = cellCentres[cellI];
722 testLevels[testI] = cellLevel[cellI];
729 features_.findHigherLevel(testCc, testLevels, maxLevel);
736 if (refineCell[cellI] == -1)
738 if (maxLevel[testI] > testLevels[testI])
740 bool reachedLimit = !markForRefine
752 Pout<<
"Stopped refining internal cells"
753 <<
" since reaching my cell limit of "
754 << mesh_.nCells()+7*nRefine <<
endl;
769 Info<<
"Reached refinement limit." <<
endl;
772 return returnReduce(nRefine-oldNRefine, sumOp<label>());
777Foam::label Foam::meshRefinement::markInternalRefinement
779 const label nAllowRefine,
785 const labelList& cellLevel = meshCutter_.cellLevel();
786 const pointField& cellCentres = mesh_.cellCentres();
788 label oldNRefine = nRefine;
792 labelList testLevels(cellLevel.size()-nRefine);
797 if (refineCell[cellI] == -1)
799 testCc[testI] = cellCentres[cellI];
800 testLevels[testI] = cellLevel[cellI];
807 shells_.findHigherLevel(testCc, testLevels, maxLevel);
814 if (refineCell[cellI] == -1)
816 if (maxLevel[testI] > testLevels[testI])
818 bool reachedLimit = !markForRefine
830 Pout<<
"Stopped refining internal cells"
831 <<
" since reaching my cell limit of "
832 << mesh_.nCells()+7*nRefine <<
endl;
847 Info<<
"Reached refinement limit." <<
endl;
850 return returnReduce(nRefine-oldNRefine, sumOp<label>());
854Foam::label Foam::meshRefinement::unmarkInternalRefinement
860 const labelList& cellLevel = meshCutter_.cellLevel();
861 const pointField& cellCentres = mesh_.cellCentres();
863 label oldNRefine = nRefine;
872 if (refineCell[cellI] >= 0)
874 testCc[testI] = cellCentres[cellI];
875 testLevels[testI] = cellLevel[cellI];
882 limitShells_.findLevel(testCc, testLevels, levelShell);
889 if (refineCell[cellI] >= 0)
891 if (levelShell[testI] != -1)
893 refineCell[cellI] = -1;
900 return returnReduce(oldNRefine-nRefine, sumOp<label>());
915 const labelList& surfIndex = surfaceIndex();
922 if (surfIndex[faceI] != -1)
924 label own = mesh_.faceOwner()[faceI];
926 if (mesh_.isInternalFace(faceI))
928 label nei = mesh_.faceNeighbour()[faceI];
930 if (refineCell[own] == -1 || refineCell[nei] == -1)
932 testFaces[nTest++] = faceI;
937 const label bFacei = faceI - mesh_.nInternalFaces();
938 if (refineCell[own] == -1 || boundaryRefineCell[bFacei] == -1)
940 testFaces[nTest++] = faceI;
945 testFaces.setSize(nTest);
952Foam::label Foam::meshRefinement::markSurfaceRefinement
954 const label nAllowRefine,
962 const labelList& cellLevel = meshCutter_.cellLevel();
964 label oldNRefine = nRefine;
972 labelList testFaces(getRefineCandidateFaces(refineCell));
997 surfaces_.findHigherIntersection
1014 label faceI = testFaces[i];
1016 label surfI = surfaceHit[i];
1026 label own = mesh_.faceOwner()[faceI];
1028 if (surfaceMinLevel[i] > cellLevel[own])
1046 if (mesh_.isInternalFace(faceI))
1048 label nei = mesh_.faceNeighbour()[faceI];
1049 if (surfaceMinLevel[i] > cellLevel[nei])
1076 Info<<
"Reached refinement limit." <<
endl;
1079 return returnReduce(nRefine-oldNRefine, sumOp<label>());
1084Foam::label Foam::meshRefinement::countMatches
1086 const List<point>& normals1,
1087 const List<point>& normals2,
1095 const vector& n1 = normals1[i];
1099 const vector& n2 = normals2[j];
1114Foam::label Foam::meshRefinement::markSurfaceCurvatureRefinement
1116 const scalar curvature,
1117 const label nAllowRefine,
1125 const labelList& cellLevel = meshCutter_.cellLevel();
1126 const pointField& cellCentres = mesh_.cellCentres();
1128 label oldNRefine = nRefine;
1142 labelList testFaces(getRefineCandidateFaces(refineCell));
1152 label faceI = testFaces[i];
1154 label own = mesh_.faceOwner()[faceI];
1156 if (mesh_.isInternalFace(faceI))
1158 label nei = mesh_.faceNeighbour()[faceI];
1160 start[i] = cellCentres[own];
1161 end[i] = cellCentres[nei];
1162 minLevel[i] =
min(cellLevel[own], cellLevel[nei]);
1166 label bFaceI = faceI - mesh_.nInternalFaces();
1168 start[i] = cellCentres[own];
1169 end[i] = neiCc[bFaceI];
1171 if (!isMasterFace[faceI])
1173 std::swap(start[i], end[i]);
1176 minLevel[i] =
min(cellLevel[own], neiLevel[bFaceI]);
1182 const vectorField smallVec(ROOTSMALL*(end-start));
1191 List<vectorList> cellSurfNormals(mesh_.nCells());
1195 List<vectorList> surfaceNormal;
1199 surfaces_.findAllIntersections
1206 labelList(surfaces_.maxLevel().size(), 0),
1207 surfaces_.maxLevel(),
1218 forAll(surfaceNormal, pointI)
1220 vectorList& pNormals = surfaceNormal[pointI];
1221 labelList& pLevel = surfaceLevel[pointI];
1223 sortedOrder(pNormals, visitOrder, normalLess(pNormals));
1225 pNormals = List<point>(pNormals, visitOrder);
1237 label faceI = testFaces[i];
1238 label own = mesh_.faceOwner()[faceI];
1240 const vectorList& fNormals = surfaceNormal[i];
1241 const labelList& fLevels = surfaceLevel[i];
1245 if (fLevels[hitI] > cellLevel[own])
1247 cellSurfLevels[own].append(fLevels[hitI]);
1248 cellSurfNormals[own].append(fNormals[hitI]);
1251 if (mesh_.isInternalFace(faceI))
1253 label nei = mesh_.faceNeighbour()[faceI];
1254 if (fLevels[hitI] > cellLevel[nei])
1256 cellSurfLevels[nei].append(fLevels[hitI]);
1257 cellSurfNormals[nei].append(fNormals[hitI]);
1271 forAll(cellSurfNormals, cellI)
1273 const vectorList& normals = cellSurfNormals[cellI];
1277 nNormals += normals.
size();
1280 reduce(nSet, sumOp<label>());
1281 reduce(nNormals, sumOp<label>());
1282 Info<<
"markSurfaceCurvatureRefinement :"
1283 <<
" cells:" << mesh_.globalData().nTotalCells()
1284 <<
" of which with normals:" << nSet
1285 <<
" ; total normals stored:" << nNormals
1291 bool reachedLimit =
false;
1300 !reachedLimit && cellI < cellSurfNormals.size();
1304 const vectorList& normals = cellSurfNormals[cellI];
1305 const labelList& levels = cellSurfLevels[cellI];
1308 for (label i = 0; !reachedLimit && i < normals.size(); i++)
1310 for (label j = i+1; !reachedLimit && j < normals.size(); j++)
1312 if ((normals[i] & normals[j]) < curvature)
1314 label maxLevel =
max(levels[i], levels[j]);
1316 if (cellLevel[cellI] < maxLevel)
1331 Pout<<
"Stopped refining since reaching my cell"
1332 <<
" limit of " << mesh_.nCells()+7*nRefine
1335 reachedLimit =
true;
1355 !reachedLimit && faceI < mesh_.nInternalFaces();
1359 label own = mesh_.faceOwner()[faceI];
1360 label nei = mesh_.faceNeighbour()[faceI];
1362 const vectorList& ownNormals = cellSurfNormals[own];
1363 const labelList& ownLevels = cellSurfLevels[own];
1364 const vectorList& neiNormals = cellSurfNormals[nei];
1365 const labelList& neiLevels = cellSurfLevels[nei];
1373 countMatches(ownNormals, neiNormals)
1374 == ownNormals.
size();
1377 countMatches(neiNormals, ownNormals)
1378 == neiNormals.size();
1381 if (!ownIsSubset && !neiIsSubset)
1384 for (label i = 0; !reachedLimit && i < ownNormals.size(); i++)
1386 for (label j = 0; !reachedLimit && j < neiNormals.size(); j++)
1389 if ((ownNormals[i] & neiNormals[j]) < curvature)
1392 if (cellLevel[own] < ownLevels[i])
1407 Pout<<
"Stopped refining since reaching"
1408 <<
" my cell limit of "
1409 << mesh_.nCells()+7*nRefine <<
endl;
1411 reachedLimit =
true;
1415 if (cellLevel[nei] < neiLevels[j])
1430 Pout<<
"Stopped refining since reaching"
1431 <<
" my cell limit of "
1432 << mesh_.nCells()+7*nRefine <<
endl;
1434 reachedLimit =
true;
1446 List<vectorList> neiSurfaceNormals;
1452 label faceI = mesh_.nInternalFaces();
1453 !reachedLimit && faceI < mesh_.nFaces();
1457 label own = mesh_.faceOwner()[faceI];
1458 label bFaceI = faceI - mesh_.nInternalFaces();
1460 const vectorList& ownNormals = cellSurfNormals[own];
1461 const labelList& ownLevels = cellSurfLevels[own];
1462 const vectorList& neiNormals = neiSurfaceNormals[bFaceI];
1469 countMatches(ownNormals, neiNormals)
1470 == ownNormals.
size();
1473 countMatches(neiNormals, ownNormals)
1474 == neiNormals.size();
1477 if (!ownIsSubset && !neiIsSubset)
1480 for (label i = 0; !reachedLimit && i < ownNormals.size(); i++)
1482 for (label j = 0; !reachedLimit && j < neiNormals.size(); j++)
1485 if ((ownNormals[i] & neiNormals[j]) < curvature)
1487 if (cellLevel[own] < ownLevels[i])
1502 Pout<<
"Stopped refining since reaching"
1503 <<
" my cell limit of "
1504 << mesh_.nCells()+7*nRefine
1507 reachedLimit =
true;
1524 Info<<
"Reached refinement limit." <<
endl;
1527 return returnReduce(nRefine-oldNRefine, sumOp<label>());
1533 const scalar planarCos,
1543 vector d = point1-point0;
1544 scalar magD =
mag(d);
1546 if (magD > mergeDistance())
1548 scalar cosAngle = (normal0 & normal1);
1551 if (cosAngle < (-1+planarCos))
1554 avg = 0.5*(normal0-normal1);
1556 else if (cosAngle > (1-planarCos))
1558 avg = 0.5*(normal0+normal1);
1566 if (
mag(avg&d) > mergeDistance())
1580 const scalar planarCos,
1594 vector d = point1-point0;
1595 scalar magD =
mag(d);
1597 if (magD > mergeDistance())
1599 scalar cosAngle = (normal0 & normal1);
1602 if (cosAngle < (-1+planarCos))
1605 avgNormal = 0.5*(normal0-normal1);
1607 else if (cosAngle > (1-planarCos))
1609 avgNormal = 0.5*(normal0+normal1);
1614 avgNormal /=
mag(avgNormal);
1636bool Foam::meshRefinement::checkProximity
1638 const scalar planarCos,
1639 const label nAllowRefine,
1641 const label surfaceLevel,
1643 const vector& surfaceNormal,
1647 label& cellMaxLevel,
1655 const labelList& cellLevel = meshCutter_.cellLevel();
1658 if (surfaceLevel > cellLevel[cellI])
1660 if (cellMaxLevel == -1)
1663 cellMaxLevel = surfaceLevel;
1665 cellMaxNormal = surfaceNormal;
1674 bool closeSurfaces = isNormalGap
1687 if (surfaceLevel > cellMaxLevel)
1689 cellMaxLevel = surfaceLevel;
1690 cellMaxLocation = surfaceLocation;
1691 cellMaxNormal = surfaceNormal;
1705 return markForRefine
1721Foam::label Foam::meshRefinement::markProximityRefinement
1723 const scalar planarCos,
1728 const label nAllowRefine,
1736 const labelList& cellLevel = meshCutter_.cellLevel();
1738 label oldNRefine = nRefine;
1749 const labelList testFaces(getRefineCandidateFaces(refineCell));
1770 labelList cellMaxLevel(mesh_.nCells(), -1);
1776 List<vectorList> surfaceLocation;
1777 List<vectorList> surfaceNormal;
1781 surfaces_.findAllIntersections
1819 label faceI = testFaces[i];
1820 label own = mesh_.faceOwner()[faceI];
1822 const labelList& fLevels = surfaceLevel[i];
1823 const vectorList& fPoints = surfaceLocation[i];
1824 const vectorList& fNormals = surfaceNormal[i];
1839 cellMaxLocation[own],
1847 if (mesh_.isInternalFace(faceI))
1849 label nei = mesh_.faceNeighbour()[faceI];
1864 cellMaxLocation[nei],
1878 labelList neiBndMaxLevel(mesh_.nBoundaryFaces());
1879 pointField neiBndMaxLocation(mesh_.nBoundaryFaces());
1880 vectorField neiBndMaxNormal(mesh_.nBoundaryFaces());
1882 for (label faceI = mesh_.nInternalFaces(); faceI < mesh_.nFaces(); faceI++)
1884 label bFaceI = faceI-mesh_.nInternalFaces();
1885 label own = mesh_.faceOwner()[faceI];
1887 neiBndMaxLevel[bFaceI] = cellMaxLevel[own];
1888 neiBndMaxLocation[bFaceI] = cellMaxLocation[own];
1889 neiBndMaxNormal[bFaceI] = cellMaxNormal[own];
1898 for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
1900 label own = mesh_.faceOwner()[faceI];
1901 label nei = mesh_.faceNeighbour()[faceI];
1903 if (cellMaxLevel[own] != -1 && cellMaxLevel[nei] != -1)
1912 cellMaxLocation[own],
1915 cellMaxLocation[nei],
1921 if (cellLevel[own] < cellMaxLevel[own])
1936 Pout<<
"Stopped refining since reaching my cell"
1937 <<
" limit of " << mesh_.nCells()+7*nRefine
1944 if (cellLevel[nei] < cellMaxLevel[nei])
1959 Pout<<
"Stopped refining since reaching my cell"
1960 <<
" limit of " << mesh_.nCells()+7*nRefine
1970 for (label faceI = mesh_.nInternalFaces(); faceI < mesh_.nFaces(); faceI++)
1972 label own = mesh_.faceOwner()[faceI];
1973 label bFaceI = faceI - mesh_.nInternalFaces();
1975 if (cellLevel[own] < cellMaxLevel[own] && neiBndMaxLevel[bFaceI] != -1)
1984 cellMaxLocation[own],
1987 neiBndMaxLocation[bFaceI],
1988 neiBndMaxNormal[bFaceI]
2005 Pout<<
"Stopped refining since reaching my cell"
2006 <<
" limit of " << mesh_.nCells()+7*nRefine
2021 Info<<
"Reached refinement limit." <<
endl;
2024 return returnReduce(nRefine-oldNRefine, sumOp<label>());
2038 const scalar curvature,
2039 const scalar planarAngle,
2041 const bool featureRefinement,
2042 const bool featureDistanceRefinement,
2043 const bool internalRefinement,
2044 const bool surfaceRefinement,
2045 const bool curvatureRefinement,
2046 const bool smallFeatureRefinement,
2047 const bool gapRefinement,
2048 const bool bigGapRefinement,
2049 const bool spreadGapSize,
2050 const label maxGlobalCells,
2051 const label maxLocalCells
2054 label totNCells = mesh_.globalData().nTotalCells();
2058 if (totNCells >= maxGlobalCells)
2060 Info<<
"No cells marked for refinement since reached limit "
2061 << maxGlobalCells <<
'.' <<
endl;
2095 labelList neiLevel(mesh_.nBoundaryFaces());
2097 calcNeighbourData(neiLevel, neiCc);
2106 if (featureRefinement)
2108 label nFeatures = markFeatureRefinement
2117 Info<<
"Marked for refinement due to explicit features "
2118 <<
": " << nFeatures <<
" cells." <<
endl;
2124 if (featureDistanceRefinement)
2126 label nShell = markInternalDistanceToFeatureRefinement
2133 Info<<
"Marked for refinement due to distance to explicit features "
2134 ": " << nShell <<
" cells." <<
endl;
2140 if (internalRefinement)
2142 label nShell = markInternalRefinement
2149 Info<<
"Marked for refinement due to refinement shells "
2150 <<
": " << nShell <<
" cells." <<
endl;
2156 if (surfaceRefinement)
2158 label nSurf = markSurfaceRefinement
2167 Info<<
"Marked for refinement due to surface intersection "
2168 <<
": " << nSurf <<
" cells." <<
endl;
2176 &&
max(shells_.maxGapLevel()) > 0
2179 label nGapSurf = markSurfaceGapRefinement
2189 Info<<
"Marked for refinement due to surface intersection"
2191 <<
": " << nGapSurf <<
" cells." <<
endl;
2201 && (curvature >= -1 && curvature <= 1)
2202 && (surfaces_.minLevel() != surfaces_.maxLevel())
2205 label nCurv = markSurfaceCurvatureRefinement
2215 Info<<
"Marked for refinement due to curvature/regions "
2216 <<
": " << nCurv <<
" cells." <<
endl;
2225 smallFeatureRefinement
2226 && (planarCos >= -1 && planarCos <= 1)
2227 &&
max(shells_.maxGapLevel()) > 0
2230 label nGap = markSmallFeatureRefinement
2240 Info<<
"Marked for refinement due to close opposite surfaces "
2241 <<
": " << nGap <<
" cells." <<
endl;
2248 const labelList& surfaceGapLevel = surfaces_.gapLevel();
2253 && (planarCos >= -1 && planarCos <= 1)
2254 && (
max(surfaceGapLevel) > -1)
2257 Info<<
"Specified gap level : " <<
max(surfaceGapLevel)
2258 <<
", planar angle " << planarAngle <<
endl;
2260 label nGap = markProximityRefinement
2274 Info<<
"Marked for refinement due to close opposite surfaces "
2275 <<
": " << nGap <<
" cells." <<
endl;
2285 && (planarCos >= -1 && planarCos <= 1)
2286 &&
max(shells_.maxGapLevel()) > 0
2293 label nGap = markInternalGapRefinement
2304 Info<<
"Marked for refinement due to opposite surfaces"
2306 <<
": " << nGap <<
" cells." <<
endl;
2313 if (limitShells_.shells().size())
2315 label nUnmarked = unmarkInternalRefinement(
refineCell, nRefine);
2318 Info<<
"Unmarked for refinement due to limit shells"
2319 <<
" : " << nUnmarked <<
" cells." <<
endl;
2328 cellsToRefine.
setSize(nRefine);
2335 cellsToRefine[nRefine++] = cellI;
2340 return cellsToRefine;
2353 meshCutter_.setRefinement(cellsToRefine, meshMod);
2357 mesh_.moving(
false);
2364 mesh_.updateMesh(map);
2381 updateMesh(map, getChangedFaces(map, cellsToRefine));
2396 const scalar maxLoadUnbalance
2400 refine(cellsToRefine);
2404 Pout<<
"Writing refined but unbalanced " << msg
2412 Pout<<
"Dumped debug data in = "
2413 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2419 Info<<
"Refined mesh in = "
2420 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2421 printMeshInfo(debug,
"After refinement " + msg);
2431 scalar nIdealCells =
2432 mesh_.globalData().nTotalCells()
2437 mag(1.0-mesh_.nCells()/nIdealCells),
2441 if (unbalance <= maxLoadUnbalance)
2443 Info<<
"Skipping balancing since max unbalance " << unbalance
2444 <<
" is less than allowable " << maxLoadUnbalance
2460 Info<<
"Balanced mesh in = "
2461 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2463 printMeshInfo(debug,
"After balancing " + msg);
2468 Pout<<
"Writing balanced " << msg
2476 Pout<<
"Dumped debug data in = "
2477 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2497 const scalar maxLoadUnbalance
2500 labelList cellsToRefine(initCellsToRefine);
2535 scalar nNewCells = scalar(mesh_.nCells() + 7*cellsToRefine.
size());
2536 scalar nIdealNewCells =
2541 mag(1.0-nNewCells/nIdealNewCells),
2545 if (unbalance <= maxLoadUnbalance)
2547 Info<<
"Skipping balancing since max unbalance " << unbalance
2548 <<
" is less than allowable " << maxLoadUnbalance
2556 cellWeights[cellsToRefine[i]] += 7;
2569 distMap().distributeCellIndices(cellsToRefine);
2571 Info<<
"Balanced mesh in = "
2572 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2587 printMeshInfo(debug,
"After balancing " + msg);
2591 Pout<<
"Writing balanced " << msg
2599 Pout<<
"Dumped debug data in = "
2600 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2611 refine(cellsToRefine);
2615 Pout<<
"Writing refined " << msg
2623 Pout<<
"Dumped debug data in = "
2624 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2630 Info<<
"Refined mesh in = "
2631 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2645 printMeshInfo(debug,
"After refinement " + msg);
2653 const label maxGlobalCells,
2654 const label maxLocalCells,
2659 const labelList& cellLevel = meshCutter_.cellLevel();
2660 const pointField& cellCentres = mesh_.cellCentres();
2662 label totNCells = mesh_.globalData().nTotalCells();
2666 if (totNCells >= maxGlobalCells)
2668 Info<<
"No cells marked for refinement since reached limit "
2669 << maxGlobalCells <<
'.' <<
endl;
2683 shells_.findDirectionalLevel
2693 forAll(insideShell, celli)
2695 if (insideShell[celli] >= 0)
2697 bool reachedLimit = !markForRefine
2709 Pout<<
"Stopped refining cells"
2710 <<
" since reaching my cell limit of "
2711 << mesh_.nCells()+nAllowRefine <<
endl;
2722 label nUnmarked = unmarkInternalRefinement(
refineCell, nRefine);
2725 Info<<
"Unmarked for refinement due to limit shells"
2726 <<
" : " << nUnmarked <<
" cells." <<
endl;
2735 cellsToRefine.
setSize(nRefine);
2742 cellsToRefine[nRefine++] = cellI;
2747 return cellsToRefine;
2764 refCells[i] =
refineCell(cellsToRefine[i], refDir);
2771 cellCuts cuts(mesh_, cellWalker, refCells);
2783 mesh_.moving(
false);
2788 mesh_.updateMesh(*morphMap);
2791 if (morphMap().hasMotionPoints())
2793 mesh_.movePoints(morphMap().preMotionPoints());
2808 updateMesh(*morphMap, getChangedFaces(*morphMap, cellsToRefine));
void setSize(const label n)
Alias for resize()
void clear()
Clear the list, i.e. set size to zero.
virtual const fileName & name() const
Get the name of the stream.
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
label nProcs() const noexcept
Number of ranks associated with PstreamBuffers.
iterator end() noexcept
Return an iterator to end traversing the UList.
void size(const label n)
Older name for setAddressableSize.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Description of cuts across cells.
Abstract base class for domain decomposition.
Sends/receives parts of mesh+fvfields to neighbouring processors. Used in load balancing.
label nTotalFaces() const noexcept
Return total number of faces in decomposed mesh. Not.
Implementation of cellLooper.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
const pointField & preMotionPoints() const
Pre-motion point positions.
bool hasMotionPoints() const
Has valid preMotionPoints?
void updateMesh(const mapPolyMesh &)
Force recalculation of locally stored data on topological change.
void setRefinement(const cellCuts &cuts, polyTopoChange &meshMod)
Do actual cutting with cut description. Inserts mesh changes.
autoPtr< mapDistributePolyMesh > balanceAndRefine(const string &msg, decompositionMethod &decomposer, fvMeshDistribute &distributor, const labelList &cellsToRefine, const scalar maxLoadUnbalance)
Balance before refining some cells.
writeType
Enumeration for what to write. Used as a bit-pattern.
labelList refineCandidates(const pointField &keepPoints, const scalar curvature, const scalar planarAngle, const bool featureRefinement, const bool featureDistanceRefinement, const bool internalRefinement, const bool surfaceRefinement, const bool curvatureRefinement, const bool smallFeatureRefinement, const bool gapRefinement, const bool bigGapRefinement, const bool spreadGapSize, const label maxGlobalCells, const label maxLocalCells) const
Calculate list of cells to refine.
debugType
Enumeration for what to debug. Used as a bit-pattern.
bool isNormalGap(const scalar planarCos, const label level0, const vector &point0, const vector &normal0, const label level1, const vector &point1, const vector &normal1) const
Is local topology a small gap normal to the test vector.
const fvMesh & mesh() const
Reference to mesh.
bool isGap(const scalar, const vector &, const vector &, const vector &, const vector &) const
Is local topology a small gap?
autoPtr< mapDistributePolyMesh > refineAndBalance(const string &msg, decompositionMethod &decomposer, fvMeshDistribute &distributor, const labelList &cellsToRefine, const scalar maxLoadUnbalance)
Refine some cells and rebalance.
autoPtr< mapPolyMesh > directionalRefine(const string &msg, const direction cmpt, const labelList &cellsToRefine)
Directionally refine in direction cmpt.
autoPtr< mapPolyMesh > refine(const labelList &cellsToRefine)
Refine some cells.
labelList directionalRefineCandidates(const label maxGlobalCells, const label maxLocalCells, const labelList ¤tLevel, const direction dir) const
Calculate list of cells to directionally refine.
normalLess(const vectorList &values)
bool operator()(const label a, const label b) const
static constexpr direction nComponents
Number of components in bool is 1.
labelList cmptType
Component type.
vectorList cmptType
Component type.
A traits class, which is primarily used for primitives.
label start() const noexcept
The start label of boundary faces in the polyMesh face list.
virtual const labelList & faceOwner() const
Return face owner.
const globalMeshData & globalData() const
Return parallel info.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const labelList & faceNeighbour() const
Return face neighbour.
Direct mesh changes based on v1.3 polyTopoChange syntax.
autoPtr< mapPolyMesh > changeMesh(polyMesh &mesh, const labelUList &patchMap, const bool inflate, const bool syncParallel=true, const bool orderCells=false, const bool orderPoints=false)
Inplace changes mesh without change of patches.
label nBoundaryFaces() const noexcept
Number of boundary faces (== nFaces - nInternalFaces)
label nInternalFaces() const noexcept
Number of internal faces.
label nFaces() const noexcept
Number of mesh faces.
const cellList & cells() const
Container with cells to refine. Refinement given as single direction.
Contains information about location on a triSurface.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
void set(List< bool > &bools, const labelUList &locations)
Set the listed locations (assign 'true').
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
List< label > labelList
A List of labels.
List< cell > cellList
A List of cells.
messageStream Info
Information stream (stdout output on master, null elsewhere)
vectorField pointField
pointField is a vectorField.
vector point
Point is a vector.
List< vector > vectorList
A List of vectors.
Ostream & endl(Ostream &os)
Add newline and flush stream.
constexpr scalar degToRad() noexcept
Multiplication factor for degrees to radians conversion.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static bool less(const vector &x, const vector &y)
To compare normals.
List< labelList > labelListList
A List of labelList.
void reduce(const List< UPstream::commsStruct > &comms, T &value, const BinaryOp &bop, const int tag, const label comm)
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Field< vector > vectorField
Specialisation of Field<T> for vector.
List< bool > boolList
A List of bools.
static constexpr const zero Zero
Global zero (0)
label checkData(const fvMesh &mesh, const instantList &timeDirs, wordList &objectNames)
Check if fields are good to use (available at all times)
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
labelList findIndices(const ListType &input, typename ListType::const_reference val, label start=0)
Linear search to find all occurrences of given element.
UIndirectList< label > labelUIndList
UIndirectList of labels.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
dimensionedScalar cos(const dimensionedScalar &ds)
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
#define forAll(list, i)
Loop across all elements in list.
A non-counting (dummy) refCount.