From 056bb4f0a0354d45d14423de1f0266d3d89a8bf9 Mon Sep 17 00:00:00 2001 From: Mattijs Janssens <ext-mjanssens@esi-group.com> Date: Thu, 30 May 2024 15:15:56 +0000 Subject: [PATCH] FaceCellWave: travel through coupled ACMI. See #3139 --- .../polyPatches/polyPatch/polyPatch.C | 18 ++++++- .../polyPatches/polyPatch/polyPatch.H | 16 ++++-- ...ureForceBaffleVelocityFvPatchVectorField.C | 30 +++++++++++- .../wallDistAddressing/wallDistAddressing.C | 37 +++++++++----- .../inverseDistanceDiffusivity.C | 1 - .../inverseVolume/inverseVolumeDiffusivity.C | 1 - .../basic/InteractionLists/InteractionLists.C | 6 +-- .../cyclicACMIPolyPatch/cyclicACMIPolyPatch.C | 29 +++++++++-- .../algorithms/MeshWave/FaceCellWave.C | 11 ++++- src/meshTools/cellDist/cellDistFuncs.C | 49 +++++++++++++++---- .../cellDist/patchWave/patchDataWave.C | 46 ++++++++++------- .../cellDist/patchWave/patchDataWave.H | 4 +- src/meshTools/cellDist/patchWave/patchWave.C | 46 ++++++++--------- src/meshTools/cellDist/patchWave/patchWave.H | 4 +- 14 files changed, 212 insertions(+), 86 deletions(-) diff --git a/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.C b/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.C index 3f2530d6da8..bc132c24abe 100644 --- a/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.C +++ b/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.C @@ -350,13 +350,15 @@ Foam::tmp<Foam::vectorField> Foam::polyPatch::faceCellCentres() const } -Foam::tmp<Foam::scalarField> Foam::polyPatch::areaFraction() const +Foam::tmp<Foam::scalarField> Foam::polyPatch::areaFraction +( + const pointField& points +) const { auto tfraction = tmp<scalarField>::New(size()); auto& fraction = tfraction.ref(); const vectorField::subField faceAreas = this->faceAreas(); - const pointField& points = this->points(); forAll(*this, facei) { @@ -369,6 +371,18 @@ Foam::tmp<Foam::scalarField> Foam::polyPatch::areaFraction() const } +const Foam::tmp<Foam::scalarField>& Foam::polyPatch::areaFraction() const +{ + return areaFraction_; +} + + +void Foam::polyPatch::areaFraction(const tmp<scalarField>& fraction) +{ + areaFraction_ = fraction; +} + + const Foam::labelUList& Foam::polyPatch::faceCells() const { if (!faceCellsPtr_) diff --git a/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.H b/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.H index c6761791975..49d2786bc07 100644 --- a/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.H +++ b/src/OpenFOAM/meshes/polyMesh/polyPatches/polyPatch/polyPatch.H @@ -88,6 +88,9 @@ class polyPatch //- Demand-driven: global edge addressing mutable std::unique_ptr<labelList> mePtr_; + //- Cached area fraction + tmp<scalarField> areaFraction_; + protected: @@ -441,9 +444,16 @@ public: //- Return face cell centres tmp<vectorField> faceCellCentres() const; - //- Return the area fraction as the ratio of the stored face area - //- and the area given by the face points - tmp<scalarField> areaFraction() const; + //- Calculate the area fraction as the ratio of the stored face + // area and the area given by the face points. + tmp<scalarField> areaFraction(const pointField& points) const; + + //- Return the cached area fraction. Usually only set for the + //- non-overlap patches on ACMI. + const tmp<scalarField>& areaFraction() const; + + //- Set cached area fraction + void areaFraction(const tmp<scalarField>&); // Addressing into mesh diff --git a/src/finiteVolume/fields/fvPatchFields/derived/activePressureForceBaffleVelocity/activePressureForceBaffleVelocityFvPatchVectorField.C b/src/finiteVolume/fields/fvPatchFields/derived/activePressureForceBaffleVelocity/activePressureForceBaffleVelocityFvPatchVectorField.C index 6f2a9897680..cb80bf364a9 100644 --- a/src/finiteVolume/fields/fvPatchFields/derived/activePressureForceBaffleVelocity/activePressureForceBaffleVelocityFvPatchVectorField.C +++ b/src/finiteVolume/fields/fvPatchFields/derived/activePressureForceBaffleVelocity/activePressureForceBaffleVelocityFvPatchVectorField.C @@ -342,6 +342,7 @@ void Foam::activePressureForceBaffleVelocityFvPatchVectorField::updateCoeffs() ).neighbFvPatch().Sf(); } + // Update this wall patch vectorField::subField Sfw = patch().patch().faceAreas(); vectorField newSfw((1 - areaFraction)*initWallSf_); @@ -350,16 +351,41 @@ void Foam::activePressureForceBaffleVelocityFvPatchVectorField::updateCoeffs() Sfw[facei] = newSfw[facei]; } const_cast<scalarField&>(patch().magSf()) = mag(patch().Sf()); + // Cache fraction + const_cast<polyPatch&>(patch().patch()).areaFraction + ( + tmp<scalarField>::New + ( + patch().patch().size(), + (1-areaFraction) + ) + ); + // Update owner side of cyclic const_cast<vectorField&>(cyclicPatch.Sf()) = areaFraction*initCyclicSf_; - const_cast<scalarField&>(cyclicPatch.magSf()) = mag(cyclicPatch.Sf()); + const_cast<polyPatch&>(cyclicPatch.patch()).areaFraction + ( + tmp<scalarField>::New + ( + cyclicPatch.patch().size(), + areaFraction + ) + ); + // Update neighbour side of cyclic const_cast<vectorField&>(nbrPatch.Sf()) = areaFraction*nbrCyclicSf_; - const_cast<scalarField&>(nbrPatch.magSf()) = mag(nbrPatch.Sf()); + const_cast<polyPatch&>(cyclicPatch.patch()).areaFraction + ( + tmp<scalarField>::New + ( + nbrPatch.patch().size(), + areaFraction + ) + ); curTimeIndex_ = this->db().time().timeIndex(); } diff --git a/src/finiteVolume/fvMesh/wallDist/wallDistAddressing/wallDistAddressing.C b/src/finiteVolume/fvMesh/wallDist/wallDistAddressing/wallDistAddressing.C index 653079e98e3..67a20853af9 100644 --- a/src/finiteVolume/fvMesh/wallDist/wallDistAddressing/wallDistAddressing.C +++ b/src/finiteVolume/fvMesh/wallDist/wallDistAddressing/wallDistAddressing.C @@ -5,7 +5,7 @@ \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- - Copyright (C) 2023 OpenCFD Ltd. + Copyright (C) 2023-2024 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -146,29 +146,40 @@ void Foam::wallDistAddressing::correct(volScalarField& y) globalWallsPtr_.reset(new globalIndex(nWalls)); const globalIndex& globalWalls = globalWallsPtr_(); - labelList seedFaces(nWalls); - List<wallPointAddressing> seedInfo(nWalls); + DynamicList<label> seedFaces(nWalls); + DynamicList<wallPointAddressing> seedInfo(nWalls); nWalls = 0; for (const label patchi : patchIDs_) { const auto& fc = C.boundaryField()[patchi]; + const auto& areaFraction = fc.patch().patch().areaFraction(); forAll(fc, patchFacei) { - seedFaces[nWalls] = fc.patch().start()+patchFacei; - seedInfo[nWalls] = wallPointAddressing + if ( - fc[patchFacei], - gt.encode + !areaFraction + || (areaFraction()[patchFacei] > 0.5) // mostly wall + ) + { + seedFaces.append(fc.patch().start()+patchFacei); + seedInfo.append ( - Pstream::myProcNo(), - nWalls, - gt.nullTransformIndex() - ), - scalar(0.0) - ); + wallPointAddressing + ( + fc[patchFacei], + gt.encode + ( + Pstream::myProcNo(), + nWalls, + gt.nullTransformIndex() + ), + scalar(0.0) + ) + ); + } nWalls++; } } diff --git a/src/fvMotionSolver/motionDiffusivity/inverseDistance/inverseDistanceDiffusivity.C b/src/fvMotionSolver/motionDiffusivity/inverseDistance/inverseDistanceDiffusivity.C index f023d3c46a7..6eea2f76757 100644 --- a/src/fvMotionSolver/motionDiffusivity/inverseDistance/inverseDistanceDiffusivity.C +++ b/src/fvMotionSolver/motionDiffusivity/inverseDistance/inverseDistanceDiffusivity.C @@ -28,7 +28,6 @@ License #include "inverseDistanceDiffusivity.H" #include "addToRunTimeSelectionTable.H" -#include "patchWave.H" #include "HashSet.H" #include "surfaceInterpolate.H" #include "zeroGradientFvPatchFields.H" diff --git a/src/fvMotionSolver/motionDiffusivity/inverseVolume/inverseVolumeDiffusivity.C b/src/fvMotionSolver/motionDiffusivity/inverseVolume/inverseVolumeDiffusivity.C index ff0582a84e4..49a6c7694f7 100644 --- a/src/fvMotionSolver/motionDiffusivity/inverseVolume/inverseVolumeDiffusivity.C +++ b/src/fvMotionSolver/motionDiffusivity/inverseVolume/inverseVolumeDiffusivity.C @@ -27,7 +27,6 @@ License #include "inverseVolumeDiffusivity.H" #include "addToRunTimeSelectionTable.H" -#include "patchWave.H" #include "HashSet.H" #include "surfaceInterpolate.H" #include "zeroGradientFvPatchFields.H" diff --git a/src/lagrangian/basic/InteractionLists/InteractionLists.C b/src/lagrangian/basic/InteractionLists/InteractionLists.C index 8bbd96b67b6..edae174aa47 100644 --- a/src/lagrangian/basic/InteractionLists/InteractionLists.C +++ b/src/lagrangian/basic/InteractionLists/InteractionLists.C @@ -299,11 +299,11 @@ void Foam::InteractionLists<ParticleType>::buildInteractionLists() { if (isA<wallPolyPatch>(patch)) { - const scalarField areaFraction(patch.areaFraction()); + const auto& areaFraction = patch.areaFraction(); - forAll(areaFraction, facei) + forAll(patch, facei) { - if (areaFraction[facei] > 0.5) + if (!areaFraction || areaFraction()[facei] > 0.5) { localWallFaces.append(facei + patch.start()); } diff --git a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C index 4a7452e93a0..404f19f6716 100644 --- a/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C +++ b/src/meshTools/AMIInterpolation/patches/cyclicACMI/cyclicACMIPolyPatch/cyclicACMIPolyPatch.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2013-2016 OpenFOAM Foundation - Copyright (C) 2017-2022 OpenCFD Ltd. + Copyright (C) 2017-2022,2024 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -217,10 +217,24 @@ void Foam::cyclicACMIPolyPatch::scalePatchFaceAreas return; } - forAll(noSf, facei) { - const scalar w = min(maxTol, max(tolerance_, mask[facei])); - noSf[facei] = noFaceArea[facei]*(scalar(1) - w); + tmp<scalarField> scale + ( + scalar(1) + - min + ( + max(mask, tolerance_), + maxTol + ) + ); + + // Scale area + forAll(noSf, facei) + { + noSf[facei] = noFaceArea[facei]*scale()[facei]; + } + // Cache scaling + const_cast<polyPatch&>(nonOverlapPatch).areaFraction(scale); } if (!createAMIFaces_) @@ -235,10 +249,15 @@ void Foam::cyclicACMIPolyPatch::scalePatchFaceAreas // Scale the coupled patch face areas vectorField::subField Sf = acmipp.faceAreas(); + tmp<scalarField> scale(max(tolerance_, mask)); + + // Scale area forAll(Sf, facei) { - Sf[facei] = faceArea[facei]*max(tolerance_, mask[facei]); + Sf[facei] = faceArea[facei]*scale()[facei]; } + // Cache scaling + const_cast<cyclicACMIPolyPatch&>(acmipp).areaFraction(scale); // Re-normalise the weights since the effect of overlap is already // accounted for in the area diff --git a/src/meshTools/algorithms/MeshWave/FaceCellWave.C b/src/meshTools/algorithms/MeshWave/FaceCellWave.C index a5487c290ab..b96792829bb 100644 --- a/src/meshTools/algorithms/MeshWave/FaceCellWave.C +++ b/src/meshTools/algorithms/MeshWave/FaceCellWave.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2017 OpenFOAM Foundation - Copyright (C) 2018-2023 OpenCFD Ltd. + Copyright (C) 2018-2024 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -821,8 +821,17 @@ void Foam::FaceCellWave<Type, TrackingData>::handleAMICyclicPatches() } // Merge into global storage + + const auto& areaFraction = patch.areaFraction(); + forAll(receiveInfo, i) { + if (areaFraction && areaFraction()[i] <= 0.5) + { + // not coupled + continue; + } + const label meshFacei = cycPatch.start()+i; const Type& newInfo = receiveInfo[i]; diff --git a/src/meshTools/cellDist/cellDistFuncs.C b/src/meshTools/cellDist/cellDistFuncs.C index 301354bc601..1b064dc391c 100644 --- a/src/meshTools/cellDist/cellDistFuncs.C +++ b/src/meshTools/cellDist/cellDistFuncs.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2020 OpenCFD Ltd. + Copyright (C) 2020,2024 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -235,6 +235,8 @@ void Foam::cellDistFuncs::correctBoundaryFaceCells Map<label>& nearestFace ) const { + const auto& pbm = mesh().boundaryMesh(); + // Size neighbours array for maximum possible (= size of largest patch) DynamicList<label> neighbours(maxPatchSize(patchIDs)); @@ -242,15 +244,22 @@ void Foam::cellDistFuncs::correctBoundaryFaceCells const vectorField& cellCentres = mesh().cellCentres(); const labelList& faceOwner = mesh().faceOwner(); - forAll(mesh().boundaryMesh(), patchi) + forAll(pbm, patchi) { if (patchIDs.found(patchi)) { - const polyPatch& patch = mesh().boundaryMesh()[patchi]; + const polyPatch& patch = pbm[patchi]; + const auto& areaFraction = patch.areaFraction(); // Check cells with face on wall forAll(patch, patchFacei) { + if (areaFraction && (areaFraction()[patchFacei] <= 0.5)) + { + // is mostly coupled + continue; + } + getPointNeighbours(patch, patchFacei, neighbours); label celli = faceOwner[patch.start() + patchFacei]; @@ -283,20 +292,42 @@ void Foam::cellDistFuncs::correctBoundaryPointCells { // Correct all (non-visited) cells with point on wall + const auto& pbm = mesh().boundaryMesh(); const vectorField& cellCentres = mesh().cellCentres(); - forAll(mesh().boundaryMesh(), patchi) + forAll(pbm, patchi) { if (patchIDs.found(patchi)) { - const polyPatch& patch = mesh().boundaryMesh()[patchi]; - + const polyPatch& patch = pbm[patchi]; + const auto& localFaces = patch.localFaces(); const labelList& meshPoints = patch.meshPoints(); const labelListList& pointFaces = patch.pointFaces(); - forAll(meshPoints, meshPointi) + bitSet isWallPoint(meshPoints.size(), true); { - const label vertI = meshPoints[meshPointi]; + const auto& areaFraction = patch.areaFraction(); + + // Check cells with face on wall + forAll(patch, patchFacei) + { + if (!areaFraction || (areaFraction()[patchFacei] <= 0.5)) + { + // face mostly coupled + isWallPoint.unset(localFaces[patchFacei]); + } + } + } + + + forAll(meshPoints, patchPointi) + { + const label vertI = meshPoints[patchPointi]; + + if (!isWallPoint[patchPointi]) + { + continue; + } const labelList& neighbours = mesh().pointCells(vertI); @@ -304,7 +335,7 @@ void Foam::cellDistFuncs::correctBoundaryPointCells { if (!nearestFace.found(celli)) { - const labelList& wallFaces = pointFaces[meshPointi]; + const labelList& wallFaces = pointFaces[patchPointi]; label minFacei = -1; diff --git a/src/meshTools/cellDist/patchWave/patchDataWave.C b/src/meshTools/cellDist/patchWave/patchDataWave.C index dfa48df725e..293018e8d91 100644 --- a/src/meshTools/cellDist/patchWave/patchDataWave.C +++ b/src/meshTools/cellDist/patchWave/patchDataWave.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2020 OpenCFD Ltd. + Copyright (C) 2020,2024 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -41,37 +41,45 @@ template<class TransferType, class TrackingData> void Foam::patchDataWave<TransferType, TrackingData>::setChangedFaces ( const labelHashSet& patchIDs, - labelList& changedFaces, - List<TransferType>& faceDist + DynamicList<label>& changedFaces, + DynamicList<TransferType>& faceDist ) const { const polyMesh& mesh = cellDistFuncs::mesh(); - label nChangedFaces = 0; - forAll(mesh.boundaryMesh(), patchi) { if (patchIDs.found(patchi)) { const polyPatch& patch = mesh.boundaryMesh()[patchi]; - const Field<Type>& patchField = initialPatchValuePtrs_[patchi]; + const auto& areaFraction = patch.areaFraction(); - forAll(patch.faceCentres(), patchFacei) - { - label meshFacei = patch.start() + patchFacei; + const auto faceCentres(patch.faceCentres()); - changedFaces[nChangedFaces] = meshFacei; + const Field<Type>& patchField = initialPatchValuePtrs_[patchi]; - faceDist[nChangedFaces] = - TransferType + forAll(faceCentres, patchFacei) + { + if + ( + !areaFraction + || (areaFraction()[patchFacei] > 0.5) // mostly wall + ) + { + label meshFacei = patch.start() + patchFacei; + changedFaces.append(meshFacei); + + faceDist.append ( - patch.faceCentres()[patchFacei], - patchField[patchFacei], - 0.0 + TransferType + ( + faceCentres[patchFacei], + patchField[patchFacei], + 0.0 + ) ); - - nChangedFaces++; + } } } } @@ -222,8 +230,8 @@ void Foam::patchDataWave<TransferType, TrackingData>::correct() // Count walls label nWalls = sumPatchSize(patchIDs_); - List<TransferType> faceDist(nWalls); - labelList changedFaces(nWalls); + DynamicList<TransferType> faceDist(nWalls); + DynamicList<label> changedFaces(nWalls); setChangedFaces(patchIDs_, changedFaces, faceDist); diff --git a/src/meshTools/cellDist/patchWave/patchDataWave.H b/src/meshTools/cellDist/patchWave/patchDataWave.H index f7f79707e57..764efcd46b8 100644 --- a/src/meshTools/cellDist/patchWave/patchDataWave.H +++ b/src/meshTools/cellDist/patchWave/patchDataWave.H @@ -115,8 +115,8 @@ private: void setChangedFaces ( const labelHashSet& patchIDs, - labelList&, - List<TransferType>& + DynamicList<label>& changedFaces, + DynamicList<TransferType>& faceDist ) const; //- Copy MeshWave values into *this diff --git a/src/meshTools/cellDist/patchWave/patchWave.C b/src/meshTools/cellDist/patchWave/patchWave.C index 6f07859038f..e8d583fb238 100644 --- a/src/meshTools/cellDist/patchWave/patchWave.C +++ b/src/meshTools/cellDist/patchWave/patchWave.C @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation + Copyright (C) 2024 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -35,50 +36,49 @@ License void Foam::patchWave::setChangedFaces ( const labelHashSet& patchIDs, - labelList& changedFaces, - List<wallPoint>& faceDist + DynamicList<label>& changedFaces, + DynamicList<wallPoint>& faceDist ) const { const polyMesh& mesh = cellDistFuncs::mesh(); - label nChangedFaces = 0; - forAll(mesh.boundaryMesh(), patchi) { if (patchIDs.found(patchi)) { const polyPatch& patch = mesh.boundaryMesh()[patchi]; - forAll(patch.faceCentres(), patchFacei) - { - label meshFacei = patch.start() + patchFacei; - - changedFaces[nChangedFaces] = meshFacei; + const auto& areaFraction = patch.areaFraction(); - faceDist[nChangedFaces] = - wallPoint - ( - patch.faceCentres()[patchFacei], - 0.0 - ); + const auto faceCentres(patch.faceCentres()); - nChangedFaces++; + forAll(faceCentres, patchFacei) + { + if + ( + !areaFraction + || (areaFraction()[patchFacei] > 0.5) // mostly wall + ) + { + changedFaces.append(patch.start() + patchFacei); + faceDist.append(wallPoint(faceCentres[patchFacei], 0.0)); + } } } } for (const label facei : sourceIDs_) { - changedFaces[nChangedFaces] = facei; + changedFaces.append(facei); - faceDist[nChangedFaces] = + faceDist.append + ( wallPoint ( mesh.faceCentres()[facei], 0.0 - ); - - nChangedFaces++; + ) + ); } } @@ -181,8 +181,8 @@ void Foam::patchWave::correct() label nPatch = sumPatchSize(patchIDs_) + sourceIDs_.size(); - List<wallPoint> faceDist(nPatch); - labelList changedFaces(nPatch); + DynamicList<wallPoint> faceDist(nPatch); + DynamicList<label> changedFaces(nPatch); // Set to faceDist information to facecentre on walls. setChangedFaces(patchIDs_, changedFaces, faceDist); diff --git a/src/meshTools/cellDist/patchWave/patchWave.H b/src/meshTools/cellDist/patchWave/patchWave.H index d7259d78a18..a99f05aa564 100644 --- a/src/meshTools/cellDist/patchWave/patchWave.H +++ b/src/meshTools/cellDist/patchWave/patchWave.H @@ -89,8 +89,8 @@ class patchWave void setChangedFaces ( const labelHashSet& patchIDs, - labelList& changedFaces, - List<wallPoint>& changedInfo + DynamicList<label>& changedFaces, + DynamicList<wallPoint>& faceDist ) const; //- Copy MeshWave cell values. Return number of illegal/unset -- GitLab