diff --git a/src/sampling/surface/distanceSurface/distanceSurface.C b/src/sampling/surface/distanceSurface/distanceSurface.C index e957d6fccc2f52cbccaf616f9786ac6470ae9765..90e934b72809371aa0ed59fac7db42f9d7ac9f71 100644 --- a/src/sampling/surface/distanceSurface/distanceSurface.C +++ b/src/sampling/surface/distanceSurface/distanceSurface.C @@ -29,7 +29,6 @@ License #include "volPointInterpolation.H" #include "addToRunTimeSelectionTable.H" #include "fvMesh.H" -#include "volumeType.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -134,8 +133,6 @@ void Foam::distanceSurface::createGeometry() const fvMesh& fvm = static_cast<const fvMesh&>(mesh_); - const labelList& own = fvm.faceOwner(); - // Distance to cell centres // ~~~~~~~~~~~~~~~~~~~~~~~~ @@ -173,35 +170,14 @@ void Foam::distanceSurface::createGeometry() if (signed_) { - List<volumeType> volType; - - surfPtr_().getVolumeType(cc, volType); + vectorField norms; + surfPtr_().getNormal(nearest, norms); - forAll(volType, i) + forAll(norms, i) { - volumeType vT = volType[i]; + const point diff(cc[i] - nearest[i].hitPoint()); - if (vT == volumeType::OUTSIDE) - { - fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint()); - } - else if (vT == volumeType::INSIDE) - { - fld[i] = -Foam::mag(cc[i] - nearest[i].hitPoint()); - } - else if (vT == volumeType::UNKNOWN) - { - // Treat as very far outside - fld[i] = GREAT; - } - else - { - FatalErrorInFunction - << "getVolumeType failure:" - << " neither INSIDE or OUTSIDE but " - << volumeType::names[vT] - << exit(FatalError); - } + fld[i] = sign(diff & norms[i]) * Foam::mag(diff); } } else @@ -223,9 +199,6 @@ void Foam::distanceSurface::createGeometry() const pointField& cc = fvm.C().boundaryField()[patchi]; fvPatchScalarField& fld = cellDistanceBf[patchi]; - const label patchStarti = fvm.boundaryMesh()[patchi].start(); - - List<pointIndexHit> nearest; surfPtr_().findNearest ( @@ -236,41 +209,14 @@ void Foam::distanceSurface::createGeometry() if (signed_) { - List<volumeType> volType; - - surfPtr_().getVolumeType(cc, volType); + vectorField norms; + surfPtr_().getNormal(nearest, norms); - forAll(volType, i) + forAll(norms, i) { - volumeType vT = volType[i]; - - if (vT == volumeType::OUTSIDE) - { - fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint()); - } - else if (vT == volumeType::INSIDE) - { - fld[i] = -Foam::mag(cc[i] - nearest[i].hitPoint()); - } - else if (vT == volumeType::UNKNOWN) - { - // Nothing known, so use the cell value. - // - this avoids spurious changes on the boundary - - // The cell value - const label meshFacei = i+patchStarti; - const scalar& cellVal = cellDistance[own[meshFacei]]; - - fld[i] = cellVal; - } - else - { - FatalErrorInFunction - << "getVolumeType failure:" - << " neither INSIDE or OUTSIDE but " - << volumeType::names[vT] - << exit(FatalError); - } + const point diff(cc[i] - nearest[i].hitPoint()); + + fld[i] = sign(diff & norms[i]) * Foam::mag(diff); } } else @@ -303,44 +249,21 @@ void Foam::distanceSurface::createGeometry() if (signed_) { - List<volumeType> volType; + vectorField norms; + surfPtr_().getNormal(nearest, norms); - surfPtr_().getVolumeType(pts, volType); - - forAll(volType, i) + forAll(norms, i) { - volumeType vT = volType[i]; + const point diff(pts[i] - nearest[i].hitPoint()); - if (vT == volumeType::OUTSIDE) - { - pointDistance_[i] = - Foam::mag(pts[i] - nearest[i].hitPoint()); - } - else if (vT == volumeType::INSIDE) - { - pointDistance_[i] = - -Foam::mag(pts[i] - nearest[i].hitPoint()); - } - else if (vT == volumeType::UNKNOWN) - { - // Treat as very far outside - pointDistance_[i] = GREAT; - } - else - { - FatalErrorInFunction - << "getVolumeType failure:" - << " neither INSIDE or OUTSIDE but " - << volumeType::names[vT] - << exit(FatalError); - } + pointDistance_[i] = sign(diff & norms[i]) * Foam::mag(diff); } } else { forAll(nearest, i) { - pointDistance_[i] = Foam::mag(pts[i]-nearest[i].hitPoint()); + pointDistance_[i] = Foam::mag(pts[i] - nearest[i].hitPoint()); } } } diff --git a/src/sampling/surface/isoSurface/isoSurface.C b/src/sampling/surface/isoSurface/isoSurface.C index adaa14655562038a9ed3a095662ceba00c3f6864..1a2f36b07ecec3b4692139cf33353d1dca2fbeae 100644 --- a/src/sampling/surface/isoSurface/isoSurface.C +++ b/src/sampling/surface/isoSurface/isoSurface.C @@ -60,14 +60,6 @@ namespace Foam } }; - - // Avoid detecting change if the cells have been marked as GREAT - // (ie, ignore them) - static inline constexpr bool ignoreValue(const scalar val) - { - return (val >= 0.5*Foam::GREAT); - } - } // End namespace Foam @@ -165,7 +157,7 @@ void Foam::isoSurface::syncUnseparatedPoints forAll(nbrPts, pointi) { - label nbrPointi = nbrPts[pointi]; + const label nbrPointi = nbrPts[pointi]; patchInfo[nbrPointi] = pointValues[meshPts[pointi]]; } @@ -314,39 +306,19 @@ bool Foam::isoSurface::isEdgeOfFaceCut const bool neiLower ) const { - // Could also count number of edges cut and return when they are > 1 - // but doesn't appear to improve anything - forAll(f, fp) { - const scalar& pt0Value = pVals[f[fp]]; - - if (ignoreValue(pt0Value)) - { - continue; - } - - const bool fpLower = (pt0Value < iso_); + const bool fpLower = (pVals[f[fp]] < iso_); - if (fpLower != ownLower || fpLower != neiLower) + if + ( + fpLower != ownLower + || fpLower != neiLower + || fpLower != (pVals[f[f.fcIndex(fp)]] < iso_) + ) { - // ++ncut; return true; } - else - { - const scalar& pt1Value = pVals[f[f.fcIndex(fp)]]; - - if (!ignoreValue(pt1Value) && (fpLower != (pt1Value < iso_))) - { - // ++ncut; - return true; - } - } - // if (ncut > 1) - // { - // return true; - // } } return false; @@ -401,17 +373,9 @@ void Foam::isoSurface::calcCutTypes faceCutType_.setSize(mesh_.nFaces()); faceCutType_ = NOTCUT; - // Avoid detecting change if the cells have been marked as GREAT - // (ie, ignore them) - for (label facei = 0; facei < mesh_.nInternalFaces(); ++facei) { const scalar& ownValue = cVals[own[facei]]; - if (ignoreValue(ownValue)) - { - continue; - } - const bool ownLower = (ownValue < iso_); scalar nbrValue; @@ -427,11 +391,6 @@ void Foam::isoSurface::calcCutTypes nbrPoint ); - if (ignoreValue(nbrValue)) - { - continue; - } - const bool neiLower = (nbrValue < iso_); if (ownLower != neiLower) @@ -503,7 +462,6 @@ void Foam::isoSurface::calcCutTypes // Propagate internal face cuts into the cells. - // For cells marked as ignore (eg, GREAT) - skip this. for (label facei = 0; facei < mesh_.nInternalFaces(); ++facei) { @@ -512,20 +470,12 @@ void Foam::isoSurface::calcCutTypes continue; } - if - ( - cellCutType_[own[facei]] == NOTCUT - && !ignoreValue(cVals[own[facei]]) - ) + if (cellCutType_[own[facei]] == NOTCUT) { cellCutType_[own[facei]] = CUT; ++nCutCells_; } - if - ( - cellCutType_[nei[facei]] == NOTCUT - && !ignoreValue(cVals[nei[facei]]) - ) + if (cellCutType_[nei[facei]] == NOTCUT) { cellCutType_[nei[facei]] = CUT; ++nCutCells_; @@ -534,8 +484,6 @@ void Foam::isoSurface::calcCutTypes // Propagate boundary face cuts into the cells. - // For cells marked as ignore (eg, GREAT) - skip this and - // also suppress the boundary face cut to prevent dangling face cuts. for (label facei = mesh_.nInternalFaces(); facei < mesh_.nFaces(); ++facei) { @@ -544,12 +492,7 @@ void Foam::isoSurface::calcCutTypes continue; } - if (ignoreValue(cVals[own[facei]])) - { - // Suppress dangling boundary face cut - faceCutType_[facei] = NOTCUT; - } - else if (cellCutType_[own[facei]] == NOTCUT) + if (cellCutType_[own[facei]] == NOTCUT) { cellCutType_[own[facei]] = CUT; ++nCutCells_; @@ -774,10 +717,8 @@ void Foam::isoSurface::calcSnappedPoint bool anyCut = false; - forAll(pFaces, i) + for (const label facei : pFaces) { - label facei = pFaces[i]; - if (faceCutType_[facei] == CUT) { anyCut = true; @@ -795,12 +736,10 @@ void Foam::isoSurface::calcSnappedPoint label nOther = 0; point otherPointSum = Zero; - forAll(pFaces, pFacei) + for (const label facei : pFaces) { // Create points for all intersections close to point // (i.e. from pyramid edges) - - label facei = pFaces[pFacei]; const face& f = mesh_.faces()[facei]; label own = mesh_.faceOwner()[facei]; diff --git a/src/sampling/surface/isoSurface/isoSurfaceCell.C b/src/sampling/surface/isoSurface/isoSurfaceCell.C index e0b1d6bbe1dfea186e11122196ea60acfbafd202..b110d1552ccca45d63ecc4620ee187ed6003aab8 100644 --- a/src/sampling/surface/isoSurface/isoSurfaceCell.C +++ b/src/sampling/surface/isoSurface/isoSurfaceCell.C @@ -44,20 +44,6 @@ namespace Foam } -// * * * * * * * * * * * * * * * Local Functions * * * * * * * * * * * * * * // - -namespace Foam -{ - // Avoid detecting change if the cells have been marked as GREAT - // (ie, ignore them) - static inline constexpr bool ignoreValue(const scalar val) - { - return (val >= 0.5*Foam::GREAT); - } - -} // End namespace Foam - - // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // Foam::scalar Foam::isoSurfaceCell::isoFraction @@ -99,11 +85,6 @@ Foam::isoSurfaceCell::cellCutType Foam::isoSurfaceCell::calcCutType const label celli ) const { - if (ignoreValue(cellValues[celli])) - { - return NOTCUT; - } - const cell& cFaces = mesh_.cells()[celli]; if (isTet.test(celli)) @@ -137,11 +118,7 @@ Foam::isoSurfaceCell::cellCutType Foam::isoSurfaceCell::calcCutType // Check pyramids cut for (const label labi : f) { - if - ( - !ignoreValue(pointValues[labi]) - && cellLower != (pointValues[labi] < iso_) - ) + if (cellLower != (pointValues[labi] < iso_)) { edgeCut = true; break; @@ -187,11 +164,7 @@ Foam::isoSurfaceCell::cellCutType Foam::isoSurfaceCell::calcCutType for (const label pointi : cPoints) { - if - ( - !ignoreValue(pointValues[pointi]) - && cellLower != (pointValues[pointi] < iso_) - ) + if (cellLower != (pointValues[pointi] < iso_)) { ++nCuts; } @@ -201,7 +174,7 @@ Foam::isoSurfaceCell::cellCutType Foam::isoSurfaceCell::calcCutType { return SPHERE; } - else if (nCuts > 1) + else { return CUT; }