From 42ff732d580768c1ee017e0af70d82107f1d2a30 Mon Sep 17 00:00:00 2001 From: Will Bainbridge Date: Fri, 8 Dec 2023 16:34:18 +0000 Subject: [PATCH] fvPatch: Added polyFaceFraction calculation This returns the fraction of the poly-face that each fv-face in the patch covers. This is useful for mapping some calculations from the polyMesh to the fvMesh in cases with non-conformal couplings. --- .../epsilonWallFunctionFvPatchScalarField.C | 25 ++++++------------- .../epsilonWallFunctionFvPatchScalarField.H | 2 +- .../omegaWallFunctionFvPatchScalarField.C | 25 ++++++------------- .../omegaWallFunctionFvPatchScalarField.H | 2 +- .../epsilonmWallFunctionFvPatchScalarField.C | 15 +++-------- .../epsilonmWallFunctionFvPatchScalarField.H | 2 +- .../fvMesh/fvPatches/fvPatch/fvPatch.C | 16 +++++++++++- .../fvMesh/fvPatches/fvPatch/fvPatch.H | 5 ++++ 8 files changed, 40 insertions(+), 52 deletions(-) diff --git a/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C b/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C index 649b12d4f2..3944e9d58a 100644 --- a/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C +++ b/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.C @@ -31,7 +31,7 @@ License // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // -Foam::scalar Foam::epsilonWallFunctionFvPatchScalarField::tolerance_ = 1e-1; +Foam::scalar Foam::epsilonWallFunctionFvPatchScalarField::tol_ = 1e-1; // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * // @@ -361,12 +361,10 @@ void Foam::epsilonWallFunctionFvPatchScalarField::updateCoeffs() internalField() ); - scalarField weights(patch().magSf()/patch().patch().magFaceAreas()); - forAll(weights, facei) - { - scalar& w = weights[facei]; - w = w <= tolerance_ ? 0 : (w - tolerance_)/(1 - tolerance_); - } + const scalarField weights + ( + max((patch().polyFaceFraction() - tol_)/(1 - tol_), scalar(0)) + ); forAll(weights, facei) { @@ -393,20 +391,11 @@ void Foam::epsilonWallFunctionFvPatchScalarField::manipulateMatrix return; } - const scalarField& epsilon0 = this->epsilon(); - - scalarField weights(patch().magSf()/patch().patch().magFaceAreas()); - forAll(weights, facei) - { - scalar& w = weights[facei]; - w = w <= tolerance_ ? 0 : (w - tolerance_)/(1 - tolerance_); - } - matrix.setValues ( patch().faceCells(), - UIndirectList(epsilon0, patch().faceCells()), - weights + UIndirectList(this->epsilon(), patch().faceCells()), + max((patch().polyFaceFraction() - tol_)/(1 - tol_), scalar(0)) ); fvPatchField::manipulateMatrix(matrix); diff --git a/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.H b/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.H index eabb47ff13..3e63896c95 100644 --- a/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.H +++ b/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonWallFunction/epsilonWallFunctionFvPatchScalarField.H @@ -86,7 +86,7 @@ protected: // Protected data //- Tolerance used in weighted calculations - static scalar tolerance_; + static scalar tol_; //- Local copy of turbulence G field scalarField G_; diff --git a/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C b/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C index 1b78609d96..4454640758 100644 --- a/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C +++ b/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.C @@ -36,7 +36,7 @@ namespace Foam // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // -scalar omegaWallFunctionFvPatchScalarField::tolerance_ = 1e-1; +scalar omegaWallFunctionFvPatchScalarField::tol_ = 1e-1; // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // @@ -402,12 +402,10 @@ void omegaWallFunctionFvPatchScalarField::updateCoeffs() internalField() ); - scalarField weights(patch().magSf()/patch().patch().magFaceAreas()); - forAll(weights, facei) - { - scalar& w = weights[facei]; - w = w <= tolerance_ ? 0 : (w - tolerance_)/(1 - tolerance_); - } + const scalarField weights + ( + max((patch().polyFaceFraction() - tol_)/(1 - tol_), scalar(0)) + ); forAll(weights, facei) { @@ -434,20 +432,11 @@ void omegaWallFunctionFvPatchScalarField::manipulateMatrix return; } - const scalarField& omega0 = this->omega(); - - scalarField weights(patch().magSf()/patch().patch().magFaceAreas()); - forAll(weights, facei) - { - scalar& w = weights[facei]; - w = w <= tolerance_ ? 0 : (w - tolerance_)/(1 - tolerance_); - } - matrix.setValues ( patch().faceCells(), - UIndirectList(omega0, patch().faceCells()), - weights + UIndirectList(this->omega(), patch().faceCells()), + max((patch().polyFaceFraction() - tol_)/(1 - tol_), scalar(0)) ); fvPatchField::manipulateMatrix(matrix); diff --git a/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.H b/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.H index f6f85aa7a6..6ccd2c612c 100644 --- a/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.H +++ b/src/MomentumTransportModels/momentumTransportModels/derivedFvPatchFields/wallFunctions/omegaWallFunctions/omegaWallFunction/omegaWallFunctionFvPatchScalarField.H @@ -86,7 +86,7 @@ protected: // Protected data //- Tolerance used in weighted calculations - static scalar tolerance_; + static scalar tol_; //- beta1 coefficient scalar beta1_; diff --git a/src/MomentumTransportModels/phaseCompressible/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonmWallFunction/epsilonmWallFunctionFvPatchScalarField.C b/src/MomentumTransportModels/phaseCompressible/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonmWallFunction/epsilonmWallFunctionFvPatchScalarField.C index b52033c44d..315e97b914 100644 --- a/src/MomentumTransportModels/phaseCompressible/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonmWallFunction/epsilonmWallFunctionFvPatchScalarField.C +++ b/src/MomentumTransportModels/phaseCompressible/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonmWallFunction/epsilonmWallFunctionFvPatchScalarField.C @@ -29,7 +29,7 @@ License // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // -Foam::scalar Foam::epsilonmWallFunctionFvPatchScalarField::tolerance_ = 1e-1; +Foam::scalar Foam::epsilonmWallFunctionFvPatchScalarField::tol_ = 1e-1; // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // @@ -85,20 +85,11 @@ void Foam::epsilonmWallFunctionFvPatchScalarField::manipulateMatrix return; } - const DimensionedField& epsilon = internalField(); - - scalarField weights(patch().magSf()/patch().patch().magFaceAreas()); - forAll(weights, facei) - { - scalar& w = weights[facei]; - w = w <= tolerance_ ? 0 : (w - tolerance_)/(1 - tolerance_); - } - matrix.setValues ( patch().faceCells(), - UIndirectList(epsilon, patch().faceCells()), - weights + UIndirectList(internalField(), patch().faceCells()), + max((patch().polyFaceFraction() - tol_)/(1 - tol_), scalar(0)) ); fvPatchField::manipulateMatrix(matrix); diff --git a/src/MomentumTransportModels/phaseCompressible/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonmWallFunction/epsilonmWallFunctionFvPatchScalarField.H b/src/MomentumTransportModels/phaseCompressible/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonmWallFunction/epsilonmWallFunctionFvPatchScalarField.H index b6150478f4..7c215e057a 100644 --- a/src/MomentumTransportModels/phaseCompressible/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonmWallFunction/epsilonmWallFunctionFvPatchScalarField.H +++ b/src/MomentumTransportModels/phaseCompressible/derivedFvPatchFields/wallFunctions/epsilonWallFunctions/epsilonmWallFunction/epsilonmWallFunctionFvPatchScalarField.H @@ -54,7 +54,7 @@ class epsilonmWallFunctionFvPatchScalarField // Private member data //- Tolerance used in weighted calculations - static scalar tolerance_; + static scalar tol_; public: diff --git a/src/finiteVolume/fvMesh/fvPatches/fvPatch/fvPatch.C b/src/finiteVolume/fvMesh/fvPatches/fvPatch/fvPatch.C index 7bc02132e3..debae5b2eb 100644 --- a/src/finiteVolume/fvMesh/fvPatches/fvPatch/fvPatch.C +++ b/src/finiteVolume/fvMesh/fvPatches/fvPatch/fvPatch.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -153,6 +153,20 @@ Foam::tmp Foam::fvPatch::delta() const } +Foam::tmp Foam::fvPatch::polyFaceFraction() const +{ + return + boundaryMesh().mesh().conformal() + ? tmp(new scalarField(size(), scalar(1))) + : magSf() + /scalarField + ( + boundaryMesh().mesh().magFaceAreas(), + boundaryMesh().mesh().polyFacesBf()[patch().index()] + ); +} + + void Foam::fvPatch::makeWeights(scalarField& w) const { w = 1.0; diff --git a/src/finiteVolume/fvMesh/fvPatches/fvPatch/fvPatch.H b/src/finiteVolume/fvMesh/fvPatches/fvPatch/fvPatch.H index 24c6a05b63..9aa666c97d 100644 --- a/src/finiteVolume/fvMesh/fvPatches/fvPatch/fvPatch.H +++ b/src/finiteVolume/fvMesh/fvPatches/fvPatch/fvPatch.H @@ -221,6 +221,11 @@ public: // to coupled-cell-centre vector is returned virtual tmp delta() const; + //- Return the fraction of the poly-face that each fv-face in this + // patch covers. Will be equal to one, unless the mesh is + // non-conformal. + tmp polyFaceFraction() const; + // Access functions for demand driven data