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())
357 "hexRef8::getLevel0EdgeLength() const" 358 ) <<
"Number of cells in mesh:" << mesh_.nCells()
359 <<
" does not equal size of cellLevel:" << cellLevel_.size()
361 <<
"This might be because of a restart with inconsistent cellLevel." 368 const scalar GREAT2 =
sqr(GREAT);
385 const label cLevel = cellLevel_[cellI];
387 const labelList& cEdges = mesh_.cellEdges(cellI);
391 label edgeI = cEdges[i];
393 if (edgeLevel[edgeI] == -1)
395 edgeLevel[edgeI] = cLevel;
397 else if (edgeLevel[edgeI] ==
labelMax)
401 else if (edgeLevel[edgeI] != cLevel)
423 const label eLevel = edgeLevel[edgeI];
425 if (eLevel >= 0 && eLevel <
labelMax)
427 const edge&
e = mesh_.edges()[edgeI];
429 scalar edgeLenSqr =
magSqr(e.
vec(mesh_.points()));
431 typEdgeLenSqr[eLevel] =
min(typEdgeLenSqr[eLevel], edgeLenSqr);
443 Pout<<
"hexRef8::getLevel0EdgeLength() :" 444 <<
" After phase1: Edgelengths (squared) per refinementlevel:" 445 << typEdgeLenSqr <<
endl;
459 const label cLevel = cellLevel_[cellI];
461 const labelList& cEdges = mesh_.cellEdges(cellI);
465 const edge&
e = mesh_.edges()[cEdges[i]];
467 scalar edgeLenSqr =
magSqr(e.
vec(mesh_.points()));
469 maxEdgeLenSqr[cLevel] =
max(maxEdgeLenSqr[cLevel], edgeLenSqr);
478 Pout<<
"hexRef8::getLevel0EdgeLength() :" 479 <<
" Crappy Edgelengths (squared) per refinementlevel:" 480 << maxEdgeLenSqr <<
endl;
487 forAll(typEdgeLenSqr, levelI)
489 if (typEdgeLenSqr[levelI] == GREAT2 && maxEdgeLenSqr[levelI] >= 0)
491 typEdgeLenSqr[levelI] = maxEdgeLenSqr[levelI];
497 Pout<<
"hexRef8::getLevel0EdgeLength() :" 498 <<
" Final Edgelengths (squared) per refinementlevel:" 499 << typEdgeLenSqr <<
endl;
503 scalar level0Size = -1;
505 forAll(typEdgeLenSqr, levelI)
507 scalar lenSqr = typEdgeLenSqr[levelI];
515 Pout<<
"hexRef8::getLevel0EdgeLength() :" 516 <<
" For level:" << levelI
518 <<
" with equivalent level0 len:" << level0Size
525 if (level0Size == -1)
546 if (cellAnchorPoints[cellI].size())
552 return cellAddedCells[cellI][index];
559 const face& f = mesh_.faces()[faceI];
567 return cellAddedCells[cellI][index];
573 Perr<<
"cell:" << cellI <<
" anchorPoints:" << cellAnchorPoints[cellI]
577 <<
"Could not find point " << pointI
578 <<
" in the anchorPoints for cell " << cellI << endl
579 <<
"Does your original mesh obey the 2:1 constraint and" 580 <<
" did you use consistentRefinement to make your cells to refine" 581 <<
" obey this constraint as well?" 594 void Foam::hexRef8::getFaceNeighbours
610 mesh_.faceOwner()[faceI],
615 if (mesh_.isInternalFace(faceI))
621 mesh_.faceNeighbour()[faceI],
641 label level = pointLevel_[f[fp]];
643 if (level < minLevel)
662 label level = pointLevel_[f[fp]];
664 if (level > maxLevel)
678 const label anchorLevel
685 if (pointLevel_[f[fp]] <= anchorLevel)
694 void Foam::hexRef8::dumpCell(
const label cellI)
const 697 Pout<<
"hexRef8 : Dumping cell as obj to " << str.
name() <<
endl;
699 const cell& cFaces = mesh_.cells()[cellI];
706 const face& f = mesh_.faces()[cFaces[i]];
710 if (pointToObjVert.
insert(f[fp], objVertI))
720 const face& f = mesh_.faces()[cFaces[i]];
724 label pointI = f[fp];
727 str <<
"l " << pointToObjVert[pointI]+1
728 <<
' ' << pointToObjVert[nexPointI]+1 <<
nl;
740 const bool searchForward,
741 const label wantedLevel
748 label pointI = f[fp];
750 if (pointLevel_[pointI] < wantedLevel)
752 dumpCell(mesh_.faceOwner()[faceI]);
753 if (mesh_.isInternalFace(faceI))
755 dumpCell(mesh_.faceNeighbour()[faceI]);
761 <<
" startFp:" << startFp
762 <<
" wantedLevel:" << wantedLevel
765 else if (pointLevel_[pointI] == wantedLevel)
780 dumpCell(mesh_.faceOwner()[faceI]);
781 if (mesh_.isInternalFace(faceI))
783 dumpCell(mesh_.faceNeighbour()[faceI]);
789 <<
" startFp:" << startFp
790 <<
" wantedLevel:" << wantedLevel
800 const face& f = mesh_.faces()[faceI];
804 return pointLevel_[f[findMaxLevel(f)]];
808 label ownLevel = cellLevel_[mesh_.faceOwner()[faceI]];
810 if (countAnchors(f, ownLevel) == 4)
814 else if (countAnchors(f, ownLevel+1) == 4)
826 void Foam::hexRef8::checkInternalOrientation
841 vector dir(neiPt - ownPt);
846 <<
"cell:" << cellI <<
" old face:" << faceI
847 <<
" newFace:" << newFace <<
endl 848 <<
" coords:" << compactPoints
849 <<
" ownPt:" << ownPt
850 <<
" neiPt:" << neiPt
854 vector fcToOwn(compactFace.
centre(compactPoints) - ownPt);
856 scalar
s = (fcToOwn&
n) / (dir&n);
858 if (s < 0.1 || s > 0.9)
861 <<
"cell:" << cellI <<
" old face:" << faceI
862 <<
" newFace:" << newFace <<
endl 863 <<
" coords:" << compactPoints
864 <<
" ownPt:" << ownPt
865 <<
" neiPt:" << neiPt
872 void Foam::hexRef8::checkBoundaryOrientation
878 const point& boundaryPt,
887 vector dir(boundaryPt - ownPt);
892 <<
"cell:" << cellI <<
" old face:" << faceI
893 <<
" newFace:" << newFace
894 <<
" coords:" << compactPoints
895 <<
" ownPt:" << ownPt
896 <<
" boundaryPt:" << boundaryPt
900 vector fcToOwn(compactFace.
centre(compactPoints) - ownPt);
902 scalar
s = (fcToOwn&dir) /
magSqr(dir);
904 if (s < 0.7 || s > 1.3)
906 WarningIn(
"checkBoundaryOrientation(..)")
907 <<
"cell:" << cellI <<
" old face:" << faceI
908 <<
" newFace:" << newFace
909 <<
" coords:" << compactPoints
910 <<
" ownPt:" << ownPt
911 <<
" boundaryPt:" << boundaryPt
921 void Foam::hexRef8::insertEdgeSplit
929 if (p0 < mesh_.nPoints() && p1 < mesh_.nPoints())
933 if (edgeI != -1 && edgeMidPoint[edgeI] != -1)
935 verts.
append(edgeMidPoint[edgeI]);
957 const bool faceOrder,
958 const label edgeMidPointI,
959 const label anchorPointI,
960 const label faceMidPointI,
969 bool changed =
false;
970 bool haveTwoAnchors =
false;
974 if (edgeMidFnd == midPointToAnchors.
end())
976 midPointToAnchors.
insert(edgeMidPointI,
edge(anchorPointI, -1));
980 edge&
e = edgeMidFnd();
982 if (anchorPointI != e[0])
991 if (e[0] != -1 && e[1] != -1)
993 haveTwoAnchors =
true;
997 bool haveTwoFaceMids =
false;
1001 if (faceMidFnd == midPointToFaceMids.
end())
1003 midPointToFaceMids.
insert(edgeMidPointI,
edge(faceMidPointI, -1));
1007 edge&
e = faceMidFnd();
1009 if (faceMidPointI != e[0])
1013 e[1] = faceMidPointI;
1018 if (e[0] != -1 && e[1] != -1)
1020 haveTwoFaceMids =
true;
1027 if (changed && haveTwoAnchors && haveTwoFaceMids)
1029 const edge& anchors = midPointToAnchors[edgeMidPointI];
1030 const edge& faceMids = midPointToFaceMids[edgeMidPointI];
1040 if (faceOrder == (mesh_.faceOwner()[faceI] == cellI))
1042 newFaceVerts.
append(faceMidPointI);
1053 newFaceVerts.
append(edgeMidPointI);
1063 newFaceVerts.
append(otherFaceMidPointI);
1064 newFaceVerts.
append(cellMidPoint[cellI]);
1068 newFaceVerts.
append(otherFaceMidPointI);
1078 newFaceVerts.
append(edgeMidPointI);
1088 newFaceVerts.
append(faceMidPointI);
1089 newFaceVerts.
append(cellMidPoint[cellI]);
1095 label anchorCell0 = getAnchorCell
1103 label anchorCell1 = getAnchorCell
1116 if (anchorCell0 < anchorCell1)
1121 ownPt = mesh_.points()[anchorPointI];
1122 neiPt = mesh_.points()[anchors.
otherVertex(anchorPointI)];
1131 ownPt = mesh_.points()[anchors.
otherVertex(anchorPointI)];
1132 neiPt = mesh_.points()[anchorPointI];
1139 if (anchorCell0 < anchorCell1)
1141 ownPt = mesh_.points()[anchorPointI];
1142 neiPt = mesh_.points()[anchors.
otherVertex(anchorPointI)];
1146 ownPt = mesh_.points()[anchors.
otherVertex(anchorPointI)];
1147 neiPt = mesh_.points()[anchorPointI];
1150 checkInternalOrientation
1161 return addInternalFace
1179 void Foam::hexRef8::createInternalFaces
1195 const cell& cFaces = mesh_.cells()[cellI];
1196 const label cLevel = cellLevel_[cellI];
1208 label nFacesAdded = 0;
1212 label faceI = cFaces[i];
1214 const face& f = mesh_.faces()[faceI];
1215 const labelList& fEdges = mesh_.faceEdges(faceI, storage);
1220 label faceMidPointI = -1;
1222 label nAnchors = countAnchors(f, cLevel);
1230 label anchorFp = -1;
1234 if (pointLevel_[f[fp]] <= cLevel)
1242 label edgeMid = findLevel
1250 label faceMid = findLevel
1259 faceMidPointI = f[faceMid];
1261 else if (nAnchors == 4)
1266 faceMidPointI = faceMidPoint[faceI];
1270 dumpCell(mesh_.faceOwner()[faceI]);
1271 if (mesh_.isInternalFace(faceI))
1273 dumpCell(mesh_.faceNeighbour()[faceI]);
1277 <<
"nAnchors:" << nAnchors
1278 <<
" faceI:" << faceI
1290 label point0 = f[fp0];
1292 if (pointLevel_[point0] <= cLevel)
1301 label edgeMidPointI = -1;
1305 if (pointLevel_[f[fp1]] <= cLevel)
1308 label edgeI = fEdges[fp0];
1310 edgeMidPointI = edgeMidPoint[edgeI];
1312 if (edgeMidPointI == -1)
1316 const labelList& cPoints = mesh_.cellPoints(cellI);
1319 <<
"cell:" << cellI <<
" cLevel:" << cLevel
1320 <<
" cell points:" << cPoints
1323 <<
" face:" << faceI
1327 <<
" faceAnchorLevel:" << faceAnchorLevel[faceI]
1328 <<
" faceMidPoint:" << faceMidPoint[faceI]
1329 <<
" faceMidPointI:" << faceMidPointI
1337 label edgeMid = findLevel(faceI, f, fp1,
true, cLevel+1);
1339 edgeMidPointI = f[edgeMid];
1342 label newFaceI = storeMidPointInfo
1365 if (nFacesAdded == 12)
1378 if (pointLevel_[f[fpMin1]] <= cLevel)
1381 label edgeI = fEdges[fpMin1];
1383 edgeMidPointI = edgeMidPoint[edgeI];
1385 if (edgeMidPointI == -1)
1389 const labelList& cPoints = mesh_.cellPoints(cellI);
1392 <<
"cell:" << cellI <<
" cLevel:" << cLevel
1393 <<
" cell points:" << cPoints
1396 <<
" face:" << faceI
1400 <<
" faceAnchorLevel:" << faceAnchorLevel[faceI]
1401 <<
" faceMidPoint:" << faceMidPoint[faceI]
1402 <<
" faceMidPointI:" << faceMidPointI
1410 label edgeMid = findLevel
1419 edgeMidPointI = f[edgeMid];
1422 newFaceI = storeMidPointInfo
1445 if (nFacesAdded == 12)
1453 if (nFacesAdded == 12)
1461 void Foam::hexRef8::walkFaceToMid
1466 const label startFp,
1470 const face& f = mesh_.faces()[faceI];
1471 const labelList& fEdges = mesh_.faceEdges(faceI);
1480 if (edgeMidPoint[fEdges[fp]] >= 0)
1482 faceVerts.
append(edgeMidPoint[fEdges[fp]]);
1487 if (pointLevel_[f[fp]] <= cLevel)
1493 else if (pointLevel_[f[fp]] == cLevel+1)
1500 else if (pointLevel_[f[fp]] == cLevel+2)
1510 void Foam::hexRef8::walkFaceFromMid
1515 const label startFp,
1519 const face& f = mesh_.faces()[faceI];
1520 const labelList& fEdges = mesh_.faceEdges(faceI);
1526 if (pointLevel_[f[fp]] <= cLevel)
1531 else if (pointLevel_[f[fp]] == cLevel+1)
1537 else if (pointLevel_[f[fp]] == cLevel+2)
1547 if (edgeMidPoint[fEdges[fp]] >= 0)
1549 faceVerts.
append(edgeMidPoint[fEdges[fp]]);
1565 Foam::label Foam::hexRef8::faceConsistentRefinement
1574 for (
label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
1576 label own = mesh_.faceOwner()[faceI];
1577 label ownLevel = cellLevel_[own] + refineCell.
get(own);
1579 label nei = mesh_.faceNeighbour()[faceI];
1580 label neiLevel = cellLevel_[nei] + refineCell.
get(nei);
1582 if (ownLevel > (neiLevel+1))
1586 refineCell.
set(nei);
1590 refineCell.
unset(own);
1594 else if (neiLevel > (ownLevel+1))
1598 refineCell.
set(own);
1602 refineCell.
unset(nei);
1611 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
1615 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
1617 neiLevel[i] = cellLevel_[own] + refineCell.
get(own);
1626 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
1627 label ownLevel = cellLevel_[own] + refineCell.
get(own);
1629 if (ownLevel > (neiLevel[i]+1))
1633 refineCell.
unset(own);
1637 else if (neiLevel[i] > (ownLevel+1))
1641 refineCell.
set(own);
1652 void Foam::hexRef8::checkWantedRefinementLevels
1660 refineCell.
set(cellsToRefine[i]);
1663 for (
label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
1665 label own = mesh_.faceOwner()[faceI];
1666 label ownLevel = cellLevel_[own] + refineCell.
get(own);
1668 label nei = mesh_.faceNeighbour()[faceI];
1669 label neiLevel = cellLevel_[nei] + refineCell.
get(nei);
1671 if (
mag(ownLevel-neiLevel) > 1)
1677 "hexRef8::checkWantedRefinementLevels(const labelList&)" 1679 <<
" current level:" << cellLevel_[own]
1680 <<
" level after refinement:" << ownLevel
1682 <<
"neighbour cell:" << nei
1683 <<
" current level:" << cellLevel_[nei]
1684 <<
" level after refinement:" << neiLevel
1686 <<
"which does not satisfy 2:1 constraints anymore." 1693 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
1697 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
1699 neiLevel[i] = cellLevel_[own] + refineCell.
get(own);
1708 label faceI = i + mesh_.nInternalFaces();
1710 label own = mesh_.faceOwner()[faceI];
1711 label ownLevel = cellLevel_[own] + refineCell.
get(own);
1713 if (
mag(ownLevel - neiLevel[i]) > 1)
1715 label patchI = mesh_.boundaryMesh().whichPatch(faceI);
1720 "hexRef8::checkWantedRefinementLevels(const labelList&)" 1721 ) <<
"Celllevel does not satisfy 2:1 constraint." 1722 <<
" On coupled face " 1724 <<
" on patch " << patchI <<
" " 1725 << mesh_.boundaryMesh()[patchI].name()
1726 <<
" owner cell " << own
1727 <<
" current level:" << cellLevel_[own]
1728 <<
" level after refinement:" << ownLevel
1730 <<
" (coupled) neighbour cell will get refinement " 1743 Pout<<
"hexRef8::setInstance(const fileName& inst) : " 1744 <<
"Resetting file instance to " << inst <<
endl;
1747 cellLevel_.instance() = inst;
1748 pointLevel_.instance() = inst;
1749 level0Edge_.instance() = inst;
1750 history_.instance() = inst;
1754 void Foam::hexRef8::collectLevelPoints
1763 if (pointLevel_[f[fp]] <= level)
1771 void Foam::hexRef8::collectLevelPoints
1781 label pointI = meshPoints[f[fp]];
1782 if (pointLevel_[pointI] <= level)
1791 bool Foam::hexRef8::matchHexShape
1794 const label cellLevel,
1798 const cell& cFaces = mesh_.cells()[cellI];
1809 label faceI = cFaces[i];
1810 const face& f = mesh_.faces()[faceI];
1813 collectLevelPoints(f, cellLevel, verts);
1814 if (verts.
size() == 4)
1816 if (mesh_.faceOwner()[faceI] != cellI)
1827 if (quads.
size() < 6)
1833 label faceI = cFaces[i];
1834 const face& f = mesh_.faces()[faceI];
1842 collectLevelPoints(f, cellLevel, verts);
1843 if (verts.
size() == 1)
1849 label pointI = f[fp];
1850 if (pointLevel_[pointI] == cellLevel+1)
1853 pointFaces.find(pointI);
1854 if (iter != pointFaces.end())
1880 if (pFaces.
size() == 4)
1886 label faceI = pFaces[pFaceI];
1887 const face& f = mesh_.faces()[faceI];
1888 if (mesh_.faceOwner()[faceI] == cellI)
1890 fourFaces[pFaceI] =
f;
1905 if (edgeLoops.size() == 1)
1911 bigFace.meshPoints(),
1912 bigFace.edgeLoops()[0],
1917 if (verts.
size() == 4)
1928 return (quads.
size() == 6);
1935 Foam::hexRef8::hexRef8(
const polyMesh& mesh,
const bool readHistory)
1943 mesh_.facesInstance(),
1956 mesh_.facesInstance(),
1969 mesh_.facesInstance(),
1981 "refinementHistory",
1982 mesh_.facesInstance(),
1988 (readHistory ? mesh_.nCells() : 0)
1990 faceRemover_(mesh_, GREAT),
1991 savedPointLevel_(0),
2013 "hexRef8::hexRef8(const polyMesh&)" 2014 ) <<
"History enabled but number of visible cells " 2017 <<
" is not equal to the number of cells in the mesh " 2024 cellLevel_.
size() != mesh_.nCells()
2025 || pointLevel_.
size() != mesh_.nPoints()
2030 "hexRef8::hexRef8(const polyMesh&)" 2031 ) <<
"Restarted from inconsistent cellLevel or pointLevel files." 2035 <<
"Number of cells in mesh:" << mesh_.nCells()
2036 <<
" does not equal size of cellLevel:" << cellLevel_.
size() <<
endl 2037 <<
"Number of points in mesh:" << mesh_.nPoints()
2038 <<
" does not equal size of pointLevel:" << pointLevel_.
size()
2057 Foam::hexRef8::hexRef8
2063 const scalar level0Edge
2072 mesh_.facesInstance(),
2085 mesh_.facesInstance(),
2098 mesh_.facesInstance(),
2108 (level0Edge >= 0 ? level0Edge : getLevel0EdgeLength())
2115 "refinementHistory",
2116 mesh_.facesInstance(),
2124 faceRemover_(mesh_, GREAT),
2125 savedPointLevel_(0),
2132 "hexRef8::hexRef8(const polyMesh&, const labelList&" 2133 ", const labelList&, const refinementHistory&)" 2134 ) <<
"History enabled but number of visible cells in it " 2136 <<
" is not equal to the number of cells in the mesh " 2142 cellLevel_.
size() != mesh_.nCells()
2143 || pointLevel_.
size() != mesh_.nPoints()
2148 "hexRef8::hexRef8(const polyMesh&, const labelList&" 2149 ", const labelList&, const refinementHistory&)" 2150 ) <<
"Incorrect cellLevel or pointLevel size." <<
endl 2151 <<
"Number of cells in mesh:" << mesh_.nCells()
2152 <<
" does not equal size of cellLevel:" << cellLevel_.
size() <<
endl 2153 <<
"Number of points in mesh:" << mesh_.nPoints()
2154 <<
" does not equal size of pointLevel:" << pointLevel_.
size()
2172 Foam::hexRef8::hexRef8
2177 const scalar level0Edge
2186 mesh_.facesInstance(),
2199 mesh_.facesInstance(),
2212 mesh_.facesInstance(),
2222 (level0Edge >= 0 ? level0Edge : getLevel0EdgeLength())
2229 "refinementHistory",
2230 mesh_.facesInstance(),
2239 faceRemover_(mesh_, GREAT),
2240 savedPointLevel_(0),
2245 cellLevel_.
size() != mesh_.nCells()
2246 || pointLevel_.
size() != mesh_.nPoints()
2251 "hexRef8::hexRef8(const polyMesh&, const labelList&" 2252 ", const labelList&)" 2253 ) <<
"Incorrect cellLevel or pointLevel size." <<
endl 2254 <<
"Number of cells in mesh:" << mesh_.nCells()
2255 <<
" does not equal size of cellLevel:" << cellLevel_.
size() <<
endl 2256 <<
"Number of points in mesh:" << mesh_.nPoints()
2257 <<
" does not equal size of pointLevel:" << pointLevel_.
size()
2290 refineCell.
set(cellsToRefine[i]);
2295 label nChanged = faceConsistentRefinement(maxSet, refineCell);
2301 Pout<<
"hexRef8::consistentRefinement : Changed " << nChanged
2302 <<
" refinement levels due to 2:1 conflicts." 2316 forAll(refineCell, cellI)
2318 if (refineCell.
get(cellI))
2327 forAll(refineCell, cellI)
2329 if (refineCell.
get(cellI))
2331 newCellsToRefine[nRefined++] = cellI;
2337 checkWantedRefinementLevels(newCellsToRefine);
2340 return newCellsToRefine;
2352 const label maxFaceDiff,
2355 const label maxPointDiff,
2359 const labelList& faceOwner = mesh_.faceOwner();
2360 const labelList& faceNeighbour = mesh_.faceNeighbour();
2363 if (maxFaceDiff <= 0)
2367 "hexRef8::consistentSlowRefinement" 2368 "(const label, const labelList&, const labelList&" 2369 ", const label, const labelList&)" 2370 ) <<
"Illegal maxFaceDiff " << maxFaceDiff <<
nl 2389 forAll(allCellInfo, cellI)
2395 maxFaceDiff*(cellLevel_[cellI]+1),
2396 maxFaceDiff*cellLevel_[cellI]
2403 label cellI = cellsToRefine[i];
2405 allCellInfo[cellI].count() = allCellInfo[cellI].refinementCount();
2415 int dummyTrackData = 0;
2423 label faceI = facesToCheck[i];
2425 if (allFaceInfo[faceI].valid(dummyTrackData))
2430 "hexRef8::consistentSlowRefinement" 2431 "(const label, const labelList&, const labelList&" 2432 ", const label, const labelList&)" 2433 ) <<
"Argument facesToCheck seems to have duplicate entries!" 2435 <<
"face:" << faceI <<
" occurs at positions " 2443 if (mesh_.isInternalFace(faceI))
2448 const refinementData& neiData = allCellInfo[faceNeighbour[faceI]];
2451 label faceRefineCount;
2454 faceCount = neiData.
count() + maxFaceDiff;
2455 faceRefineCount = faceCount + maxFaceDiff;
2459 faceCount = ownData.
count() + maxFaceDiff;
2460 faceRefineCount = faceCount + maxFaceDiff;
2463 seedFaces.append(faceI);
2464 seedFacesInfo.append
2472 allFaceInfo[faceI] = seedFacesInfo.last();
2476 label faceCount = ownData.
count() + maxFaceDiff;
2477 label faceRefineCount = faceCount + maxFaceDiff;
2479 seedFaces.append(faceI);
2480 seedFacesInfo.append
2488 allFaceInfo[faceI] = seedFacesInfo.last();
2496 forAll(faceNeighbour, faceI)
2499 if (!allFaceInfo[faceI].valid(dummyTrackData))
2501 label own = faceOwner[faceI];
2502 label nei = faceNeighbour[faceI];
2506 if (allCellInfo[own].count() > allCellInfo[nei].count())
2508 allFaceInfo[faceI].updateFace
2518 seedFacesInfo.append(allFaceInfo[faceI]);
2520 else if (allCellInfo[own].count() < allCellInfo[nei].count())
2522 allFaceInfo[faceI].updateFace
2531 seedFaces.append(faceI);
2532 seedFacesInfo.append(allFaceInfo[faceI]);
2540 for (
label faceI = mesh_.nInternalFaces(); faceI < mesh_.nFaces(); faceI++)
2543 if (!allFaceInfo[faceI].valid(dummyTrackData))
2545 label own = faceOwner[faceI];
2558 seedFaces.append(faceI);
2559 seedFacesInfo.append(faceData);
2577 Pout<<
"hexRef8::consistentSlowRefinement : Seeded " 2578 << seedFaces.size() <<
" faces between cells with different" 2579 <<
" refinement level." <<
endl;
2583 levelCalc.
setFaceInfo(seedFaces.shrink(), seedFacesInfo.shrink());
2585 seedFacesInfo.clear();
2588 levelCalc.
iterate(mesh_.globalData().nTotalFaces()+1);
2596 if (maxPointDiff == -1)
2604 labelList maxPointCount(mesh_.nPoints(), 0);
2606 forAll(maxPointCount, pointI)
2608 label& pLevel = maxPointCount[pointI];
2610 const labelList& pCells = mesh_.pointCells(pointI);
2614 pLevel =
max(pLevel, allCellInfo[pCells[i]].count());
2635 label pointI = pointsToCheck[i];
2640 const labelList& pCells = mesh_.pointCells(pointI);
2644 label cellI = pCells[pCellI];
2652 maxPointCount[pointI]
2653 > cellInfo.
count() + maxFaceDiff*maxPointDiff
2661 const cell& cFaces = mesh_.cells()[cellI];
2665 label faceI = cFaces[cFaceI];
2678 if (faceData.
count() > allFaceInfo[faceI].count())
2680 changedFacesInfo.insert(faceI, faceData);
2687 label nChanged = changedFacesInfo.size();
2697 seedFaces.setCapacity(changedFacesInfo.size());
2698 seedFacesInfo.setCapacity(changedFacesInfo.size());
2702 seedFaces.append(iter.key());
2703 seedFacesInfo.append(iter());
2710 for (
label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
2712 label own = mesh_.faceOwner()[faceI];
2715 + (allCellInfo[own].isRefined() ? 1 : 0);
2717 label nei = mesh_.faceNeighbour()[faceI];
2720 + (allCellInfo[nei].isRefined() ? 1 : 0);
2722 if (
mag(ownLevel-neiLevel) > 1)
2728 "hexRef8::consistentSlowRefinement" 2729 "(const label, const labelList&, const labelList&" 2730 ", const label, const labelList&)" 2732 <<
" current level:" << cellLevel_[own]
2733 <<
" current refData:" << allCellInfo[own]
2734 <<
" level after refinement:" << ownLevel
2736 <<
"neighbour cell:" << nei
2737 <<
" current level:" << cellLevel_[nei]
2738 <<
" current refData:" << allCellInfo[nei]
2739 <<
" level after refinement:" << neiLevel
2741 <<
"which does not satisfy 2:1 constraints anymore." <<
nl 2742 <<
"face:" << faceI <<
" faceRefData:" << allFaceInfo[faceI]
2751 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
2752 labelList neiCount(mesh_.nFaces()-mesh_.nInternalFaces());
2753 labelList neiRefCount(mesh_.nFaces()-mesh_.nInternalFaces());
2757 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
2758 neiLevel[i] = cellLevel_[own];
2759 neiCount[i] = allCellInfo[own].count();
2760 neiRefCount[i] = allCellInfo[own].refinementCount();
2771 label faceI = i+mesh_.nInternalFaces();
2773 label own = mesh_.faceOwner()[faceI];
2776 + (allCellInfo[own].isRefined() ? 1 : 0);
2780 + ((neiCount[i] >= neiRefCount[i]) ? 1 : 0);
2782 if (
mag(ownLevel - nbrLevel) > 1)
2785 label patchI = mesh_.boundaryMesh().whichPatch(faceI);
2789 "hexRef8::consistentSlowRefinement" 2790 "(const label, const labelList&, const labelList&" 2791 ", const label, const labelList&)" 2792 ) <<
"Celllevel does not satisfy 2:1 constraint." 2793 <<
" On coupled face " 2795 <<
" refData:" << allFaceInfo[faceI]
2796 <<
" on patch " << patchI <<
" " 2797 << mesh_.boundaryMesh()[patchI].name() <<
nl 2798 <<
"owner cell " << own
2799 <<
" current level:" << cellLevel_[own]
2800 <<
" current count:" << allCellInfo[own].count()
2801 <<
" current refCount:" 2802 << allCellInfo[own].refinementCount()
2803 <<
" level after refinement:" << ownLevel
2805 <<
"(coupled) neighbour cell" 2806 <<
" has current level:" << neiLevel[i]
2807 <<
" current count:" << neiCount[i]
2808 <<
" current refCount:" << neiRefCount[i]
2809 <<
" level after refinement:" << nbrLevel
2819 forAll(allCellInfo, cellI)
2821 if (allCellInfo[cellI].isRefined())
2831 forAll(allCellInfo, cellI)
2833 if (allCellInfo[cellI].isRefined())
2835 newCellsToRefine[nRefined++] = cellI;
2841 Pout<<
"hexRef8::consistentSlowRefinement : From " 2842 << cellsToRefine.
size() <<
" to " << newCellsToRefine.
size()
2843 <<
" cells to refine." <<
endl;
2846 return newCellsToRefine;
2852 const label maxFaceDiff,
2857 const labelList& faceOwner = mesh_.faceOwner();
2858 const labelList& faceNeighbour = mesh_.faceNeighbour();
2860 if (maxFaceDiff <= 0)
2864 "hexRef8::consistentSlowRefinement2" 2865 "(const label, const labelList&, const labelList&)" 2866 ) <<
"Illegal maxFaceDiff " << maxFaceDiff <<
nl 2885 int dummyTrackData = 0;
2891 label cellI = cellsToRefine[i];
2896 mesh_.cellCentres()[cellI],
2901 forAll(allCellInfo, cellI)
2903 if (!allCellInfo[cellI].valid(dummyTrackData))
2908 mesh_.cellCentres()[cellI],
2924 label faceI = facesToCheck[i];
2926 if (allFaceInfo[faceI].valid(dummyTrackData))
2931 "hexRef8::consistentSlowRefinement2" 2932 "(const label, const labelList&, const labelList&)" 2933 ) <<
"Argument facesToCheck seems to have duplicate entries!" 2935 <<
"face:" << faceI <<
" occurs at positions " 2940 label own = faceOwner[faceI];
2944 allCellInfo[own].valid(dummyTrackData)
2945 ? allCellInfo[own].originLevel()
2949 if (!mesh_.isInternalFace(faceI))
2953 const point& fc = mesh_.faceCentres()[faceI];
2962 allFaceInfo[faceI].updateFace
2974 label nei = faceNeighbour[faceI];
2978 allCellInfo[nei].valid(dummyTrackData)
2979 ? allCellInfo[nei].originLevel()
2983 if (ownLevel == neiLevel)
2986 allFaceInfo[faceI].updateFace
2995 allFaceInfo[faceI].updateFace
3008 allFaceInfo[faceI].updateFace
3017 allFaceInfo[faceI].updateFace
3029 seedFacesInfo.append(allFaceInfo[faceI]);
3036 forAll(faceNeighbour, faceI)
3039 if (!allFaceInfo[faceI].valid(dummyTrackData))
3041 label own = faceOwner[faceI];
3045 allCellInfo[own].valid(dummyTrackData)
3046 ? allCellInfo[own].originLevel()
3050 label nei = faceNeighbour[faceI];
3054 allCellInfo[nei].valid(dummyTrackData)
3055 ? allCellInfo[nei].originLevel()
3059 if (ownLevel > neiLevel)
3063 allFaceInfo[faceI].updateFace
3072 seedFacesInfo.append(allFaceInfo[faceI]);
3074 else if (neiLevel > ownLevel)
3076 seedFaces.append(faceI);
3077 allFaceInfo[faceI].updateFace
3086 seedFacesInfo.append(allFaceInfo[faceI]);
3092 seedFacesInfo.shrink();
3102 mesh_.globalData().nTotalCells()+1,
3143 label cellI = cellsToRefine[i];
3145 allCellInfo[cellI].originLevel() =
sizeof(
label)*8-2;
3146 allCellInfo[cellI].origin() = cc[cellI];
3153 forAll(allCellInfo, cellI)
3155 label wanted = allCellInfo[cellI].wantedLevel(cc[cellI]);
3157 if (wanted > cellLevel_[cellI]+1)
3159 refineCell.
set(cellI);
3162 faceConsistentRefinement(
true, refineCell);
3166 label nChanged = faceConsistentRefinement(
true, refineCell);
3172 Pout<<
"hexRef8::consistentSlowRefinement2 : Changed " << nChanged
3173 <<
" refinement levels due to 2:1 conflicts." 3186 forAll(refineCell, cellI)
3189 if (refineCell[cellI])
3198 forAll(refineCell, cellI)
3201 if (refineCell[cellI])
3203 newCellsToRefine[nRefined++] = cellI;
3209 Pout<<
"hexRef8::consistentSlowRefinement2 : From " 3210 << cellsToRefine.
size() <<
" to " << newCellsToRefine.
size()
3211 <<
" cells to refine." <<
endl;
3216 cellSet cellsIn(mesh_,
"cellsToRefineIn", cellsToRefine);
3217 Pout<<
"hexRef8::consistentSlowRefinement2 : writing " 3218 << cellsIn.
size() <<
" to cellSet " 3223 cellSet cellsOut(mesh_,
"cellsToRefineOut", newCellsToRefine);
3224 Pout<<
"hexRef8::consistentSlowRefinement2 : writing " 3225 << cellsOut.
size() <<
" to cellSet " 3232 forAll(newCellsToRefine, i)
3234 refineCell.
set(newCellsToRefine[i]);
3238 label nChanged = faceConsistentRefinement(
true, refineCell);
3243 mesh_,
"cellsToRefineOut2", newCellsToRefine.
size()
3245 forAll(refineCell, cellI)
3247 if (refineCell.
get(cellI))
3249 cellsOut2.insert(cellI);
3252 Pout<<
"hexRef8::consistentSlowRefinement2 : writing " 3253 << cellsOut2.size() <<
" to cellSet " 3254 << cellsOut2.objectPath() <<
endl;
3260 forAll(refineCell, cellI)
3262 if (refineCell.
get(cellI) && !savedRefineCell.
get(cellI))
3267 "hexRef8::consistentSlowRefinement2" 3268 "(const label, const labelList&, const labelList&)" 3269 ) <<
"Cell:" << cellI <<
" cc:" 3270 << mesh_.cellCentres()[cellI]
3271 <<
" was not marked for refinement but does not obey" 3272 <<
" 2:1 constraints." 3279 return newCellsToRefine;
3292 Pout<<
"hexRef8::setRefinement :" 3293 <<
" Checking initial mesh just to make sure" <<
endl;
3302 savedPointLevel_.
clear();
3303 savedCellLevel_.
clear();
3308 forAll(cellLevel_, cellI)
3310 newCellLevel.append(cellLevel_[cellI]);
3313 forAll(pointLevel_, pointI)
3315 newPointLevel.append(pointLevel_[pointI]);
3321 Pout<<
"hexRef8::setRefinement :" 3322 <<
" Allocating " << cellLabels.
size() <<
" cell midpoints." 3330 labelList cellMidPoint(mesh_.nCells(), -1);
3334 label cellI = cellLabels[i];
3336 label anchorPointI = mesh_.faces()[mesh_.cells()[cellI][0]][0];
3342 mesh_.cellCentres()[cellI],
3349 newPointLevel(cellMidPoint[cellI]) = cellLevel_[cellI]+1;
3355 cellSet splitCells(mesh_,
"splitCells", cellLabels.
size());
3357 forAll(cellMidPoint, cellI)
3359 if (cellMidPoint[cellI] >= 0)
3361 splitCells.insert(cellI);
3365 Pout<<
"hexRef8::setRefinement : Dumping " << splitCells.size()
3366 <<
" cells to split to cellSet " << splitCells.objectPath()
3379 Pout<<
"hexRef8::setRefinement :" 3380 <<
" Allocating edge midpoints." 3389 labelList edgeMidPoint(mesh_.nEdges(), -1);
3392 forAll(cellMidPoint, cellI)
3394 if (cellMidPoint[cellI] >= 0)
3396 const labelList& cEdges = mesh_.cellEdges(cellI);
3400 label edgeI = cEdges[i];
3402 const edge&
e = mesh_.edges()[edgeI];
3406 pointLevel_[e[0]] <= cellLevel_[cellI]
3407 && pointLevel_[e[1]] <= cellLevel_[cellI]
3410 edgeMidPoint[edgeI] = 12345;
3437 forAll(edgeMidPoint, edgeI)
3439 if (edgeMidPoint[edgeI] >= 0)
3442 edgeMids[edgeI] = mesh_.edges()[edgeI].centre(mesh_.points());
3450 point(-GREAT, -GREAT, -GREAT)
3455 forAll(edgeMidPoint, edgeI)
3457 if (edgeMidPoint[edgeI] >= 0)
3462 const edge&
e = mesh_.edges()[edgeI];
3475 newPointLevel(edgeMidPoint[edgeI]) =
3488 OFstream str(mesh_.time().path()/
"edgeMidPoint.obj");
3490 forAll(edgeMidPoint, edgeI)
3492 if (edgeMidPoint[edgeI] >= 0)
3494 const edge&
e = mesh_.edges()[edgeI];
3500 Pout<<
"hexRef8::setRefinement :" 3501 <<
" Dumping edge centres to split to file " << str.
name() <<
endl;
3511 Pout<<
"hexRef8::setRefinement :" 3512 <<
" Allocating face midpoints." 3518 labelList faceAnchorLevel(mesh_.nFaces());
3520 for (
label faceI = 0; faceI < mesh_.nFaces(); faceI++)
3527 labelList faceMidPoint(mesh_.nFaces(), -1);
3532 for (
label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
3534 if (faceAnchorLevel[faceI] >= 0)
3536 label own = mesh_.faceOwner()[faceI];
3537 label ownLevel = cellLevel_[own];
3538 label newOwnLevel = ownLevel + (cellMidPoint[own] >= 0 ? 1 : 0);
3540 label nei = mesh_.faceNeighbour()[faceI];
3541 label neiLevel = cellLevel_[nei];
3542 label newNeiLevel = neiLevel + (cellMidPoint[nei] >= 0 ? 1 : 0);
3546 newOwnLevel > faceAnchorLevel[faceI]
3547 || newNeiLevel > faceAnchorLevel[faceI]
3550 faceMidPoint[faceI] = 12345;
3563 labelList newNeiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
3567 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
3568 label ownLevel = cellLevel_[own];
3569 label newOwnLevel = ownLevel + (cellMidPoint[own] >= 0 ? 1 : 0);
3571 newNeiLevel[i] = newOwnLevel;
3581 label faceI = i+mesh_.nInternalFaces();
3583 if (faceAnchorLevel[faceI] >= 0)
3585 label own = mesh_.faceOwner()[faceI];
3586 label ownLevel = cellLevel_[own];
3587 label newOwnLevel = ownLevel + (cellMidPoint[own] >= 0 ? 1 : 0);
3591 newOwnLevel > faceAnchorLevel[faceI]
3592 || newNeiLevel[i] > faceAnchorLevel[faceI]
3595 faceMidPoint[faceI] = 12345;
3620 mesh_.nFaces()-mesh_.nInternalFaces(),
3621 point(-GREAT, -GREAT, -GREAT)
3626 label faceI = i+mesh_.nInternalFaces();
3628 if (faceMidPoint[faceI] >= 0)
3630 bFaceMids[i] = mesh_.faceCentres()[faceI];
3640 forAll(faceMidPoint, faceI)
3642 if (faceMidPoint[faceI] >= 0)
3647 const face& f = mesh_.faces()[faceI];
3654 faceI < mesh_.nInternalFaces()
3655 ? mesh_.faceCentres()[faceI]
3656 : bFaceMids[faceI-mesh_.nInternalFaces()]
3666 newPointLevel(faceMidPoint[faceI]) = faceAnchorLevel[faceI]+1;
3673 faceSet splitFaces(mesh_,
"splitFaces", cellLabels.
size());
3675 forAll(faceMidPoint, faceI)
3677 if (faceMidPoint[faceI] >= 0)
3679 splitFaces.insert(faceI);
3683 Pout<<
"hexRef8::setRefinement : Dumping " << splitFaces.size()
3684 <<
" faces to split to faceSet " << splitFaces.objectPath() <<
endl;
3705 Pout<<
"hexRef8::setRefinement :" 3706 <<
" Finding cell anchorPoints (8 per cell)" 3717 labelList nAnchorPoints(mesh_.nCells(), 0);
3719 forAll(cellMidPoint, cellI)
3721 if (cellMidPoint[cellI] >= 0)
3723 cellAnchorPoints[cellI].
setSize(8);
3727 forAll(pointLevel_, pointI)
3729 const labelList& pCells = mesh_.pointCells(pointI);
3733 label cellI = pCells[pCellI];
3737 cellMidPoint[cellI] >= 0
3738 && pointLevel_[pointI] <= cellLevel_[cellI]
3741 if (nAnchorPoints[cellI] == 8)
3746 "hexRef8::setRefinement(const labelList&" 3747 ", polyTopoChange&)" 3748 ) <<
"cell " << cellI
3749 <<
" of level " << cellLevel_[cellI]
3750 <<
" uses more than 8 points of equal or" 3751 <<
" lower level" <<
nl 3752 <<
"Points so far:" << cellAnchorPoints[cellI]
3756 cellAnchorPoints[cellI][nAnchorPoints[cellI]++] = pointI;
3761 forAll(cellMidPoint, cellI)
3763 if (cellMidPoint[cellI] >= 0)
3765 if (nAnchorPoints[cellI] != 8)
3769 const labelList& cPoints = mesh_.cellPoints(cellI);
3773 "hexRef8::setRefinement(const labelList&" 3774 ", polyTopoChange&)" 3775 ) <<
"cell " << cellI
3776 <<
" of level " << cellLevel_[cellI]
3777 <<
" does not seem to have 8 points of equal or" 3778 <<
" lower level" <<
endl 3779 <<
"cellPoints:" << cPoints <<
endl 3794 Pout<<
"hexRef8::setRefinement :" 3795 <<
" Adding cells (1 per anchorPoint)" 3802 forAll(cellAnchorPoints, cellI)
3804 const labelList& cAnchors = cellAnchorPoints[cellI];
3806 if (cAnchors.
size() == 8)
3808 labelList& cAdded = cellAddedCells[cellI];
3815 newCellLevel[cellI] = cellLevel_[cellI]+1;
3818 for (
label i = 1; i < 8; i++)
3828 mesh_.cellZones().whichZone(cellI)
3832 newCellLevel(cAdded[i]) = cellLevel_[cellI]+1;
3847 Pout<<
"hexRef8::setRefinement :" 3848 <<
" Marking faces to be handled" 3856 forAll(cellMidPoint, cellI)
3858 if (cellMidPoint[cellI] >= 0)
3860 const cell& cFaces = mesh_.cells()[cellI];
3864 affectedFace.set(cFaces[i]);
3869 forAll(faceMidPoint, faceI)
3871 if (faceMidPoint[faceI] >= 0)
3873 affectedFace.set(faceI);
3877 forAll(edgeMidPoint, edgeI)
3879 if (edgeMidPoint[edgeI] >= 0)
3881 const labelList& eFaces = mesh_.edgeFaces(edgeI);
3885 affectedFace.set(eFaces[i]);
3897 Pout<<
"hexRef8::setRefinement : Splitting faces" <<
endl;
3900 forAll(faceMidPoint, faceI)
3902 if (faceMidPoint[faceI] >= 0 && affectedFace.get(faceI))
3908 const face& f = mesh_.faces()[faceI];
3912 bool modifiedFace =
false;
3914 label anchorLevel = faceAnchorLevel[faceI];
3920 label pointI = f[fp];
3922 if (pointLevel_[pointI] <= anchorLevel)
3928 faceVerts.
append(pointI);
3944 faceVerts.
append(faceMidPoint[faceI]);
3977 if (mesh_.isInternalFace(faceI))
3979 label oldOwn = mesh_.faceOwner()[faceI];
3980 label oldNei = mesh_.faceNeighbour()[faceI];
3982 checkInternalOrientation
3987 mesh_.cellCentres()[oldOwn],
3988 mesh_.cellCentres()[oldNei],
3994 label oldOwn = mesh_.faceOwner()[faceI];
3996 checkBoundaryOrientation
4001 mesh_.cellCentres()[oldOwn],
4002 mesh_.faceCentres()[faceI],
4011 modifiedFace =
true;
4013 modFace(meshMod, faceI, newFace, own, nei);
4017 addFace(meshMod, faceI, newFace, own, nei);
4023 affectedFace.unset(faceI);
4033 Pout<<
"hexRef8::setRefinement :" 4034 <<
" Adding edge splits to unsplit faces" 4041 forAll(edgeMidPoint, edgeI)
4043 if (edgeMidPoint[edgeI] >= 0)
4047 const labelList& eFaces = mesh_.edgeFaces(edgeI, eFacesStorage);
4051 label faceI = eFaces[i];
4053 if (faceMidPoint[faceI] < 0 && affectedFace.get(faceI))
4057 const face& f = mesh_.faces()[faceI];
4058 const labelList& fEdges = mesh_.faceEdges
4068 newFaceVerts.append(f[fp]);
4070 label edgeI = fEdges[fp];
4072 if (edgeMidPoint[edgeI] >= 0)
4074 newFaceVerts.
append(edgeMidPoint[edgeI]);
4083 label anchorFp = findMinLevel(f);
4100 if (mesh_.isInternalFace(faceI))
4102 label oldOwn = mesh_.faceOwner()[faceI];
4103 label oldNei = mesh_.faceNeighbour()[faceI];
4105 checkInternalOrientation
4110 mesh_.cellCentres()[oldOwn],
4111 mesh_.cellCentres()[oldNei],
4117 label oldOwn = mesh_.faceOwner()[faceI];
4119 checkBoundaryOrientation
4124 mesh_.cellCentres()[oldOwn],
4125 mesh_.faceCentres()[faceI],
4131 modFace(meshMod, faceI, newFace, own, nei);
4134 affectedFace.unset(faceI);
4146 Pout<<
"hexRef8::setRefinement :" 4147 <<
" Changing owner/neighbour for otherwise unaffected faces" 4151 forAll(affectedFace, faceI)
4153 if (affectedFace.get(faceI))
4155 const face& f = mesh_.faces()[faceI];
4159 label anchorFp = findMinLevel(f);
4173 modFace(meshMod, faceI, f, own, nei);
4176 affectedFace.unset(faceI);
4191 Pout<<
"hexRef8::setRefinement :" 4192 <<
" Create new internal faces for split cells" 4196 forAll(cellMidPoint, cellI)
4198 if (cellMidPoint[cellI] >= 0)
4223 forAll(cellMidPoint, cellI)
4225 if (cellMidPoint[cellI] >= 0)
4227 minPointI =
min(minPointI, cellMidPoint[cellI]);
4228 maxPointI =
max(maxPointI, cellMidPoint[cellI]);
4231 forAll(faceMidPoint, faceI)
4233 if (faceMidPoint[faceI] >= 0)
4235 minPointI =
min(minPointI, faceMidPoint[faceI]);
4236 maxPointI =
max(maxPointI, faceMidPoint[faceI]);
4239 forAll(edgeMidPoint, edgeI)
4241 if (edgeMidPoint[edgeI] >= 0)
4243 minPointI =
min(minPointI, edgeMidPoint[edgeI]);
4244 maxPointI =
max(maxPointI, edgeMidPoint[edgeI]);
4248 if (minPointI !=
labelMax && minPointI != mesh_.nPoints())
4251 <<
"Added point labels not consecutive to existing mesh points." 4253 <<
"mesh_.nPoints():" << mesh_.nPoints()
4254 <<
" minPointI:" << minPointI
4255 <<
" maxPointI:" << maxPointI
4260 pointLevel_.
transfer(newPointLevel);
4275 Pout<<
"hexRef8::setRefinement :" 4276 <<
" Updating refinement history to " << cellLevel_.
size()
4277 <<
" cells" <<
endl;
4283 forAll(cellAddedCells, cellI)
4285 const labelList& addedCells = cellAddedCells[cellI];
4287 if (addedCells.
size())
4301 label cellI = cellLabels[i];
4303 refinedCells[i].
transfer(cellAddedCells[cellI]);
4306 return refinedCells;
4317 savedPointLevel_.
resize(2*pointsToStore.
size());
4320 label pointI = pointsToStore[i];
4321 savedPointLevel_.
insert(pointI, pointLevel_[pointI]);
4324 savedCellLevel_.
resize(2*cellsToStore.
size());
4327 label cellI = cellsToStore[i];
4328 savedCellLevel_.
insert(cellI, cellLevel_[cellI]);
4340 updateMesh(map, dummyMap, dummyMap, dummyMap);
4358 Pout<<
"hexRef8::updateMesh :" 4359 <<
" Updating various lists" 4368 Pout<<
"hexRef8::updateMesh :" 4371 <<
" nCells:" << mesh_.nCells()
4373 <<
" cellLevel_:" << cellLevel_.
size()
4376 <<
" nPoints:" << mesh_.nPoints()
4378 <<
" pointLevel_:" << pointLevel_.
size()
4382 if (reverseCellMap.
size() == cellLevel_.
size())
4388 reorder(reverseCellMap, mesh_.nCells(), -1, cellLevel_);
4396 forAll(cellMap, newCellI)
4398 label oldCellI = cellMap[newCellI];
4402 newCellLevel[newCellI] = -1;
4406 newCellLevel[newCellI] = cellLevel_[oldCellI];
4416 label newCellI = iter.key();
4417 label storedCellI = iter();
4421 if (fnd == savedCellLevel_.
end())
4423 FatalErrorIn(
"hexRef8::updateMesh(const mapPolyMesh&)")
4424 <<
"Problem : trying to restore old value for new cell " 4425 << newCellI <<
" but cannot find old cell " << storedCellI
4426 <<
" in map of stored values " << savedCellLevel_
4429 cellLevel_[newCellI] = fnd();
4450 if (reversePointMap.
size() == pointLevel_.
size())
4453 reorder(reversePointMap, mesh_.nPoints(), -1, pointLevel_);
4462 forAll(pointMap, newPointI)
4464 label oldPointI = pointMap[newPointI];
4466 if (oldPointI == -1)
4475 newPointLevel[newPointI] = -1;
4479 newPointLevel[newPointI] = pointLevel_[oldPointI];
4482 pointLevel_.
transfer(newPointLevel);
4489 label newPointI = iter.key();
4490 label storedPointI = iter();
4494 if (fnd == savedPointLevel_.
end())
4496 FatalErrorIn(
"hexRef8::updateMesh(const mapPolyMesh&)")
4497 <<
"Problem : trying to restore old value for new point " 4498 << newPointI <<
" but cannot find old point " 4500 <<
" in map of stored values " << savedPointLevel_
4503 pointLevel_[newPointI] = fnd();
4532 cellShapesPtr_.clear();
4547 Pout<<
"hexRef8::subset :" 4548 <<
" Updating various lists" 4556 "hexRef8::subset(const labelList&, const labelList&" 4557 ", const labelList&)" 4558 ) <<
"Subsetting will not work in combination with unrefinement." 4560 <<
"Proceed at your own risk." <<
endl;
4568 forAll(cellMap, newCellI)
4570 newCellLevel[newCellI] = cellLevel_[cellMap[newCellI]];
4579 <<
"cellLevel_ contains illegal value -1 after mapping:" 4589 forAll(pointMap, newPointI)
4591 newPointLevel[newPointI] = pointLevel_[pointMap[newPointI]];
4594 pointLevel_.
transfer(newPointLevel);
4600 <<
"pointLevel_ contains illegal value -1 after mapping:" 4609 history_.
subset(pointMap, faceMap, cellMap);
4619 cellShapesPtr_.clear();
4628 Pout<<
"hexRef8::distribute :" 4629 <<
" Updating various lists" 4648 cellShapesPtr_.clear();
4654 const scalar smallDim = 1
e-6 * mesh_.bounds().mag();
4658 Pout<<
"hexRef8::checkMesh : Using matching tolerance smallDim:" 4659 << smallDim <<
endl;
4669 labelList nei(mesh_.nFaces()-mesh_.nInternalFaces());
4672 nei[i] = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
4695 label own = mesh_.faceOwner()[faceI];
4696 label bFaceI = faceI-mesh_.nInternalFaces();
4702 <<
"Faces do not seem to be correct across coupled" 4703 <<
" boundaries" <<
endl 4704 <<
"Coupled face " << faceI
4705 <<
" between owner " << own
4706 <<
" on patch " << pp.
name()
4707 <<
" and coupled neighbour " << nei[bFaceI]
4708 <<
" has two faces connected to it:" 4724 scalarField neiFaceAreas(mesh_.nFaces()-mesh_.nInternalFaces());
4727 neiFaceAreas[i] =
mag(mesh_.faceAreas()[i+mesh_.nInternalFaces()]);
4735 label faceI = i+mesh_.nInternalFaces();
4737 const scalar magArea =
mag(mesh_.faceAreas()[faceI]);
4739 if (
mag(magArea - neiFaceAreas[i]) > smallDim)
4741 const face& f = mesh_.faces()[faceI];
4742 label patchI = mesh_.boundaryMesh().whichPatch(faceI);
4744 dumpCell(mesh_.faceOwner()[faceI]);
4747 <<
"Faces do not seem to be correct across coupled" 4748 <<
" boundaries" <<
endl 4749 <<
"Coupled face " << faceI
4750 <<
" on patch " << patchI
4751 <<
" " << mesh_.boundaryMesh()[patchI].
name()
4753 <<
" has face area:" << magArea
4754 <<
" (coupled) neighbour face area differs:" 4756 <<
" to within tolerance " << smallDim
4766 labelList nVerts(mesh_.nFaces()-mesh_.nInternalFaces());
4770 nVerts[i] = mesh_.faces()[i+mesh_.nInternalFaces()].size();
4778 label faceI = i+mesh_.nInternalFaces();
4780 const face& f = mesh_.faces()[faceI];
4782 if (f.
size() != nVerts[i])
4784 dumpCell(mesh_.faceOwner()[faceI]);
4786 label patchI = mesh_.boundaryMesh().whichPatch(faceI);
4789 <<
"Faces do not seem to be correct across coupled" 4790 <<
" boundaries" <<
endl 4791 <<
"Coupled face " << faceI
4792 <<
" on patch " << patchI
4793 <<
" " << mesh_.boundaryMesh()[patchI].
name()
4795 <<
" has size:" << f.
size()
4796 <<
" (coupled) neighbour face has size:" 4808 pointField anchorPoints(mesh_.nFaces()-mesh_.nInternalFaces());
4812 label faceI = i+mesh_.nInternalFaces();
4813 const point& fc = mesh_.faceCentres()[faceI];
4814 const face& f = mesh_.faces()[faceI];
4815 const vector anchorVec(mesh_.points()[f[0]] - fc);
4817 anchorPoints[i] = anchorVec;
4827 label faceI = i+mesh_.nInternalFaces();
4828 const point& fc = mesh_.faceCentres()[faceI];
4829 const face& f = mesh_.faces()[faceI];
4830 const vector anchorVec(mesh_.points()[f[0]] - fc);
4832 if (
mag(anchorVec - anchorPoints[i]) > smallDim)
4834 dumpCell(mesh_.faceOwner()[faceI]);
4836 label patchI = mesh_.boundaryMesh().whichPatch(faceI);
4839 <<
"Faces do not seem to be correct across coupled" 4840 <<
" boundaries" <<
endl 4841 <<
"Coupled face " << faceI
4842 <<
" on patch " << patchI
4843 <<
" " << mesh_.boundaryMesh()[patchI].
name()
4845 <<
" has anchor vector:" << anchorVec
4846 <<
" (coupled) neighbour face anchor vector differs:" 4848 <<
" to within tolerance " << smallDim
4856 Pout<<
"hexRef8::checkMesh : Returning" <<
endl;
4863 const label maxPointDiff,
4869 Pout<<
"hexRef8::checkRefinementLevels :" 4870 <<
" Checking 2:1 refinement level" <<
endl;
4875 cellLevel_.
size() != mesh_.nCells()
4876 || pointLevel_.
size() != mesh_.nPoints()
4879 FatalErrorIn(
"hexRef8::checkRefinementLevels(const label)")
4880 <<
"cellLevel size should be number of cells" 4881 <<
" and pointLevel size should be number of points."<<
nl 4882 <<
"cellLevel:" << cellLevel_.
size()
4883 <<
" mesh.nCells():" << mesh_.nCells() <<
nl 4884 <<
"pointLevel:" << pointLevel_.
size()
4885 <<
" mesh.nPoints():" << mesh_.nPoints()
4895 for (
label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
4897 label own = mesh_.faceOwner()[faceI];
4898 label nei = mesh_.faceNeighbour()[faceI];
4900 if (
mag(cellLevel_[own] - cellLevel_[nei]) > 1)
4907 "hexRef8::checkRefinementLevels(const label)" 4908 ) <<
"Celllevel does not satisfy 2:1 constraint." <<
nl 4909 <<
"On face " << faceI <<
" owner cell " << own
4910 <<
" has refinement " << cellLevel_[own]
4911 <<
" neighbour cell " << nei <<
" has refinement " 4918 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
4922 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
4924 neiLevel[i] = cellLevel_[own];
4932 label faceI = i+mesh_.nInternalFaces();
4934 label own = mesh_.faceOwner()[faceI];
4936 if (
mag(cellLevel_[own] - neiLevel[i]) > 1)
4940 label patchI = mesh_.boundaryMesh().whichPatch(faceI);
4944 "hexRef8::checkRefinementLevels(const label)" 4945 ) <<
"Celllevel does not satisfy 2:1 constraint." 4946 <<
" On coupled face " << faceI
4947 <<
" on patch " << patchI <<
" " 4948 << mesh_.boundaryMesh()[patchI].name()
4949 <<
" owner cell " << own <<
" has refinement " 4951 <<
" (coupled) neighbour cell has refinement " 4974 forAll(syncPointLevel, pointI)
4976 if (pointLevel_[pointI] != syncPointLevel[pointI])
4980 "hexRef8::checkRefinementLevels(const label)" 4981 ) <<
"PointLevel is not consistent across coupled patches." 4983 <<
"point:" << pointI <<
" coord:" << mesh_.points()[pointI]
4984 <<
" has level " << pointLevel_[pointI]
4985 <<
" whereas the coupled point has level " 4986 << syncPointLevel[pointI]
4994 if (maxPointDiff != -1)
4997 labelList maxPointLevel(mesh_.nPoints(), 0);
4999 forAll(maxPointLevel, pointI)
5001 const labelList& pCells = mesh_.pointCells(pointI);
5003 label& pLevel = maxPointLevel[pointI];
5007 pLevel =
max(pLevel, cellLevel_[pCells[i]]);
5023 label pointI = pointsToCheck[i];
5025 const labelList& pCells = mesh_.pointCells(pointI);
5029 label cellI = pCells[i];
5033 mag(cellLevel_[cellI]-maxPointLevel[pointI])
5041 "hexRef8::checkRefinementLevels(const label)" 5042 ) <<
"Too big a difference between" 5043 <<
" point-connected cells." <<
nl 5045 <<
" cellLevel:" << cellLevel_[cellI]
5046 <<
" uses point:" << pointI
5047 <<
" coord:" << mesh_.points()[pointI]
5048 <<
" which is also used by a cell with level:" 5049 << maxPointLevel[pointI]
5126 if (cellShapesPtr_.empty())
5130 Pout<<
"hexRef8::cellShapes() : calculating splitHex cellShapes." 5131 <<
" cellLevel:" << cellLevel_.
size()
5132 <<
" pointLevel:" << pointLevel_.
size()
5139 label nSplitHex = 0;
5140 label nUnrecognised = 0;
5142 forAll(cellLevel_, cellI)
5144 if (meshShapes[cellI].model().index() == 0)
5146 label level = cellLevel_[cellI];
5150 bool haveQuads = matchHexShape
5160 cellShapesPtr_()[cellI] = degenerateMatcher::match(faces);
5171 Pout<<
"hexRef8::cellShapes() :" 5172 <<
" nCells:" << mesh_.nCells() <<
" of which" <<
nl 5173 <<
" primitive:" << (mesh_.nCells()-nSplitHex-nUnrecognised)
5175 <<
" split-hex:" << nSplitHex <<
nl 5176 <<
" poly :" << nUnrecognised <<
nl 5180 return cellShapesPtr_();
5200 Pout<<
"hexRef8::getSplitPoints :" 5201 <<
" Calculating unrefineable points" <<
endl;
5208 <<
"Only call if constructed with history capability" 5216 labelList splitMaster(mesh_.nPoints(), -1);
5217 labelList splitMasterLevel(mesh_.nPoints(), 0);
5222 for (
label pointI = 0; pointI < mesh_.nPoints(); pointI++)
5224 const labelList& pCells = mesh_.pointCells(pointI);
5226 if (pCells.
size() != 8)
5228 splitMaster[pointI] = -2;
5235 forAll(visibleCells, cellI)
5237 const labelList& cPoints = mesh_.cellPoints(cellI);
5239 if (visibleCells[cellI] != -1 && history_.
parentIndex(cellI) >= 0)
5246 label pointI = cPoints[i];
5248 label masterCellI = splitMaster[pointI];
5250 if (masterCellI == -1)
5257 splitMaster[pointI] = parentIndex;
5258 splitMasterLevel[pointI] = cellLevel_[cellI] - 1;
5260 else if (masterCellI == -2)
5266 (masterCellI != parentIndex)
5267 || (splitMasterLevel[pointI] != cellLevel_[cellI] - 1)
5272 splitMaster[pointI] = -2;
5281 label pointI = cPoints[i];
5283 splitMaster[pointI] = -2;
5291 label faceI = mesh_.nInternalFaces();
5292 faceI < mesh_.nFaces();
5296 const face& f = mesh_.faces()[faceI];
5300 splitMaster[f[fp]] = -2;
5307 label nSplitPoints = 0;
5309 forAll(splitMaster, pointI)
5311 if (splitMaster[pointI] >= 0)
5320 forAll(splitMaster, pointI)
5322 if (splitMaster[pointI] >= 0)
5324 splitPoints[nSplitPoints++] = pointI;
5402 Pout<<
"hexRef8::consistentUnrefinement :" 5403 <<
" Determining 2:1 consistent unrefinement" <<
endl;
5410 "hexRef8::consistentUnrefinement(const labelList&, const bool" 5411 ) <<
"maxSet not implemented yet." 5423 forAll(pointsToUnrefine, i)
5425 label pointI = pointsToUnrefine[i];
5427 unrefinePoint.set(pointI);
5438 forAll(unrefinePoint, pointI)
5440 if (unrefinePoint.get(pointI))
5442 const labelList& pCells = mesh_.pointCells(pointI);
5446 unrefineCell.set(pCells[j]);
5459 for (
label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
5461 label own = mesh_.faceOwner()[faceI];
5462 label ownLevel = cellLevel_[own] - unrefineCell.get(own);
5464 label nei = mesh_.faceNeighbour()[faceI];
5465 label neiLevel = cellLevel_[nei] - unrefineCell.get(nei);
5467 if (ownLevel < (neiLevel-1))
5474 unrefineCell.set(nei);
5485 if (unrefineCell.get(own) == 0)
5491 unrefineCell.unset(own);
5495 else if (neiLevel < (ownLevel-1))
5499 unrefineCell.set(own);
5503 if (unrefineCell.get(nei) == 0)
5509 unrefineCell.unset(nei);
5517 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
5521 label own = mesh_.faceOwner()[i+mesh_.nInternalFaces()];
5523 neiLevel[i] = cellLevel_[own] - unrefineCell.get(own);
5531 label faceI = i+mesh_.nInternalFaces();
5532 label own = mesh_.faceOwner()[faceI];
5533 label ownLevel = cellLevel_[own] - unrefineCell.get(own);
5535 if (ownLevel < (neiLevel[i]-1))
5539 if (unrefineCell.get(own) == 0)
5545 unrefineCell.unset(own);
5549 else if (neiLevel[i] < (ownLevel-1))
5553 if (unrefineCell.get(own) == 1)
5559 unrefineCell.set(own);
5569 Pout<<
"hexRef8::consistentUnrefinement :" 5570 <<
" Changed " << nChanged
5571 <<
" refinement levels due to 2:1 conflicts." 5585 forAll(unrefinePoint, pointI)
5587 if (unrefinePoint.get(pointI))
5589 const labelList& pCells = mesh_.pointCells(pointI);
5593 if (!unrefineCell.get(pCells[j]))
5595 unrefinePoint.unset(pointI);
5607 forAll(unrefinePoint, pointI)
5609 if (unrefinePoint.get(pointI))
5618 forAll(unrefinePoint, pointI)
5620 if (unrefinePoint.get(pointI))
5622 newPointsToUnrefine[nSet++] = pointI;
5626 return newPointsToUnrefine;
5640 "hexRef8::setUnrefinement(const labelList&, polyTopoChange&)" 5641 ) <<
"Only call if constructed with history capability" 5647 Pout<<
"hexRef8::setUnrefinement :" 5648 <<
" Checking initial mesh just to make sure" <<
endl;
5652 forAll(cellLevel_, cellI)
5654 if (cellLevel_[cellI] < 0)
5658 "hexRef8::setUnrefinement" 5660 "const labelList&, " 5663 ) <<
"Illegal cell level " << cellLevel_[cellI]
5664 <<
" for cell " << cellI
5671 pointSet pSet(mesh_,
"splitPoints", splitPointLabels);
5674 cellSet cSet(mesh_,
"splitPointCells", splitPointLabels.
size());
5676 forAll(splitPointLabels, i)
5678 const labelList& pCells = mesh_.pointCells(splitPointLabels[i]);
5682 cSet.insert(pCells[j]);
5687 Pout<<
"hexRef8::setRefinement : Dumping " << pSet.
size()
5689 << cSet.size() <<
" cells for unrefinement to" <<
nl 5691 <<
" cellSet " << cSet.objectPath()
5703 forAll(splitPointLabels, i)
5709 splitFaces.insert(pFaces[j]);
5723 if (facesToRemove.
size() != splitFaces.size())
5727 "hexRef8::setUnrefinement(const labelList&, polyTopoChange&)" 5728 ) <<
"Ininitial set of split points to unrefine does not" 5729 <<
" seem to be consistent or not mid points of refined cells" 5738 forAll(splitPointLabels, i)
5740 label pointI = splitPointLabels[i];
5744 const labelList& pCells = mesh_.pointCells(pointI);
5747 if (pCells.
size() != 8)
5751 "hexRef8::setUnrefinement(const labelList&, polyTopoChange&)" 5752 ) <<
"splitPoint " << pointI
5753 <<
" should have 8 cells using it. It has " << pCells
5766 label cellI = pCells[j];
5768 label region = cellRegion[cellI];
5773 <<
"Ininitial set of split points to unrefine does not" 5774 <<
" seem to be consistent or not mid points" 5775 <<
" of refined cells" <<
nl 5776 <<
"cell:" << cellI <<
" on splitPoint " << pointI
5777 <<
" has no region to be merged into" 5781 if (masterCellI != cellRegionMaster[region])
5784 <<
"cell:" << cellI <<
" on splitPoint:" << pointI
5785 <<
" in region " << region
5786 <<
" has master:" << cellRegionMaster[region]
5787 <<
" which is not the lowest numbered cell" 5788 <<
" among the pointCells:" << pCells
5808 forAll(splitPointLabels, i)
5810 label pointI = splitPointLabels[i];
5812 const labelList& pCells = mesh_.pointCells(pointI);
5818 cellLevel_[pCells[j]]--;
5836 && pointLevel_.
write()
5837 && level0Edge_.
write();
5841 writeOk = writeOk && history_.
write();
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
dimensionedScalar sqrt(const dimensionedScalar &ds)
void checkRefinementLevels(const label maxPointDiff, const labelList &pointsToCheck) const
Debug: Check 2:1 consistency across faces.
const boolList & flipMap() const
Return face flip map.
label setAction(const topoAction &action)
For compatibility with polyTopoChange: set topological action.
A face addition data class. A face can be inflated either from a point or from another face and can e...
readOption readOpt() const
virtual bool coupled() const
Return true if this patch is geometrically coupled (i.e. faces and.
vector point
Point is a vector.
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
An ordered pair of two objects of type <T> with first() and second() elements.
void setFaceInfo(const labelList &changedFaces, const List< Type > &changedFacesInfo)
Set initial changed faces.
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 ))
scalar level0EdgeLength() const
Typical edge length between unrefined points.
const word & name() const
Return name.
void unset(const PackedList< 1 > &)
Unset specified bits.
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
dimensioned< scalar > mag(const dimensioned< Type > &)
labelList consistentSlowRefinement(const label maxFaceDiff, const labelList &cellsToRefine, const labelList &facesToCheck, const label maxPointDiff, const labelList &pointsToCheck) const
Like consistentRefinement but slower:
face reverseFace() const
Return face with reverse direction.
word name(const complex &)
Return a string representation of a complex.
A cell is defined as a list of faces with extra functionality.
Reduction class. If x and y are not equal assign value.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
const labelList & reversePointMap() const
Reverse point map.
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
An STL-conforming iterator.
T & last()
Return the last element of the list.
An STL-conforming hash table.
void setInstance(const fileName &inst)
void set(const PackedList< 1 > &)
Set specified bits.
A subset of mesh faces organised as a primitive patch.
void storeSplit(const label cellI, const labelList &addedCells)
Store splitting of cell into 8.
A collection of cell labels.
void clear()
Clear the addressed list, i.e. set the size to zero.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
void resize(const label newSize)
Resize the hash table for efficiency.
const labelList & visibleCells() const
Per cell in the current mesh (i.e. visible) either -1 (unrefined)
void size(const label)
Override size to be inconsistent with allocated storage.
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
errorManipArg< error, int > exit(error &err, const int errNo=1)
label fcIndex(const label i) const
Return the forward circular index, i.e. the next index.
label getAnchorLevel(const label faceI) const
Gets level such that the face has four points <= level.
static void syncEdgeList(const polyMesh &, List< T > &, const CombineOp &cop, const T &nullValue, const TransformOp &top)
Synchronize values on all mesh edges.
prefixOSstream Perr(cerr,"Perr")
All refinement history. Used in unrefinement.
void storeData(const labelList &pointsToStore, const labelList &facesToStore, const labelList &cellsToStore)
Signal points/face/cells for which to store data.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
static void syncFaceList(const polyMesh &mesh, UList< T > &l, const CombineOp &cop)
Synchronize values on all mesh faces.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
point centre(const pointField &) const
Centre point of face.
const cellShapeList & cellShapes() const
Utility: get hexes as cell shapes.
const DynamicList< point > & points() const
Points. Shrunk after constructing mesh (or calling of compact())
label compatibleRemoves(const labelList &inPiercedFaces, labelList &cellRegion, labelList &cellRegionMaster, labelList &outPiercedFaces) const
Find faces including those with cells which have the same mastercell.
void distribute(const mapDistributePolyMesh &)
Force recalculation of locally stored data for mesh distribution.
A patch is a list of labels that address the faces in the global face list.
void updateMesh(const mapPolyMesh &)
Force recalculation of locally stored data on topological change.
label whichFace(const label globalCellID) const
Helper function to re-direct to zone::localID(...)
void resize(const label nCells)
Extend/shrink storage. additional visibleCells_ elements get.
labelList consistentRefinement(const labelList &cellsToRefine, const bool maxSet) const
Given valid mesh and current cell level and proposed.
Class containing data for cell addition.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
A topoSetSource to select a faceSet from cells.
virtual bool write() const
Write using setting from DB.
A face is a list of labels corresponding to mesh vertices.
void distribute(const mapDistributePolyMesh &)
Update local numbering for mesh redistribution.
const double e
Elementary charge.
bool updateFace(const polyMesh &, const label thisFaceI, const label neighbourCellI, const refinementData &neighbourInfo, const scalar tol, TrackingData &td)
Influence of neighbouring cell.
void setSize(const label)
Reset size of List.
Ostream & endl(Ostream &os)
Add newline and flush stream.
labelList consistentUnrefinement(const labelList &pointsToUnrefine, const bool maxSet) const
Given proposed.
Wave propagation of information through grid. Every iteration information goes through one layer of c...
labelList getSplitPoints() const
Return the points at the centre of top-level split cells.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
#define WarningIn(functionName)
Report a warning using Foam::Warning.
void clear()
Clear all entries from table.
static void syncEdgePositions(const polyMesh &mesh, List< point > &l, const CombineOp &cop, const point &nullValue)
Synchronize locations on all mesh edges.
bool write() const
Force writing refinement+history to polyMesh directory.
A List with indirect addressing.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
Transfers refinement levels such that slow transition between levels is maintained. Used in FaceCellWave.
Xfer< List< T > > xfer()
Transfer contents to the Xfer container as a plain List.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
virtual Ostream & write(const token &)=0
Write next token to stream.
label otherVertex(const label a) const
Given one vertex, return the other.
fileName objectPath() const
Return complete path + object name.
virtual bool read()
Read object.
static void swapBoundaryFaceList(const polyMesh &mesh, UList< T > &l)
Swap coupled boundary face values.
void updateMesh(const mapPolyMesh &)
Update local numbering for changed mesh.
static const label labelMin
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Class describing modification of a face.
label size() const
Return number of elements in table.
vector vec(const pointField &) const
Return the vector (end - start)
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
errorManip< error > abort(error &err)
void updateMesh(const mapPolyMesh &)
Update numbering for mesh changes.
static void syncPointList(const polyMesh &, List< T > &, const CombineOp &cop, const T &nullValue, const TransformOp &top)
Synchronize values on all mesh points.
virtual const fileName & name() const
Return the name of the stream.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
A list of faces which address into the list of points.
void checkMesh() const
Debug: Check coupled mesh for correctness.
label start() const
Return start label of this patch in the polyMesh face list.
void operator()(label &x, const label y) const
void distribute(const mapDistributePolyMesh &)
Update local numbering for mesh redistribution.
void setRefinement(const labelList &piercedFaces, const labelList &cellRegion, const labelList &cellRegionMaster, polyTopoChange &) const
Play commands into polyTopoChange to remove faces.
void flip()
Flip the face in-place.
virtual const fileName & name() const
Return the name of the stream.
Holds information regarding type of cell. Used in inside/outside determination in cellClassification...
void combineCells(const label masterCellI, const labelList &combinedCells)
Store combining 8 cells into master.
Mesh consisting of general polyhedral cells.
void subset(const labelList &pointMap, const labelList &faceMap, const labelList &cellMap)
Update local numbering for subsetted mesh.
vector normal(const pointField &) const
Vector normal; magnitude is equal to area of face.
labelListList setRefinement(const labelList &cells, polyTopoChange &)
Insert refinement. All selected cells will be split into 8.
List< label > labelList
A List of labels.
A class for handling file names.
const labelList & pointMap() const
Old point map.
Pair< label > labelPair
Label pair.
static fileCheckTypes fileModificationChecking
void reverse(UList< T > &, const label n)
const labelList & reverseCellMap() const
Reverse cell map.
Direct mesh changes based on v1.3 polyTopoChange syntax.
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
label rcIndex(const label i) const
Return the reverse circular index, i.e. the previous index.
A List obtained as a section of another List.
List< cellShape > cellShapeList
List of cellShapes and PtrList of List of cellShape.
void append(const T &)
Append an element at the end of the list.
Class containing data for point addition.
bool headerOk()
Read and check header info.
Transfers refinement levels such that slow transition between levels is maintained. Used in FaceCellWave.
iterator find(const label &)
Find and return an iterator set at the hashedEntry.
label iterate(const label maxIter)
Iterate until no changes or maxIter reached. Returns actual.
bool active() const
Is there unrefinement history. Note that this will fall over if.
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]
static void syncBoundaryFacePositions(const polyMesh &mesh, UList< point > &l, const CombineOp &cop)
Synchronize locations on boundary faces only.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
fileCheckTypes
Types of communications.
label findEdge(const edgeList &edges, const labelList &candidates, const label v0, const label v1)
Return edge among candidates that uses the two vertices.
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
void distributeCellData(List< T > &lst) const
Distribute list of cell data.
label nOldPoints() const
Number of old points.
Type gMax(const FieldField< Field, Type > &f)
labelList consistentSlowRefinement2(const label maxFaceDiff, const labelList &cellsToRefine, const labelList &facesToCheck) const
Like consistentSlowRefinement but uses different meshWave.
labelList findIndices(const ListType &, typename ListType::const_reference, const label start=0)
Find all occurences of given element. Linear search.
void setUnrefinement(const labelList &splitPointLabels, polyTopoChange &)
Remove some refinement. Needs to be supplied output of.
void subset(const labelList &pointMap, const labelList &faceMap, const labelList &cellMap)
Update numbering for subsetting.
defineTypeNameAndDebug(combustionModel, 0)
prefixOSstream Pout(cout,"Pout")
unsigned int get(const label) const
Get value at index I.
label nOldCells() const
Number of old cells.
label refinementCount() const
point centre(const pointField &) const
Return centre (centroid)
static const label labelMax
const labelList & cellMap() const
Old cell map.
label parentIndex(const label cellI) const
Get parent of cell.
void distributePointData(List< T > &lst) const
Distribute list of point data.