diff --git a/src/finiteVolume/fvMesh/wallDist/wallDistAddressing/wallDistAddressing.C b/src/finiteVolume/fvMesh/wallDist/wallDistAddressing/wallDistAddressing.C index dad43a7d9fb11ec9f56a3fb37eccdfaf53e536a4..4f2189723cc630c868cec3b6bb57a5584a6986f2 100644 --- a/src/finiteVolume/fvMesh/wallDist/wallDistAddressing/wallDistAddressing.C +++ b/src/finiteVolume/fvMesh/wallDist/wallDistAddressing/wallDistAddressing.C @@ -234,16 +234,25 @@ void Foam::wallDistAddressing::correct(volScalarField& y) if (correctWalls_) { cellToWallFace.reserve(nWalls); - correctBoundaryFaceCells - ( - patchSet, - y, - cellToWallFace - ); - correctBoundaryPointCells + //correctBoundaryFaceCells + //( + // patchSet, + // y, + // cellToWallFace + //); + //correctBoundaryPointCells + //( + // patchSet, + // y, + // cellToWallFace + //); + + // Correct across multiple patches + correctBoundaryCells ( - patchSet, + patchIDs_, + true, // do point-connected cells as well y, cellToWallFace ); diff --git a/src/meshTools/cellDist/cellDistFuncs.C b/src/meshTools/cellDist/cellDistFuncs.C index 9a3a1b00d60a92038e68264e576d569fa3f32ecc..15f276ce846514cb3332d579c6a113c361a2aa31 100644 --- a/src/meshTools/cellDist/cellDistFuncs.C +++ b/src/meshTools/cellDist/cellDistFuncs.C @@ -29,6 +29,7 @@ License #include "cellDistFuncs.H" #include "polyMesh.H" #include "polyBoundaryMesh.H" +#include "uindirectPrimitivePatch.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -56,134 +57,6 @@ Foam::labelHashSet Foam::cellDistFuncs::getPatchIDs } -// Return smallest true distance from p to any of wallFaces. -// Note that even if normal hits face we still check other faces. -// Note that wallFaces is untruncated and we explicitly pass in size. -Foam::scalar Foam::cellDistFuncs::smallestDist -( - const point& p, - const polyPatch& patch, - const labelUList& wallFaces, - label& minFacei -) const -{ - const pointField& points = patch.points(); - - scalar minDist = GREAT; - minFacei = -1; - - for (const label patchFacei : wallFaces) - { - const pointHit curHit = patch[patchFacei].nearestPoint(p, points); - - if (curHit.distance() < minDist) - { - minDist = curHit.distance(); - minFacei = patch.start() + patchFacei; - } - } - - return minDist; -} - - -// Get point neighbours of facei (including facei). Returns number of faces. -// Note: does not allocate storage but does use linear search to determine -// uniqueness. For polygonal faces this might be quite inefficient. -void Foam::cellDistFuncs::getPointNeighbours -( - const primitivePatch& patch, - const label patchFacei, - DynamicList<label>& neighbours -) const -{ - neighbours.clear(); - - // Add myself - neighbours.append(patchFacei); - - // Add all face neighbours - const labelList& faceNeighbours = patch.faceFaces()[patchFacei]; - - for (const label nbr : faceNeighbours) - { - neighbours.push_uniq(nbr); - } - - // Add all point-only neighbours by linear searching in edge neighbours. - // Assumes that point-only neighbours are not using multiple points on - // face. - - const face& f = patch.localFaces()[patchFacei]; - - forAll(f, fp) - { - label pointi = f[fp]; - - const labelList& pointNbs = patch.pointFaces()[pointi]; - - for (const label facei : pointNbs) - { - // Check for facei in edge-neighbours part of neighbours - neighbours.push_uniq(facei); - } - } - - - if (debug) - { - // Check for duplicates - - // Use hashSet to determine nbs. - labelHashSet nbs(4*f.size()); - - forAll(f, fp) - { - const labelList& pointNbs = patch.pointFaces()[f[fp]]; - nbs.insert(pointNbs); - } - - // Subtract ours. - for (const label nb : neighbours) - { - if (!nbs.found(nb)) - { - SeriousErrorInFunction - << "getPointNeighbours : patchFacei:" << patchFacei - << " verts:" << f << endl; - - forAll(f, fp) - { - SeriousErrorInFunction - << "point:" << f[fp] << " pointFaces:" - << patch.pointFaces()[f[fp]] << endl; - } - - for (const label facei : neighbours) - { - SeriousErrorInFunction - << "fast nbr:" << facei - << endl; - } - - FatalErrorInFunction - << "Problem: fast pointNeighbours routine included " << nb - << " which is not in proper neighbour list " << nbs.toc() - << abort(FatalError); - } - nbs.erase(nb); - } - - if (nbs.size()) - { - FatalErrorInFunction - << "Problem: fast pointNeighbours routine did not find " - << nbs.toc() << abort(FatalError); - } - } -} - - // size of largest patch (out of supplied subset of patches) Foam::label Foam::cellDistFuncs::maxPatchSize ( @@ -275,7 +148,7 @@ void Foam::cellDistFuncs::correctBoundaryFaceCells ); // Store wallCell and its nearest neighbour - nearestFace.insert(celli, minFacei); + nearestFace.insert(celli, patch.start()+minFacei); } } } @@ -348,7 +221,7 @@ void Foam::cellDistFuncs::correctBoundaryPointCells ); // Store wallCell and its nearest neighbour - nearestFace.insert(celli, minFacei); + nearestFace.insert(celli, patch.start()+minFacei); } } } @@ -357,4 +230,142 @@ void Foam::cellDistFuncs::correctBoundaryPointCells } +void Foam::cellDistFuncs::correctBoundaryCells +( + const labelList& patchIDs, + const bool doPointCells, + scalarField& wallDistCorrected, + Map<label>& nearestFace +) const +{ + label nWalls = 0; + { + for (const label patchi : patchIDs) + { + nWalls += mesh().boundaryMesh()[patchi].size(); + } + } + + + DynamicList<label> faceLabels(nWalls); + { + for (const label patchi : patchIDs) + { + const auto& patch = mesh().boundaryMesh()[patchi]; + forAll(patch, i) + { + faceLabels.append(patch.start()+i); + } + } + } + + const uindirectPrimitivePatch wallPatch + ( + UIndirectList<face>(mesh().faces(), faceLabels), + mesh_.points() + ); + + + // Correct all cells with face on wall + const vectorField& cellCentres = mesh().cellCentres(); + + DynamicList<label> neighbours; + + nWalls = 0; + for (const label patchi : patchIDs) + { + const auto& patch = mesh().boundaryMesh()[patchi]; + const auto areaFraction(patch.areaFraction()); + const labelUList& faceCells = patch.faceCells(); + + // Check cells with face on wall + forAll(patch, patchFacei) + { + if (areaFraction && (areaFraction()[patchFacei] <= 0.5)) + { + // For cyclicACMI: more cyclic than wall + } + else + { + getPointNeighbours(wallPatch, nWalls, neighbours); + + const label celli = faceCells[patchFacei]; + + label minFacei = -1; + wallDistCorrected[celli] = smallestDist + ( + cellCentres[celli], + wallPatch, + neighbours, + minFacei + ); + + // Store wallCell and its nearest neighbour + nearestFace.insert(celli, nWalls+minFacei); + } + + nWalls++; + } + } + + // Correct all cells with a point on the wall + if (doPointCells) + { + const auto& meshPoints = wallPatch.meshPoints(); + const auto& localFaces = wallPatch.localFaces(); + + bitSet isWallPoint(meshPoints.size(), true); + + nWalls = 0; + for (const label patchi : patchIDs) + { + const auto& patch = mesh().boundaryMesh()[patchi]; + const auto areaFraction(patch.areaFraction()); + + // Check cells with face on wall + forAll(patch, patchFacei) + { + if (areaFraction && (areaFraction()[patchFacei] <= 0.5)) + { + // For cyclicACMI: more cyclic than wall + isWallPoint.unset(localFaces[nWalls]); + } + + nWalls++; + } + } + + const auto& pointFaces = wallPatch.pointFaces(); + + for (const label patchPointi : isWallPoint) + { + const label verti = meshPoints[patchPointi]; + + const labelList& neighbours = mesh().pointCells(verti); + + for (const label celli : neighbours) + { + if (!nearestFace.found(celli)) + { + const labelList& wallFaces = pointFaces[patchPointi]; + + label minFacei = -1; + + wallDistCorrected[celli] = smallestDist + ( + cellCentres[celli], + wallPatch, + wallFaces, + minFacei + ); + + // Store wallCell and its nearest neighbour + nearestFace.insert(celli, wallPatch.addressing()[minFacei]); + } + } + } + } +} + + // ************************************************************************* // diff --git a/src/meshTools/cellDist/cellDistFuncs.H b/src/meshTools/cellDist/cellDistFuncs.H index 6935a1c93d31dfd15e118ed3f1d7997b1c2b1572..2c9955a88b26ef0d618082267b6a40802bd949b9 100644 --- a/src/meshTools/cellDist/cellDistFuncs.H +++ b/src/meshTools/cellDist/cellDistFuncs.H @@ -103,22 +103,22 @@ public: template<class Type> labelHashSet getPatchIDs() const; - //- Calculate smallest true distance (and face index) + //- Calculate smallest true distance (and patch face index) // from pt to faces wallFaces. - // For efficiency reasons we still pass in patch instead of extracting - // it from mesh_ + template<class PatchType> scalar smallestDist ( const point& p, - const polyPatch& patch, + const PatchType& patch, const labelUList& wallFaces, - label& meshFacei + label& patchFacei ) const; //- Get faces sharing point with face on patch + template<class PatchType> void getPointNeighbours ( - const primitivePatch&, + const PatchType&, const label patchFacei, DynamicList<label>& ) const; @@ -147,6 +147,17 @@ public: scalarField& wallDistCorrected, Map<label>& nearestFace ) const; + + //- Correct all cells connected to any of the patches in patchIDs. Sets + // - cell values in wallDistCorrected + // - (mesh) face that contains the nearest point + void correctBoundaryCells + ( + const labelList& patchIDs, + const bool doPointCells, + scalarField& wallDistCorrected, + Map<label>& nearestFace + ) const; }; diff --git a/src/meshTools/cellDist/cellDistFuncsTemplates.C b/src/meshTools/cellDist/cellDistFuncsTemplates.C index 1b74f101540becfcc4d6a69fd5f63ecd4efdce24..86e225077713f68032a5e075a90574f9f6fd13a0 100644 --- a/src/meshTools/cellDist/cellDistFuncsTemplates.C +++ b/src/meshTools/cellDist/cellDistFuncsTemplates.C @@ -5,7 +5,7 @@ \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- - Copyright (C) 2011-2016 OpenFOAM Foundation + Copyright (C) 2011-2016,2024 OpenFOAM Foundation ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -50,4 +50,135 @@ Foam::labelHashSet Foam::cellDistFuncs::getPatchIDs() const } +template<class PatchType> +Foam::scalar Foam::cellDistFuncs::smallestDist +( + const point& p, + const PatchType& patch, + const labelUList& wallFaces, + label& minFacei +) const +{ + // Return smallest true distance from p to any of wallFaces. + // Note that even if normal hits face we still check other faces. + + const pointField& points = patch.points(); + + scalar minDist = GREAT; + minFacei = -1; + + for (const label patchFacei : wallFaces) + { + const pointHit curHit = patch[patchFacei].nearestPoint(p, points); + + if (curHit.distance() < minDist) + { + minDist = curHit.distance(); + minFacei = patchFacei; + } + } + + return minDist; +} + + +template<class PatchType> +void Foam::cellDistFuncs::getPointNeighbours +( + const PatchType& patch, + const label patchFacei, + DynamicList<label>& neighbours +) const +{ + // Get point neighbours of facei (including facei). Returns number of faces. + // Note: does not allocate storage but does use linear search to determine + // uniqueness. For polygonal faces this might be quite inefficient. + + neighbours.clear(); + + // Add myself + neighbours.append(patchFacei); + + // Add all face neighbours + const labelList& faceNeighbours = patch.faceFaces()[patchFacei]; + + for (const label nbr : faceNeighbours) + { + neighbours.push_uniq(nbr); + } + + // Add all point-only neighbours by linear searching in edge neighbours. + // Assumes that point-only neighbours are not using multiple points on + // face. + + const face& f = patch.localFaces()[patchFacei]; + + forAll(f, fp) + { + label pointi = f[fp]; + + const labelList& pointNbs = patch.pointFaces()[pointi]; + + for (const label facei : pointNbs) + { + // Check for facei in edge-neighbours part of neighbours + neighbours.push_uniq(facei); + } + } + + + if (debug) + { + // Check for duplicates + + // Use hashSet to determine nbs. + labelHashSet nbs(4*f.size()); + + forAll(f, fp) + { + const labelList& pointNbs = patch.pointFaces()[f[fp]]; + nbs.insert(pointNbs); + } + + // Subtract ours. + for (const label nb : neighbours) + { + if (!nbs.found(nb)) + { + SeriousErrorInFunction + << "getPointNeighbours : patchFacei:" << patchFacei + << " verts:" << f << endl; + + forAll(f, fp) + { + SeriousErrorInFunction + << "point:" << f[fp] << " pointFaces:" + << patch.pointFaces()[f[fp]] << endl; + } + + for (const label facei : neighbours) + { + SeriousErrorInFunction + << "fast nbr:" << facei + << endl; + } + + FatalErrorInFunction + << "Problem: fast pointNeighbours routine included " << nb + << " which is not in proper neighbour list " << nbs.toc() + << abort(FatalError); + } + nbs.erase(nb); + } + + if (nbs.size()) + { + FatalErrorInFunction + << "Problem: fast pointNeighbours routine did not find " + << nbs.toc() << abort(FatalError); + } + } +} + + // ************************************************************************* // diff --git a/src/meshTools/cellDist/patchWave/patchDataWave.C b/src/meshTools/cellDist/patchWave/patchDataWave.C index 9f943225c3e90476870f2db550b9813bca889549..2c9070995cacbcf2c0fee17eaefe7b57d2d37149 100644 --- a/src/meshTools/cellDist/patchWave/patchDataWave.C +++ b/src/meshTools/cellDist/patchWave/patchDataWave.C @@ -264,16 +264,25 @@ void Foam::patchDataWave<TransferType, TrackingData>::correct() Map<label> nearestFace(2 * nWalls); // Get distance and indices of nearest face - correctBoundaryFaceCells + //correctBoundaryFaceCells + //( + // patchIDs_, + // distance_, + // nearestFace + //); + + //correctBoundaryPointCells + //( + // patchIDs_, + // distance_, + // nearestFace + //); + + // Correct across multiple patches + correctBoundaryCells ( - patchIDs_, - distance_, - nearestFace - ); - - correctBoundaryPointCells - ( - patchIDs_, + patchIDs_.sortedToc(), + true, // do point-connected cells as well distance_, nearestFace ); diff --git a/src/meshTools/cellDist/patchWave/patchWave.C b/src/meshTools/cellDist/patchWave/patchWave.C index deeb062d40757a08fa4876427a7411017b68c8db..c5b5e89b109e7dc397503e09217c1eae0c036678 100644 --- a/src/meshTools/cellDist/patchWave/patchWave.C +++ b/src/meshTools/cellDist/patchWave/patchWave.C @@ -204,16 +204,25 @@ void Foam::patchWave::correct() { Map<label> nearestFace(2*nPatch); - correctBoundaryFaceCells + //correctBoundaryFaceCells + //( + // patchIDs_, + // distance_, + // nearestFace + //); + + //correctBoundaryPointCells + //( + // patchIDs_, + // distance_, + // nearestFace + //); + + // Correct across multiple patches + correctBoundaryCells ( - patchIDs_, - distance_, - nearestFace - ); - - correctBoundaryPointCells - ( - patchIDs_, + patchIDs_.sortedToc(), + true, // do point-connected cells as well distance_, nearestFace );