diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriver.H b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriver.H index 298ca934d19206356ee53c3c17851aac1d9bed40..ee94b0b0156265b28f7c57b89975407c41f49188 100644 --- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriver.H +++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriver.H @@ -510,7 +510,7 @@ private: const motionSmoother& meshMover, const label nSmoothNormals, const label nSmoothSurfaceNormals, - const scalar minMedianAxisAngleCos, + const scalar minMedialAxisAngleCos, const scalar featureAngle, pointVectorField& dispVec, diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C index a2d7a39a3a6eb5fb1680c47636489fa90271c847..e8c3d38cd6c30778ca58b80a12459cc797d48ef0 100644 --- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C +++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoLayerDriverShrink.C @@ -864,7 +864,7 @@ void Foam::autoLayerDriver::medialAxisSmoothingInfo const motionSmoother& meshMover, const label nSmoothNormals, const label nSmoothSurfaceNormals, - const scalar minMedianAxisAngleCos, + const scalar minMedialAxisAngleCos, const scalar featureAngle, pointVectorField& dispVec, @@ -1049,7 +1049,7 @@ void Foam::autoLayerDriver::medialAxisSmoothingInfo { // Unvisited point. See above about nUnvisit warning } - else if (isMaxEdge(pointWallDist, edgeI, minMedianAxisAngleCos)) + else if (isMaxEdge(pointWallDist, edgeI, minMedialAxisAngleCos)) { // Both end points of edge have very different nearest wall // point. Mark both points as medial axis points. diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C index 8ac5eadea544b715175d1cc8db3f6b12886e6a0c..b466e17de2fcb627b72eed5c0dc293f3de3ad8f8 100644 --- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C +++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C @@ -154,10 +154,6 @@ Foam::layerParameters::layerParameters ( readScalar(dict.lookup("maxThicknessToMedialRatio")) ), - minMedianAxisAngleCos_ - ( - Foam::cos(degToRad(readScalar(dict.lookup("minMedianAxisAngle")))) - ), nBufferCellsNoExtrude_ ( readLabel(dict.lookup("nBufferCellsNoExtrude")) @@ -175,6 +171,18 @@ Foam::layerParameters::layerParameters ) ) { + word angleKey = "minMedialAxisAngle"; + if (!dict.found(angleKey)) + { + // Backwards compatibility + angleKey = "minMedianAxisAngle"; + } + minMedialAxisAngleCos_ = Foam::cos + ( + degToRad(readScalar(dict.lookup(angleKey))) + ); + + // Detect layer specification mode label nSpec = 0; diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.H b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.H index 3b4d09f76f362e745428ac681ed765b272b70065..d9d2860238fff705520dd5ae45e92d14bc050bb5 100644 --- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.H +++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.H @@ -134,7 +134,7 @@ private: scalar maxThicknessToMedialRatio_; - scalar minMedianAxisAngleCos_; + scalar minMedialAxisAngleCos_; label nBufferCellsNoExtrude_; @@ -348,9 +348,9 @@ public: } //- Angle used to pick up medial axis points - scalar minMedianAxisAngleCos() const + scalar minMedialAxisAngleCos() const { - return minMedianAxisAngleCos_; + return minMedialAxisAngleCos_; } label nSnap() const diff --git a/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.C b/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.C index db2dedf89ff6928130fe3744e684306aa4dcc7e9..eb223c88414a6c294290464d790e15b4629d74a1 100644 --- a/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.C +++ b/src/mesh/autoMesh/autoHexMesh/externalDisplacementMeshMover/medialAxisMeshMover.C @@ -402,7 +402,7 @@ bool Foam::medialAxisMeshMover::isMaxEdge void Foam::medialAxisMeshMover::update(const dictionary& coeffDict) { Info<< typeName - << " : Calculate distance to Medial Axis ..." << endl; + << " : Calculating distance to Medial Axis ..." << endl; const pointField& points = mesh().points(); @@ -420,9 +420,15 @@ void Foam::medialAxisMeshMover::update(const dictionary& coeffDict) ); //- When is medial axis - const scalar minMedianAxisAngleCos = Foam::cos + word angleKey = "minMedialAxisAngle"; + if (!coeffDict.found(angleKey)) + { + // Backwards compatibility + angleKey = "minMedianAxisAngle"; + } + scalar minMedialAxisAngleCos = Foam::cos ( - degToRad(readScalar(coeffDict.lookup("minMedianAxisAngle"))) + degToRad(readScalar(coeffDict.lookup(angleKey))) ); //- Feature angle when to stop adding layers @@ -442,6 +448,13 @@ void Foam::medialAxisMeshMover::update(const dictionary& coeffDict) coeffDict.lookup("nSmoothNormals") ); + //- Number of edges walking out + const label nMedialAxisIter = coeffDict.lookupOrDefault<label> + ( + "nMedialAxisIter", + mesh().globalData().nTotalPoints() + ); + // Predetermine mesh edges // ~~~~~~~~~~~~~~~~~~~~~~~ @@ -512,10 +525,10 @@ void Foam::medialAxisMeshMover::update(const dictionary& coeffDict) wallInfo, pointWallDist, edgeWallDist, - mesh().globalData().nTotalPoints(), // max iterations + 0, // max iterations dummyTrackData ); - + wallDistCalc.iterate(nMedialAxisIter); label nUnvisit = returnReduce ( @@ -525,16 +538,27 @@ void Foam::medialAxisMeshMover::update(const dictionary& coeffDict) if (nUnvisit > 0) { - WarningIn - ( - "medialAxisMeshMover::" - "medialAxisSmoothingInfo(..)" - ) << "Walking did not visit all points." << nl - << " Did not visit " << nUnvisit - << " out of " << mesh().globalData().nTotalPoints() - << " points. This is not necessarily a problem" << nl - << " and might be due to faceZones splitting of part" - << " of the domain." << nl << endl; + if (nMedialAxisIter > 0) + { + Info<< typeName + << " : Limited walk to " << nMedialAxisIter + << " steps. Not visited " << nUnvisit + << " out of " << mesh().globalData().nTotalPoints() + << " points" << endl; + } + else + { + WarningIn + ( + "medialAxisMeshMover::" + "medialAxisSmoothingInfo(..)" + ) << "Walking did not visit all points." << nl + << " Did not visit " << nUnvisit + << " out of " << mesh().globalData().nTotalPoints() + << " points. This is not necessarily a problem" << nl + << " and might be due to faceZones splitting of part" + << " of the domain." << nl << endl; + } } } @@ -564,8 +588,29 @@ void Foam::medialAxisMeshMover::update(const dictionary& coeffDict) ) { // Unvisited point. See above about nUnvisit warning + forAll(e, ep) + { + label pointI = e[ep]; + + if (!pointMedialDist[pointI].valid(dummyTrackData)) + { + maxPoints.append(pointI); + maxInfo.append + ( + pointData + ( + points[pointI], + 0.0, + pointI, // passive data + vector::zero // passive data + ) + ); + pointMedialDist[pointI] = maxInfo.last(); + } + } + } - else if (isMaxEdge(pointWallDist, edgeI, minMedianAxisAngleCos)) + else if (isMaxEdge(pointWallDist, edgeI, minMedialAxisAngleCos)) { // Both end points of edge have very different nearest wall // point. Mark both points as medial axis points. @@ -653,7 +698,8 @@ void Foam::medialAxisMeshMover::update(const dictionary& coeffDict) if (pvf.fixesValue()) { // Disable all movement on fixedValue patchFields - Info<< "Inserting all points on patch " << pp.name() + Info<< typeName + << " : Inserting all points on patch " << pp.name() << endl; forAll(meshPoints, i) @@ -683,7 +729,8 @@ void Foam::medialAxisMeshMover::update(const dictionary& coeffDict) // normal as the passive vector. Note that this points // out of the originating wall so inside of the domain // on this patch. - Info<< "Inserting points on patch " << pp.name() + Info<< typeName + << " : Inserting points on patch " << pp.name() << " if angle to nearest layer patch > " << slipFeatureAngle << " degrees." << endl; @@ -753,15 +800,29 @@ void Foam::medialAxisMeshMover::update(const dictionary& coeffDict) pointMedialDist, edgeMedialDist, - mesh().globalData().nTotalPoints(), // max iterations + 0, // max iterations dummyTrackData ); + medialDistCalc.iterate(2*nMedialAxisIter); + // Extract medial axis distance as pointScalarField forAll(pointMedialDist, pointI) { - medialDist_[pointI] = Foam::sqrt(pointMedialDist[pointI].distSqr()); - medialVec_[pointI] = pointMedialDist[pointI].origin(); + if (pointMedialDist[pointI].valid(dummyTrackData)) + { + medialDist_[pointI] = Foam::sqrt + ( + pointMedialDist[pointI].distSqr() + ); + medialVec_[pointI] = pointMedialDist[pointI].origin(); + } + else + { + // Unvisited. Do as if on medial axis so unmoving + medialDist_[pointI] = 0.0; + medialVec_[pointI] = point(1, 0, 0); + } } } @@ -1380,7 +1441,8 @@ void Foam::medialAxisMeshMover::findIsolatedRegions } reduce(nPointCounter, sumOp<label>()); - Info<< "Number isolated points extrusion stopped : "<< nPointCounter + Info<< typeName + << " : Number of isolated points extrusion stopped : "<< nPointCounter << endl; } @@ -1632,6 +1694,13 @@ void Foam::medialAxisMeshMover::calculateDisplacement coeffDict.lookup("nSmoothThickness") ); + //- Number of edges walking out + const label nMedialAxisIter = coeffDict.lookupOrDefault<label> + ( + "nMedialAxisIter", + mesh().globalData().nTotalPoints() + ); + // Precalulate master points/edge (only relevant for shared points/edges) const PackedBoolList isMasterPoint(syncTools::getMasterPoints(mesh())); @@ -1677,7 +1746,8 @@ void Foam::medialAxisMeshMover::calculateDisplacement + ".obj" ) ); - Info<< "Writing points with too large an extrusion distance to " + Info<< typeName + << " : Writing points with too large an extrusion distance to " << str().name() << endl; } @@ -1694,8 +1764,9 @@ void Foam::medialAxisMeshMover::calculateDisplacement + ".obj" ) ); - Info<< "Writing points with too large an extrusion distance to " - << medialVecStr().name() << endl; + Info<< typeName + << " : Writing medial axis vectors on points with too large" + << " an extrusion distance to " << medialVecStr().name() << endl; } forAll(meshPoints, patchPointI) @@ -1721,7 +1792,7 @@ void Foam::medialAxisMeshMover::calculateDisplacement if (thicknessRatio > maxThicknessToMedialRatio) { // Truncate thickness. - if (debug) + if (debug&2) { Pout<< "truncating displacement at " << mesh().points()[pointI] @@ -1770,7 +1841,7 @@ void Foam::medialAxisMeshMover::calculateDisplacement } reduce(numThicknessRatioExclude, sumOp<label>()); - Info<< typeName << " : Reduce layer thickness at " + Info<< typeName << " : Reducing layer thickness at " << numThicknessRatioExclude << " nodes where thickness to medial axis distance is large " << endl; @@ -1849,9 +1920,10 @@ void Foam::medialAxisMeshMover::calculateDisplacement wallInfo, pointWallDist, edgeWallDist, - mesh().globalData().nTotalPoints(), // max iterations + 0, // max iterations dummyTrackData ); + wallDistCalc.iterate(nMedialAxisIter); } @@ -1917,10 +1989,13 @@ bool Foam::medialAxisMeshMover::shrinkMesh for (label iter = 0; iter < 2*nSnap ; iter++) { - Info<< "Iteration " << iter << endl; + Info<< typeName + << " : Iteration " << iter << endl; if (iter == nSnap) { - Info<< "Displacement scaling for error reduction set to 0." << endl; + Info<< typeName + << " : Displacement scaling for error reduction set to 0." + << endl; oldErrorReduction = meshMover_.setErrorReduction(0.0); } @@ -1937,7 +2012,7 @@ bool Foam::medialAxisMeshMover::shrinkMesh ) ) { - Info<< typeName << " :" << " Successfully moved mesh" << endl; + Info<< typeName << " : Successfully moved mesh" << endl; meshOk = true; break; }