52 Foam::label Foam::meshRefinement::createBaffle
57 polyTopoChange& meshMod
60 const face&
f = mesh_.
faces()[faceI];
62 bool zoneFlip =
false;
66 const faceZone& fZone = mesh_.
faceZones()[zoneID];
67 zoneFlip = fZone.flipMap()[fZone.whichFace(faceI)];
95 "meshRefinement::createBaffle"
96 "(const label, const label, const label, polyTopoChange&)"
97 ) <<
"No neighbour patch for internal face " << faceI
102 bool reverseFlip =
false;
105 reverseFlip = !zoneFlip;
108 dupFaceI = meshMod.setAction
130 Foam::label Foam::meshRefinement::getBafflePatch
136 const polyBoundaryMesh&
patches = mesh_.boundaryMesh();
145 forAll(mesh_.faces()[faceI], fp)
147 label pointI = mesh_.faces()[faceI][fp];
149 forAll(mesh_.pointFaces()[pointI], pf)
151 label pFaceI = mesh_.pointFaces()[pointI][pf];
153 label patchI = patches.whichPatch(pFaceI);
155 if (patchI != -1 && !patches[patchI].coupled())
159 else if (facePatch[pFaceI] != -1)
161 return facePatch[pFaceI];
168 const cell& ownFaces = mesh_.cells()[mesh_.faceOwner()[faceI]];
172 label cFaceI = ownFaces[i];
174 label patchI = patches.whichPatch(cFaceI);
176 if (patchI != -1 && !patches[patchI].coupled())
180 else if (facePatch[cFaceI] != -1)
182 return facePatch[cFaceI];
186 if (mesh_.isInternalFace(faceI))
188 const cell& neiFaces = mesh_.cells()[mesh_.faceNeighbour()[faceI]];
192 label cFaceI = neiFaces[i];
194 label patchI = patches.whichPatch(cFaceI);
196 if (patchI != -1 && !patches[patchI].coupled())
200 else if (facePatch[cFaceI] != -1)
202 return facePatch[cFaceI];
209 "meshRefinement::getBafflePatch(const labelList&, const label)"
210 ) <<
"Could not find boundary face neighbouring internal face "
211 << faceI <<
" with face centre " << mesh_.faceCentres()[faceI]
213 <<
"Using arbitrary patch " << 0 <<
" instead." <<
endl;
220 void Foam::meshRefinement::getBafflePatches
230 autoPtr<OFstream> str;
232 if (debug&OBJINTERSECTIONS)
238 mesh_.time().path()/
timeName()/
"intersections.obj"
242 Pout<<
"getBafflePatches : Writing surface intersections to file "
246 const pointField& cellCentres = mesh_.cellCentres();
249 const wordList& cellZoneNames = surfaces_.cellZoneNames();
251 labelList surfacesToBaffle(cellZoneNames.size());
253 forAll(cellZoneNames, surfI)
255 if (cellZoneNames[surfI].size())
259 Pout<<
"getBafflePatches : Not baffling surface "
260 << surfaces_.names()[surfI] <<
endl;
265 surfacesToBaffle[baffleI++] = surfI;
268 surfacesToBaffle.setSize(baffleI);
270 ownPatch.setSize(mesh_.nFaces());
272 neiPatch.setSize(mesh_.nFaces());
289 label faceI = testFaces[i];
291 label own = mesh_.faceOwner()[faceI];
293 if (mesh_.isInternalFace(faceI))
295 start[i] = cellCentres[own];
296 end[i] = cellCentres[mesh_.faceNeighbour()[faceI]];
300 start[i] = cellCentres[own];
301 end[i] = neiCc[faceI-mesh_.nInternalFaces()];
309 List<pointIndexHit> hit1;
312 List<pointIndexHit> hit2;
314 surfaces_.findNearestIntersection
331 label faceI = testFaces[i];
333 if (hit1[i].hit() && hit2[i].hit())
345 str()<<
"l " << vertI-3 <<
' ' << vertI-2 <<
nl;
346 str()<<
"l " << vertI-2 <<
' ' << vertI-1 <<
nl;
347 str()<<
"l " << vertI-1 <<
' ' << vertI <<
nl;
351 ownPatch[faceI] = globalToPatch
353 surfaces_.globalRegion(surface1[i], region1[i])
355 neiPatch[faceI] = globalToPatch
357 surfaces_.globalRegion(surface2[i], region2[i])
360 if (ownPatch[faceI] == -1 || neiPatch[faceI] == -1)
390 ownPatch.
size() != mesh_.nFaces()
391 || neiPatch.
size() != mesh_.nFaces()
396 "meshRefinement::createBaffles"
397 "(const labelList&, const labelList&)"
398 ) <<
"Illegal size :"
399 <<
" ownPatch:" << ownPatch.
size()
400 <<
" neiPatch:" << neiPatch.
size()
401 <<
". Should be number of faces:" << mesh_.nFaces()
412 forAll(syncedOwnPatch, faceI)
416 (ownPatch[faceI] == -1 && syncedOwnPatch[faceI] != -1)
417 || (neiPatch[faceI] == -1 && syncedNeiPatch[faceI] != -1)
422 "meshRefinement::createBaffles"
423 "(const labelList&, const labelList&)"
424 ) <<
"Non synchronised at face:" << faceI
425 <<
" on patch:" << mesh_.boundaryMesh().whichPatch(faceI)
426 <<
" fc:" << mesh_.faceCentres()[faceI] <<
endl
427 <<
"ownPatch:" << ownPatch[faceI]
428 <<
" syncedOwnPatch:" << syncedOwnPatch[faceI]
429 <<
" neiPatch:" << neiPatch[faceI]
430 <<
" syncedNeiPatch:" << syncedNeiPatch[faceI]
443 if (ownPatch[faceI] != -1)
462 mesh_.updateMesh(map);
465 if (map().hasMotionPoints())
467 mesh_.movePoints(map().preMotionPoints());
477 mesh_.setInstance(oldInstance());
482 faceSet baffledFacesSet(mesh_,
"baffledFacesSet", 2*nBaffles);
484 const labelList& reverseFaceMap = map().reverseFaceMap();
485 const labelList& faceMap = map().faceMap();
488 forAll(ownPatch, oldFaceI)
490 label faceI = reverseFaceMap[oldFaceI];
492 if (ownPatch[oldFaceI] != -1 && faceI >= 0)
494 const cell& ownFaces = mesh_.cells()[mesh_.faceOwner()[faceI]];
498 baffledFacesSet.
insert(ownFaces[i]);
505 label oldFaceI = faceMap[faceI];
507 if (oldFaceI >= 0 && reverseFaceMap[oldFaceI] != faceI)
509 const cell& ownFaces = mesh_.cells()[mesh_.faceOwner()[faceI]];
513 baffledFacesSet.
insert(ownFaces[i]);
517 baffledFacesSet.
sync(mesh_);
519 updateMesh(map, baffledFacesSet.
toc());
550 label otherFaceI = duplicateFace[i];
552 if (otherFaceI != -1 && i < otherFaceI)
554 label meshFace0 = testFaces[i];
556 label meshFace1 = testFaces[otherFaceI];
561 (patch0 != -1 && isA<processorPolyPatch>(patches[patch0]))
562 || (patch1 != -1 && isA<processorPolyPatch>(patches[patch1]))
567 "meshRefinement::getDuplicateFaces"
568 "(const bool, const labelList&)"
569 ) <<
"One of two duplicate faces is on"
570 <<
" processorPolyPatch."
571 <<
"This is not allowed." <<
nl
572 <<
"Face:" << meshFace0
573 <<
" is on patch:" << patches[patch0].
name()
575 <<
"Face:" << meshFace1
576 <<
" is on patch:" << patches[patch1].
name()
580 duplicateFaces[dupI++] =
labelPair(meshFace0, meshFace1);
583 duplicateFaces.setSize(dupI);
586 <<
" pairs of duplicate faces." <<
nl <<
endl;
591 faceSet duplicateFaceSet(mesh_,
"duplicateFaces", 2*dupI);
595 duplicateFaceSet.
insert(duplicateFaces[i][0]);
596 duplicateFaceSet.
insert(duplicateFaces[i][1]);
598 Pout<<
"Writing duplicate faces (baffles) to faceSet "
600 duplicateFaceSet.
write();
603 return duplicateFaces;
628 labelList nBafflesPerEdge(mesh_.nEdges(), 0);
644 label faceI = pp.
start();
648 const labelList& fEdges = mesh_.faceEdges(faceI);
652 nBafflesPerEdge[fEdges[fEdgeI]]++;
660 DynamicList<label> fe0;
661 DynamicList<label> fe1;
669 const labelList& fEdges0 = mesh_.faceEdges(couples[i].first(), fe0);
673 nBafflesPerEdge[fEdges0[fEdgeI]] += 1000000;
676 const labelList& fEdges1 = mesh_.faceEdges(couples[i].second(), fe1);
680 nBafflesPerEdge[fEdges1[fEdgeI]] += 1000000;
698 List<labelPair> filteredCouples(couples.
size());
707 patches.whichPatch(couple.first())
708 == patches.whichPatch(couple.second())
711 const labelList& fEdges = mesh_.faceEdges(couples[i].first());
715 label edgeI = fEdges[fEdgeI];
717 if (nBafflesPerEdge[edgeI] == 2*1000000+2*1)
719 filteredCouples[filterI++] = couples[i];
725 filteredCouples.
setSize(filterI);
733 return filteredCouples;
746 const faceList& faces = mesh_.faces();
747 const labelList& faceOwner = mesh_.faceOwner();
752 label face0 = couples[i].first();
753 label face1 = couples[i].second();
757 label own0 = faceOwner[face0];
758 label own1 = faceOwner[face1];
760 if (face1 < 0 || own0 < own1)
763 label zoneID = faceZones.
whichZone(face0);
764 bool zoneFlip =
false;
768 const faceZone& fZone = faceZones[zoneID];
772 label nei = (face1 < 0 ? -1 : own1);
794 label zoneID = faceZones.
whichZone(face1);
795 bool zoneFlip =
false;
799 const faceZone& fZone = faceZones[zoneID];
826 mesh_.updateMesh(map);
829 if (map().hasMotionPoints())
831 mesh_.movePoints(map().preMotionPoints());
841 mesh_.setInstance(oldInstance());
852 label newFace0 = map().reverseFaceMap()[couples[i].first()];
855 newExposedFaces[newI++] = newFace0;
858 label newFace1 = map().reverseFaceMap()[couples[i].second()];
861 newExposedFaces[newI++] = newFace1;
864 newExposedFaces.setSize(newI);
865 updateMesh(map, newExposedFaces);
872 void Foam::meshRefinement::findCellZoneGeometric
881 const pointField& cellCentres = mesh_.cellCentres();
882 const labelList& faceOwner = mesh_.faceOwner();
883 const labelList& faceNeighbour = mesh_.faceNeighbour();
894 forAll(insideSurfaces, cellI)
896 if (cellToZone[cellI] == -2)
898 label surfI = insideSurfaces[cellI];
902 cellToZone[cellI] = surfaceToCellZone[surfI];
915 label nCandidates = 0;
916 forAll(namedSurfaceIndex, faceI)
918 label surfI = namedSurfaceIndex[faceI];
922 if (mesh_.isInternalFace(faceI))
936 forAll(namedSurfaceIndex, faceI)
938 label surfI = namedSurfaceIndex[faceI];
942 label own = faceOwner[faceI];
943 const point& ownCc = cellCentres[own];
945 if (mesh_.isInternalFace(faceI))
947 label nei = faceNeighbour[faceI];
948 const point& neiCc = cellCentres[nei];
950 const vector d = 1
E-4*(neiCc - ownCc);
951 candidatePoints[nCandidates++] = ownCc-
d;
952 candidatePoints[nCandidates++] = neiCc+
d;
956 const point& neiFc = mesh_.faceCentres()[faceI];
958 const vector d = 1
E-4*(neiFc - ownCc);
959 candidatePoints[nCandidates++] = ownCc-
d;
978 forAll(namedSurfaceIndex, faceI)
980 label surfI = namedSurfaceIndex[faceI];
984 label own = faceOwner[faceI];
986 if (mesh_.isInternalFace(faceI))
988 label ownSurfI = insideSurfaces[nCandidates++];
991 cellToZone[own] = surfaceToCellZone[ownSurfI];
994 label neiSurfI = insideSurfaces[nCandidates++];
997 label nei = faceNeighbour[faceI];
999 cellToZone[nei] = surfaceToCellZone[neiSurfI];
1004 label ownSurfI = insideSurfaces[nCandidates++];
1007 cellToZone[own] = surfaceToCellZone[ownSurfI];
1018 for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
1020 label ownZone = cellToZone[mesh_.faceOwner()[faceI]];
1021 label neiZone = cellToZone[mesh_.faceNeighbour()[faceI]];
1023 if (namedSurfaceIndex[faceI] == -1 && (ownZone != neiZone))
1029 max(ownZone, neiZone)
1034 labelList neiCellZone(mesh_.nFaces()-mesh_.nInternalFaces());
1035 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
1039 const polyPatch& pp = patches[patchI];
1045 label faceI = pp.start()+i;
1046 label ownZone = cellToZone[mesh_.faceOwner()[faceI]];
1047 neiCellZone[faceI-mesh_.nInternalFaces()] = ownZone;
1055 const polyPatch& pp = patches[patchI];
1061 label faceI = pp.start()+i;
1062 label ownZone = cellToZone[mesh_.faceOwner()[faceI]];
1063 label neiZone = neiCellZone[faceI-mesh_.nInternalFaces()];
1065 if (namedSurfaceIndex[faceI] == -1 && (ownZone != neiZone))
1071 max(ownZone, neiZone)
1089 bool Foam::meshRefinement::calcRegionToZone
1091 const label surfZoneI,
1092 const label ownRegion,
1093 const label neiRegion,
1098 bool changed =
false;
1101 if (ownRegion != neiRegion)
1108 if (regionToCellZone[ownRegion] == -2)
1110 if (regionToCellZone[neiRegion] == surfZoneI)
1114 regionToCellZone[ownRegion] = -1;
1117 else if (regionToCellZone[neiRegion] != -2)
1121 regionToCellZone[ownRegion] = surfZoneI;
1125 else if (regionToCellZone[neiRegion] == -2)
1127 if (regionToCellZone[ownRegion] == surfZoneI)
1131 regionToCellZone[neiRegion] = -1;
1134 else if (regionToCellZone[ownRegion] != -2)
1138 regionToCellZone[neiRegion] = surfZoneI;
1151 void Foam::meshRefinement::findCellZoneTopo
1153 const point& keepPoint,
1160 boolList blockedFace(mesh_.nFaces());
1162 forAll(namedSurfaceIndex, faceI)
1164 if (namedSurfaceIndex[faceI] == -1)
1166 blockedFace[faceI] =
false;
1170 blockedFace[faceI] =
true;
1176 regionSplit cellRegion(mesh_, blockedFace);
1177 blockedFace.clear();
1183 labelList regionToCellZone(cellRegion.nRegions(), -2);
1188 forAll(cellToZone, cellI)
1190 if (cellToZone[cellI] != -2)
1192 regionToCellZone[cellRegion[cellI]] = cellToZone[cellI];
1199 label keepRegionI = -1;
1201 label cellI = mesh_.findCell(keepPoint);
1205 keepRegionI = cellRegion[cellI];
1207 reduce(keepRegionI, maxOp<label>());
1209 Info<<
"Found point " << keepPoint <<
" in cell " << cellI
1210 <<
" in global region " << keepRegionI
1211 <<
" out of " << cellRegion.nRegions() <<
" regions." <<
endl;
1213 if (keepRegionI == -1)
1217 "meshRefinement::findCellZoneTopo"
1218 "(const point&, const labelList&, const labelList&, labelList&)"
1219 ) <<
"Point " << keepPoint
1220 <<
" is not inside the mesh." <<
nl
1221 <<
"Bounding box of the mesh:" << mesh_.globalData().bb()
1226 if (regionToCellZone[keepRegionI] == -2)
1228 regionToCellZone[keepRegionI] = -1;
1247 bool changed =
false;
1251 for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
1253 label surfI = namedSurfaceIndex[faceI];
1259 bool changedCell = calcRegionToZone
1261 surfaceToCellZone[surfI],
1262 cellRegion[mesh_.faceOwner()[faceI]],
1263 cellRegion[mesh_.faceNeighbour()[faceI]],
1267 changed = changed | changedCell;
1273 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
1276 labelList neiCellRegion(mesh_.nFaces()-mesh_.nInternalFaces());
1279 const polyPatch& pp = patches[patchI];
1285 label faceI = pp.start()+i;
1286 neiCellRegion[faceI-mesh_.nInternalFaces()] =
1287 cellRegion[mesh_.faceOwner()[faceI]];
1297 const polyPatch& pp = patches[patchI];
1303 label faceI = pp.start()+i;
1305 label surfI = namedSurfaceIndex[faceI];
1309 bool changedCell = calcRegionToZone
1311 surfaceToCellZone[surfI],
1312 cellRegion[mesh_.faceOwner()[faceI]],
1313 neiCellRegion[faceI-mesh_.nInternalFaces()],
1317 changed = changed | changedCell;
1330 forAll(regionToCellZone, regionI)
1332 label zoneI = regionToCellZone[regionI];
1338 "meshRefinement::findCellZoneTopo"
1339 "(const point&, const labelList&, const labelList&, labelList&)"
1340 ) <<
"For region " << regionI <<
" haven't set cell zone."
1347 forAll(regionToCellZone, regionI)
1349 Pout<<
"Region " << regionI
1350 <<
" becomes cellZone:" << regionToCellZone[regionI]
1356 forAll(cellToZone, cellI)
1358 cellToZone[cellI] = regionToCellZone[cellRegion[cellI]];
1365 void Foam::meshRefinement::makeConsistentFaceIndex
1371 const labelList& faceOwner = mesh_.faceOwner();
1372 const labelList& faceNeighbour = mesh_.faceNeighbour();
1374 for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
1376 label ownZone = cellToZone[faceOwner[faceI]];
1377 label neiZone = cellToZone[faceNeighbour[faceI]];
1379 if (ownZone == neiZone && namedSurfaceIndex[faceI] != -1)
1381 namedSurfaceIndex[faceI] = -1;
1383 else if (ownZone != neiZone && namedSurfaceIndex[faceI] == -1)
1386 <<
"Different cell zones on either side of face " << faceI
1387 <<
" at " << mesh_.faceCentres()[faceI]
1388 <<
" but face not marked with a surface."
1393 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
1396 labelList neiCellZone(mesh_.nFaces()-mesh_.nInternalFaces());
1399 const polyPatch& pp = patches[patchI];
1405 label faceI = pp.start()+i;
1406 neiCellZone[faceI-mesh_.nInternalFaces()] =
1407 cellToZone[mesh_.faceOwner()[faceI]];
1416 const polyPatch& pp = patches[patchI];
1422 label faceI = pp.start()+i;
1424 label ownZone = cellToZone[faceOwner[faceI]];
1425 label neiZone = neiCellZone[faceI-mesh_.nInternalFaces()];
1427 if (ownZone == neiZone && namedSurfaceIndex[faceI] != -1)
1429 namedSurfaceIndex[faceI] = -1;
1431 else if (ownZone != neiZone && namedSurfaceIndex[faceI] == -1)
1434 <<
"Different cell zones on either side of face "
1435 << faceI <<
" at " << mesh_.faceCentres()[faceI]
1436 <<
" but face not marked with a surface."
1450 const bool handleSnapProblems,
1451 const bool removeEdgeConnectedCells,
1453 const bool mergeFreeStanding,
1457 const point& keepPoint
1466 Info<<
"Introducing baffles for "
1468 <<
" faces that are intersected by the surface." <<
nl <<
endl;
1471 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
1472 pointField neiCc(mesh_.nFaces()-mesh_.nInternalFaces());
1473 calcNeighbourData(neiLevel, neiCc);
1491 createBaffles(ownPatch, neiPatch);
1499 Info<<
"Created baffles in = "
1502 printMeshInfo(debug,
"After introducing baffles");
1508 write(debug, runTime.
path()/
"baffles");
1509 Pout<<
"Dumped debug data in = "
1519 if (handleSnapProblems)
1522 <<
"Introducing baffles to block off problem cells" <<
nl
1523 <<
"----------------------------------------------" <<
nl
1528 markFacesOnProblemCells
1531 removeEdgeConnectedCells,
1537 Info<<
"Analyzed problem cells in = "
1542 faceSet problemTopo(mesh_,
"problemFacesTopo", 100);
1546 markFacesOnProblemCells
1549 removeEdgeConnectedCells,
1554 forAll(facePatchTopo, faceI)
1556 if (facePatchTopo[faceI] != -1)
1558 problemTopo.
insert(faceI);
1561 Pout<<
"Dumping " << problemTopo.
size()
1563 problemTopo.
write();
1566 Info<<
"Introducing baffles to delete problem cells." <<
nl <<
endl;
1574 createBaffles(facePatch, facePatch);
1581 Info<<
"Created baffles in = "
1584 printMeshInfo(debug,
"After introducing baffles");
1588 Pout<<
"Writing extra baffled mesh to time "
1590 write(debug, runTime.
path()/
"extraBaffles");
1591 Pout<<
"Dumped debug data in = "
1601 <<
"Remove unreachable sections of mesh" <<
nl
1602 <<
"-----------------------------------" <<
nl
1610 splitMeshRegions(keepPoint);
1617 Info<<
"Split mesh in = "
1620 printMeshInfo(debug,
"After subsetting");
1627 Pout<<
"Dumped debug data in = "
1635 if (mergeFreeStanding)
1638 <<
"Merge free-standing baffles" <<
nl
1639 <<
"---------------------------" <<
nl
1650 filterDuplicateFaces
1654 identity(mesh_.nFaces()-mesh_.nInternalFaces())
1655 +mesh_.nInternalFaces()
1660 label nCouples = couples.
size();
1663 Info<<
"Detected free-standing baffles : " << nCouples <<
endl;
1670 mergeBaffles(couples);
1678 Info<<
"Merged free-standing baffles in = "
1687 const label nBufferLayers,
1689 const point& keepPoint
1696 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
1697 pointField neiCc(mesh_.nFaces()-mesh_.nInternalFaces());
1698 calcNeighbourData(neiLevel, neiCc);
1712 boolList blockedFace(mesh_.nFaces(),
false);
1716 if (ownPatch[faceI] != -1 || neiPatch[faceI] != -1)
1718 blockedFace[faceI] =
true;
1725 blockedFace.clear();
1728 label keepRegionI = -1;
1730 label cellI = mesh_.findCell(keepPoint);
1734 keepRegionI = cellRegion[cellI];
1738 Info<<
"Found point " << keepPoint <<
" in cell " << cellI
1739 <<
" in global region " << keepRegionI
1740 <<
" out of " << cellRegion.
nRegions() <<
" regions." <<
endl;
1742 if (keepRegionI == -1)
1746 "meshRefinement::splitMesh"
1747 "(const label, const labelList&, const point&)"
1748 ) <<
"Point " << keepPoint
1749 <<
" is not inside the mesh." <<
nl
1750 <<
"Bounding box of the mesh:" << mesh_.globalData().bb()
1761 const labelList& faceOwner = mesh_.faceOwner();
1762 const labelList& faceNeighbour = mesh_.faceNeighbour();
1765 label defaultPatch = 0;
1766 if (globalToPatch.
size())
1768 defaultPatch = globalToPatch[0];
1771 for (label i = 0; i < nBufferLayers; i++)
1775 labelList pointBaffle(mesh_.nPoints(), -1);
1777 forAll(faceNeighbour, faceI)
1779 const face& f = mesh_.faces()[faceI];
1781 label ownRegion = cellRegion[faceOwner[faceI]];
1782 label neiRegion = cellRegion[faceNeighbour[faceI]];
1784 if (ownRegion == keepRegionI && neiRegion != keepRegionI)
1791 pointBaffle[f[fp]] =
max(defaultPatch, ownPatch[faceI]);
1794 else if (ownRegion != keepRegionI && neiRegion == keepRegionI)
1796 label newPatchI = neiPatch[faceI];
1797 if (newPatchI == -1)
1799 newPatchI =
max(defaultPatch, ownPatch[faceI]);
1803 pointBaffle[f[fp]] = newPatchI;
1809 label faceI = mesh_.nInternalFaces();
1810 faceI < mesh_.nFaces();
1814 const face& f = mesh_.faces()[faceI];
1816 label ownRegion = cellRegion[faceOwner[faceI]];
1818 if (ownRegion == keepRegionI)
1822 pointBaffle[f[fp]] =
max(defaultPatch, ownPatch[faceI]);
1842 forAll(pointFaces, pointI)
1844 if (pointBaffle[pointI] != -1)
1850 label faceI = pFaces[pFaceI];
1852 if (ownPatch[faceI] == -1)
1854 ownPatch[faceI] = pointBaffle[pointI];
1868 if (ownPatch[faceI] != -1)
1870 label own = faceOwner[faceI];
1872 if (cellRegion[own] != keepRegionI)
1874 cellRegion[own] = keepRegionI;
1876 const cell& ownFaces = mesh_.cells()[own];
1879 if (ownPatch[ownFaces[j]] == -1)
1881 newOwnPatch[ownFaces[j]] = ownPatch[faceI];
1885 if (mesh_.isInternalFace(faceI))
1887 label nei = faceNeighbour[faceI];
1889 if (cellRegion[nei] != keepRegionI)
1891 cellRegion[nei] = keepRegionI;
1893 const cell& neiFaces = mesh_.cells()[nei];
1896 if (ownPatch[neiFaces[j]] == -1)
1898 newOwnPatch[neiFaces[j]] = ownPatch[faceI];
1918 forAll(cellRegion, cellI)
1920 if (cellRegion[cellI] != keepRegionI)
1922 cellsToRemove.append(cellI);
1925 cellsToRemove.shrink();
1927 label nCellsToKeep = mesh_.nCells() - cellsToRemove.size();
1930 Info<<
"Keeping all cells in region " << keepRegionI
1931 <<
" containing point " << keepPoint <<
endl
1932 <<
"Selected for keeping : " << nCellsToKeep
1933 <<
" cells." <<
endl;
1941 labelList exposedPatches(exposedFaces.size());
1945 label faceI = exposedFaces[i];
1947 if (ownPatch[faceI] != -1)
1949 exposedPatches[i] = ownPatch[faceI];
1953 WarningIn(
"meshRefinement::splitMesh(..)")
1954 <<
"For exposed face " << faceI
1955 <<
" fc:" << mesh_.faceCentres()[faceI]
1956 <<
" found no patch." << endl
1957 <<
" Taking patch " << defaultPatch
1958 <<
" instead." <<
endl;
1959 exposedPatches[i] = defaultPatch;
1963 return doRemoveCells
1990 Info<<
"dupNonManifoldPoints : Found : " << nNonManifPoints
1991 <<
" non-manifold points (out of "
1992 << mesh_.globalData().nTotalPoints()
2005 mesh_.updateMesh(map);
2008 if (map().hasMotionPoints())
2010 mesh_.movePoints(map().preMotionPoints());
2020 mesh_.setInstance(oldInstance());
2034 const point& keepPoint,
2035 const bool allowFreeStandingZoneFaces
2038 const wordList& cellZoneNames = surfaces_.cellZoneNames();
2039 const wordList& faceZoneNames = surfaces_.faceZoneNames();
2041 labelList namedSurfaces(surfaces_.getNamedSurfaces());
2047 label surfI = namedSurfaces[i];
2049 isNamedSurface[surfI] =
true;
2051 Info<<
"Surface : " << surfaces_.names()[surfI] <<
nl
2052 <<
" faceZone : " << faceZoneNames[surfI] <<
nl
2053 <<
" cellZone : " << cellZoneNames[surfI] <<
endl;
2065 label surfI = namedSurfaces[i];
2067 label zoneI = faceZones.
findZoneID(faceZoneNames[surfI]);
2071 zoneI = faceZones.
size();
2078 faceZoneNames[surfI],
2089 Pout<<
"Faces on " << surfaces_.names()[surfI]
2090 <<
" will go into faceZone " << zoneI <<
endl;
2092 surfaceToFaceZone[surfI] = zoneI;
2101 for (label procI = 1; procI < allFaceZones.
size(); procI++)
2103 if (allFaceZones[procI] != allFaceZones[0])
2107 "meshRefinement::zonify"
2108 "(const label, const point&)"
2109 ) <<
"Zones not synchronised among processors." <<
nl
2110 <<
" Processor0 has faceZones:" << allFaceZones[0]
2111 <<
" , processor" << procI
2112 <<
" has faceZones:" << allFaceZones[procI]
2124 label surfI = namedSurfaces[i];
2126 label zoneI = cellZones.
findZoneID(cellZoneNames[surfI]);
2130 zoneI = cellZones.
size();
2137 cellZoneNames[surfI],
2147 Pout<<
"Cells inside " << surfaces_.names()[surfI]
2148 <<
" will go into cellZone " << zoneI <<
endl;
2150 surfaceToCellZone[surfI] = zoneI;
2159 for (label procI = 1; procI < allCellZones.
size(); procI++)
2161 if (allCellZones[procI] != allCellZones[0])
2165 "meshRefinement::zonify"
2166 "(const label, const point&)"
2167 ) <<
"Zones not synchronised among processors." <<
nl
2168 <<
" Processor0 has cellZones:" << allCellZones[0]
2169 <<
" , processor" << procI
2170 <<
" has cellZones:" << allCellZones[procI]
2178 const pointField& cellCentres = mesh_.cellCentres();
2179 const labelList& faceOwner = mesh_.faceOwner();
2180 const labelList& faceNeighbour = mesh_.faceNeighbour();
2189 labelList namedSurfaceIndex(mesh_.nFaces(), -1);
2196 labelList neiLevel(mesh_.nFaces()-mesh_.nInternalFaces());
2197 pointField neiCc(mesh_.nFaces()-mesh_.nInternalFaces());
2198 calcNeighbourData(neiLevel, neiCc);
2208 labelList testFaces(intersectedFaces());
2218 label faceI = testFaces[i];
2220 if (mesh_.isInternalFace(faceI))
2222 start[i] = cellCentres[faceOwner[faceI]];
2223 end[i] = cellCentres[faceNeighbour[faceI]];
2227 start[i] = cellCentres[faceOwner[faceI]];
2228 end[i] = neiCc[faceI-mesh_.nInternalFaces()];
2245 surfaces_.findNearestIntersection
2262 label faceI = testFaces[i];
2264 if (surface1[i] != -1)
2267 namedSurfaceIndex[faceI] = surface1[i];
2268 nSurfFaces[surface1[i]]++;
2270 else if (surface2[i] != -1)
2272 namedSurfaceIndex[faceI] = surface2[i];
2273 nSurfFaces[surface2[i]]++;
2293 forAll(nSurfFaces, surfI)
2296 << surfaces_.names()[surfI]
2297 <<
" nZoneFaces:" << nSurfFaces[surfI] <<
nl;
2308 labelList closedNamedSurfaces(surfaces_.getClosedNamedSurfaces());
2314 labelList cellToZone(mesh_.nCells(), -2);
2320 if (closedNamedSurfaces.
size())
2322 findCellZoneGeometric
2324 closedNamedSurfaces,
2348 if (!allowFreeStandingZoneFaces)
2350 makeConsistentFaceIndex(cellToZone, namedSurfaceIndex);
2360 for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
2362 label surfI = namedSurfaceIndex[faceI];
2367 label ownZone = cellToZone[faceOwner[faceI]];
2368 label neiZone = cellToZone[faceNeighbour[faceI]];
2371 if (ownZone ==
max(ownZone, neiZone))
2384 mesh_.faces()[faceI],
2387 faceNeighbour[faceI],
2391 surfaceToFaceZone[surfI],
2399 labelList neiCellZone(mesh_.nFaces()-mesh_.nInternalFaces(), -1);
2408 label faceI = pp.
start()+i;
2409 neiCellZone[faceI-mesh_.nInternalFaces()] =
2410 cellToZone[mesh_.faceOwner()[faceI]];
2424 label faceI = pp.
start();
2428 label surfI = namedSurfaceIndex[faceI];
2432 label ownZone = cellToZone[faceOwner[faceI]];
2433 label neiZone = neiCellZone[faceI-mesh_.nInternalFaces()];
2437 label maxZone =
max(ownZone, neiZone);
2443 else if (ownZone == neiZone)
2447 flip = !isMasterFace[faceI];
2449 else if (neiZone == maxZone)
2462 mesh_.faces()[faceI],
2469 surfaceToFaceZone[surfI],
2482 forAll(cellToZone, cellI)
2484 label zoneI = cellToZone[cellI];
2504 mesh_.updateMesh(map);
2507 if (map().hasMotionPoints())
2509 mesh_.movePoints(map().preMotionPoints());
2519 mesh_.setInstance(oldInstance());