From 07a9ee86f323a4491ff2fb281d6de6bfd2d4befc Mon Sep 17 00:00:00 2001 From: Andrew Heather <> Date: Wed, 10 Aug 2022 10:59:39 +0100 Subject: [PATCH] ENH: Code refactoring --- .../SpalartAllmarasDDES/SpalartAllmarasDDES.C | 51 ++++++++--------- .../SpalartAllmarasDDES/SpalartAllmarasDDES.H | 6 +- .../SpalartAllmarasDES/SpalartAllmarasDES.C | 50 +++++++++++++++++ .../SpalartAllmarasDES/SpalartAllmarasDES.H | 14 ++++- .../DES/kOmegaSSTDDES/kOmegaSSTDDES.C | 39 ++++--------- .../DES/kOmegaSSTDDES/kOmegaSSTDDES.H | 4 +- .../DES/kOmegaSSTDES/kOmegaSSTDES.C | 55 +++++++++++++++++++ .../DES/kOmegaSSTDES/kOmegaSSTDES.H | 20 +++++++ 8 files changed, 173 insertions(+), 66 deletions(-) diff --git a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.C b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.C index b744117066..4c75b3c72d 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.C +++ b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.C @@ -61,7 +61,27 @@ tmp SpalartAllmarasDDES::Stilda const volScalarField& dTilda ) const { - tmp St = + if (this->useSigma_) + { + const volScalarField& lRAS(this->y_); + const volScalarField fv2(this->fv2(chi, fv1)); + const volScalarField lLES(this->lengthScaleLES(chi, fv1)); + const volScalarField Omega(this->Omega(gradU)); + const volScalarField Ssigma(this->Ssigma(gradU)); + const volScalarField SsigmaDES + ( + Omega - fd(mag(gradU))*pos(lRAS - lLES)*(Omega - Ssigma) + ); + + return + max + ( + SsigmaDES + fv2*this->nuTilda_/sqr(this->kappa_*dTilda), + this->Cs_*SsigmaDES + ); + } + + return SpalartAllmarasBase>::Stilda ( chi, @@ -69,23 +89,6 @@ tmp SpalartAllmarasDDES::Stilda gradU, dTilda ); - - if (useSigma_) - { - const volScalarField& lRAS(this->y_); - const volScalarField lLES(this->lengthScaleLES(chi, fv1)); - const volScalarField Omega(this->Omega(gradU)); - const volScalarField Ssigma(this->Ssigma(gradU)); - - return - max - ( - St - fd(mag(gradU))*pos(lRAS - lLES)*(Omega - Ssigma), - this->Cs_*Omega - ); - } - - return St; } @@ -136,18 +139,9 @@ SpalartAllmarasDDES::SpalartAllmarasDDES type ), - useSigma_ - ( - Switch::getOrAddToDict - ( - "useSigma", - this->coeffDict_, - false - ) - ), Cd1_ ( - useSigma_ ? + this->useSigma_ ? dimensioned::getOrAddToDict ( "Cd1Sigma", @@ -185,7 +179,6 @@ bool SpalartAllmarasDDES::read() { if (SpalartAllmarasDES::read()) { - useSigma_.readIfPresent("useSigma", this->coeffDict()); Cd1_.readIfPresent(this->coeffDict()); Cd2_.readIfPresent(this->coeffDict()); diff --git a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.H b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.H index 5f6c652add..17d3ede0ee 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.H +++ b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.H @@ -72,7 +72,7 @@ class SpalartAllmarasDDES { // Private Member Functions - //- Delay function + //- Shielding function tmp fd(const volScalarField& magGradU) const; //- No copy construct @@ -86,9 +86,6 @@ protected: // Protected data - //- Switch to activate grey-area enhanced sigma-DDES - Switch useSigma_; - // Model coefficients dimensionedScalar Cd1_; @@ -97,6 +94,7 @@ protected: // Protected Member Functions + //- Production term virtual tmp Stilda ( const volScalarField& chi, diff --git a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.C b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.C index 3de6413d39..3d9b3d6704 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.C +++ b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.C @@ -95,6 +95,46 @@ tmp SpalartAllmarasDES::lengthScaleLES } +template +tmp SpalartAllmarasDES::Stilda +( + const volScalarField& chi, + const volScalarField& fv1, + const volTensorField& gradU, + const volScalarField& dTilda +) const +{ + if (useSigma_) + { + const volScalarField& lRAS = this->y_; + const volScalarField fv2(this->fv2(chi, fv1)); + const volScalarField lLES(this->lengthScaleLES(chi, fv1)); + const volScalarField Omega(this->Omega(gradU)); + const volScalarField Ssigma(this->Ssigma(gradU)); + const volScalarField SsigmaDES + ( + pos(dTilda - lRAS)*Omega + (scalar(1) - pos(dTilda - lRAS))*Ssigma + ); + + return + max + ( + SsigmaDES + fv2*this->nuTilda_/sqr(this->kappa_*dTilda), + this->Cs_*SsigmaDES + ); + } + + return + SpalartAllmarasBase>::Stilda + ( + chi, + fv1, + gradU, + dTilda + ); +} + + template tmp SpalartAllmarasDES::dTilda ( @@ -139,6 +179,15 @@ SpalartAllmarasDES::SpalartAllmarasDES propertiesName ), + useSigma_ + ( + Switch::getOrAddToDict + ( + "useSigma", + this->coeffDict_, + false + ) + ), CDES_ ( dimensioned::getOrAddToDict @@ -190,6 +239,7 @@ bool SpalartAllmarasDES::read() { if (SpalartAllmarasBase>::read()) { + useSigma_.readIfPresent("useSigma", this->coeffDict()); CDES_.readIfPresent(this->coeffDict()); lowReCorrection_.readIfPresent("lowReCorrection", this->coeffDict()); fwStar_.readIfPresent(this->coeffDict()); diff --git a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.H b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.H index 8069f66ec7..6936c59f03 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.H +++ b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.H @@ -96,6 +96,9 @@ protected: // Protected data + //- Switch to activate grey-area enhanced sigma-(D)DES + Switch useSigma_; + // Model constants // DES coefficient @@ -108,7 +111,7 @@ protected: // Protected Member Functions - //- Shielding function + //- Low Reynolds number correction function virtual tmp psi ( const volScalarField& chi, @@ -122,6 +125,15 @@ protected: const volScalarField& fv1 ) const; + //- Production term + virtual tmp Stilda + ( + const volScalarField& chi, + const volScalarField& fv1, + const volTensorField& gradU, + const volScalarField& dTilda + ) const; + //- Length scale virtual tmp dTilda ( diff --git a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.C b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.C index 445c4de7fe..2eff1bb0d5 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.C +++ b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.C @@ -58,20 +58,18 @@ tmp kOmegaSSTDDES::S2 ) const { tmp tS2 = - this->kOmegaSSTDES::S2(F1, gradU); + kOmegaSSTBase>::S2(F1, gradU); - if (useSigma_) + if (this->useSigma_) { volScalarField& S2 = tS2.ref(); const volScalarField CDES(this->CDES(F1)); const volScalarField Ssigma(this->Ssigma(gradU)); S2 -= - ( - fd(mag(gradU)) - *pos(this->lengthScaleRAS() - this->lengthScaleLES(CDES)) - *(S2 - sqr(Ssigma)) - ); + fd(mag(gradU)) + *pos(this->lengthScaleRAS() - this->lengthScaleLES(CDES)) + *(S2 - sqr(Ssigma)); } return tS2; @@ -104,21 +102,13 @@ tmp kOmegaSSTDDES::GbyNu0 const volScalarField& S2 ) const { - tmp tGbyNu0 = - this->kOmegaSSTDES::GbyNu0(gradU, F1, S2); - - if (useSigma_) + if (this->useSigma_) { - volScalarField::Internal& GbyNu0 = tGbyNu0.ref(); - const volScalarField CDES(this->CDES(F1)); - - GbyNu0 -= - fd(mag(gradU))()() - *pos(this->lengthScaleRAS()()() - this->lengthScaleLES(CDES)()()) - *(GbyNu0 - S2()); + return S2(); } - return tGbyNu0; + return + kOmegaSSTBase>::GbyNu0(gradU, F1, S2); } @@ -161,18 +151,9 @@ kOmegaSSTDDES::kOmegaSSTDDES type ), - useSigma_ - ( - Switch::getOrAddToDict - ( - "useSigma", - this->coeffDict_, - false - ) - ), Cd1_ ( - useSigma_ ? + this->useSigma_ ? dimensioned::getOrAddToDict ( "Cd1Sigma", diff --git a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.H b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.H index 6a6369121f..f572216a3a 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.H +++ b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.H @@ -70,6 +70,7 @@ class kOmegaSSTDDES { // Private Member Functions + //- Shielding function tmp fd(const volScalarField& magGradU) const; //- No copy construct @@ -83,9 +84,6 @@ protected: // Protected data - //- Switch to activate grey-area enhanced sigma-DDES - Switch useSigma_; - // Model coefficients dimensionedScalar Cd1_; diff --git a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.C b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.C index 7518509a30..edae6c7b57 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.C +++ b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.C @@ -75,6 +75,33 @@ tmp kOmegaSSTDES::r } +template +tmp kOmegaSSTDES::S2 +( + const volScalarField& F1, + const volTensorField& gradU +) const +{ + tmp tS2 = + kOmegaSSTBase>::S2(F1, gradU); + + if (this->useSigma_) + { + volScalarField& S2 = tS2.ref(); + const volScalarField CDES(this->CDES(F1)); + const volScalarField dTilda(this->dTilda(mag(gradU), CDES)); + const volScalarField lengthScaleRAS(this->lengthScaleRAS()); + const volScalarField Ssigma(this->Ssigma(gradU)); + + S2 = + pos(dTilda - lengthScaleRAS)*S2 + + (scalar(1) - pos(dTilda - lengthScaleRAS))*sqr(Ssigma); + } + + return tS2; +} + + template tmp kOmegaSSTDES::dTilda ( @@ -98,6 +125,24 @@ tmp kOmegaSSTDES::epsilonByk } +template +tmp kOmegaSSTDES::GbyNu0 +( + const volTensorField& gradU, + const volScalarField& F1, + const volScalarField& S2 +) const +{ + if (this->useSigma_) + { + return S2(); + } + + return + kOmegaSSTBase>::GbyNu0(gradU, F1, S2); +} + + template tmp kOmegaSSTDES::GbyNu ( @@ -137,6 +182,15 @@ kOmegaSSTDES::kOmegaSSTDES propertiesName ), + useSigma_ + ( + Switch::getOrAddToDict + ( + "useSigma", + this->coeffDict_, + false + ) + ), kappa_ ( dimensioned::getOrAddToDict @@ -188,6 +242,7 @@ bool kOmegaSSTDES::read() { if (kOmegaSSTBase>::read()) { + useSigma_.readIfPresent("useSigma", this->coeffDict()); kappa_.readIfPresent(this->coeffDict()); CDESkom_.readIfPresent(this->coeffDict()); CDESkeps_.readIfPresent(this->coeffDict()); diff --git a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.H b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.H index fd7d9b7587..114633c919 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.H +++ b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.H @@ -86,9 +86,14 @@ protected: // Protected data + //- Switch to activate grey-area enhanced sigma-(D)DES + Switch useSigma_; + // Model coefficients dimensionedScalar kappa_; + + //- DES coefficients dimensionedScalar CDESkom_; dimensionedScalar CDESkeps_; @@ -110,6 +115,13 @@ protected: const volScalarField& magGradU ) const; + //- Return square of strain rate + virtual tmp S2 + ( + const volScalarField& F1, + const volTensorField& gradU + ) const; + //- Length scale virtual tmp dTilda ( @@ -124,6 +136,14 @@ protected: const volTensorField& gradU ) const; + //- Return (G/nu)_0 + virtual tmp GbyNu0 + ( + const volTensorField& gradU, + const volScalarField& F1, + const volScalarField& S2 + ) const; + //- Return G/nu virtual tmp GbyNu (