diff --git a/applications/utilities/mesh/advanced/removeFaces/removeFaces.C b/applications/utilities/mesh/advanced/removeFaces/removeFaces.C index 4133fbb9e..c3ea4ed34 100644 --- a/applications/utilities/mesh/advanced/removeFaces/removeFaces.C +++ b/applications/utilities/mesh/advanced/removeFaces/removeFaces.C @@ -147,11 +147,15 @@ int main(int argc, char *argv[]) // Topo changes container directTopoChange meshMod(mesh); + // Dummy point region master + labelList pointRegionMaster(cellRegionMaster.size(), -1); + // Insert mesh refinement into directTopoChange. faceRemover.setRefinement ( facesToRemove, cellRegion, + pointRegionMaster, cellRegionMaster, meshMod ); diff --git a/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/hexRef8.C b/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/hexRef8.C index d26289f85..c91ffdd08 100644 --- a/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/hexRef8.C +++ b/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/hexRef8.C @@ -5345,10 +5345,17 @@ void Foam::hexRef8::setUnrefinement // Insert all commands to combine cells. Never fails so don't have to // test for success. + + // Added dummy pointRegionMaster, which will be ignored. + // Probably, splitPointLabels should be used + // HJ, 6/Sep/2019 + labelList pointRegionMaster(cellRegionMaster.size(), label(-1)); + faceRemover_.setRefinement ( facesToRemove, cellRegion, + pointRegionMaster, cellRegionMaster, meshMod ); diff --git a/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFaces.C b/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFaces.C index 863ed6e49..d93ffad97 100644 --- a/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFaces.C +++ b/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFaces.C @@ -149,10 +149,14 @@ Foam::Xfer Foam::removeFaces::affectedFaces { const label& region = cellRegion[cellI]; - if (region != -1 && (cellI != cellRegionMaster[region])) + // Bug fix: map coarse cell from all fine neighbours + // All cells will be removed and replacew with a new one from + // the central point. HJ, 9/Sep/2019 + if (region != -1) { // Get this cell (list of cell faces) and mark all of its faces const labelList& cFaces = meshCells[cellI]; + forAll(cFaces, cFaceI) { affectedFace[cFaces[cFaceI]] = true; @@ -284,7 +288,7 @@ Foam::label Foam::removeFaces::compatibleRemoves ( const labelList& facesToRemove, labelList& cellRegion, - labelList& regionMaster, + labelList& cellRegionMaster, labelList& newFacesToRemove ) const { @@ -298,8 +302,8 @@ Foam::label Foam::removeFaces::compatibleRemoves cellRegion.setSize(nCells); cellRegion = -1; - regionMaster.setSize(nCells); - regionMaster = -1; + cellRegionMaster.setSize(nCells); + cellRegionMaster = -1; // Count regions label nRegions = 0; @@ -342,7 +346,7 @@ Foam::label Foam::removeFaces::compatibleRemoves // Make owner (lowest numbered!) the master of the region and // increment the number of regions - regionMaster[nRegions] = own; + cellRegionMaster[nRegions] = own; ++nRegions; } else @@ -352,7 +356,7 @@ Foam::label Foam::removeFaces::compatibleRemoves // See if owner becomes the master of the region (if its index // is lower than the current master of the region) - regionMaster[region1] = min(own, regionMaster[region1]); + cellRegionMaster[region1] = min(own, cellRegionMaster[region1]); } } else @@ -400,24 +404,29 @@ Foam::label Foam::removeFaces::compatibleRemoves freedRegion = region0; } - label master0 = regionMaster[region0]; - label master1 = regionMaster[region1]; + label master0 = cellRegionMaster[region0]; + label master1 = cellRegionMaster[region1]; - regionMaster[freedRegion] = -1; - regionMaster[keptRegion] = min(master0, master1); + cellRegionMaster[freedRegion] = -1; + cellRegionMaster[keptRegion] = min(master0, master1); } } } // Set size - regionMaster.setSize(nRegions); + cellRegionMaster.setSize(nRegions); + // Note: + // It is legal to have cellRegionMaster = -1 if the region has been created + // and then abandoned because it has been merged with another region + // HJ, 6/Sep/2019 + // Various checks, additional scope for clarity and memory management // - master is lowest numbered in any region // - regions have more than 1 cell { // Number of cells per region - labelList nCells(regionMaster.size(), 0); + labelList nCells(cellRegionMaster.size(), 0); // Loop through cell regions forAll(cellRegion, cellI) @@ -430,15 +439,12 @@ Foam::label Foam::removeFaces::compatibleRemoves // Region found, increment number of cells per region ++nCells[r]; - if (cellI < regionMaster[r]) + if (cellI < cellRegionMaster[r]) { - FatalErrorIn - ( - "removeFaces::compatibleRemoves(const labelList&" - ", labelList&, labelList&, labelList&)" - ) << "Not lowest numbered! Cell: " << cellI + FatalErrorInFunction + << "Not lowest numbered! Cell: " << cellI << " region: " << r - << " region master: " << regionMaster[r] + << " region master: " << cellRegionMaster[r] << abort(FatalError); } } @@ -449,11 +455,8 @@ Foam::label Foam::removeFaces::compatibleRemoves { if (nCells[regionI] == 1) { - FatalErrorIn - ( - "removeFaces::compatibleRemoves(const labelList&" - ", labelList&, labelList&, labelList&)" - ) << "Region " << regionI + FatalErrorInFunction + << "Region " << regionI << " has only " << nCells[regionI] << " cell in it." << abort(FatalError); } @@ -464,9 +467,9 @@ Foam::label Foam::removeFaces::compatibleRemoves // Count number of used regions label nUsedRegions = 0; - forAll(regionMaster, i) + forAll(cellRegionMaster, i) { - if (regionMaster[i] != -1) + if (cellRegionMaster[i] != -1) { ++nUsedRegions; } diff --git a/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFaces.H b/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFaces.H index 650d075f2..4a462b556 100644 --- a/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFaces.H +++ b/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFaces.H @@ -28,8 +28,6 @@ Description Given list of faces to remove insert all the topology changes. Contains helper function to get consistent set of faces to remove. - Not very well tested in parallel. - SourceFiles removeFaces.C @@ -185,6 +183,7 @@ public: // and always merge (if on same patch) removeFaces(const polyMesh&, const scalar minCos); + // Member Functions //- Given a set of faces to pierce calculates: @@ -210,7 +209,9 @@ public: ( const labelList& piercedFaces, const labelList& cellRegion, - const labelList& cellRegionMaster, + const labelList& pointRegionMaster, + labelList& cellRegionMaster, + TopoChangeEngine& ref ) const; diff --git a/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFacesTemplates.C b/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFacesTemplates.C index f9790d5cc..9004c37b7 100644 --- a/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFacesTemplates.C +++ b/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/removeFacesTemplates.C @@ -67,18 +67,9 @@ void Foam::removeFaces::mergeFaces if (fp.edgeLoops().size() != 1) { writeOBJ(fp, mesh_.time().path()/"facesToBeMerged.obj"); - FatalErrorIn - ( - "template" - "void removeFaces::mergeFaces" - "\n(" - "\n const labelList& cellRegion," - "\n const labelList& cellRegionMaster," - "\n const labelHashSet& pointsToRemove," - "\n const labelList& faceLabels," - "\n const TopoChangeEngine& ref" - "\n) const" - ) << "Cannot merge faces " << faceLabels + + FatalErrorInFunction + << "Cannot merge faces " << faceLabels << " into single face since outside vertices " << fp.edgeLoops() << " do not form single loop but form " << fp.edgeLoops().size() << " loops instead." << abort(FatalError); @@ -131,18 +122,9 @@ void Foam::removeFaces::mergeFaces if (masterIndex == -1) { writeOBJ(fp, mesh_.time().path()/"facesToBeMerged.obj"); - FatalErrorIn - ( - "template" - "void removeFaces::mergeFaces" - "\n(" - "\n const labelList& cellRegion," - "\n const labelList& cellRegionMaster," - "\n const labelHashSet& pointsToRemove," - "\n const labelList& faceLabels," - "\n const TopoChangeEngine& ref" - "\n) const" - ) << "Problem." << abort(FatalError); + + FatalErrorInFunction + << "Problem." << abort(FatalError); } @@ -312,7 +294,8 @@ void Foam::removeFaces::setRefinement ( const labelList& faceLabels, const labelList& cellRegion, - const labelList& cellRegionMaster, + const labelList& pointRegionMaster, + labelList& cellRegionMaster, TopoChangeEngine& ref ) const @@ -343,18 +326,8 @@ void Foam::removeFaces::setRefinement if (!mesh_.isInternalFace(faceI)) { - FatalErrorIn - ( - "template" - "\n" - "void removeFaces::setRefinement" - "\n(" - "\n const labelList& faceLabels," - "\n const labelList& cellRegion," - "\n const labelList& cellRegionMaster," - "\n TopoChangeEngine& ref" - "\n)" - ) << "Face " << faceI << " is not an internal faces, therefore" + FatalErrorInFunction + << "Face " << faceI << " is not an internal faces, therefore" << " it cannot be removed. Check faceLabels argument." << abort(FatalError); } @@ -446,18 +419,8 @@ void Foam::removeFaces::setRefinement mesh_.write(); // Write data for debugging - FatalErrorIn - ( - "template" - "\n" - "void removeFaces::setRefinement" - "\n(" - "\n const labelList& faceLabels," - "\n const labelList& cellRegion," - "\n const labelList& cellRegionMaster," - "\n TopoChangeEngine& ref" - "\n)" - ) << "Problem : edge has too few face neighbours:" + FatalErrorInFunction + << "Problem : edge has too few face neighbours:" << eFaces << endl << "edge:" << edgeI << " vertices:" << e @@ -586,18 +549,8 @@ void Foam::removeFaces::setRefinement const edge& e = mesh_.edges()[edgeI]; // Only found one boundary face. Problem - FatalErrorIn - ( - "template" - "\n" - "void removeFaces::setRefinement" - "\n(" - "\n const labelList& faceLabels," - "\n const labelList& cellRegion," - "\n const labelList& cellRegionMaster," - "\n TopoChangeEngine& ref" - "\n)" - ) << "Problem : edge would have one boundary face" + FatalErrorInFunction + << "Problem : edge would have one boundary face" << " and one internal face using it." << endl << "Your remove pattern is probably incorrect." << endl << "edge:" << edgeI @@ -620,18 +573,7 @@ void Foam::removeFaces::setRefinement { const edge& e = mesh_.edges()[edgeI]; - FatalErrorIn - ( - "template" - "\n" - "void removeFaces::setRefinement" - "\n(" - "\n const labelList& faceLabels," - "\n const labelList& cellRegion," - "\n const labelList& cellRegionMaster," - "\n TopoChangeEngine& ref" - "\n)" - ) + FatalErrorInFunction << "Problem : edge would get 1 face using it only" << " edge:" << edgeI << " nFaces:" << nFacesPerEdge[edgeI] @@ -779,18 +721,8 @@ void Foam::removeFaces::setRefinement if (nRegion < 1) { - FatalErrorIn - ( - "template" - "\n" - "void removeFaces::setRefinement" - "\n(" - "\n const labelList& faceLabels," - "\n const labelList& cellRegion," - "\n const labelList& cellRegionMaster," - "\n TopoChangeEngine& ref" - "\n)" - ) << "Problem with region number." << abort(FatalError); + FatalErrorInFunction + << "Problem with region number." << abort(FatalError); } else if (nRegion == 1) { @@ -837,18 +769,8 @@ void Foam::removeFaces::setRefinement { if (nbrRegion != myRegion) { - FatalErrorIn - ( - "template" - "\n" - "void removeFaces::setRefinement" - "\n(" - "\n const labelList& faceLabels," - "\n const labelList& cellRegion," - "\n const labelList& cellRegionMaster," - "\n TopoChangeEngine& ref" - "\n)" - ) << "Inconsistent face region across coupled patches." + FatalErrorInFunction + << "Inconsistent face region across coupled patches." << endl << "This side has for faceI:" << faceI << " region:" << myRegion << endl @@ -868,18 +790,8 @@ void Foam::removeFaces::setRefinement // Second visit of this region if (toNbrRegion[myRegion] != nbrRegion) { - FatalErrorIn - ( - "template" - "\n" - "void removeFaces::setRefinement" - "\n(" - "\n const labelList& faceLabels," - "\n const labelList& cellRegion," - "\n const labelList& cellRegionMaster," - "\n TopoChangeEngine& ref" - "\n)" - ) << "Inconsistent face region across coupled patches." + FatalErrorInFunction + << "Inconsistent face region across coupled patches." << endl << "This side has for faceI:" << faceI << " region:" << myRegion @@ -938,18 +850,8 @@ void Foam::removeFaces::setRefinement { if (nEdgesPerPoint[pointI] == 1) { - FatalErrorIn - ( - "template" - "\n" - "void removeFaces::setRefinement" - "\n(" - "\n const labelList& faceLabels," - "\n const labelList& cellRegion," - "\n const labelList& cellRegionMaster," - "\n TopoChangeEngine& ref" - "\n)" - ) << "Problem : point would get 1 edge using it only." + FatalErrorInFunction + << "Problem : point would get 1 edge using it only." << " pointI:" << pointI << " coord:" << mesh_.points()[pointI] << abort(FatalError); @@ -1060,11 +962,43 @@ void Foam::removeFaces::setRefinement ref.setAction(polyRemovePoint(pointI, -1)); } + // Add master cells for correct mapping + forAll (cellRegionMaster, regionI) + { + // Note: it is legal to have cellRegionMaster = -1 if the region + // has been created and them marged into another region. + // Such masters will also have pointRegionMaster = -1 and should + // be ignored. HJ, 6/Sep/2019 + + // Additionally protect for old directTopoChangers which do not + // identify points for mapping. Non-existent pointRegionMaster + // is rejected + if (cellRegionMaster[regionI] > -1 && pointRegionMaster[regionI] > -1) + { + // Add master cell from master point for correct mapping + cellRegionMaster[regionI] = + ref.setAction + ( + polyAddCell + ( + pointRegionMaster[regionI], // masterPointID + -1, // masterEdgeID + -1, // masterFaceID + -1, // masterCellID + mesh_.cellZones().whichZone(cellRegionMaster[regionI]) + ) + ); + } + } + // Remove cells forAll(cellRegion, cellI) { label region = cellRegion[cellI]; + // Old check is acceptable: for mapping from point, the cellRegionMaster + // has been replaced in polyAddPoint + // HJ, 6/Sep/2019 if (region != -1 && (cellI != cellRegionMaster[region])) { ref.setAction(polyRemoveCell(cellI, cellRegionMaster[region])); @@ -1087,18 +1021,8 @@ void Foam::removeFaces::setRefinement if (rFaces.size() <= 1) { - FatalErrorIn - ( - "template" - "\n" - "void removeFaces::setRefinement" - "\n(" - "\n const labelList& faceLabels," - "\n const labelList& cellRegion," - "\n const labelList& cellRegionMaster," - "\n TopoChangeEngine& ref" - "\n)" - ) << "Region: " << regionI + FatalErrorInFunction + << "Region: " << regionI << " contains only these faces: " << rFaces << abort(FatalError); } diff --git a/src/dynamicMesh/dynamicMesh/directTopoChange/meshCut/modifiers/undoableMeshCutter/undoableMeshCutter.C b/src/dynamicMesh/dynamicMesh/directTopoChange/meshCut/modifiers/undoableMeshCutter/undoableMeshCutter.C index 1be3b18e1..8df525248 100644 --- a/src/dynamicMesh/dynamicMesh/directTopoChange/meshCut/modifiers/undoableMeshCutter/undoableMeshCutter.C +++ b/src/dynamicMesh/dynamicMesh/directTopoChange/meshCut/modifiers/undoableMeshCutter/undoableMeshCutter.C @@ -634,10 +634,16 @@ Foam::labelList Foam::undoableMeshCutter::removeSplitFaces // Insert all commands to combine cells. Never fails so don't have to // test for success. + + // Added dummy pointRegionMaster, which will be ignored. + // HJ, 6/Sep/2019 + labelList pointRegionMaster(cellRegionMaster.size(), label(-1)); + faceRemover().setRefinement ( facesToRemove, cellRegion, + pointRegionMaster, cellRegionMaster, meshMod ); diff --git a/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/polyhedralRefinement/polyhedralRefinement.C b/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/polyhedralRefinement/polyhedralRefinement.C index 64862038f..4443c612b 100644 --- a/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/polyhedralRefinement/polyhedralRefinement.C +++ b/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/polyhedralRefinement/polyhedralRefinement.C @@ -1260,10 +1260,8 @@ void Foam::polyhedralRefinement::setRefinementInstruction // cell points const labelList& cPoints = meshCellPoints[cellI]; - FatalErrorIn - ( - "polyhedralRefinement::setRefinementInstruction(...)" - ) << "Cell " << cellI + FatalErrorInFunction + << "Cell " << cellI << " of level " << cellLevel_[cellI] << " does not seem to have enough points of " << " lower level" << endl @@ -1967,11 +1965,73 @@ void Foam::polyhedralRefinement::setUnrefinementInstruction } } + // Find point region master for every cell region. This is the central point + // from which the coarse cell will be made + // The property of the point region master is that all cells that touch it + // have the same cell region index + // HJ, 6/Sep/2019 + labelList pointRegionMaster(cellRegionMaster.size(), label(-1)); + + // Get point-cell addressing + const labelListList& pc = mesh_.pointCells(); + + forAll (splitPointsToUnrefine_, i) + { + const labelList& curPc = pc[splitPointsToUnrefine_[i]]; + + label curRegion = -1; + + forAll (curPc, curPcI) + { + if (curRegion == -1) + { + // First region found. Grab it + curRegion = cellRegion[curPc[curPcI]]; + } + else + { + // Region already found. Check that all other cells that + // touch this point have the same region + if (curRegion != cellRegion[curPc[curPcI]]) + { + // Error: different region cells touching in split point + // This is not a valid unrefinement pattern + FatalErrorIn + ( + "polyhedralRefinement::setUnrefinementInstruction" + "(polyTopoChange& ref)" + ) << "Different region cells touching in split point." + << abort(FatalError); + } + } + } + + // Record point region master + if (curRegion > -1) + { + pointRegionMaster[curRegion] = splitPointsToUnrefine_[i]; + } + else + { + // Error: Cannot find region for point + FatalErrorIn + ( + "polyhedralRefinement::setUnrefinementInstruction" + "(polyTopoChange& ref)" + ) << "Different region cells touching in split point." + << abort(FatalError); + } + } + + Info<< "pointRegionMaster: " << pointRegionMaster << nl + << " cellRegionMaster: " << cellRegionMaster << endl; + // Insert all commands to combine cells faceRemover_.setRefinement ( facesToRemove, cellRegion, + pointRegionMaster, cellRegionMaster, ref ); diff --git a/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/prismatic2DRefinement/prismatic2DRefinement.C b/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/prismatic2DRefinement/prismatic2DRefinement.C index e6e3677d5..fb928849a 100644 --- a/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/prismatic2DRefinement/prismatic2DRefinement.C +++ b/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/prismatic2DRefinement/prismatic2DRefinement.C @@ -2329,7 +2329,7 @@ void Foam::prismatic2DRefinement::setRefinementInstruction if (minPointI != labelMax && minPointI != mesh_.nPoints()) { - FatalErrorIn("prismatic2DRefinement::setRefinementInstruction(...)") + FatalErrorInFunction << "Added point labels not consecutive to existing mesh points." << nl << "mesh_.nPoints():" << mesh_.nPoints() @@ -2352,10 +2352,8 @@ void Foam::prismatic2DRefinement::setUnrefinementInstruction // Check whether the refinementLevelIndicator is valid if (refinementLevelIndicator_.size() != mesh_.nCells()) { - FatalErrorIn - ( - "prismatic2DRefinement::setUnrefinementInstruction(...)" - ) << "Refinement level indicator list has invalid size: " + FatalErrorInFunction + << "Refinement level indicator list has invalid size: " << refinementLevelIndicator_.size() << ", number of cells: " << mesh_.nCells() << nl @@ -2379,11 +2377,8 @@ void Foam::prismatic2DRefinement::setUnrefinementInstruction { if (cellLevel_[cellI] < 0) { - FatalErrorIn - ( - "prismatic2DRefinement::setUnrefinementInstruction" - "(polyTopoChange& ref)" - ) << "Illegal cell level " << cellLevel_[cellI] + FatalErrorInFunction + << "Illegal cell level " << cellLevel_[cellI] << " for cell " << cellI << abort(FatalError); } @@ -2526,11 +2521,70 @@ void Foam::prismatic2DRefinement::setUnrefinementInstruction } } + // Find point region master for every cell region. This is the central point + // from which the coarse cell will be made + // The property of the point region master is that all cells that touch it + // have the same cell region index + // HJ, 6/Sep/2019 + labelList pointRegionMaster(cellRegionMaster.size(), label(-1)); + + // Get point-cell addressing + const labelListList& pc = mesh_.pointCells(); + + forAll (splitPointsToUnrefine_, i) + { + const labelList& curPc = pc[splitPointsToUnrefine_[i]]; + + label curRegion = -1; + + forAll (curPc, curPcI) + { + if (curRegion == -1) + { + // First region found. Grab it + curRegion = cellRegion[curPc[curPcI]]; + } + else + { + // Region already found. Check that all other cells that + // touch this point have the same region + if (curRegion != cellRegion[curPc[curPcI]]) + { + // Error: different region cells touching in split point + // This is not a valid unrefinement pattern + FatalErrorIn + ( + "polyhedralRefinement::setUnrefinementInstruction" + "(polyTopoChange& ref)" + ) << "Different region cells touching in split point." + << abort(FatalError); + } + } + } + + // Record point region master + if (curRegion > -1) + { + pointRegionMaster[curRegion] = splitPointsToUnrefine_[i]; + } + else + { + // Error: Cannot find region for point + FatalErrorIn + ( + "polyhedralRefinement::setUnrefinementInstruction" + "(polyTopoChange& ref)" + ) << "Different region cells touching in split point." + << abort(FatalError); + } + } + // Insert all commands to combine cells faceRemover_.setRefinement ( facesToRemove, cellRegion, + pointRegionMaster, cellRegionMaster, ref ); diff --git a/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/refinement/refinement.C b/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/refinement/refinement.C index 43f709ef5..742417278 100644 --- a/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/refinement/refinement.C +++ b/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/refinement/refinement.C @@ -23,6 +23,7 @@ License Author Vuko Vukcevic, Wikki Ltd. All rights reserved. + Hrvoje Jasak, Wikki Ltd. Notes Generalisation of hexRef8 for polyhedral cells and refactorisation into mesh @@ -423,6 +424,54 @@ Foam::label Foam::refinement::countAnchors return nAnchors; } +void Foam::refinement::adjustRefLevel +( + label& curNewCellLevel, + const label oldCellI +) +{ + if (oldCellI == -1) + { + // This cell is inflated (does not originate from other cell), set + // cell level to -1 + curNewCellLevel = -1; + } + else + { + // This cell has either been added based on another cell or it + // hasn't changed. Update new cell level according to refinement + // level indicator and old cell level + + // Get refinement status of the old cell + const label& refStatus = refinementLevelIndicator_[oldCellI]; + + if (refStatus == UNREFINED) + { + // New cell has been obtained by unrefining other cells - this + // is the remaining "master" cell. Decrement cell level + curNewCellLevel = cellLevel_[oldCellI] - 1; + } + else if (refStatus == UNCHANGED) + { + // Cell hasn't been changed during this refinement, copy old + // cell level + curNewCellLevel = cellLevel_[oldCellI]; + } + else if (refStatus == REFINED) + { + // Cell has been refined, increment cell level + curNewCellLevel = cellLevel_[oldCellI] + 1; + } + else + { + FatalErrorInFunction + << "Invalid refinement status detected: " + << refStatus << nl + << "Old cell index: " << oldCellI << abort(FatalError); + } + } +} + void Foam::refinement::checkInternalOrientation ( @@ -1386,51 +1435,21 @@ void Foam::refinement::updateMesh(const mapPolyMesh& map) // Loop through all new cells forAll (cellMap, newCellI) { - // Get index of the corresponding old cell - const label& oldCellI = cellMap[newCellI]; + adjustRefLevel(newCellLevel[newCellI], cellMap[newCellI]); + } - if (oldCellI == -1) - { - // This cell is inflated (does not originate from other cell), set - // cell level to -1 - newCellLevel[newCellI] = -1; - } - else - { - // This cell has either been added based on another cell or it - // hasn't changed. Update new cell level according to refinement - // level indicator and old cell level + // Do cells from points: unrefinement + // Note: should other types be done as well? + // HJ, 9/Sep/2019 + const List& cellsFromPoints = map.cellsFromPointsMap(); - // Get refinement status of the old cell - const label& refStatus = refinementLevelIndicator_[oldCellI]; - - if (refStatus == UNREFINED) - { - // New cell has been obtained by unrefining other cells - this - // is the remaining "master" cell. Decrement cell level - newCellLevel[newCellI] = cellLevel_[oldCellI] - 1; - } - else if (refStatus == UNCHANGED) - { - // Cell hasn't been changed during this refinement, copy old - // cell level - newCellLevel[newCellI] = cellLevel_[oldCellI]; - } - else if (refStatus == REFINED) - { - // Cell has been refined, increment cell level - newCellLevel[newCellI] = cellLevel_[oldCellI] + 1; - } - else - { - FatalErrorIn - ( - "refinement::updateMesh(const mapPolyMesh& map)" - ) << "Invalid refinement status detected: " << refStatus << nl - << "Old cell index: " << oldCellI << nl - << "New cell index: " << newCellI << abort(FatalError); - } - } + forAll (cellsFromPoints, cfpI) + { + adjustRefLevel + ( + newCellLevel[cellsFromPoints[cfpI].index()], + cellsFromPoints[cfpI].masterObjects()[0] + ); } // Transfer the new cell level into the data member diff --git a/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/refinement/refinement.H b/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/refinement/refinement.H index f361ca707..32e769820 100644 --- a/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/refinement/refinement.H +++ b/src/dynamicMesh/dynamicMesh/polyMeshModifiers/refinement/refinement/refinement.H @@ -232,7 +232,14 @@ protected: const label anchorLevel ) const; + //- Adjust cell refinement level after topo change + void adjustRefLevel + ( + label& curNewCellLevel, + const label oldCellI + ); + // Debug functions //- Check orientation of added internal face