60 x = (x ==
y) ? x : value;
68 void Foam::hexRef8::reorder
89 newElems[newI] = elems[i];
97 void Foam::hexRef8::getFaceInfo
107 if (!mesh_.isInternalFace(facei))
109 patchID = mesh_.boundaryMesh().whichPatch(facei);
112 zoneID = mesh_.faceZones().whichZone(facei);
118 const faceZone& fZone = mesh_.faceZones()[zoneID];
135 label patchID, zoneID, zoneFlip;
137 getFaceInfo(facei, patchID, zoneID, zoneFlip);
141 if ((nei == -1) || (own < nei))
194 const label meshFacei,
195 const label meshPointi,
201 if (mesh_.isInternalFace(meshFacei))
287 void Foam::hexRef8::modFace
296 label patchID, zoneID, zoneFlip;
298 getFaceInfo(facei, patchID, zoneID, zoneFlip);
302 (own != mesh_.faceOwner()[facei])
304 mesh_.isInternalFace(facei)
305 && (nei != mesh_.faceNeighbour()[facei])
307 || (newFace != mesh_.faces()[facei])
310 if ((nei == -1) || (own < nei))
351 Foam::scalar Foam::hexRef8::getLevel0EdgeLength()
const 353 if (cellLevel_.size() != mesh_.nCells())
356 <<
"Number of cells in mesh:" << mesh_.nCells()
357 <<
" does not equal size of cellLevel:" << cellLevel_.size()
359 <<
"This might be because of a restart with inconsistent cellLevel." 366 const scalar great2 =
sqr(great);
383 const label cLevel = cellLevel_[celli];
385 const labelList& cEdges = mesh_.cellEdges(celli);
389 label edgeI = cEdges[i];
391 if (edgeLevel[edgeI] == -1)
393 edgeLevel[edgeI] = cLevel;
395 else if (edgeLevel[edgeI] ==
labelMax)
399 else if (edgeLevel[edgeI] != cLevel)
421 const label eLevel = edgeLevel[edgeI];
423 if (eLevel >= 0 && eLevel <
labelMax)
425 const edge&
e = mesh_.edges()[edgeI];
427 scalar edgeLenSqr =
magSqr(e.
vec(mesh_.points()));
429 typEdgeLenSqr[eLevel] =
min(typEdgeLenSqr[eLevel], edgeLenSqr);
441 Pout<<
"hexRef8::getLevel0EdgeLength() :" 442 <<
" After phase1: Edgelengths (squared) per refinementlevel:" 443 << typEdgeLenSqr <<
endl;
457 const label cLevel = cellLevel_[celli];
459 const labelList& cEdges = mesh_.cellEdges(celli);
463 const edge&
e = mesh_.edges()[cEdges[i]];
465 scalar edgeLenSqr =
magSqr(e.
vec(mesh_.points()));
467 maxEdgeLenSqr[cLevel] =
max(maxEdgeLenSqr[cLevel], edgeLenSqr);
476 Pout<<
"hexRef8::getLevel0EdgeLength() :" 477 <<
" Crappy Edgelengths (squared) per refinementlevel:" 478 << maxEdgeLenSqr <<
endl;
485 forAll(typEdgeLenSqr, levelI)
487 if (typEdgeLenSqr[levelI] == great2 && maxEdgeLenSqr[levelI] >= 0)
489 typEdgeLenSqr[levelI] = maxEdgeLenSqr[levelI];
495 Pout<<
"hexRef8::getLevel0EdgeLength() :" 496 <<
" Final Edgelengths (squared) per refinementlevel:" 497 << typEdgeLenSqr <<
endl;
501 scalar level0Size = -1;
503 forAll(typEdgeLenSqr, levelI)
505 scalar lenSqr = typEdgeLenSqr[levelI];
513 Pout<<
"hexRef8::getLevel0EdgeLength() :" 514 <<
" For level:" << levelI
516 <<
" with equivalent level0 len:" << level0Size
523 if (level0Size == -1)
544 if (cellAnchorPoints[celli].size())
550 return cellAddedCells[celli][index];
557 const face& f = mesh_.faces()[facei];
565 return cellAddedCells[celli][index];
571 Perr<<
"cell:" << celli <<
" anchorPoints:" << cellAnchorPoints[celli]
575 <<
"Could not find point " << pointi
576 <<
" in the anchorPoints for cell " << celli << endl
577 <<
"Does your original mesh obey the 2:1 constraint and" 578 <<
" did you use consistentRefinement to make your cells to refine" 579 <<
" obey this constraint as well?" 592 void Foam::hexRef8::getFaceNeighbours
608 mesh_.faceOwner()[facei],
613 if (mesh_.isInternalFace(facei))
619 mesh_.faceNeighbour()[facei],
639 label level = pointLevel_[f[fp]];
641 if (level < minLevel)
660 label level = pointLevel_[f[fp]];
662 if (level > maxLevel)
676 const label anchorLevel
683 if (pointLevel_[f[fp]] <= anchorLevel)
692 void Foam::hexRef8::dumpCell(
const label celli)
const 695 Pout<<
"hexRef8 : Dumping cell as obj to " << str.
name() <<
endl;
697 const cell& cFaces = mesh_.cells()[celli];
704 const face& f = mesh_.faces()[cFaces[i]];
708 if (pointToObjVert.
insert(f[fp], objVertI))
718 const face& f = mesh_.faces()[cFaces[i]];
722 label pointi = f[fp];
725 str <<
"l " << pointToObjVert[pointi]+1
726 <<
' ' << pointToObjVert[nexPointi]+1 <<
nl;
738 const bool searchForward,
739 const label wantedLevel
746 label pointi = f[fp];
748 if (pointLevel_[pointi] < wantedLevel)
750 dumpCell(mesh_.faceOwner()[facei]);
751 if (mesh_.isInternalFace(facei))
753 dumpCell(mesh_.faceNeighbour()[facei]);
759 <<
" startFp:" << startFp
760 <<
" wantedLevel:" << wantedLevel
763 else if (pointLevel_[pointi] == wantedLevel)
778 dumpCell(mesh_.faceOwner()[facei]);
779 if (mesh_.isInternalFace(facei))
781 dumpCell(mesh_.faceNeighbour()[facei]);
787 <<
" startFp:" << startFp
788 <<
" wantedLevel:" << wantedLevel
798 const face& f = mesh_.faces()[facei];
802 return pointLevel_[f[findMaxLevel(f)]];
806 label ownLevel = cellLevel_[mesh_.faceOwner()[facei]];
808 if (countAnchors(f, ownLevel) == 4)
812 else if (countAnchors(f, ownLevel+1) == 4)
824 void Foam::hexRef8::checkInternalOrientation
837 const vector a(compactFace.
area(compactPoints));
839 const vector dir(neiPt - ownPt);
844 <<
"cell:" << celli <<
" old face:" << facei
845 <<
" newFace:" << newFace <<
endl 846 <<
" coords:" << compactPoints
847 <<
" ownPt:" << ownPt
848 <<
" neiPt:" << neiPt
852 const vector fcToOwn(compactFace.
centre(compactPoints) - ownPt);
854 const scalar
s = (fcToOwn & a) / (dir & a);
856 if (s < 0.1 || s > 0.9)
859 <<
"cell:" << celli <<
" old face:" << facei
860 <<
" newFace:" << newFace <<
endl 861 <<
" coords:" << compactPoints
862 <<
" ownPt:" << ownPt
863 <<
" neiPt:" << neiPt
870 void Foam::hexRef8::checkBoundaryOrientation
876 const point& boundaryPt,
883 const vector a(compactFace.
area(compactPoints));
885 const vector dir(boundaryPt - ownPt);
890 <<
"cell:" << celli <<
" old face:" << facei
891 <<
" newFace:" << newFace
892 <<
" coords:" << compactPoints
893 <<
" ownPt:" << ownPt
894 <<
" boundaryPt:" << boundaryPt
898 const vector fcToOwn(compactFace.
centre(compactPoints) - ownPt);
900 const scalar
s = (fcToOwn&dir) /
magSqr(dir);
902 if (s < 0.7 || s > 1.3)
905 <<
"cell:" << celli <<
" old face:" << facei
906 <<
" newFace:" << newFace
907 <<
" coords:" << compactPoints
908 <<
" ownPt:" << ownPt
909 <<
" boundaryPt:" << boundaryPt
918 void Foam::hexRef8::insertEdgeSplit
926 if (p0 < mesh_.nPoints() && p1 < mesh_.nPoints())
930 if (edgeI != -1 && edgeMidPoint[edgeI] != -1)
932 verts.
append(edgeMidPoint[edgeI]);
954 const bool faceOrder,
955 const label edgeMidPointi,
956 const label anchorPointi,
957 const label faceMidPointi,
966 bool changed =
false;
967 bool haveTwoAnchors =
false;
971 if (edgeMidFnd == midPointToAnchors.
end())
973 midPointToAnchors.
insert(edgeMidPointi,
edge(anchorPointi, -1));
977 edge&
e = edgeMidFnd();
979 if (anchorPointi != e[0])
988 if (e[0] != -1 && e[1] != -1)
990 haveTwoAnchors =
true;
994 bool haveTwoFaceMids =
false;
998 if (faceMidFnd == midPointToFaceMids.
end())
1000 midPointToFaceMids.
insert(edgeMidPointi,
edge(faceMidPointi, -1));
1004 edge&
e = faceMidFnd();
1006 if (faceMidPointi != e[0])
1010 e[1] = faceMidPointi;
1015 if (e[0] != -1 && e[1] != -1)
1017 haveTwoFaceMids =
true;
1024 if (changed && haveTwoAnchors && haveTwoFaceMids)
1026 const edge& anchors = midPointToAnchors[edgeMidPointi];
1027 const edge& faceMids = midPointToFaceMids[edgeMidPointi];
1037 if (faceOrder == (mesh_.faceOwner()[facei] == celli))
1039 newFaceVerts.
append(faceMidPointi);
1050 newFaceVerts.
append(edgeMidPointi);
1060 newFaceVerts.
append(otherFaceMidPointi);
1061 newFaceVerts.
append(cellMidPoint[celli]);
1065 newFaceVerts.
append(otherFaceMidPointi);
1075 newFaceVerts.
append(edgeMidPointi);
1085 newFaceVerts.
append(faceMidPointi);
1086 newFaceVerts.
append(cellMidPoint[celli]);
1092 label anchorCell0 = getAnchorCell
1100 label anchorCell1 = getAnchorCell
1113 if (anchorCell0 < anchorCell1)
1118 ownPt = mesh_.points()[anchorPointi];
1119 neiPt = mesh_.points()[anchors.
otherVertex(anchorPointi)];
1128 ownPt = mesh_.points()[anchors.
otherVertex(anchorPointi)];
1129 neiPt = mesh_.points()[anchorPointi];
1136 if (anchorCell0 < anchorCell1)
1138 ownPt = mesh_.points()[anchorPointi];
1139 neiPt = mesh_.points()[anchors.
otherVertex(anchorPointi)];
1143 ownPt = mesh_.points()[anchors.
otherVertex(anchorPointi)];
1144 neiPt = mesh_.points()[anchorPointi];
1147 checkInternalOrientation
1158 return addInternalFace
1176 void Foam::hexRef8::createInternalFaces
1192 const cell& cFaces = mesh_.cells()[celli];
1193 const label cLevel = cellLevel_[celli];
1205 label nFacesAdded = 0;
1209 label facei = cFaces[i];
1211 const face& f = mesh_.faces()[facei];
1212 const labelList& fEdges = mesh_.faceEdges(facei, storage);
1217 label faceMidPointi = -1;
1219 label nAnchors = countAnchors(f, cLevel);
1227 label anchorFp = -1;
1231 if (pointLevel_[f[fp]] <= cLevel)
1239 label edgeMid = findLevel
1247 label faceMid = findLevel
1256 faceMidPointi = f[faceMid];
1258 else if (nAnchors == 4)
1263 faceMidPointi = faceMidPoint[facei];
1267 dumpCell(mesh_.faceOwner()[facei]);
1268 if (mesh_.isInternalFace(facei))
1270 dumpCell(mesh_.faceNeighbour()[facei]);
1274 <<
"nAnchors:" << nAnchors
1275 <<
" facei:" << facei
1287 label point0 = f[fp0];
1289 if (pointLevel_[point0] <= cLevel)
1298 label edgeMidPointi = -1;
1302 if (pointLevel_[f[fp1]] <= cLevel)
1305 label edgeI = fEdges[fp0];
1307 edgeMidPointi = edgeMidPoint[edgeI];
1309 if (edgeMidPointi == -1)
1313 const labelList& cPoints = mesh_.cellPoints(celli);
1316 <<
"cell:" << celli <<
" cLevel:" << cLevel
1317 <<
" cell points:" << cPoints
1320 <<
" face:" << facei
1324 <<
" faceAnchorLevel:" << faceAnchorLevel[facei]
1325 <<
" faceMidPoint:" << faceMidPoint[facei]
1326 <<
" faceMidPointi:" << faceMidPointi
1334 label edgeMid = findLevel(facei, f, fp1,
true, cLevel+1);
1336 edgeMidPointi = f[edgeMid];
1339 label newFacei = storeMidPointInfo
1362 if (nFacesAdded == 12)
1375 if (pointLevel_[f[fpMin1]] <= cLevel)
1378 label edgeI = fEdges[fpMin1];
1380 edgeMidPointi = edgeMidPoint[edgeI];
1382 if (edgeMidPointi == -1)
1386 const labelList& cPoints = mesh_.cellPoints(celli);
1389 <<
"cell:" << celli <<
" cLevel:" << cLevel
1390 <<
" cell points:" << cPoints
1393 <<
" face:" << facei
1397 <<
" faceAnchorLevel:" << faceAnchorLevel[facei]
1398 <<
" faceMidPoint:" << faceMidPoint[facei]
1399 <<
" faceMidPointi:" << faceMidPointi
1407 label edgeMid = findLevel
1416 edgeMidPointi = f[edgeMid];
1419 newFacei = storeMidPointInfo
1442 if (nFacesAdded == 12)
1450 if (nFacesAdded == 12)
1458 void Foam::hexRef8::walkFaceToMid
1463 const label startFp,
1467 const face& f = mesh_.faces()[facei];
1468 const labelList& fEdges = mesh_.faceEdges(facei);
1477 if (edgeMidPoint[fEdges[fp]] >= 0)
1479 faceVerts.
append(edgeMidPoint[fEdges[fp]]);
1484 if (pointLevel_[f[fp]] <= cLevel)
1490 else if (pointLevel_[f[fp]] == cLevel+1)
1497 else if (pointLevel_[f[fp]] == cLevel+2)
1507 void Foam::hexRef8::walkFaceFromMid
1512 const label startFp,
1516 const face& f = mesh_.faces()[facei];
1517 const labelList& fEdges = mesh_.faceEdges(facei);
1523 if (pointLevel_[f[fp]] <= cLevel)
1528 else if (pointLevel_[f[fp]] == cLevel+1)
1534 else if (pointLevel_[f[fp]] == cLevel+2)
1544 if (edgeMidPoint[fEdges[fp]] >= 0)
1546 faceVerts.
append(edgeMidPoint[fEdges[fp]]);
1562 Foam::label Foam::hexRef8::faceConsistentRefinement
1571 for (
label facei = 0; facei < mesh_.nInternalFaces(); facei++)
1573 label own = mesh_.faceOwner()[facei];
1574 label ownLevel = cellLevel_[own] + refineCell.
get(own);
1576 label nei = mesh_.faceNeighbour()[facei];
1577 label neiLevel = cellLevel_[nei] + refineCell.
get(nei);
1579 if (ownLevel > (neiLevel+1))
1583 refineCell.
set(nei);
1587 refineCell.
unset(own);
1591 else if (neiLevel > (ownLevel+1))
1595 refineCell.
set(own);
1599 refineCell.
unset(nei);
1608 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
1612 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
1614 neiLevel[i] = cellLevel_[own] + refineCell.
get(own);
1623 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
1624 label ownLevel = cellLevel_[own] + refineCell.
get(own);
1626 if (ownLevel > (neiLevel[i]+1))
1630 refineCell.
unset(own);
1634 else if (neiLevel[i] > (ownLevel+1))
1638 refineCell.
set(own);
1649 void Foam::hexRef8::checkWantedRefinementLevels
1657 refineCell.
set(cellsToRefine[i]);
1660 for (
label facei = 0; facei < mesh_.nInternalFaces(); facei++)
1662 label own = mesh_.faceOwner()[facei];
1663 label ownLevel = cellLevel_[own] + refineCell.
get(own);
1665 label nei = mesh_.faceNeighbour()[facei];
1666 label neiLevel = cellLevel_[nei] + refineCell.
get(nei);
1668 if (
mag(ownLevel-neiLevel) > 1)
1674 <<
" current level:" << cellLevel_[own]
1675 <<
" level after refinement:" << ownLevel
1677 <<
"neighbour cell:" << nei
1678 <<
" current level:" << cellLevel_[nei]
1679 <<
" level after refinement:" << neiLevel
1681 <<
"which does not satisfy 2:1 constraints anymore." 1688 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
1692 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
1694 neiLevel[i] = cellLevel_[own] + refineCell.
get(own);
1703 label facei = i + mesh_.nInternalFaces();
1705 label own = mesh_.faceOwner()[facei];
1706 label ownLevel = cellLevel_[own] + refineCell.
get(own);
1708 if (
mag(ownLevel - neiLevel[i]) > 1)
1710 label patchi = mesh_.boundaryMesh().whichPatch(facei);
1714 <<
"Celllevel does not satisfy 2:1 constraint." 1715 <<
" On coupled face " 1717 <<
" on patch " << patchi <<
" " 1718 << mesh_.boundaryMesh()[
patchi].name()
1719 <<
" owner cell " << own
1720 <<
" current level:" << cellLevel_[own]
1721 <<
" level after refinement:" << ownLevel
1723 <<
" (coupled) neighbour cell will get refinement " 1736 Pout<<
"hexRef8::setInstance(const fileName& inst) : " 1737 <<
"Resetting file instance to " << inst <<
endl;
1740 cellLevel_.instance() = inst;
1741 pointLevel_.instance() = inst;
1742 level0Edge_.instance() = inst;
1743 history_.instance() = inst;
1747 void Foam::hexRef8::collectLevelPoints
1756 if (pointLevel_[f[fp]] <= level)
1764 void Foam::hexRef8::collectLevelPoints
1774 label pointi = meshPoints[f[fp]];
1775 if (pointLevel_[pointi] <= level)
1784 bool Foam::hexRef8::matchHexShape
1787 const label cellLevel,
1791 const cell& cFaces = mesh_.cells()[celli];
1802 label facei = cFaces[i];
1803 const face& f = mesh_.faces()[facei];
1806 collectLevelPoints(f, cellLevel, verts);
1807 if (verts.
size() == 4)
1809 if (mesh_.faceOwner()[facei] != celli)
1820 if (quads.
size() < 6)
1826 label facei = cFaces[i];
1827 const face& f = mesh_.faces()[facei];
1835 collectLevelPoints(f, cellLevel, verts);
1836 if (verts.
size() == 1)
1842 label pointi = f[fp];
1843 if (pointLevel_[pointi] == cellLevel+1)
1846 pointFaces.find(pointi);
1847 if (iter != pointFaces.end())
1873 if (pFaces.
size() == 4)
1879 label facei = pFaces[pFacei];
1880 const face& f = mesh_.faces()[facei];
1881 if (mesh_.faceOwner()[facei] == celli)
1883 fourFaces[pFacei] =
f;
1898 if (edgeLoops.size() == 1)
1904 bigFace.meshPoints(),
1905 bigFace.edgeLoops()[0],
1910 if (verts.
size() == 4)
1921 return (quads.
size() == 6);
1928 Foam::hexRef8::hexRef8(
const polyMesh& mesh,
const bool readHistory)
1936 mesh_.facesInstance(),
1949 mesh_.facesInstance(),
1962 mesh_.facesInstance(),
1974 "refinementHistory",
1975 mesh_.facesInstance(),
1982 (readHistory ? mesh_.nCells() : 0)
1984 faceRemover_(mesh_, great),
1985 savedPointLevel_(0),
2006 <<
"History enabled but number of visible cells " 2009 <<
" is not equal to the number of cells in the mesh " 2016 cellLevel_.
size() != mesh_.nCells()
2017 || pointLevel_.
size() != mesh_.nPoints()
2021 <<
"Restarted from inconsistent cellLevel or pointLevel files." 2025 <<
"Number of cells in mesh:" << mesh_.nCells()
2026 <<
" does not equal size of cellLevel:" << cellLevel_.
size() <<
endl 2027 <<
"Number of points in mesh:" << mesh_.nPoints()
2028 <<
" does not equal size of pointLevel:" << pointLevel_.
size()
2047 Foam::hexRef8::hexRef8
2053 const scalar level0Edge
2062 mesh_.facesInstance(),
2075 mesh_.facesInstance(),
2088 mesh_.facesInstance(),
2098 (level0Edge >= 0 ? level0Edge : getLevel0EdgeLength())
2105 "refinementHistory",
2106 mesh_.facesInstance(),
2114 faceRemover_(mesh_, great),
2115 savedPointLevel_(0),
2121 <<
"History enabled but number of visible cells in it " 2123 <<
" is not equal to the number of cells in the mesh " 2129 cellLevel_.
size() != mesh_.nCells()
2130 || pointLevel_.
size() != mesh_.nPoints()
2134 <<
"Incorrect cellLevel or pointLevel size." <<
endl 2135 <<
"Number of cells in mesh:" << mesh_.nCells()
2136 <<
" does not equal size of cellLevel:" << cellLevel_.
size() <<
endl 2137 <<
"Number of points in mesh:" << mesh_.nPoints()
2138 <<
" does not equal size of pointLevel:" << pointLevel_.
size()
2156 Foam::hexRef8::hexRef8
2161 const scalar level0Edge
2170 mesh_.facesInstance(),
2183 mesh_.facesInstance(),
2196 mesh_.facesInstance(),
2206 (level0Edge >= 0 ? level0Edge : getLevel0EdgeLength())
2213 "refinementHistory",
2214 mesh_.facesInstance(),
2224 faceRemover_(mesh_, great),
2225 savedPointLevel_(0),
2230 cellLevel_.
size() != mesh_.nCells()
2231 || pointLevel_.
size() != mesh_.nPoints()
2235 <<
"Incorrect cellLevel or pointLevel size." <<
endl 2236 <<
"Number of cells in mesh:" << mesh_.nCells()
2237 <<
" does not equal size of cellLevel:" << cellLevel_.
size() <<
endl 2238 <<
"Number of points in mesh:" << mesh_.nPoints()
2239 <<
" does not equal size of pointLevel:" << pointLevel_.
size()
2272 refineCell.
set(cellsToRefine[i]);
2277 label nChanged = faceConsistentRefinement(maxSet, refineCell);
2283 Pout<<
"hexRef8::consistentRefinement : Changed " << nChanged
2284 <<
" refinement levels due to 2:1 conflicts." 2298 forAll(refineCell, celli)
2300 if (refineCell.
get(celli))
2309 forAll(refineCell, celli)
2311 if (refineCell.
get(celli))
2313 newCellsToRefine[nRefined++] = celli;
2319 checkWantedRefinementLevels(newCellsToRefine);
2322 return newCellsToRefine;
2334 const label maxFaceDiff,
2337 const label maxPointDiff,
2341 const labelList& faceOwner = mesh_.faceOwner();
2342 const labelList& faceNeighbour = mesh_.faceNeighbour();
2345 if (maxFaceDiff <= 0)
2348 <<
"Illegal maxFaceDiff " << maxFaceDiff <<
nl 2367 forAll(allCellInfo, celli)
2373 maxFaceDiff*(cellLevel_[celli]+1),
2374 maxFaceDiff*cellLevel_[celli]
2381 label celli = cellsToRefine[i];
2383 allCellInfo[celli].count() = allCellInfo[celli].refinementCount();
2393 int dummyTrackData = 0;
2401 label facei = facesToCheck[i];
2403 if (allFaceInfo[facei].valid(dummyTrackData))
2407 <<
"Argument facesToCheck seems to have duplicate entries!" 2409 <<
"face:" << facei <<
" occurs at positions " 2417 if (mesh_.isInternalFace(facei))
2422 const refinementData& neiData = allCellInfo[faceNeighbour[facei]];
2425 label faceRefineCount;
2428 faceCount = neiData.
count() + maxFaceDiff;
2429 faceRefineCount = faceCount + maxFaceDiff;
2433 faceCount = ownData.
count() + maxFaceDiff;
2434 faceRefineCount = faceCount + maxFaceDiff;
2437 seedFaces.append(facei);
2438 seedFacesInfo.append
2446 allFaceInfo[facei] = seedFacesInfo.last();
2450 label faceCount = ownData.
count() + maxFaceDiff;
2451 label faceRefineCount = faceCount + maxFaceDiff;
2453 seedFaces.append(facei);
2454 seedFacesInfo.append
2462 allFaceInfo[facei] = seedFacesInfo.last();
2470 forAll(faceNeighbour, facei)
2473 if (!allFaceInfo[facei].valid(dummyTrackData))
2475 label own = faceOwner[facei];
2476 label nei = faceNeighbour[facei];
2480 if (allCellInfo[own].count() > allCellInfo[nei].count())
2482 allFaceInfo[facei].updateFace
2492 seedFacesInfo.append(allFaceInfo[facei]);
2494 else if (allCellInfo[own].count() < allCellInfo[nei].count())
2496 allFaceInfo[facei].updateFace
2505 seedFaces.append(facei);
2506 seedFacesInfo.append(allFaceInfo[facei]);
2514 for (
label facei = mesh_.nInternalFaces(); facei < mesh_.nFaces(); facei++)
2517 if (!allFaceInfo[facei].valid(dummyTrackData))
2519 label own = faceOwner[facei];
2532 seedFaces.append(facei);
2533 seedFacesInfo.append(faceData);
2551 Pout<<
"hexRef8::consistentSlowRefinement : Seeded " 2552 << seedFaces.size() <<
" faces between cells with different" 2553 <<
" refinement level." <<
endl;
2557 levelCalc.
setFaceInfo(seedFaces.shrink(), seedFacesInfo.shrink());
2559 seedFacesInfo.clear();
2562 levelCalc.
iterate(mesh_.globalData().nTotalFaces()+1);
2570 if (maxPointDiff == -1)
2578 labelList maxPointCount(mesh_.nPoints(), 0);
2580 forAll(maxPointCount, pointi)
2582 label& pLevel = maxPointCount[pointi];
2584 const labelList& pCells = mesh_.pointCells(pointi);
2588 pLevel =
max(pLevel, allCellInfo[pCells[i]].count());
2609 label pointi = pointsToCheck[i];
2614 const labelList& pCells = mesh_.pointCells(pointi);
2618 label celli = pCells[pCelli];
2626 maxPointCount[pointi]
2627 > cellInfo.
count() + maxFaceDiff*maxPointDiff
2635 const cell& cFaces = mesh_.cells()[celli];
2639 label facei = cFaces[cFacei];
2652 if (faceData.
count() > allFaceInfo[facei].count())
2654 changedFacesInfo.insert(facei, faceData);
2661 label nChanged = changedFacesInfo.size();
2671 seedFaces.setCapacity(changedFacesInfo.size());
2672 seedFacesInfo.setCapacity(changedFacesInfo.size());
2676 seedFaces.append(iter.key());
2677 seedFacesInfo.append(iter());
2684 for (
label facei = 0; facei < mesh_.nInternalFaces(); facei++)
2686 label own = mesh_.faceOwner()[facei];
2689 + (allCellInfo[own].isRefined() ? 1 : 0);
2691 label nei = mesh_.faceNeighbour()[facei];
2694 + (allCellInfo[nei].isRefined() ? 1 : 0);
2696 if (
mag(ownLevel-neiLevel) > 1)
2702 <<
" current level:" << cellLevel_[own]
2703 <<
" current refData:" << allCellInfo[own]
2704 <<
" level after refinement:" << ownLevel
2706 <<
"neighbour cell:" << nei
2707 <<
" current level:" << cellLevel_[nei]
2708 <<
" current refData:" << allCellInfo[nei]
2709 <<
" level after refinement:" << neiLevel
2711 <<
"which does not satisfy 2:1 constraints anymore." <<
nl 2712 <<
"face:" << facei <<
" faceRefData:" << allFaceInfo[facei]
2721 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
2722 labelList neiCount(mesh_.nFaces()-mesh_.nInternalFaces());
2723 labelList neiRefCount(mesh_.nFaces()-mesh_.nInternalFaces());
2727 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
2728 neiLevel[i] = cellLevel_[own];
2729 neiCount[i] = allCellInfo[own].count();
2730 neiRefCount[i] = allCellInfo[own].refinementCount();
2741 label facei = i+mesh_.nInternalFaces();
2743 label own = mesh_.faceOwner()[facei];
2746 + (allCellInfo[own].isRefined() ? 1 : 0);
2750 + ((neiCount[i] >= neiRefCount[i]) ? 1 : 0);
2752 if (
mag(ownLevel - nbrLevel) > 1)
2755 label patchi = mesh_.boundaryMesh().whichPatch(facei);
2758 <<
"Celllevel does not satisfy 2:1 constraint." 2759 <<
" On coupled face " 2761 <<
" refData:" << allFaceInfo[facei]
2762 <<
" on patch " << patchi <<
" " 2763 << mesh_.boundaryMesh()[
patchi].name() <<
nl 2764 <<
"owner cell " << own
2765 <<
" current level:" << cellLevel_[own]
2766 <<
" current count:" << allCellInfo[own].count()
2767 <<
" current refCount:" 2768 << allCellInfo[own].refinementCount()
2769 <<
" level after refinement:" << ownLevel
2771 <<
"(coupled) neighbour cell" 2772 <<
" has current level:" << neiLevel[i]
2773 <<
" current count:" << neiCount[i]
2774 <<
" current refCount:" << neiRefCount[i]
2775 <<
" level after refinement:" << nbrLevel
2785 forAll(allCellInfo, celli)
2787 if (allCellInfo[celli].isRefined())
2797 forAll(allCellInfo, celli)
2799 if (allCellInfo[celli].isRefined())
2801 newCellsToRefine[nRefined++] = celli;
2807 Pout<<
"hexRef8::consistentSlowRefinement : From " 2808 << cellsToRefine.
size() <<
" to " << newCellsToRefine.
size()
2809 <<
" cells to refine." <<
endl;
2812 return newCellsToRefine;
2818 const label maxFaceDiff,
2823 const labelList& faceOwner = mesh_.faceOwner();
2824 const labelList& faceNeighbour = mesh_.faceNeighbour();
2826 if (maxFaceDiff <= 0)
2829 <<
"Illegal maxFaceDiff " << maxFaceDiff <<
nl 2848 int dummyTrackData = 0;
2854 label celli = cellsToRefine[i];
2859 mesh_.cellCentres()[celli],
2864 forAll(allCellInfo, celli)
2866 if (!allCellInfo[celli].valid(dummyTrackData))
2871 mesh_.cellCentres()[celli],
2887 label facei = facesToCheck[i];
2889 if (allFaceInfo[facei].valid(dummyTrackData))
2893 <<
"Argument facesToCheck seems to have duplicate entries!" 2895 <<
"face:" << facei <<
" occurs at positions " 2900 label own = faceOwner[facei];
2904 allCellInfo[own].valid(dummyTrackData)
2905 ? allCellInfo[own].originLevel()
2909 if (!mesh_.isInternalFace(facei))
2913 const point& fc = mesh_.faceCentres()[facei];
2922 allFaceInfo[facei].updateFace
2934 label nei = faceNeighbour[facei];
2938 allCellInfo[nei].valid(dummyTrackData)
2939 ? allCellInfo[nei].originLevel()
2943 if (ownLevel == neiLevel)
2946 allFaceInfo[facei].updateFace
2955 allFaceInfo[facei].updateFace
2968 allFaceInfo[facei].updateFace
2977 allFaceInfo[facei].updateFace
2989 seedFacesInfo.append(allFaceInfo[facei]);
2996 forAll(faceNeighbour, facei)
2999 if (!allFaceInfo[facei].valid(dummyTrackData))
3001 label own = faceOwner[facei];
3005 allCellInfo[own].valid(dummyTrackData)
3006 ? allCellInfo[own].originLevel()
3010 label nei = faceNeighbour[facei];
3014 allCellInfo[nei].valid(dummyTrackData)
3015 ? allCellInfo[nei].originLevel()
3019 if (ownLevel > neiLevel)
3023 allFaceInfo[facei].updateFace
3032 seedFacesInfo.append(allFaceInfo[facei]);
3034 else if (neiLevel > ownLevel)
3036 seedFaces.append(facei);
3037 allFaceInfo[facei].updateFace
3046 seedFacesInfo.append(allFaceInfo[facei]);
3052 seedFacesInfo.shrink();
3062 mesh_.globalData().nTotalCells()+1,
3103 label celli = cellsToRefine[i];
3105 allCellInfo[celli].originLevel() =
sizeof(
label)*8-2;
3106 allCellInfo[celli].origin() = cc[celli];
3113 forAll(allCellInfo, celli)
3115 label wanted = allCellInfo[celli].wantedLevel(cc[celli]);
3117 if (wanted > cellLevel_[celli]+1)
3119 refineCell.
set(celli);
3122 faceConsistentRefinement(
true, refineCell);
3126 label nChanged = faceConsistentRefinement(
true, refineCell);
3132 Pout<<
"hexRef8::consistentSlowRefinement2 : Changed " << nChanged
3133 <<
" refinement levels due to 2:1 conflicts." 3146 forAll(refineCell, celli)
3149 if (refineCell[celli])
3158 forAll(refineCell, celli)
3161 if (refineCell[celli])
3163 newCellsToRefine[nRefined++] = celli;
3169 Pout<<
"hexRef8::consistentSlowRefinement2 : From " 3170 << cellsToRefine.
size() <<
" to " << newCellsToRefine.
size()
3171 <<
" cells to refine." <<
endl;
3176 cellSet cellsIn(mesh_,
"cellsToRefineIn", cellsToRefine);
3177 Pout<<
"hexRef8::consistentSlowRefinement2 : writing " 3178 << cellsIn.
size() <<
" to cellSet " 3183 cellSet cellsOut(mesh_,
"cellsToRefineOut", newCellsToRefine);
3184 Pout<<
"hexRef8::consistentSlowRefinement2 : writing " 3185 << cellsOut.
size() <<
" to cellSet " 3192 forAll(newCellsToRefine, i)
3194 refineCell.
set(newCellsToRefine[i]);
3198 label nChanged = faceConsistentRefinement(
true, refineCell);
3203 mesh_,
"cellsToRefineOut2", newCellsToRefine.
size()
3205 forAll(refineCell, celli)
3207 if (refineCell.
get(celli))
3209 cellsOut2.insert(celli);
3212 Pout<<
"hexRef8::consistentSlowRefinement2 : writing " 3213 << cellsOut2.size() <<
" to cellSet " 3214 << cellsOut2.objectPath() <<
endl;
3220 forAll(refineCell, celli)
3222 if (refineCell.
get(celli) && !savedRefineCell.
get(celli))
3226 <<
"Cell:" << celli <<
" cc:" 3227 << mesh_.cellCentres()[celli]
3228 <<
" was not marked for refinement but does not obey" 3229 <<
" 2:1 constraints." 3236 return newCellsToRefine;
3249 Pout<<
"hexRef8::setRefinement :" 3250 <<
" Checking initial mesh just to make sure" <<
endl;
3259 savedPointLevel_.
clear();
3260 savedCellLevel_.
clear();
3265 forAll(cellLevel_, celli)
3267 newCellLevel.append(cellLevel_[celli]);
3270 forAll(pointLevel_, pointi)
3272 newPointLevel.append(pointLevel_[pointi]);
3278 Pout<<
"hexRef8::setRefinement :" 3279 <<
" Allocating " << cellLabels.
size() <<
" cell midpoints." 3287 labelList cellMidPoint(mesh_.nCells(), -1);
3291 label celli = cellLabels[i];
3293 label anchorPointi = mesh_.faces()[mesh_.cells()[celli][0]][0];
3299 mesh_.cellCentres()[celli],
3306 newPointLevel(cellMidPoint[celli]) = cellLevel_[celli]+1;
3312 cellSet splitCells(mesh_,
"splitCells", cellLabels.
size());
3314 forAll(cellMidPoint, celli)
3316 if (cellMidPoint[celli] >= 0)
3318 splitCells.insert(celli);
3322 Pout<<
"hexRef8::setRefinement : Dumping " << splitCells.size()
3323 <<
" cells to split to cellSet " << splitCells.objectPath()
3336 Pout<<
"hexRef8::setRefinement :" 3337 <<
" Allocating edge midpoints." 3346 labelList edgeMidPoint(mesh_.nEdges(), -1);
3349 forAll(cellMidPoint, celli)
3351 if (cellMidPoint[celli] >= 0)
3353 const labelList& cEdges = mesh_.cellEdges(celli);
3357 label edgeI = cEdges[i];
3359 const edge&
e = mesh_.edges()[edgeI];
3363 pointLevel_[e[0]] <= cellLevel_[celli]
3364 && pointLevel_[e[1]] <= cellLevel_[celli]
3367 edgeMidPoint[edgeI] = 12345;
3394 forAll(edgeMidPoint, edgeI)
3396 if (edgeMidPoint[edgeI] >= 0)
3399 edgeMids[edgeI] = mesh_.edges()[edgeI].centre(mesh_.points());
3407 point(-great, -great, -great)
3412 forAll(edgeMidPoint, edgeI)
3414 if (edgeMidPoint[edgeI] >= 0)
3419 const edge&
e = mesh_.edges()[edgeI];
3432 newPointLevel(edgeMidPoint[edgeI]) =
3445 OFstream str(mesh_.time().path()/
"edgeMidPoint.obj");
3447 forAll(edgeMidPoint, edgeI)
3449 if (edgeMidPoint[edgeI] >= 0)
3451 const edge&
e = mesh_.edges()[edgeI];
3457 Pout<<
"hexRef8::setRefinement :" 3458 <<
" Dumping edge centres to split to file " << str.
name() <<
endl;
3468 Pout<<
"hexRef8::setRefinement :" 3469 <<
" Allocating face midpoints." 3475 labelList faceAnchorLevel(mesh_.nFaces());
3477 for (
label facei = 0; facei < mesh_.nFaces(); facei++)
3479 faceAnchorLevel[facei] =
faceLevel(facei);
3484 labelList faceMidPoint(mesh_.nFaces(), -1);
3489 for (
label facei = 0; facei < mesh_.nInternalFaces(); facei++)
3491 if (faceAnchorLevel[facei] >= 0)
3493 label own = mesh_.faceOwner()[facei];
3494 label ownLevel = cellLevel_[own];
3495 label newOwnLevel = ownLevel + (cellMidPoint[own] >= 0 ? 1 : 0);
3497 label nei = mesh_.faceNeighbour()[facei];
3498 label neiLevel = cellLevel_[nei];
3499 label newNeiLevel = neiLevel + (cellMidPoint[nei] >= 0 ? 1 : 0);
3503 newOwnLevel > faceAnchorLevel[facei]
3504 || newNeiLevel > faceAnchorLevel[facei]
3507 faceMidPoint[facei] = 12345;
3520 labelList newNeiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
3524 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
3525 label ownLevel = cellLevel_[own];
3526 label newOwnLevel = ownLevel + (cellMidPoint[own] >= 0 ? 1 : 0);
3528 newNeiLevel[i] = newOwnLevel;
3538 label facei = i+mesh_.nInternalFaces();
3540 if (faceAnchorLevel[facei] >= 0)
3542 label own = mesh_.faceOwner()[facei];
3543 label ownLevel = cellLevel_[own];
3544 label newOwnLevel = ownLevel + (cellMidPoint[own] >= 0 ? 1 : 0);
3548 newOwnLevel > faceAnchorLevel[facei]
3549 || newNeiLevel[i] > faceAnchorLevel[facei]
3552 faceMidPoint[facei] = 12345;
3577 mesh_.nFaces()-mesh_.nInternalFaces(),
3578 point(-great, -great, -great)
3583 label facei = i+mesh_.nInternalFaces();
3585 if (faceMidPoint[facei] >= 0)
3587 bFaceMids[i] = mesh_.faceCentres()[facei];
3597 forAll(faceMidPoint, facei)
3599 if (faceMidPoint[facei] >= 0)
3604 const face& f = mesh_.faces()[facei];
3611 facei < mesh_.nInternalFaces()
3612 ? mesh_.faceCentres()[facei]
3613 : bFaceMids[facei-mesh_.nInternalFaces()]
3623 newPointLevel(faceMidPoint[facei]) = faceAnchorLevel[facei]+1;
3630 faceSet splitFaces(mesh_,
"splitFaces", cellLabels.
size());
3632 forAll(faceMidPoint, facei)
3634 if (faceMidPoint[facei] >= 0)
3636 splitFaces.insert(facei);
3640 Pout<<
"hexRef8::setRefinement : Dumping " << splitFaces.size()
3641 <<
" faces to split to faceSet " << splitFaces.objectPath() <<
endl;
3662 Pout<<
"hexRef8::setRefinement :" 3663 <<
" Finding cell anchorPoints (8 per cell)" 3674 labelList nAnchorPoints(mesh_.nCells(), 0);
3676 forAll(cellMidPoint, celli)
3678 if (cellMidPoint[celli] >= 0)
3680 cellAnchorPoints[celli].
setSize(8);
3684 forAll(pointLevel_, pointi)
3686 const labelList& pCells = mesh_.pointCells(pointi);
3690 label celli = pCells[pCelli];
3694 cellMidPoint[celli] >= 0
3695 && pointLevel_[pointi] <= cellLevel_[celli]
3698 if (nAnchorPoints[celli] == 8)
3703 <<
" of level " << cellLevel_[celli]
3704 <<
" uses more than 8 points of equal or" 3705 <<
" lower level" <<
nl 3706 <<
"Points so far:" << cellAnchorPoints[celli]
3710 cellAnchorPoints[celli][nAnchorPoints[celli]++] = pointi;
3715 forAll(cellMidPoint, celli)
3717 if (cellMidPoint[celli] >= 0)
3719 if (nAnchorPoints[celli] != 8)
3723 const labelList& cPoints = mesh_.cellPoints(celli);
3727 <<
" of level " << cellLevel_[celli]
3728 <<
" does not seem to have 8 points of equal or" 3729 <<
" lower level" <<
endl 3730 <<
"cellPoints:" << cPoints <<
endl 3745 Pout<<
"hexRef8::setRefinement :" 3746 <<
" Adding cells (1 per anchorPoint)" 3753 forAll(cellAnchorPoints, celli)
3755 const labelList& cAnchors = cellAnchorPoints[celli];
3757 if (cAnchors.
size() == 8)
3759 labelList& cAdded = cellAddedCells[celli];
3766 newCellLevel[celli] = cellLevel_[celli]+1;
3769 for (
label i = 1; i < 8; i++)
3779 mesh_.cellZones().whichZone(celli)
3783 newCellLevel(cAdded[i]) = cellLevel_[celli]+1;
3798 Pout<<
"hexRef8::setRefinement :" 3799 <<
" Marking faces to be handled" 3807 forAll(cellMidPoint, celli)
3809 if (cellMidPoint[celli] >= 0)
3811 const cell& cFaces = mesh_.cells()[celli];
3815 affectedFace.set(cFaces[i]);
3820 forAll(faceMidPoint, facei)
3822 if (faceMidPoint[facei] >= 0)
3824 affectedFace.set(facei);
3828 forAll(edgeMidPoint, edgeI)
3830 if (edgeMidPoint[edgeI] >= 0)
3832 const labelList& eFaces = mesh_.edgeFaces(edgeI);
3836 affectedFace.set(eFaces[i]);
3848 Pout<<
"hexRef8::setRefinement : Splitting faces" <<
endl;
3851 forAll(faceMidPoint, facei)
3853 if (faceMidPoint[facei] >= 0 && affectedFace.get(facei))
3859 const face& f = mesh_.faces()[facei];
3863 bool modifiedFace =
false;
3865 label anchorLevel = faceAnchorLevel[facei];
3871 label pointi = f[fp];
3873 if (pointLevel_[pointi] <= anchorLevel)
3879 faceVerts.
append(pointi);
3895 faceVerts.
append(faceMidPoint[facei]);
3928 if (mesh_.isInternalFace(facei))
3930 label oldOwn = mesh_.faceOwner()[facei];
3931 label oldNei = mesh_.faceNeighbour()[facei];
3933 checkInternalOrientation
3938 mesh_.cellCentres()[oldOwn],
3939 mesh_.cellCentres()[oldNei],
3945 label oldOwn = mesh_.faceOwner()[facei];
3947 checkBoundaryOrientation
3952 mesh_.cellCentres()[oldOwn],
3953 mesh_.faceCentres()[facei],
3962 modifiedFace =
true;
3964 modFace(meshMod, facei, newFace, own, nei);
3968 addFace(meshMod, facei, newFace, own, nei);
3974 affectedFace.unset(facei);
3984 Pout<<
"hexRef8::setRefinement :" 3985 <<
" Adding edge splits to unsplit faces" 3992 forAll(edgeMidPoint, edgeI)
3994 if (edgeMidPoint[edgeI] >= 0)
3998 const labelList& eFaces = mesh_.edgeFaces(edgeI, eFacesStorage);
4002 label facei = eFaces[i];
4004 if (faceMidPoint[facei] < 0 && affectedFace.get(facei))
4008 const face& f = mesh_.faces()[facei];
4009 const labelList& fEdges = mesh_.faceEdges
4019 newFaceVerts.append(f[fp]);
4021 label edgeI = fEdges[fp];
4023 if (edgeMidPoint[edgeI] >= 0)
4025 newFaceVerts.
append(edgeMidPoint[edgeI]);
4034 label anchorFp = findMinLevel(f);
4051 if (mesh_.isInternalFace(facei))
4053 label oldOwn = mesh_.faceOwner()[facei];
4054 label oldNei = mesh_.faceNeighbour()[facei];
4056 checkInternalOrientation
4061 mesh_.cellCentres()[oldOwn],
4062 mesh_.cellCentres()[oldNei],
4068 label oldOwn = mesh_.faceOwner()[facei];
4070 checkBoundaryOrientation
4075 mesh_.cellCentres()[oldOwn],
4076 mesh_.faceCentres()[facei],
4082 modFace(meshMod, facei, newFace, own, nei);
4085 affectedFace.unset(facei);
4097 Pout<<
"hexRef8::setRefinement :" 4098 <<
" Changing owner/neighbour for otherwise unaffected faces" 4102 forAll(affectedFace, facei)
4104 if (affectedFace.get(facei))
4106 const face& f = mesh_.faces()[facei];
4110 label anchorFp = findMinLevel(f);
4124 modFace(meshMod, facei, f, own, nei);
4127 affectedFace.unset(facei);
4142 Pout<<
"hexRef8::setRefinement :" 4143 <<
" Create new internal faces for split cells" 4147 forAll(cellMidPoint, celli)
4149 if (cellMidPoint[celli] >= 0)
4174 forAll(cellMidPoint, celli)
4176 if (cellMidPoint[celli] >= 0)
4178 minPointi =
min(minPointi, cellMidPoint[celli]);
4179 maxPointi =
max(maxPointi, cellMidPoint[celli]);
4182 forAll(faceMidPoint, facei)
4184 if (faceMidPoint[facei] >= 0)
4186 minPointi =
min(minPointi, faceMidPoint[facei]);
4187 maxPointi =
max(maxPointi, faceMidPoint[facei]);
4190 forAll(edgeMidPoint, edgeI)
4192 if (edgeMidPoint[edgeI] >= 0)
4194 minPointi =
min(minPointi, edgeMidPoint[edgeI]);
4195 maxPointi =
max(maxPointi, edgeMidPoint[edgeI]);
4199 if (minPointi !=
labelMax && minPointi != mesh_.nPoints())
4202 <<
"Added point labels not consecutive to existing mesh points." 4204 <<
"mesh_.nPoints():" << mesh_.nPoints()
4205 <<
" minPointi:" << minPointi
4206 <<
" maxPointi:" << maxPointi
4211 pointLevel_.
transfer(newPointLevel);
4226 Pout<<
"hexRef8::setRefinement :" 4227 <<
" Updating refinement history to " << cellLevel_.
size()
4228 <<
" cells" <<
endl;
4234 forAll(cellAddedCells, celli)
4236 const labelList& addedCells = cellAddedCells[celli];
4238 if (addedCells.
size())
4252 label celli = cellLabels[i];
4254 refinedCells[i].
transfer(cellAddedCells[celli]);
4257 return refinedCells;
4268 savedPointLevel_.
resize(2*pointsToStore.
size());
4271 label pointi = pointsToStore[i];
4272 savedPointLevel_.
insert(pointi, pointLevel_[pointi]);
4275 savedCellLevel_.
resize(2*cellsToStore.
size());
4278 label celli = cellsToStore[i];
4279 savedCellLevel_.
insert(celli, cellLevel_[celli]);
4291 updateMesh(map, dummyMap, dummyMap, dummyMap);
4309 Pout<<
"hexRef8::updateMesh :" 4310 <<
" Updating various lists" 4319 Pout<<
"hexRef8::updateMesh :" 4322 <<
" nCells:" << mesh_.nCells()
4324 <<
" cellLevel_:" << cellLevel_.
size()
4327 <<
" nPoints:" << mesh_.nPoints()
4329 <<
" pointLevel_:" << pointLevel_.
size()
4333 if (reverseCellMap.
size() == cellLevel_.
size())
4339 reorder(reverseCellMap, mesh_.nCells(), -1, cellLevel_);
4347 forAll(cellMap, newCelli)
4349 label oldCelli = cellMap[newCelli];
4353 newCellLevel[newCelli] = -1;
4357 newCellLevel[newCelli] = cellLevel_[oldCelli];
4367 label newCelli = iter.key();
4368 label storedCelli = iter();
4372 if (fnd == savedCellLevel_.
end())
4375 <<
"Problem : trying to restore old value for new cell " 4376 << newCelli <<
" but cannot find old cell " << storedCelli
4377 <<
" in map of stored values " << savedCellLevel_
4380 cellLevel_[newCelli] = fnd();
4401 if (reversePointMap.
size() == pointLevel_.
size())
4404 reorder(reversePointMap, mesh_.nPoints(), -1, pointLevel_);
4417 if (oldPointi == -1)
4430 newPointLevel[
newPointi] = pointLevel_[oldPointi];
4433 pointLevel_.
transfer(newPointLevel);
4441 label storedPointi = iter();
4445 if (fnd == savedPointLevel_.
end())
4448 <<
"Problem : trying to restore old value for new point " 4449 << newPointi <<
" but cannot find old point " 4451 <<
" in map of stored values " << savedPointLevel_
4483 cellShapesPtr_.clear();
4498 Pout<<
"hexRef8::subset :" 4499 <<
" Updating various lists" 4506 <<
"Subsetting will not work in combination with unrefinement." 4508 <<
"Proceed at your own risk." <<
endl;
4516 forAll(cellMap, newCelli)
4518 newCellLevel[newCelli] = cellLevel_[cellMap[newCelli]];
4527 <<
"cellLevel_ contains illegal value -1 after mapping:" 4542 pointLevel_.
transfer(newPointLevel);
4548 <<
"pointLevel_ contains illegal value -1 after mapping:" 4557 history_.
subset(pointMap, faceMap, cellMap);
4567 cellShapesPtr_.clear();
4576 Pout<<
"hexRef8::distribute :" 4577 <<
" Updating various lists" 4596 cellShapesPtr_.clear();
4602 const scalar smallDim = 1
e-6 * mesh_.bounds().mag();
4606 Pout<<
"hexRef8::checkMesh : Using matching tolerance smallDim:" 4607 << smallDim <<
endl;
4617 labelList nei(mesh_.nFaces()-mesh_.nInternalFaces());
4620 nei[i] = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
4643 label own = mesh_.faceOwner()[facei];
4644 label bFacei = facei-mesh_.nInternalFaces();
4650 <<
"Faces do not seem to be correct across coupled" 4651 <<
" boundaries" <<
endl 4652 <<
"Coupled face " << facei
4653 <<
" between owner " << own
4654 <<
" on patch " << pp.
name()
4655 <<
" and coupled neighbour " << nei[bFacei]
4656 <<
" has two faces connected to it:" 4672 scalarField neiFaceAreas(mesh_.nFaces()-mesh_.nInternalFaces());
4675 neiFaceAreas[i] =
mag(mesh_.faceAreas()[i+mesh_.nInternalFaces()]);
4683 label facei = i+mesh_.nInternalFaces();
4685 const scalar magArea =
mag(mesh_.faceAreas()[facei]);
4687 if (
mag(magArea - neiFaceAreas[i]) > smallDim)
4689 const face& f = mesh_.faces()[facei];
4690 label patchi = mesh_.boundaryMesh().whichPatch(facei);
4692 dumpCell(mesh_.faceOwner()[facei]);
4695 <<
"Faces do not seem to be correct across coupled" 4696 <<
" boundaries" <<
endl 4697 <<
"Coupled face " << facei
4698 <<
" on patch " << patchi
4699 <<
" " << mesh_.boundaryMesh()[
patchi].
name()
4701 <<
" has face area:" << magArea
4702 <<
" (coupled) neighbour face area differs:" 4704 <<
" to within tolerance " << smallDim
4714 labelList nVerts(mesh_.nFaces()-mesh_.nInternalFaces());
4718 nVerts[i] = mesh_.faces()[i+mesh_.nInternalFaces()].size();
4726 label facei = i+mesh_.nInternalFaces();
4728 const face& f = mesh_.faces()[facei];
4730 if (f.
size() != nVerts[i])
4732 dumpCell(mesh_.faceOwner()[facei]);
4734 label patchi = mesh_.boundaryMesh().whichPatch(facei);
4737 <<
"Faces do not seem to be correct across coupled" 4738 <<
" boundaries" <<
endl 4739 <<
"Coupled face " << facei
4740 <<
" on patch " << patchi
4741 <<
" " << mesh_.boundaryMesh()[
patchi].
name()
4743 <<
" has size:" << f.
size()
4744 <<
" (coupled) neighbour face has size:" 4756 pointField anchorPoints(mesh_.nFaces()-mesh_.nInternalFaces());
4760 label facei = i+mesh_.nInternalFaces();
4761 const point& fc = mesh_.faceCentres()[facei];
4762 const face& f = mesh_.faces()[facei];
4763 const vector anchorVec(mesh_.points()[f[0]] - fc);
4765 anchorPoints[i] = anchorVec;
4775 label facei = i+mesh_.nInternalFaces();
4776 const point& fc = mesh_.faceCentres()[facei];
4777 const face& f = mesh_.faces()[facei];
4778 const vector anchorVec(mesh_.points()[f[0]] - fc);
4780 if (
mag(anchorVec - anchorPoints[i]) > smallDim)
4782 dumpCell(mesh_.faceOwner()[facei]);
4784 label patchi = mesh_.boundaryMesh().whichPatch(facei);
4787 <<
"Faces do not seem to be correct across coupled" 4788 <<
" boundaries" <<
endl 4789 <<
"Coupled face " << facei
4790 <<
" on patch " << patchi
4791 <<
" " << mesh_.boundaryMesh()[
patchi].
name()
4793 <<
" has anchor vector:" << anchorVec
4794 <<
" (coupled) neighbour face anchor vector differs:" 4796 <<
" to within tolerance " << smallDim
4804 Pout<<
"hexRef8::checkMesh : Returning" <<
endl;
4811 const label maxPointDiff,
4817 Pout<<
"hexRef8::checkRefinementLevels :" 4818 <<
" Checking 2:1 refinement level" <<
endl;
4823 cellLevel_.
size() != mesh_.nCells()
4824 || pointLevel_.
size() != mesh_.nPoints()
4828 <<
"cellLevel size should be number of cells" 4829 <<
" and pointLevel size should be number of points."<<
nl 4830 <<
"cellLevel:" << cellLevel_.
size()
4831 <<
" mesh.nCells():" << mesh_.nCells() <<
nl 4832 <<
"pointLevel:" << pointLevel_.
size()
4833 <<
" mesh.nPoints():" << mesh_.nPoints()
4843 for (
label facei = 0; facei < mesh_.nInternalFaces(); facei++)
4845 label own = mesh_.faceOwner()[facei];
4846 label nei = mesh_.faceNeighbour()[facei];
4848 if (
mag(cellLevel_[own] - cellLevel_[nei]) > 1)
4854 <<
"Celllevel does not satisfy 2:1 constraint." <<
nl 4855 <<
"On face " << facei <<
" owner cell " << own
4856 <<
" has refinement " << cellLevel_[own]
4857 <<
" neighbour cell " << nei <<
" has refinement " 4864 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
4868 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
4870 neiLevel[i] = cellLevel_[own];
4878 label facei = i+mesh_.nInternalFaces();
4880 label own = mesh_.faceOwner()[facei];
4882 if (
mag(cellLevel_[own] - neiLevel[i]) > 1)
4886 label patchi = mesh_.boundaryMesh().whichPatch(facei);
4889 <<
"Celllevel does not satisfy 2:1 constraint." 4890 <<
" On coupled face " << facei
4891 <<
" on patch " << patchi <<
" " 4892 << mesh_.boundaryMesh()[
patchi].name()
4893 <<
" owner cell " << own <<
" has refinement " 4895 <<
" (coupled) neighbour cell has refinement " 4918 forAll(syncPointLevel, pointi)
4920 if (pointLevel_[pointi] != syncPointLevel[pointi])
4923 <<
"PointLevel is not consistent across coupled patches." 4925 <<
"point:" << pointi <<
" coord:" << mesh_.points()[pointi]
4926 <<
" has level " << pointLevel_[pointi]
4927 <<
" whereas the coupled point has level " 4928 << syncPointLevel[pointi]
4936 if (maxPointDiff != -1)
4939 labelList maxPointLevel(mesh_.nPoints(), 0);
4941 forAll(maxPointLevel, pointi)
4943 const labelList& pCells = mesh_.pointCells(pointi);
4945 label& pLevel = maxPointLevel[pointi];
4949 pLevel =
max(pLevel, cellLevel_[pCells[i]]);
4965 label pointi = pointsToCheck[i];
4967 const labelList& pCells = mesh_.pointCells(pointi);
4971 label celli = pCells[i];
4975 mag(cellLevel_[celli]-maxPointLevel[pointi])
4982 <<
"Too big a difference between" 4983 <<
" point-connected cells." <<
nl 4985 <<
" cellLevel:" << cellLevel_[celli]
4986 <<
" uses point:" << pointi
4987 <<
" coord:" << mesh_.points()[pointi]
4988 <<
" which is also used by a cell with level:" 4989 << maxPointLevel[pointi]
5064 if (cellShapesPtr_.empty())
5068 Pout<<
"hexRef8::cellShapes() : calculating splitHex cellShapes." 5069 <<
" cellLevel:" << cellLevel_.
size()
5070 <<
" pointLevel:" << pointLevel_.
size()
5077 label nSplitHex = 0;
5078 label nUnrecognised = 0;
5080 forAll(cellLevel_, celli)
5082 if (meshShapes[celli].model().index() == 0)
5084 label level = cellLevel_[celli];
5088 bool haveQuads = matchHexShape
5098 cellShapesPtr_()[celli] = degenerateMatcher::match(faces);
5109 Pout<<
"hexRef8::cellShapes() :" 5110 <<
" nCells:" << mesh_.nCells() <<
" of which" <<
nl 5111 <<
" primitive:" << (mesh_.nCells()-nSplitHex-nUnrecognised)
5113 <<
" split-hex:" << nSplitHex <<
nl 5114 <<
" poly :" << nUnrecognised <<
nl 5118 return cellShapesPtr_();
5131 Pout<<
"hexRef8::getSplitPoints :" 5132 <<
" Calculating unrefineable points" <<
endl;
5139 <<
"Only call if constructed with history capability" 5147 labelList splitMaster(mesh_.nPoints(), -1);
5148 labelList splitMasterLevel(mesh_.nPoints(), 0);
5153 for (
label pointi = 0; pointi < mesh_.nPoints(); pointi++)
5155 const labelList& pCells = mesh_.pointCells(pointi);
5157 if (pCells.
size() != 8)
5159 splitMaster[pointi] = -2;
5166 forAll(visibleCells, celli)
5168 const labelList& cPoints = mesh_.cellPoints(celli);
5170 if (visibleCells[celli] != -1 && history_.
parentIndex(celli) >= 0)
5177 label pointi = cPoints[i];
5179 label masterCelli = splitMaster[pointi];
5181 if (masterCelli == -1)
5188 splitMaster[pointi] = parentIndex;
5189 splitMasterLevel[pointi] = cellLevel_[celli] - 1;
5191 else if (masterCelli == -2)
5197 (masterCelli != parentIndex)
5198 || (splitMasterLevel[pointi] != cellLevel_[celli] - 1)
5203 splitMaster[pointi] = -2;
5212 label pointi = cPoints[i];
5214 splitMaster[pointi] = -2;
5222 label facei = mesh_.nInternalFaces();
5223 facei < mesh_.nFaces();
5227 const face& f = mesh_.faces()[facei];
5231 splitMaster[f[fp]] = -2;
5238 label nSplitPoints = 0;
5240 forAll(splitMaster, pointi)
5242 if (splitMaster[pointi] >= 0)
5251 forAll(splitMaster, pointi)
5253 if (splitMaster[pointi] >= 0)
5255 splitPoints[nSplitPoints++] = pointi;
5333 Pout<<
"hexRef8::consistentUnrefinement :" 5334 <<
" Determining 2:1 consistent unrefinement" <<
endl;
5340 <<
"maxSet not implemented yet." 5352 forAll(pointsToUnrefine, i)
5354 label pointi = pointsToUnrefine[i];
5356 unrefinePoint.set(pointi);
5367 forAll(unrefinePoint, pointi)
5369 if (unrefinePoint.get(pointi))
5371 const labelList& pCells = mesh_.pointCells(pointi);
5375 unrefineCell.set(pCells[j]);
5388 for (
label facei = 0; facei < mesh_.nInternalFaces(); facei++)
5390 label own = mesh_.faceOwner()[facei];
5391 label ownLevel = cellLevel_[own] - unrefineCell.get(own);
5393 label nei = mesh_.faceNeighbour()[facei];
5394 label neiLevel = cellLevel_[nei] - unrefineCell.get(nei);
5396 if (ownLevel < (neiLevel-1))
5403 unrefineCell.set(nei);
5414 if (unrefineCell.get(own) == 0)
5420 unrefineCell.unset(own);
5424 else if (neiLevel < (ownLevel-1))
5428 unrefineCell.set(own);
5432 if (unrefineCell.get(nei) == 0)
5438 unrefineCell.unset(nei);
5446 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
5450 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
5452 neiLevel[i] = cellLevel_[own] - unrefineCell.get(own);
5460 label facei = i+mesh_.nInternalFaces();
5461 label own = mesh_.faceOwner()[facei];
5462 label ownLevel = cellLevel_[own] - unrefineCell.get(own);
5464 if (ownLevel < (neiLevel[i]-1))
5468 if (unrefineCell.get(own) == 0)
5474 unrefineCell.unset(own);
5478 else if (neiLevel[i] < (ownLevel-1))
5482 if (unrefineCell.get(own) == 1)
5488 unrefineCell.set(own);
5498 Pout<<
"hexRef8::consistentUnrefinement :" 5499 <<
" Changed " << nChanged
5500 <<
" refinement levels due to 2:1 conflicts." 5514 forAll(unrefinePoint, pointi)
5516 if (unrefinePoint.get(pointi))
5518 const labelList& pCells = mesh_.pointCells(pointi);
5522 if (!unrefineCell.get(pCells[j]))
5524 unrefinePoint.unset(pointi);
5536 forAll(unrefinePoint, pointi)
5538 if (unrefinePoint.get(pointi))
5547 forAll(unrefinePoint, pointi)
5549 if (unrefinePoint.get(pointi))
5551 newPointsToUnrefine[nSet++] = pointi;
5555 return newPointsToUnrefine;
5568 <<
"Only call if constructed with history capability" 5574 Pout<<
"hexRef8::setUnrefinement :" 5575 <<
" Checking initial mesh just to make sure" <<
endl;
5579 forAll(cellLevel_, celli)
5581 if (cellLevel_[celli] < 0)
5584 <<
"Illegal cell level " << cellLevel_[celli]
5585 <<
" for cell " << celli
5592 pointSet pSet(mesh_,
"splitPoints", splitPointLabels);
5595 cellSet cSet(mesh_,
"splitPointCells", splitPointLabels.
size());
5597 forAll(splitPointLabels, i)
5599 const labelList& pCells = mesh_.pointCells(splitPointLabels[i]);
5603 cSet.insert(pCells[j]);
5608 Pout<<
"hexRef8::setRefinement : Dumping " << pSet.
size()
5610 << cSet.size() <<
" cells for unrefinement to" <<
nl 5612 <<
" cellSet " << cSet.objectPath()
5624 forAll(splitPointLabels, i)
5630 splitFaces.insert(pFaces[j]);
5644 if (facesToRemove.
size() != splitFaces.size())
5647 <<
"Ininitial set of split points to unrefine does not" 5648 <<
" seem to be consistent or not mid points of refined cells" 5657 forAll(splitPointLabels, i)
5659 label pointi = splitPointLabels[i];
5663 const labelList& pCells = mesh_.pointCells(pointi);
5666 if (pCells.
size() != 8)
5669 <<
"splitPoint " << pointi
5670 <<
" should have 8 cells using it. It has " << pCells
5683 label celli = pCells[j];
5685 label region = cellRegion[celli];
5690 <<
"Ininitial set of split points to unrefine does not" 5691 <<
" seem to be consistent or not mid points" 5692 <<
" of refined cells" <<
nl 5693 <<
"cell:" << celli <<
" on splitPoint " << pointi
5694 <<
" has no region to be merged into" 5698 if (masterCelli != cellRegionMaster[region])
5701 <<
"cell:" << celli <<
" on splitPoint:" << pointi
5702 <<
" in region " << region
5703 <<
" has master:" << cellRegionMaster[region]
5704 <<
" which is not the lowest numbered cell" 5705 <<
" among the pointCells:" << pCells
5725 forAll(splitPointLabels, i)
5727 label pointi = splitPointLabels[i];
5729 const labelList& pCells = mesh_.pointCells(pointi);
5735 cellLevel_[pCells[j]]--;
5752 cellLevel_.
write(valid)
5753 && pointLevel_.
write(valid)
5754 && level0Edge_.
write(valid);
5758 writeOk = writeOk && history_.
write(valid);
fileCheckTypes
Enumeration defining the file checking options.
void distribute(const mapDistributePolyMesh &)
Update local numbering for mesh redistribution.
virtual const fileName & name() const
Return the name of the stream.
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
#define forAll(list, i)
Loop across all elements in list.
const labelList & reversePointMap() const
Reverse point map.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const word & name() const
Return name.
A class for handling file names.
label nOldCells() const
Number of old cells.
const labelList & cellMap() const
Old cell map.
virtual label iterate(const label maxIter)
Iterate until no changes or maxIter reached. Returns actual.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Class describing modification of a face.
unsigned int get(const label) const
Get value at index I.
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
A face is a list of labels corresponding to mesh vertices.
void distributePointData(List< T > &lst) const
Distribute list of point data.
void setFaceInfo(const labelList &changedFaces, const List< Type > &changedFacesInfo)
Set initial changed faces.
void operator()(label &x, const label y) const
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
void setUnrefinement(const labelList &splitPointLabels, polyTopoChange &)
Remove some refinement. Needs to be supplied output of.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
label rcIndex(const label i) const
Return the reverse circular index, i.e. the previous index.
bool active() const
Is there unrefinement history?
dimensionedSymmTensor sqr(const dimensionedVector &dv)
bool updateFace(const polyMesh &, const label thisFacei, const label neighbourCelli, const refinementData &neighbourInfo, const scalar tol, TrackingData &td)
Influence of neighbouring cell.
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
void size(const label)
Override size to be inconsistent with allocated storage.
const boolList & flipMap() const
Return face flip map.
label fcIndex(const label i) const
Return the forward circular index, i.e. the next index.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Wave propagation of information through grid. Every iteration information goes through one layer of c...
bool headerOk()
Read and check header info.
void checkMesh() const
Debug: Check coupled mesh for correctness.
point centre(const pointField &) const
Centre point of face.
void updateMesh(const mapPolyMesh &)
Force recalculation of locally stored data on topological change.
void storeData(const labelList &pointsToStore, const labelList &facesToStore, const labelList &cellsToStore)
Signal points/face/cells for which to store data.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
labelList consistentSlowRefinement2(const label maxFaceDiff, const labelList &cellsToRefine, const labelList &facesToCheck) const
Like consistentSlowRefinement but uses different meshWave.
Class containing data for point addition.
label whichFace(const label globalCellID) const
Helper function to re-direct to zone::localID(...)
labelListList setRefinement(const labelList &cells, polyTopoChange &)
Insert refinement. All selected cells will be split into 8.
label otherVertex(const label a) const
Given one vertex, return the other.
Holds information regarding type of cell. Used in inside/outside determination in cellClassification...
labelList findIndices(const ListType &, typename ListType::const_reference, const label start=0)
Find all occurrences of given element. Linear search.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
labelList consistentUnrefinement(const labelList &pointsToUnrefine, const bool maxSet) const
Given proposed.
void distributeCellData(List< T > &lst) const
Distribute list of cell data.
label size() const
Return number of elements in table.
scalar level0EdgeLength() const
Typical edge length between unrefined points.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
A face addition data class. A face can be inflated either from a point or from another face and can e...
point centre(const pointField &) const
Return centre (centroid)
void subset(const labelList &pointMap, const labelList &faceMap, const labelList &cellMap)
Update local numbering for subsetted mesh.
label refinementCount() const
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
iterator find(const label &)
Find and return an iterator set at the hashedEntry.
void setInstance(const fileName &inst)
void distribute(const mapDistributePolyMesh &)
Force recalculation of locally stored data for mesh distribution.
A topoSetSource to select a faceSet from cells.
An STL-conforming iterator.
Reduction class. If x and y are not equal assign value.
Class containing data for cell addition.
A list of faces which address into the list of points.
A List obtained as a section of another List.
face reverseFace() const
Return face with reverse direction.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
List< cellShape > cellShapeList
List of cellShapes and PtrList of List of cellShape.
An ordered pair of two objects of type <T> with first() and second() elements.
virtual bool coupled() const
Return true if this patch is geometrically coupled (i.e. faces and.
void combineCells(const label masterCelli, const labelList &combinedCells)
Store combining 8 cells into master.
void set(const PackedList< 1 > &)
Set specified bits.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
labelList consistentRefinement(const labelList &cellsToRefine, const bool maxSet) const
Given valid mesh and current cell level and proposed.
virtual const fileName & name() const
Return the name of the stream.
const DynamicList< point > & points() const
Points. Shrunk after constructing mesh (or calling of compact())
void clear()
Clear all entries from table.
void append(const T &)
Append an element at the end of the list.
void subset(const labelList &pointMap, const labelList &faceMap, const labelList &cellMap)
Update numbering for subsetting.
Xfer< List< T > > xfer()
Transfer contents to the Xfer container as a plain List.
labelList consistentSlowRefinement(const label maxFaceDiff, const labelList &cellsToRefine, const labelList &facesToCheck, const label maxPointDiff, const labelList &pointsToCheck) const
Like consistentRefinement but slower:
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Pair< label > labelPair
Label pair.
static const label labelMax
List< label > labelList
A List of labels.
An STL-conforming hash table.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
errorManip< error > abort(error &err)
void updateMesh(const mapPolyMesh &)
Update local numbering for changed mesh.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
vector area(const pointField &) const
Return vector area.
vector vec(const pointField &) const
Return the vector (end - start)
void reverse(UList< T > &, const label n)
virtual bool read()
Read object. If global number of visible cells > 0 becomes active.
const labelList & reverseCellMap() const
Reverse cell map.
Type gMax(const FieldField< Field, Type > &f)
defineTypeNameAndDebug(combustionModel, 0)
static fileCheckTypes fileModificationChecking
Type of file modification checking.
prefixOSstream Perr(cerr, "Perr")
const labelList & pointMap() const
Old point map.
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
void storeSplit(const label celli, const labelList &addedCells)
Store splitting of cell into 8.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
void updateMesh(const mapPolyMesh &)
Update numbering for mesh changes.
void distribute(const mapDistributePolyMesh &)
Update local numbering for mesh redistribution.
const labelList & visibleCells() const
Per cell in the current mesh (i.e. visible) either -1 (unrefined)
Transfers refinement levels such that slow transition between levels is maintained. Used in FaceCellWave.
void setSize(const label)
Reset size of List.
vector point
Point is a vector.
label faceLevel(const label facei) const
Gets level such that the face has four points <= level.
#define WarningInFunction
Report a warning using Foam::Warning.
bool write(const bool valid=true) const
Force writing refinement+history to polyMesh directory.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
A cell is defined as a list of faces with extra functionality.
A collection of cell labels.
prefixOSstream Pout(cout, "Pout")
label parentIndex(const label celli) const
Get parent of cell.
A List with indirect addressing.
Direct mesh changes based on v1.3 polyTopoChange syntax.
label start() const
Return start label of this patch in the polyMesh face list.
void flip()
Flip the face in-place.
labelList getSplitPoints() const
Return the points at the centre of top-level split cells.
void setRefinement(const labelList &piercedFaces, const labelList &cellRegion, const labelList &cellRegionMaster, polyTopoChange &) const
Play commands into polyTopoChange to remove faces.
label compatibleRemoves(const labelList &inPiercedFaces, labelList &cellRegion, labelList &cellRegionMaster, labelList &outPiercedFaces) const
Find faces including those with cells which have the same mastercell.
void unset(const PackedList< 1 > &)
Unset specified bits.
dimensioned< scalar > mag(const dimensioned< Type > &)
void resize(const label newSize)
Resize the hash table for efficiency.
virtual Ostream & write(const token &)=0
Write next token to stream.
const doubleScalar e
Elementary charge.
Mesh consisting of general polyhedral cells.
void checkRefinementLevels(const label maxPointDiff, const labelList &pointsToCheck) const
Debug: Check 2:1 consistency across faces.
A subset of mesh faces organised as a primitive patch.
All refinement history. Used in unrefinement.
A patch is a list of labels that address the faces in the global face list.
T & last()
Return the last element of the list.
virtual bool write(const bool valid=true) const
Write using setting from DB.
void clear()
Clear the addressed list, i.e. set the size to zero.
readOption readOpt() const
Info<< "Finished reading KIVA file"<< endl;cellShapeList cellShapes(nPoints);labelList cellZoning(nPoints, -1);const cellModel &hex=*(cellModeller::lookup("hex"));labelList hexLabels(8);label activeCells=0;labelList pointMap(nPoints);forAll(pointMap, i){ pointMap[i]=i;}for(label i=0;i< nPoints;i++){ if(f[i] > 0.0) { hexLabels[0]=i;hexLabels[1]=i1tab[i];hexLabels[2]=i3tab[i1tab[i]];hexLabels[3]=i3tab[i];hexLabels[4]=i8tab[i];hexLabels[5]=i1tab[i8tab[i]];hexLabels[6]=i3tab[i1tab[i8tab[i]]];hexLabels[7]=i3tab[i8tab[i]];cellShapes[activeCells]=cellShape(hex, hexLabels);edgeList edges=cellShapes[activeCells].edges();forAll(edges, ei) { if(edges[ei].mag(points)< small) { label start=pointMap[edges[ei].start()];while(start !=pointMap[start]) { start=pointMap[start];} label end=pointMap[edges[ei].end()];while(end !=pointMap[end]) { end=pointMap[end];} label minLabel=min(start, end);pointMap[start]=pointMap[end]=minLabel;} } cellZoning[activeCells]=idreg[i];activeCells++;}}cellShapes.setSize(activeCells);cellZoning.setSize(activeCells);forAll(cellShapes, celli){ cellShape &cs=cellShapes[celli];forAll(cs, i) { cs[i]=pointMap[cs[i]];} cs.collapse();}label bcIDs[11]={-1, 0, 2, 4, -1, 5, -1, 6, 7, 8, 9};const label nBCs=12;const word *kivaPatchTypes[nBCs]={ &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &symmetryPolyPatch::typeName, &wedgePolyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &symmetryPolyPatch::typeName, &oldCyclicPolyPatch::typeName};enum patchTypeNames{ PISTON, VALVE, LINER, CYLINDERHEAD, AXIS, WEDGE, INFLOW, OUTFLOW, PRESIN, PRESOUT, SYMMETRYPLANE, CYCLIC};const char *kivaPatchNames[nBCs]={ "piston", "valve", "liner", "cylinderHead", "axis", "wedge", "inflow", "outflow", "presin", "presout", "symmetryPlane", "cyclic"};List< SLList< face > > pFaces[nBCs]
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
fileName objectPath() const
Return complete path + object name.
void resize(const label nCells)
Extend/shrink storage. additional visibleCells_ elements get.
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
const cellShapeList & cellShapes() const
Utility: get hexes as cell shapes.
static const label labelMin
label nOldPoints() const
Number of old points.
Transfers refinement levels such that slow transition between levels is maintained. Used in FaceCellWave.
label setAction(const topoAction &action)
For compatibility with polyTopoChange: set topological action.