80 return less(values_[a], values_[
b]);
117 const polyTopoChangeMap& map,
121 const polyMesh&
mesh = map.mesh();
126 label nMasterChanged = 0;
144 PackedBoolList oldRefineCell(map.nOldCells());
146 forAll(oldCellsToRefine, i)
148 oldRefineCell.set(oldCellsToRefine[i], 1u);
152 PackedBoolList refinedInternalFace(nInternalFaces);
156 for (
label facei = 0; facei < nInternalFaces; facei++)
158 const label oldOwn = map.cellMap()[faceOwner[facei]];
159 const label oldNei = map.cellMap()[faceNeighbour[facei]];
164 && oldRefineCell.get(oldOwn) == 0u
166 && oldRefineCell.get(oldNei) == 0u
173 refinedInternalFace.set(facei, 1u);
186 label facei = pp.start();
190 label oldOwn = map.cellMap()[faceOwner[facei]];
192 if (oldOwn >= 0 && oldRefineCell.get(oldOwn) == 0u)
198 refinedBoundaryFace[facei - nInternalFaces] =
true;
218 forAll(refinedInternalFace, facei)
220 if (refinedInternalFace.get(facei) == 1u)
222 const cell& ownFaces =
cells[faceOwner[facei]];
225 changedFace[ownFaces[owni]] =
true;
227 const cell& neiFaces =
cells[faceNeighbour[facei]];
230 changedFace[neiFaces[neii]] =
true;
235 forAll(refinedBoundaryFace, i)
237 if (refinedBoundaryFace[i])
239 const cell& ownFaces =
cells[faceOwner[i+nInternalFaces]];
242 changedFace[ownFaces[owni]] =
true;
263 forAll(changedFace, facei)
265 if (changedFace[facei] && isMasterFace[facei])
275 Pout<<
"getChangedFaces : Detected "
276 <<
" local:" << changedFaces.size()
277 <<
" global:" <<
returnReduce(nMasterChanged, sumOp<label>())
281 faceSet changedFacesSet(
mesh,
"changedFaces", changedFaces);
282 Pout<<
"getChangedFaces : Writing " << changedFaces.size()
283 <<
" changed faces to faceSet " << changedFacesSet.
name()
285 changedFacesSet.
write();
295 bool Foam::meshRefinement::markForRefine
297 const label markValue,
298 const label nAllowRefine,
306 cellValue = markValue;
310 return nRefine <= nAllowRefine;
314 void Foam::meshRefinement::markFeatureCellLevel
338 Cloud<trackedParticle> startPointCloud
342 IDLList<trackedParticle>()
354 const label celli = mesh_.cellTree().findInside(insidePoint);
362 const edgeMesh& featureMesh = features_[feati];
363 const label featureLevel = features_.levels()[feati][0];
368 const label nRegions = featureMesh.regions(edgeRegion);
371 PackedBoolList regionVisited(nRegions);
377 forAll(pointEdges, pointi)
379 if (pointEdges[pointi].size() != 2)
383 Pout<<
"Adding particle from point:" << pointi
384 <<
" coord:" << featureMesh.points()[pointi]
385 <<
" since number of emanating edges:"
386 << pointEdges[pointi].size()
391 startPointCloud.addParticle
398 featureMesh.points()[pointi],
407 if (pointEdges[pointi].size() > 0)
409 label e0 = pointEdges[pointi][0];
410 label regioni = edgeRegion[e0];
411 regionVisited[regioni] = 1u;
419 forAll(featureMesh.edges(), edgei)
421 if (regionVisited.set(edgeRegion[edgei], 1u))
423 const edge&
e = featureMesh.edges()[edgei];
424 const label pointi =
e.start();
427 Pout<<
"Adding particle from point:" << pointi
428 <<
" coord:" << featureMesh.points()[pointi]
429 <<
" on circular region:" << edgeRegion[edgei]
434 startPointCloud.addParticle
441 featureMesh.points()[pointi],
456 maxFeatureLevel =
labelList(mesh_.nCells(), -1);
459 List<PackedBoolList> featureEdgeVisited(features_.size());
463 featureEdgeVisited[feati].setSize(features_[feati].edges().size());
464 featureEdgeVisited[feati] = 0u;
468 trackedParticle::trackingData td
471 2.0*mesh_.bounds().mag(),
482 Pout<<
"Tracking " << startPointCloud.size()
483 <<
" particles over distance " << td.maxTrackLen_
484 <<
" to find the starting cell" <<
endl;
486 startPointCloud.move(startPointCloud, td);
490 maxFeatureLevel = -1;
493 featureEdgeVisited[feati] = 0u;
497 Cloud<trackedParticle> cloud
501 IDLList<trackedParticle>()
506 Pout<<
"Constructing cloud for cell marking" <<
endl;
509 forAllIter(Cloud<trackedParticle>, startPointCloud, iter)
511 const trackedParticle& startTp = iter();
513 const label feati = startTp.i();
514 const label pointi = startTp.j();
516 const edgeMesh& featureMesh = features_[feati];
517 const labelList& pEdges = featureMesh.pointEdges()[pointi];
522 const label edgei = pEdges[pEdgei];
524 if (featureEdgeVisited[feati].set(edgei, 1u))
529 const edge&
e = featureMesh.edges()[edgei];
530 label otherPointi =
e.otherVertex(pointi);
532 trackedParticle* tp(
new trackedParticle(startTp));
533 tp->start() = tp->position(mesh_);
534 tp->end() = featureMesh.points()[otherPointi];
535 tp->j() = otherPointi;
540 Pout<<
"Adding particle for point:" << pointi
541 <<
" coord:" << tp->position(mesh_)
542 <<
" feature:" << feati
543 <<
" to track to:" << tp->end()
547 cloud.addParticle(tp);
552 startPointCloud.clear();
560 Pout<<
"Tracking " << cloud.size()
561 <<
" particles over distance " << td.maxTrackLen_
562 <<
" to mark cells" <<
endl;
564 cloud.move(cloud, td);
567 forAllIter(Cloud<trackedParticle>, cloud, iter)
569 trackedParticle& tp = iter();
571 const label feati = tp.i();
572 const label pointi = tp.j();
574 const edgeMesh& featureMesh = features_[feati];
575 const labelList& pEdges = featureMesh.pointEdges()[pointi];
580 bool keepParticle =
false;
584 const label edgei = pEdges[i];
586 if (featureEdgeVisited[feati].set(edgei, 1u))
591 const edge&
e = featureMesh.edges()[edgei];
592 const label otherPointi =
e.otherVertex(pointi);
594 tp.start() = tp.position(mesh_);
595 tp.end() = featureMesh.points()[otherPointi];
596 tp.j() = otherPointi;
607 cloud.deleteParticle(tp);
614 Pout<<
"Remaining particles " << cloud.size() <<
endl;
642 Foam::label Foam::meshRefinement::markFeatureRefinement
645 const label nAllowRefine,
658 const labelList& cellLevel = meshCutter_.cellLevel();
660 const label oldNRefine = nRefine;
662 forAll(maxFeatureLevel, celli)
664 if (maxFeatureLevel[celli] > cellLevel[celli])
690 Info<<
"Reached refinement limit." <<
endl;
693 return returnReduce(nRefine-oldNRefine, sumOp<label>());
698 Foam::label Foam::meshRefinement::markInternalDistanceToFeatureRefinement
700 const label nAllowRefine,
706 const labelList& cellLevel = meshCutter_.cellLevel();
707 const pointField& cellCentres = mesh_.cellCentres();
710 if (features_.maxDistance() <= 0.0)
715 label oldNRefine = nRefine;
719 labelList testLevels(cellLevel.size()-nRefine);
724 if (refineCell[celli] == -1)
726 testCc[testi] = cellCentres[celli];
727 testLevels[testi] = cellLevel[celli];
734 features_.findHigherLevel(testCc, testLevels, maxLevel);
741 if (refineCell[celli] == -1)
743 if (maxLevel[testi] > testLevels[testi])
745 const bool reachedLimit = !markForRefine
757 Pout<<
"Stopped refining internal cells"
758 <<
" since reaching my cell limit of "
759 << mesh_.nCells()+7*nRefine <<
endl;
774 Info<<
"Reached refinement limit." <<
endl;
777 return returnReduce(nRefine-oldNRefine, sumOp<label>());
781 Foam::label Foam::meshRefinement::markInternalRefinement
783 const label nAllowRefine,
789 const labelList& cellLevel = meshCutter_.cellLevel();
790 const pointField& cellCentres = mesh_.cellCentres();
792 const label oldNRefine = nRefine;
796 labelList testLevels(cellLevel.size()-nRefine);
801 if (refineCell[celli] == -1)
803 testCc[testi] = cellCentres[celli];
804 testLevels[testi] = cellLevel[celli];
811 shells_.findHigherLevel
815 meshCutter_.level0EdgeLength(),
824 if (refineCell[celli] == -1)
826 if (maxLevel[testi] > testLevels[testi])
828 bool reachedLimit = !markForRefine
840 Pout<<
"Stopped refining internal cells"
841 <<
" since reaching my cell limit of "
842 << mesh_.nCells()+7*nRefine <<
endl;
857 Info<<
"Reached refinement limit." <<
endl;
860 return returnReduce(nRefine-oldNRefine, sumOp<label>());
873 forAll(surfaceIndex_, facei)
875 if (surfaceIndex_[facei] != -1)
877 label own = mesh_.faceOwner()[facei];
879 if (mesh_.isInternalFace(facei))
881 const label nei = mesh_.faceNeighbour()[facei];
883 if (refineCell[own] == -1 || refineCell[nei] == -1)
885 testFaces[nTest++] = facei;
890 if (refineCell[own] == -1)
892 testFaces[nTest++] = facei;
897 testFaces.setSize(nTest);
903 Foam::label Foam::meshRefinement::markSurfaceRefinement
905 const label nAllowRefine,
913 const labelList& cellLevel = meshCutter_.cellLevel();
914 const pointField& cellCentres = mesh_.cellCentres();
916 const label oldNRefine = nRefine;
924 labelList testFaces(getRefineCandidateFaces(refineCell));
935 const label facei = testFaces[i];
937 const label own = mesh_.faceOwner()[facei];
939 if (mesh_.isInternalFace(facei))
941 label nei = mesh_.faceNeighbour()[facei];
943 start[i] = cellCentres[own];
944 end[i] = cellCentres[nei];
945 minLevel[i] =
min(cellLevel[own], cellLevel[nei]);
949 const label bFacei = facei - mesh_.nInternalFaces();
951 start[i] = cellCentres[own];
952 end[i] = neiCc[bFacei];
953 minLevel[i] =
min(cellLevel[own], neiLevel[bFacei]);
970 surfaces_.findHigherIntersection
986 const label facei = testFaces[i];
988 const label surfi = surfaceHit[i];
998 const label own = mesh_.faceOwner()[facei];
1000 if (surfaceMinLevel[i] > cellLevel[own])
1018 if (mesh_.isInternalFace(facei))
1020 const label nei = mesh_.faceNeighbour()[facei];
1022 if (surfaceMinLevel[i] > cellLevel[nei])
1049 Info<<
"Reached refinement limit." <<
endl;
1052 return returnReduce(nRefine-oldNRefine, sumOp<label>());
1059 const List<point>& normals1,
1060 const List<point>& normals2,
1068 const vector& n1 = normals1[i];
1072 const vector& n2 = normals2[j];
1087 Foam::label Foam::meshRefinement::markSurfaceCurvatureRefinement
1089 const scalar curvature,
1090 const label nAllowRefine,
1098 const labelList& cellLevel = meshCutter_.cellLevel();
1099 const pointField& cellCentres = mesh_.cellCentres();
1101 const label oldNRefine = nRefine;
1115 labelList testFaces(getRefineCandidateFaces(refineCell));
1124 const label facei = testFaces[i];
1125 const label own = mesh_.faceOwner()[facei];
1127 if (mesh_.isInternalFace(facei))
1129 const label nei = mesh_.faceNeighbour()[facei];
1131 start[i] = cellCentres[own];
1132 end[i] = cellCentres[nei];
1133 minLevel[i] =
min(cellLevel[own], cellLevel[nei]);
1137 const label bFacei = facei - mesh_.nInternalFaces();
1139 start[i] = cellCentres[own];
1140 end[i] = neiCc[bFacei];
1142 if (!isMasterFace[facei])
1144 Swap(start[i], end[i]);
1147 minLevel[i] =
min(cellLevel[own], neiLevel[bFacei]);
1153 const vectorField smallVec(rootSmall*(end-start));
1162 List<vectorList> cellSurfNormals(mesh_.nCells());
1166 List<vectorList> surfaceNormal;
1171 surfaces_.findAllHigherIntersections
1178 surfaces_.maxLevel(),
1189 forAll(surfaceNormal, pointi)
1191 vectorList& pNormals = surfaceNormal[pointi];
1192 labelList& pLevel = surfaceLevel[pointi];
1194 sortedOrder(pNormals, visitOrder, normalLess(pNormals));
1196 pNormals = List<point>(pNormals, visitOrder);
1197 pLevel = UIndirectList<label>(pLevel, visitOrder);
1208 const label facei = testFaces[i];
1209 const label own = mesh_.faceOwner()[facei];
1211 const vectorList& fNormals = surfaceNormal[i];
1212 const labelList& fLevels = surfaceLevel[i];
1216 if (fLevels[hiti] > cellLevel[own])
1218 cellSurfLevels[own].append(fLevels[hiti]);
1219 cellSurfNormals[own].append(fNormals[hiti]);
1222 if (mesh_.isInternalFace(facei))
1224 label nei = mesh_.faceNeighbour()[facei];
1225 if (fLevels[hiti] > cellLevel[nei])
1227 cellSurfLevels[nei].append(fLevels[hiti]);
1228 cellSurfNormals[nei].append(fNormals[hiti]);
1242 forAll(cellSurfNormals, celli)
1244 const vectorList& normals = cellSurfNormals[celli];
1248 nNormals += normals.size();
1251 reduce(nSet, sumOp<label>());
1252 reduce(nNormals, sumOp<label>());
1253 Info<<
"markSurfaceCurvatureRefinement :"
1254 <<
" cells:" << mesh_.globalData().nTotalCells()
1255 <<
" of which with normals:" << nSet
1256 <<
" ; total normals stored:" << nNormals
1262 bool reachedLimit =
false;
1271 !reachedLimit && celli < cellSurfNormals.size();
1275 const vectorList& normals = cellSurfNormals[celli];
1276 const labelList& levels = cellSurfLevels[celli];
1279 for (
label i = 0; !reachedLimit && i < normals.size(); i++)
1281 for (
label j = i+1; !reachedLimit && j < normals.size(); j++)
1283 if ((normals[i] & normals[j]) < curvature)
1285 const label maxLevel =
max(levels[i], levels[j]);
1287 if (cellLevel[celli] < maxLevel)
1302 Pout<<
"Stopped refining since reaching my cell"
1303 <<
" limit of " << mesh_.nCells()+7*nRefine
1306 reachedLimit =
true;
1326 !reachedLimit && facei < mesh_.nInternalFaces();
1330 const label own = mesh_.faceOwner()[facei];
1331 const label nei = mesh_.faceNeighbour()[facei];
1333 const vectorList& ownNormals = cellSurfNormals[own];
1334 const labelList& ownLevels = cellSurfLevels[own];
1335 const vectorList& neiNormals = cellSurfNormals[nei];
1336 const labelList& neiLevels = cellSurfLevels[nei];
1344 countMatches(ownNormals, neiNormals)
1345 == ownNormals.
size();
1348 countMatches(neiNormals, ownNormals)
1349 == neiNormals.size();
1352 if (!ownisSubset && !neiisSubset)
1355 for (
label i = 0; !reachedLimit && i < ownNormals.size(); i++)
1357 for (
label j = 0; !reachedLimit && j < neiNormals.size(); j++)
1360 if ((ownNormals[i] & neiNormals[j]) < curvature)
1363 if (cellLevel[own] < ownLevels[i])
1378 Pout<<
"Stopped refining since reaching"
1379 <<
" my cell limit of "
1380 << mesh_.nCells()+7*nRefine <<
endl;
1382 reachedLimit =
true;
1386 if (cellLevel[nei] < neiLevels[j])
1401 Pout<<
"Stopped refining since reaching"
1402 <<
" my cell limit of "
1403 << mesh_.nCells()+7*nRefine <<
endl;
1405 reachedLimit =
true;
1417 List<vectorList> neiSurfaceNormals;
1423 label facei = mesh_.nInternalFaces();
1424 !reachedLimit && facei < mesh_.nFaces();
1428 label own = mesh_.faceOwner()[facei];
1429 label bFacei = facei - mesh_.nInternalFaces();
1431 const vectorList& ownNormals = cellSurfNormals[own];
1432 const labelList& ownLevels = cellSurfLevels[own];
1433 const vectorList& neiNormals = neiSurfaceNormals[bFacei];
1440 countMatches(ownNormals, neiNormals)
1441 == ownNormals.size();
1444 countMatches(neiNormals, ownNormals)
1445 == neiNormals.size();
1448 if (!ownisSubset && !neiisSubset)
1451 for (
label i = 0; !reachedLimit && i < ownNormals.size(); i++)
1453 for (
label j = 0; !reachedLimit && j < neiNormals.size(); j++)
1456 if ((ownNormals[i] & neiNormals[j]) < curvature)
1458 if (cellLevel[own] < ownLevels[i])
1473 Pout<<
"Stopped refining since reaching"
1474 <<
" my cell limit of "
1475 << mesh_.nCells()+7*nRefine
1478 reachedLimit =
true;
1495 Info<<
"Reached refinement limit." <<
endl;
1498 return returnReduce(nRefine-oldNRefine, sumOp<label>());
1504 const scalar planarCos,
1514 const vector d = point1 - point0;
1515 const scalar magD =
mag(d);
1517 if (magD > mergeDistance())
1519 scalar cosAngle = (normal0 & normal1);
1522 if (cosAngle < (-1+planarCos))
1525 avg = 0.5*(normal0-normal1);
1527 else if (cosAngle > (1-planarCos))
1529 avg = 0.5*(normal0+normal1);
1537 if (
mag(avg&d) > mergeDistance())
1561 const scalar planarCos,
1571 vector d = point1 - point0;
1572 const scalar magD =
mag(d);
1574 if (magD > mergeDistance())
1576 scalar cosAngle = (normal0 & normal1);
1579 if (cosAngle < (-1+planarCos))
1582 avg = 0.5*(normal0-normal1);
1584 else if (cosAngle > (1-planarCos))
1586 avg = 0.5*(normal0+normal1);
1616 bool Foam::meshRefinement::checkProximity
1618 const scalar planarCos,
1619 const label nAllowRefine,
1621 const label surfaceLevel,
1623 const vector& surfaceNormal,
1627 label& cellMaxLevel,
1635 const labelList& cellLevel = meshCutter_.cellLevel();
1638 if (surfaceLevel > cellLevel[celli])
1640 if (cellMaxLevel == -1)
1643 cellMaxLevel = surfaceLevel;
1645 cellMaxNormal = surfaceNormal;
1654 bool closeSurfaces = isNormalGap
1665 if (surfaceLevel > cellMaxLevel)
1667 cellMaxLevel = surfaceLevel;
1668 cellMaxLocation = surfaceLocation;
1669 cellMaxNormal = surfaceNormal;
1683 return markForRefine
1699 Foam::label Foam::meshRefinement::markProximityRefinement
1701 const scalar planarCos,
1702 const label nAllowRefine,
1710 const labelList& cellLevel = meshCutter_.cellLevel();
1711 const pointField& cellCentres = mesh_.cellCentres();
1713 const label oldNRefine = nRefine;
1724 labelList testFaces(getRefineCandidateFaces(refineCell));
1733 const label facei = testFaces[i];
1734 const label own = mesh_.faceOwner()[facei];
1736 if (mesh_.isInternalFace(facei))
1738 const label nei = mesh_.faceNeighbour()[facei];
1740 start[i] = cellCentres[own];
1741 end[i] = cellCentres[nei];
1742 minLevel[i] =
min(cellLevel[own], cellLevel[nei]);
1746 const label bFacei = facei - mesh_.nInternalFaces();
1748 start[i] = cellCentres[own];
1749 end[i] = neiCc[bFacei];
1750 minLevel[i] =
min(cellLevel[own], neiLevel[bFacei]);
1756 const vectorField smallVec(rootSmall*(end-start));
1765 labelList cellMaxLevel(mesh_.nCells(), -1);
1771 List<vectorList> surfaceLocation;
1772 List<vectorList> surfaceNormal;
1776 surfaces_.findAllHigherIntersections
1783 surfaces_.gapLevel(),
1813 label facei = testFaces[i];
1814 label own = mesh_.faceOwner()[facei];
1816 const labelList& fLevels = surfaceLevel[i];
1817 const vectorList& fPoints = surfaceLocation[i];
1818 const vectorList& fNormals = surfaceNormal[i];
1833 cellMaxLocation[own],
1841 if (mesh_.isInternalFace(facei))
1843 label nei = mesh_.faceNeighbour()[facei];
1858 cellMaxLocation[nei],
1872 labelList neiBndMaxLevel(mesh_.nFaces()-mesh_.nInternalFaces());
1873 pointField neiBndMaxLocation(mesh_.nFaces()-mesh_.nInternalFaces());
1874 vectorField neiBndMaxNormal(mesh_.nFaces()-mesh_.nInternalFaces());
1876 for (
label facei = mesh_.nInternalFaces(); facei < mesh_.nFaces(); facei++)
1878 const label bFacei = facei-mesh_.nInternalFaces();
1879 const label own = mesh_.faceOwner()[facei];
1881 neiBndMaxLevel[bFacei] = cellMaxLevel[own];
1882 neiBndMaxLocation[bFacei] = cellMaxLocation[own];
1883 neiBndMaxNormal[bFacei] = cellMaxNormal[own];
1892 for (
label facei = 0; facei < mesh_.nInternalFaces(); facei++)
1894 const label own = mesh_.faceOwner()[facei];
1895 const label nei = mesh_.faceNeighbour()[facei];
1897 if (cellMaxLevel[own] != -1 && cellMaxLevel[nei] != -1)
1905 cellMaxLocation[own],
1907 cellMaxLocation[nei],
1913 if (cellLevel[own] < cellMaxLevel[own])
1928 Pout<<
"Stopped refining since reaching my cell"
1929 <<
" limit of " << mesh_.nCells()+7*nRefine
1936 if (cellLevel[nei] < cellMaxLevel[nei])
1951 Pout<<
"Stopped refining since reaching my cell"
1952 <<
" limit of " << mesh_.nCells()+7*nRefine
1962 for (
label facei = mesh_.nInternalFaces(); facei < mesh_.nFaces(); facei++)
1964 const label own = mesh_.faceOwner()[facei];
1965 const label bFacei = facei - mesh_.nInternalFaces();
1967 if (cellLevel[own] < cellMaxLevel[own] && neiBndMaxLevel[bFacei] != -1)
1975 cellMaxLocation[own],
1977 neiBndMaxLocation[bFacei],
1978 neiBndMaxNormal[bFacei]
1995 Pout<<
"Stopped refining since reaching my cell"
1996 <<
" limit of " << mesh_.nCells()+7*nRefine
2011 Info<<
"Reached refinement limit." <<
endl;
2014 return returnReduce(nRefine-oldNRefine, sumOp<label>());
2028 const scalar curvature,
2029 const scalar planarAngle,
2031 const bool featureRefinement,
2032 const bool featureDistanceRefinement,
2033 const bool internalRefinement,
2034 const bool surfaceRefinement,
2035 const bool curvatureRefinement,
2036 const bool gapRefinement,
2037 const label maxGlobalCells,
2038 const label maxLocalCells
2041 const label totNCells = mesh_.globalData().nTotalCells();
2045 if (totNCells >= maxGlobalCells)
2047 Info<<
"No cells marked for refinement since reached limit "
2048 << maxGlobalCells <<
'.' <<
endl;
2082 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
2083 pointField neiCc(mesh_.nFaces()-mesh_.nInternalFaces());
2084 calcNeighbourData(neiLevel, neiCc);
2091 if (featureRefinement)
2093 const label nFeatures = markFeatureRefinement
2102 Info<<
"Marked for refinement due to explicit features "
2103 <<
": " << nFeatures <<
" cells." <<
endl;
2109 if (featureDistanceRefinement)
2111 const label nCellsFeat = markInternalDistanceToFeatureRefinement
2118 Info<<
"Marked for refinement due to distance to explicit features "
2119 ": " << nCellsFeat <<
" cells." <<
endl;
2125 if (internalRefinement)
2127 const label nShell = markInternalRefinement
2134 Info<<
"Marked for refinement due to refinement shells "
2135 <<
": " << nShell <<
" cells." <<
endl;
2141 if (surfaceRefinement)
2143 const label nSurf = markSurfaceRefinement
2152 Info<<
"Marked for refinement due to surface intersection "
2153 <<
": " << nSurf <<
" cells." <<
endl;
2162 && (curvature >= -1 && curvature <= 1)
2163 && (surfaces_.minLevel() != surfaces_.maxLevel())
2166 const label nCurv = markSurfaceCurvatureRefinement
2176 Info<<
"Marked for refinement due to curvature/regions "
2177 <<
": " << nCurv <<
" cells." <<
endl;
2186 && (planarCos >= -1 && planarCos <= 1)
2187 && (
max(surfaces_.gapLevel()) > -1)
2190 Info<<
"Specified gap level : " <<
max(surfaces_.gapLevel())
2191 <<
", planar angle " << planarAngle <<
endl;
2193 const label nGap = markProximityRefinement
2203 Info<<
"Marked for refinement due to close opposite surfaces "
2204 <<
": " << nGap <<
" cells." <<
endl;
2212 cellsToRefine.
setSize(nRefine);
2219 cellsToRefine[nRefine++] = celli;
2224 return cellsToRefine;
2237 meshCutter_.setRefinement(cellsToRefine, meshMod);
2243 mesh_.topoChange(map);
2246 if (map().hasMotionPoints())
2248 mesh_.movePoints(map().preMotionPoints());
2257 mesh_.setInstance(
name());
2273 const scalar maxLoadUnbalance
2277 refine(cellsToRefine);
2281 Pout<<
"Writing refined but unbalanced " << msg
2282 <<
" mesh to time " <<
name() <<
endl;
2287 mesh_.time().path()/
name()
2289 Pout<<
"Dumped debug data in = "
2290 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2296 Info<<
"Refined mesh in = "
2297 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2298 printMeshInfo(debug,
"After refinement " + msg);
2308 const scalar nIdealCells =
2309 mesh_.globalData().nTotalCells()
2314 mag(1.0-mesh_.nCells()/nIdealCells),
2318 if (unbalance <= maxLoadUnbalance)
2320 Info<<
"Skipping balancing since max unbalance " << unbalance
2321 <<
" is less than allowable " << maxLoadUnbalance
2337 Info<<
"Balanced mesh in = "
2338 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2340 printMeshInfo(debug,
"After balancing " + msg);
2345 Pout<<
"Writing balanced " << msg
2346 <<
" mesh to time " <<
name() <<
endl;
2351 mesh_.time().path()/
name()
2353 Pout<<
"Dumped debug data in = "
2354 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2374 const scalar maxLoadUnbalance
2377 labelList cellsToRefine(initCellsToRefine);
2412 const scalar nNewCells =
2413 scalar(mesh_.nCells() + 7*cellsToRefine.
size());
2414 const scalar nIdealNewCells =
2418 mag(1.0-nNewCells/nIdealNewCells),
2422 if (unbalance <= maxLoadUnbalance)
2424 Info<<
"Skipping balancing since max unbalance " << unbalance
2425 <<
" is less than allowable " << maxLoadUnbalance
2433 cellWeights[cellsToRefine[i]] += 7;
2446 distMap().distributeCellIndices(cellsToRefine);
2448 Info<<
"Balanced mesh in = "
2449 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2464 printMeshInfo(debug,
"After balancing " + msg);
2468 Pout<<
"Writing balanced " << msg
2469 <<
" mesh to time " <<
name() <<
endl;
2474 mesh_.time().path()/
name()
2476 Pout<<
"Dumped debug data in = "
2477 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2488 refine(cellsToRefine);
2492 Pout<<
"Writing refined " << msg
2493 <<
" mesh to time " <<
name() <<
endl;
2498 mesh_.time().path()/
name()
2500 Pout<<
"Dumped debug data in = "
2501 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2507 Info<<
"Refined mesh in = "
2508 << mesh_.time().cpuTimeIncrement() <<
" s" <<
endl;
2522 printMeshInfo(debug,
"After refinement " + msg);
#define forAll(list, i)
Loop across all elements in list.
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
void size(const label)
Override size to be inconsistent with allocated storage.
void clear()
Clear the list, i.e. set size to zero.
void setSize(const label)
Reset size of List.
virtual const fileName & name() const
Return the name of the stream.
virtual Ostream & write(const char)=0
Write character.
static label nProcs(const label communicator=0)
Number of processes in parallel run.
static const direction nComponents
Number of components in this vector space.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Abstract base class for decomposition.
Sends/receives parts of mesh+fvfields to neighbouring processors. Used in load balancing.
label nTotalFaces() const
Return total number of faces in decomposed mesh. Not.
autoPtr< polyDistributionMap > refineAndBalance(const string &msg, decompositionMethod &decomposer, fvMeshDistribute &distributor, const labelList &cellsToRefine, const scalar maxLoadUnbalance)
Refine some cells and rebalance.
labelList refineCandidates(const List< point > &insidePoints, const scalar curvature, const scalar planarAngle, const bool featureRefinement, const bool featureDistanceRefinement, const bool internalRefinement, const bool surfaceRefinement, const bool curvatureRefinement, const bool gapRefinement, const label maxGlobalCells, const label maxLocalCells) const
Calculate list of cells to refine.
const fvMesh & mesh() const
Reference to mesh.
autoPtr< polyTopoChangeMap > refine(const labelList &cellsToRefine)
Refine some cells.
bool isGap(const scalar, const vector &, const vector &, const vector &, const vector &) const
Is local topology a small gap?
autoPtr< polyDistributionMap > balanceAndRefine(const string &msg, decompositionMethod &decomposer, fvMeshDistribute &distributor, const labelList &cellsToRefine, const scalar maxLoadUnbalance)
Balance before refining some cells.
bool isNormalGap(const scalar, const vector &, const vector &, const vector &, const vector &) const
Is local topology a small gap normal to the test vector.
normalLess(const vectorList &values)
bool operator()(const label a, const label b) const
labelList cmptType
Component type.
vectorList cmptType
Component type.
Traits class for primitives.
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< polyTopoChangeMap > changeMesh(polyMesh &mesh, const bool inflate, const bool syncParallel=true, const bool orderCells=false, const bool orderPoints=false)
Inplace changes mesh without change of patches.
label nInternalFaces() const
const cellList & cells() const
Container with cells to refine. Refinement given as single direction.
Contains information about location on a triSurface.
labelList getChangedFaces(const polyMesh &mesh, const labelHashSet &patchIDs)
Get initial set of changed faces.
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
List< label > labelList
A List of labels.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
List< cell > cellList
list of cells
Ostream & endl(Ostream &os)
Add newline and flush stream.
labelList findIndices(const ListType &, typename ListType::const_reference, const label start=0)
Find all occurrences of given element. Linear search.
vectorField pointField
pointField is a vectorField.
List< bool > boolList
Bool container classes.
vector point
Point is a vector.
List< vector > vectorList
A List of vectors.
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
Vector< scalar > vector
A scalar version of the templated Vector.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
List< labelList > labelListList
A List of labelList.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
Field< vector > vectorField
Specialisation of Field<T> for vector.
static bool less(const vector &x, const vector &y)
To compare normals.
prefixOSstream Pout(cout, "Pout")
static const label labelMax
void sortedOrder(const UList< T > &, labelList &order)
Generate the (stable) sort order for the list.
word name(const complex &)
Return a string representation of a complex.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dimensionedScalar cos(const dimensionedScalar &ds)
scalar degToRad(const scalar deg)
Conversion from degrees to radians.