diff --git a/src/TurbulenceModels/turbulenceModels/Base/SpalartAllmaras/SpalartAllmarasBase.C b/src/TurbulenceModels/turbulenceModels/Base/SpalartAllmaras/SpalartAllmarasBase.C new file mode 100644 index 0000000000..6400f85e27 --- /dev/null +++ b/src/TurbulenceModels/turbulenceModels/Base/SpalartAllmaras/SpalartAllmarasBase.C @@ -0,0 +1,471 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2011-2016 OpenFOAM Foundation + Copyright (C) 2016-2022 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "SpalartAllmarasBase.H" +#include "wallDist.H" +#include "bound.H" +#include "fvOptions.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +template +tmp SpalartAllmarasBase::chi() const +{ + return nuTilda_/this->nu(); +} + + +template +tmp SpalartAllmarasBase::fv1 +( + const volScalarField& chi +) const +{ + const volScalarField chi3("chi3", pow3(chi)); + return chi3/(chi3 + pow3(Cv1_)); +} + + +template +tmp SpalartAllmarasBase::fv2 +( + const volScalarField& chi, + const volScalarField& fv1 +) const +{ + return scalar(1) - chi/(scalar(1) + chi*fv1); +} + + +template +tmp SpalartAllmarasBase::ft2 +( + const volScalarField& chi +) const +{ + return Ct3_*exp(-Ct4_*sqr(chi)); +} + + +template +tmp SpalartAllmarasBase::Omega +( + const volTensorField& gradU +) const +{ + return sqrt(2.0)*mag(skew(gradU)); +} + + +template +tmp SpalartAllmarasBase::r +( + const volScalarField& nur, + const volScalarField& Stilda, + const volScalarField& dTilda +) const +{ + const dimensionedScalar eps("SMALL", Stilda.dimensions(), SMALL); + + tmp tr = + min(nur/(max(Stilda, eps)*sqr(kappa_*dTilda)), scalar(10)); + + tr.ref().boundaryFieldRef() == 0; + + return tr; +} + + +template +tmp SpalartAllmarasBase::fw +( + const volScalarField& Stilda, + const volScalarField& dTilda +) const +{ + const volScalarField::Internal r(this->r(nuTilda_, Stilda, dTilda)()()); + const volScalarField::Internal g(r + Cw2_*(pow6(r) - r)); + + return g*pow((1 + pow6(Cw3_))/(pow6(g) + pow6(Cw3_)), 1.0/6.0); +} + + +template +tmp SpalartAllmarasBase::Stilda +( + const volScalarField& chi, + const volScalarField& fv1, + const volTensorField& gradU, + const volScalarField& dTilda +) const +{ + const volScalarField Omega(this->Omega(gradU)); + + return + max + ( + Omega + fv2(chi, fv1)*nuTilda_/sqr(kappa_*dTilda), + Cs_*Omega + ); +} + + +template +void SpalartAllmarasBase::correctNut +( + const volScalarField& fv1 +) +{ + this->nut_ = nuTilda_*fv1; + this->nut_.correctBoundaryConditions(); + fv::options::New(this->mesh_).correct(this->nut_); +} + + +template +void SpalartAllmarasBase::correctNut() +{ + correctNut(fv1(this->chi())); +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template +SpalartAllmarasBase::SpalartAllmarasBase +( + const word& type, + const alphaField& alpha, + const rhoField& rho, + const volVectorField& U, + const surfaceScalarField& alphaRhoPhi, + const surfaceScalarField& phi, + const transportModel& transport, + const word& propertiesName +) +: + BasicEddyViscosityModel + ( + type, + alpha, + rho, + U, + alphaRhoPhi, + phi, + transport, + propertiesName + ), + + sigmaNut_ + ( + dimensioned::getOrAddToDict + ( + "sigmaNut", + this->coeffDict_, + 0.66666 + ) + ), + kappa_ + ( + dimensioned::getOrAddToDict + ( + "kappa", + this->coeffDict_, + 0.41 + ) + ), + Cb1_ + ( + dimensioned::getOrAddToDict + ( + "Cb1", + this->coeffDict_, + 0.1355 + ) + ), + Cb2_ + ( + dimensioned::getOrAddToDict + ( + "Cb2", + this->coeffDict_, + 0.622 + ) + ), + Cw1_(Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_), + Cw2_ + ( + dimensioned::getOrAddToDict + ( + "Cw2", + this->coeffDict_, + 0.3 + ) + ), + Cw3_ + ( + dimensioned::getOrAddToDict + ( + "Cw3", + this->coeffDict_, + 2.0 + ) + ), + Cv1_ + ( + dimensioned::getOrAddToDict + ( + "Cv1", + this->coeffDict_, + 7.1 + ) + ), + Cs_ + ( + dimensioned::getOrAddToDict + ( + "Cs", + this->coeffDict_, + 0.3 + ) + ), + ck_ + ( + dimensioned::getOrAddToDict + ( + "ck", + this->coeffDict_, + 0.07 + ) + ), + Ct3_ + ( + dimensioned::getOrAddToDict + ( + "Ct3", + this->coeffDict_, + 1.2 + ) + ), + Ct4_ + ( + dimensioned::getOrAddToDict + ( + "Ct4", + this->coeffDict_, + 0.5 + ) + ), + + nuTilda_ + ( + IOobject + ( + "nuTilda", + this->runTime_.timeName(), + this->mesh_, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + this->mesh_ + ), + + y_(wallDist::New(this->mesh_).y()) +{ + if (mag(Ct3_.value()) > SMALL) + { + Info<< " ft2 term: active" << nl; + } + else + { + Info<< " ft2 term: inactive" << nl; + } +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template +bool SpalartAllmarasBase::read() +{ + if (BasicEddyViscosityModel::read()) + { + sigmaNut_.readIfPresent(this->coeffDict()); + kappa_.readIfPresent(*this); + + Cb1_.readIfPresent(this->coeffDict()); + Cb2_.readIfPresent(this->coeffDict()); + Cw1_ = Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_; + Cw2_.readIfPresent(this->coeffDict()); + Cw3_.readIfPresent(this->coeffDict()); + Cv1_.readIfPresent(this->coeffDict()); + Cs_.readIfPresent(this->coeffDict()); + + ck_.readIfPresent(this->coeffDict()); + + Ct3_.readIfPresent(this->coeffDict()); + Ct4_.readIfPresent(this->coeffDict()); + + if (mag(Ct3_.value()) > SMALL) + { + Info<< " ft2 term: active" << nl; + } + else + { + Info<< " ft2 term: inactive" << nl; + } + + return true; + } + + return false; +} + + +template +tmp +SpalartAllmarasBase::DnuTildaEff() const +{ + return tmp::New + ( + IOobject::groupName("DnuTildaEff", this->alphaRhoPhi_.group()), + (nuTilda_ + this->nu())/sigmaNut_ + ); +} + + +template +tmp SpalartAllmarasBase::k() const +{ + // (B:Eq. 4.50) + const scalar Cmu = 0.09; + const auto fv1 = this->fv1(chi()); + + return tmp::New + ( + IOobject::groupName("k", this->alphaRhoPhi_.group()), + cbrt(fv1)*nuTilda_*::sqrt(scalar(2)/Cmu)*mag(symm(fvc::grad(this->U_))) + ); +} + +template +tmp +SpalartAllmarasBase::epsilon() const +{ + // (B:Eq. 4.50) + const scalar Cmu = 0.09; + const auto fv1 = this->fv1(chi()); + const dimensionedScalar nutSMALL(nuTilda_.dimensions(), SMALL); + + return tmp::New + ( + IOobject::groupName("epsilon", this->alphaRhoPhi_.group()), + sqrt(fv1)*sqr(::sqrt(Cmu)*this->k())/(nuTilda_ + this->nut_ + nutSMALL) + ); +} + + +template +tmp SpalartAllmarasBase::omega() const +{ + // (P:p. 384) + const scalar betaStar = 0.09; + const dimensionedScalar k0(sqr(dimLength/dimTime), SMALL); + + return tmp::New + ( + IOobject::groupName("omega", this->alphaRhoPhi_.group()), + this->epsilon()/(betaStar*(this->k() + k0)) + ); +} + + +template +void SpalartAllmarasBase::correct() +{ + if (!this->turbulence_) + { + return; + } + + // Local references + const alphaField& alpha = this->alpha_; + const rhoField& rho = this->rho_; + const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_; + const volVectorField& U = this->U_; + fv::options& fvOptions(fv::options::New(this->mesh_)); + + BasicEddyViscosityModel::correct(); + + const volScalarField chi(this->chi()); + const volScalarField fv1(this->fv1(chi)); + const volScalarField ft2(this->ft2(chi)); + + tmp tgradU = fvc::grad(U); + volScalarField dTilda(this->dTilda(chi, fv1, tgradU())); + volScalarField Stilda(this->Stilda(chi, fv1, tgradU(), dTilda)); + tgradU.clear(); + + tmp nuTildaEqn + ( + fvm::ddt(alpha, rho, nuTilda_) + + fvm::div(alphaRhoPhi, nuTilda_) + - fvm::laplacian(alpha*rho*DnuTildaEff(), nuTilda_) + - Cb2_/sigmaNut_*alpha()*rho()*magSqr(fvc::grad(nuTilda_)()()) + == + Cb1_*alpha()*rho()*Stilda()*nuTilda_()*(scalar(1) - ft2()) + - fvm::Sp + ( + (Cw1_*fw(Stilda, dTilda) - Cb1_/sqr(kappa_)*ft2()) + *alpha()*rho()*nuTilda_()/sqr(dTilda()), + nuTilda_ + ) + + fvOptions(alpha, rho, nuTilda_) + ); + + nuTildaEqn.ref().relax(); + fvOptions.constrain(nuTildaEqn.ref()); + solve(nuTildaEqn); + fvOptions.correct(nuTilda_); + bound(nuTilda_, dimensionedScalar(nuTilda_.dimensions(), Zero)); + nuTilda_.correctBoundaryConditions(); + + correctNut(); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/TurbulenceModels/turbulenceModels/Base/SpalartAllmaras/SpalartAllmarasBase.H b/src/TurbulenceModels/turbulenceModels/Base/SpalartAllmaras/SpalartAllmarasBase.H new file mode 100644 index 0000000000..9aca849ea4 --- /dev/null +++ b/src/TurbulenceModels/turbulenceModels/Base/SpalartAllmaras/SpalartAllmarasBase.H @@ -0,0 +1,226 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | www.openfoam.com + \\/ M anipulation | +------------------------------------------------------------------------------- + Copyright (C) 2011-2016 OpenFOAM Foundation + Copyright (C) 2015-2022 OpenCFD Ltd. +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Class + Foam::LESModels::SpalartAllmarasBase + +Group + grpDESTurbulence + +Description + SpalartAllmarasBase DES turbulence model for incompressible and + compressible flows + + Reference: + \verbatim + Spalart, P. R., Jou, W. H., Strelets, M., & Allmaras, S. R. (1997). + Comments on the feasibility of LES for wings, and on a hybrid + RANS/LES approach. + Advances in DNS/LES, 1, 4-8. + \endverbatim + +SourceFiles + SpalartAllmarasBase.C + +\*---------------------------------------------------------------------------*/ + +#ifndef Foam_SpalartAllmarasBase_H +#define Foam_SpalartAllmarasBase_H + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class SpalartAllmarasBase Declaration +\*---------------------------------------------------------------------------*/ + +template +class SpalartAllmarasBase +: + public BasicEddyViscosityModel +{ + // Private Member Functions + + //- No copy construct + SpalartAllmarasBase(const SpalartAllmarasBase&) = delete; + + //- No copy assignment + void operator=(const SpalartAllmarasBase&) = delete; + + +protected: + + // Protected data + + // Model constants + + dimensionedScalar sigmaNut_; + dimensionedScalar kappa_; + + dimensionedScalar Cb1_; + dimensionedScalar Cb2_; + dimensionedScalar Cw1_; + dimensionedScalar Cw2_; + dimensionedScalar Cw3_; + dimensionedScalar Cv1_; + dimensionedScalar Cs_; + dimensionedScalar ck_; + + + dimensionedScalar Ct3_; + dimensionedScalar Ct4_; + + + // Fields + + volScalarField nuTilda_; + + //- Wall distance + // Note: different to wall distance in parent RASModel + // which is for near-wall cells only + const volScalarField& y_; + + + // Protected Member Functions + + tmp chi() const; + + tmp fv1(const volScalarField& chi) const; + + tmp fv2 + ( + const volScalarField& chi, + const volScalarField& fv1 + ) const; + + tmp ft2(const volScalarField& chi) const; + + tmp Omega(const volTensorField& gradU) const; + + tmp r + ( + const volScalarField& nur, + const volScalarField& Stilda, + const volScalarField& dTilda + ) const; + + tmp fw + ( + const volScalarField& Stilda, + const volScalarField& dTilda + ) const; + + virtual tmp Stilda + ( + const volScalarField& chi, + const volScalarField& fv1, + const volTensorField& gradU, + const volScalarField& dTilda + ) const; + + //- Length scale + virtual tmp dTilda + ( + const volScalarField& chi, + const volScalarField& fv1, + const volTensorField& gradU + ) const = 0; + + void correctNut(const volScalarField& fv1); + virtual void correctNut(); + + +public: + + typedef typename BasicEddyViscosityModel::alphaField alphaField; + typedef typename BasicEddyViscosityModel::rhoField rhoField; + typedef typename BasicEddyViscosityModel::transportModel transportModel; + + + // Constructors + + //- Construct from components + SpalartAllmarasBase + ( + const word& type, + const alphaField& alpha, + const rhoField& rho, + const volVectorField& U, + const surfaceScalarField& alphaRhoPhi, + const surfaceScalarField& phi, + const transportModel& transport, + const word& propertiesName = turbulenceModel::propertiesName + ); + + + //- Destructor + virtual ~SpalartAllmarasBase() = default; + + + // Member Functions + + //- Re-read model coefficients if they have changed + virtual bool read(); + + //- Return the effective diffusivity for nuTilda + tmp DnuTildaEff() const; + + //- Return the (estimated) turbulent kinetic energy + virtual tmp k() const; + + //- Return the (estimated) turbulent kinetic energy dissipation rate + virtual tmp epsilon() const; + + //- Return the (estimated) specific dissipation rate + virtual tmp omega() const; + + tmp nuTilda() const + { + return nuTilda_; + } + + //- Correct nuTilda and related properties + virtual void correct(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository + #include "SpalartAllmarasBase.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/TurbulenceModels/turbulenceModels/DES/DESModel/DESModel.C b/src/TurbulenceModels/turbulenceModels/DES/DESModel/DESModel.C index ae5a815b7d..94ed812c67 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/DESModel/DESModel.C +++ b/src/TurbulenceModels/turbulenceModels/DES/DESModel/DESModel.C @@ -63,7 +63,7 @@ DESModel::DESModel {} -// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // template DESModel::~DESModel() diff --git a/src/TurbulenceModels/turbulenceModels/DES/DESModel/DESModel.H b/src/TurbulenceModels/turbulenceModels/DES/DESModel/DESModel.H index 4422bf3274..fd1a452518 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/DESModel/DESModel.H +++ b/src/TurbulenceModels/turbulenceModels/DES/DESModel/DESModel.H @@ -92,7 +92,7 @@ public: //- Destructor - virtual ~DESModel(); + virtual ~DESModel() = default; // Public Member Functions diff --git a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.C b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.C index aa955595b5..f0af3c94c0 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.C +++ b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2015-2020 OpenCFD Ltd. + Copyright (C) 2015-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -37,41 +37,15 @@ namespace LESModels // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * // -template -tmp SpalartAllmarasDDES::rd -( - const volScalarField& magGradU -) const -{ - tmp tr - ( - min - ( - this->nuEff() - /( - max - ( - magGradU, - dimensionedScalar("SMALL", magGradU.dimensions(), SMALL) - ) - *sqr(this->kappa_*this->y_) - ), - scalar(10) - ) - ); - tr.ref().boundaryFieldRef() == 0.0; - - return tr; -} - - template tmp SpalartAllmarasDDES::fd ( const volScalarField& magGradU ) const { - return 1 - tanh(pow(Cd1_*rd(magGradU), Cd2_)); + return + 1 + - tanh(pow(this->Cd1_*this->r(this->nuEff(), magGradU, this->y_), Cd2_)); } @@ -86,7 +60,7 @@ tmp SpalartAllmarasDDES::dTilda ) const { const volScalarField& lRAS(this->y_); - const volScalarField lLES(this->psi(chi, fv1)*this->CDES_*this->delta()); + const volScalarField lLES(this->lengthScaleLES(chi, fv1)); return max ( diff --git a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.H b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.H index 414d278759..29f279d353 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.H +++ b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDDES/SpalartAllmarasDDES.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2019 OpenCFD Ltd. + Copyright (C) 2019-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -48,10 +48,10 @@ SourceFiles \*---------------------------------------------------------------------------*/ -#ifndef SpalartAllmarasDDES_H -#define SpalartAllmarasDDES_H +#ifndef Foam_SpalartAllmarasDDES_H +#define Foam_SpalartAllmarasDDES_H -#include "SpalartAllmarasDES.H" +#include "SpalartAllmarasBase.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -71,10 +71,9 @@ class SpalartAllmarasDDES { // Private Member Functions + //- Delay function tmp fd(const volScalarField& magGradU) const; - tmp rd(const volScalarField& magGradU) const; - //- No copy construct SpalartAllmarasDDES(const SpalartAllmarasDDES&) = delete; diff --git a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.C b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.C index e474d3065f..1b17db41f6 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.C +++ b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2016-2020 OpenCFD Ltd. + Copyright (C) 2016-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -38,120 +38,6 @@ namespace LESModels // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // -template -tmp SpalartAllmarasDES::chi() const -{ - return nuTilda_/this->nu(); -} - - -template -tmp SpalartAllmarasDES::fv1 -( - const volScalarField& chi -) const -{ - const volScalarField chi3("chi3", pow3(chi)); - return chi3/(chi3 + pow3(Cv1_)); -} - - -template -tmp SpalartAllmarasDES::fv2 -( - const volScalarField& chi, - const volScalarField& fv1 -) const -{ - return 1.0 - chi/(1.0 + chi*fv1); -} - - -template -tmp SpalartAllmarasDES::ft2 -( - const volScalarField& chi -) const -{ - return Ct3_*exp(-Ct4_*sqr(chi)); -} - - -template -tmp SpalartAllmarasDES::Omega -( - const volTensorField& gradU -) const -{ - return sqrt(2.0)*mag(skew(gradU)); -} - - -template -tmp SpalartAllmarasDES::Stilda -( - const volScalarField& chi, - const volScalarField& fv1, - const volScalarField& Omega, - const volScalarField& dTilda -) const -{ - return - ( - max - ( - Omega - + fv2(chi, fv1)*nuTilda_/sqr(kappa_*dTilda), - Cs_*Omega - ) - ); -} - - -template -tmp SpalartAllmarasDES::r -( - const volScalarField& nur, - const volScalarField& Omega, - const volScalarField& dTilda -) const -{ - tmp tr - ( - min - ( - nur - /( - max - ( - Omega, - dimensionedScalar("SMALL", Omega.dimensions(), SMALL) - ) - *sqr(kappa_*dTilda) - ), - scalar(10) - ) - ); - tr.ref().boundaryFieldRef() == 0.0; - - return tr; -} - - -template -tmp SpalartAllmarasDES::fw -( - const volScalarField& Omega, - const volScalarField& dTilda -) const -{ - const volScalarField r(this->r(nuTilda_, Omega, dTilda)); - const volScalarField g(r + Cw2_*(pow6(r) - r)); - - return g*pow((1 + pow6(Cw3_))/(pow6(g) + pow6(Cw3_)), 1.0/6.0); -} - - template tmp SpalartAllmarasDES::psi ( @@ -159,26 +45,23 @@ tmp SpalartAllmarasDES::psi const volScalarField& fv1 ) const { - tmp tpsi + auto tpsi = tmp::New ( - new volScalarField + IOobject ( - IOobject - ( - type() + ":psi", - this->time().timeName(), - this->mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE - ), + IOobject::scopedName(type(), "psi"), + this->time().timeName(), this->mesh(), - dimensionedScalar("one", dimless, 1) - ) + IOobject::NO_READ, + IOobject::NO_WRITE + ), + this->mesh(), + dimensionedScalar(dimless, Zero) ); if (lowReCorrection_) { - volScalarField& psi = tpsi.ref(); + auto& psi = tpsi.ref(); const volScalarField fv2(this->fv2(chi, fv1)); const volScalarField ft2(this->ft2(chi)); @@ -189,7 +72,8 @@ tmp SpalartAllmarasDES::psi min ( scalar(100), - (1 - Cb1_/(Cw1_*sqr(kappa_)*fwStar_)*(ft2 + (1 - ft2)*fv2)) + (1 - this->Cb1_/(this->Cw1_*sqr(this->kappa_)*fwStar_) + *(ft2 + (1 - ft2)*fv2)) /max(SMALL, (fv1*max(scalar(1e-10), 1 - ft2))) ) ); @@ -199,6 +83,17 @@ tmp SpalartAllmarasDES::psi } +template +tmp SpalartAllmarasDES::lengthScaleLES +( + const volScalarField& chi, + const volScalarField& fv1 +) const +{ + return psi(chi, fv1)*CDES_*this->delta(); +} + + template tmp SpalartAllmarasDES::dTilda ( @@ -207,33 +102,15 @@ tmp SpalartAllmarasDES::dTilda const volTensorField& gradU ) const { - tmp tdTilda(psi(chi, fv1)*CDES_*this->delta()); - min(tdTilda.ref().ref(), tdTilda(), y_); + // Initialise with LES length scale + tmp tdTilda = lengthScaleLES(chi, fv1); + + // Take min vs wall distance + min(tdTilda.ref(), tdTilda(), this->y_); return tdTilda; } -template -void SpalartAllmarasDES::correctNut -( - const volScalarField& fv1 -) -{ - this->nut_ = nuTilda_*fv1; - this->nut_.correctBoundaryConditions(); - fv::options::New(this->mesh_).correct(this->nut_); - - BasicTurbulenceModel::correctNut(); -} - - -template -void SpalartAllmarasDES::correctNut() -{ - correctNut(fv1(this->chi())); -} - - // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // template @@ -249,7 +126,7 @@ SpalartAllmarasDES::SpalartAllmarasDES const word& type ) : - DESModel + SpalartAllmarasBase> ( type, alpha, @@ -261,79 +138,6 @@ SpalartAllmarasDES::SpalartAllmarasDES propertiesName ), - sigmaNut_ - ( - dimensioned::getOrAddToDict - ( - "sigmaNut", - this->coeffDict_, - 0.66666 - ) - ), - kappa_ - ( - dimensioned::getOrAddToDict - ( - "kappa", - this->coeffDict_, - 0.41 - ) - ), - Cb1_ - ( - dimensioned::getOrAddToDict - ( - "Cb1", - this->coeffDict_, - 0.1355 - ) - ), - Cb2_ - ( - dimensioned::getOrAddToDict - ( - "Cb2", - this->coeffDict_, - 0.622 - ) - ), - Cw1_(Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_), - Cw2_ - ( - dimensioned::getOrAddToDict - ( - "Cw2", - this->coeffDict_, - 0.3 - ) - ), - Cw3_ - ( - dimensioned::getOrAddToDict - ( - "Cw3", - this->coeffDict_, - 2.0 - ) - ), - Cv1_ - ( - dimensioned::getOrAddToDict - ( - "Cv1", - this->coeffDict_, - 7.1 - ) - ), - Cs_ - ( - dimensioned::getOrAddToDict - ( - "Cs", - this->coeffDict_, - 0.3 - ) - ), CDES_ ( dimensioned::getOrAddToDict @@ -343,15 +147,6 @@ SpalartAllmarasDES::SpalartAllmarasDES 0.65 ) ), - ck_ - ( - dimensioned::getOrAddToDict - ( - "ck", - this->coeffDict_, - 0.07 - ) - ), lowReCorrection_ ( Switch::getOrAddToDict @@ -361,24 +156,6 @@ SpalartAllmarasDES::SpalartAllmarasDES true ) ), - Ct3_ - ( - dimensioned::getOrAddToDict - ( - "Ct3", - this->coeffDict_, - 1.2 - ) - ), - Ct4_ - ( - dimensioned::getOrAddToDict - ( - "Ct4", - this->coeffDict_, - 0.5 - ) - ), fwStar_ ( dimensioned::getOrAddToDict @@ -387,22 +164,7 @@ SpalartAllmarasDES::SpalartAllmarasDES this->coeffDict_, 0.424 ) - ), - - nuTilda_ - ( - IOobject - ( - "nuTilda", - this->runTime_.timeName(), - this->mesh_, - IOobject::MUST_READ, - IOobject::AUTO_WRITE - ), - this->mesh_ - ), - - y_(wallDist::New(this->mesh_).y()) + ) { if (type == typeName) { @@ -416,25 +178,10 @@ SpalartAllmarasDES::SpalartAllmarasDES template bool SpalartAllmarasDES::read() { - if (DESModel::read()) + if (SpalartAllmarasBase>::read()) { - sigmaNut_.readIfPresent(this->coeffDict()); - kappa_.readIfPresent(*this); - - Cb1_.readIfPresent(this->coeffDict()); - Cb2_.readIfPresent(this->coeffDict()); - Cw1_ = Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_; - Cw2_.readIfPresent(this->coeffDict()); - Cw3_.readIfPresent(this->coeffDict()); - Cv1_.readIfPresent(this->coeffDict()); - Cs_.readIfPresent(this->coeffDict()); - CDES_.readIfPresent(this->coeffDict()); - ck_.readIfPresent(this->coeffDict()); - lowReCorrection_.readIfPresent("lowReCorrection", this->coeffDict()); - Ct3_.readIfPresent(this->coeffDict()); - Ct4_.readIfPresent(this->coeffDict()); fwStar_.readIfPresent(this->coeffDict()); return true; @@ -445,124 +192,24 @@ bool SpalartAllmarasDES::read() template -tmp SpalartAllmarasDES:: -DnuTildaEff() const -{ - return tmp - ( - new volScalarField("DnuTildaEff", (nuTilda_ + this->nu())/sigmaNut_) - ); -} - - -template -tmp SpalartAllmarasDES::k() const +tmp +SpalartAllmarasDES::LESRegion() const { const volScalarField chi(this->chi()); const volScalarField fv1(this->fv1(chi)); - tmp tdTilda + return tmp::New ( - new volScalarField + IOobject ( - IOobject - ( - "dTilda", - this->mesh_.time().timeName(), - this->mesh_, - IOobject::NO_READ, - IOobject::NO_WRITE - ), + IOobject::scopedName("DES", "LESRegion"), + this->mesh_.time().timeName(), this->mesh_, - dimensionedScalar(dimLength, Zero), - zeroGradientFvPatchField::typeName - ) + IOobject::NO_READ, + IOobject::NO_WRITE + ), + neg(dTilda(chi, fv1, fvc::grad(this->U_)) - this->y_) ); - volScalarField& dTildaF = tdTilda.ref(); - dTildaF = dTilda(chi, fv1, fvc::grad(this->U_)); - dTildaF.correctBoundaryConditions(); - - return sqr(this->nut()/ck_/dTildaF); -} - - -template -tmp SpalartAllmarasDES::LESRegion() const -{ - const volScalarField chi(this->chi()); - const volScalarField fv1(this->fv1(chi)); - - tmp tLESRegion - ( - new volScalarField - ( - IOobject - ( - "DES::LESRegion", - this->mesh_.time().timeName(), - this->mesh_, - IOobject::NO_READ, - IOobject::NO_WRITE - ), - neg(dTilda(chi, fv1, fvc::grad(this->U_)) - y_) - ) - ); - - return tLESRegion; -} - - -template -void SpalartAllmarasDES::correct() -{ - if (!this->turbulence_) - { - return; - } - - // Local references - const alphaField& alpha = this->alpha_; - const rhoField& rho = this->rho_; - const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_; - const volVectorField& U = this->U_; - fv::options& fvOptions(fv::options::New(this->mesh_)); - - DESModel::correct(); - - const volScalarField chi(this->chi()); - const volScalarField fv1(this->fv1(chi)); - const volScalarField ft2(this->ft2(chi)); - - tmp tgradU = fvc::grad(U); - const volScalarField Omega(this->Omega(tgradU())); - const volScalarField dTilda(this->dTilda(chi, fv1, tgradU())); - const volScalarField Stilda(this->Stilda(chi, fv1, Omega, dTilda)); - - tmp nuTildaEqn - ( - fvm::ddt(alpha, rho, nuTilda_) - + fvm::div(alphaRhoPhi, nuTilda_) - - fvm::laplacian(alpha*rho*DnuTildaEff(), nuTilda_) - - Cb2_/sigmaNut_*alpha*rho*magSqr(fvc::grad(nuTilda_)) - == - Cb1_*alpha*rho*Stilda*nuTilda_*(scalar(1) - ft2) - - fvm::Sp - ( - (Cw1_*fw(Stilda, dTilda) - Cb1_/sqr(kappa_)*ft2) - *alpha*rho*nuTilda_/sqr(dTilda), - nuTilda_ - ) - + fvOptions(alpha, rho, nuTilda_) - ); - - nuTildaEqn.ref().relax(); - fvOptions.constrain(nuTildaEqn.ref()); - solve(nuTildaEqn); - fvOptions.correct(nuTilda_); - bound(nuTilda_, dimensionedScalar(nuTilda_.dimensions(), Zero)); - nuTilda_.correctBoundaryConditions(); - - correctNut(); } diff --git a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.H b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.H index 094ccf99b8..8069f66ec7 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.H +++ b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasDES/SpalartAllmarasDES.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2015-2019 OpenCFD Ltd. + Copyright (C) 2015-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -61,9 +61,10 @@ SourceFiles \*---------------------------------------------------------------------------*/ -#ifndef SpalartAllmarasDES_H -#define SpalartAllmarasDES_H +#ifndef Foam_SpalartAllmarasDES_H +#define Foam_SpalartAllmarasDES_H +#include "SpalartAllmarasBase.H" #include "DESModel.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -80,7 +81,7 @@ namespace LESModels template class SpalartAllmarasDES : - public DESModel + public SpalartAllmarasBase> { // Private Member Functions @@ -97,76 +98,25 @@ protected: // Model constants - dimensionedScalar sigmaNut_; - dimensionedScalar kappa_; - - dimensionedScalar Cb1_; - dimensionedScalar Cb2_; - dimensionedScalar Cw1_; - dimensionedScalar Cw2_; - dimensionedScalar Cw3_; - dimensionedScalar Cv1_; - dimensionedScalar Cs_; + // DES coefficient dimensionedScalar CDES_; - dimensionedScalar ck_; - // Low Reynolds number correction - - Switch lowReCorrection_; - dimensionedScalar Ct3_; - dimensionedScalar Ct4_; - dimensionedScalar fwStar_; - - - // Fields - - volScalarField nuTilda_; - - //- Wall distance - // Note: different to wall distance in parent RASModel - // which is for near-wall cells only - const volScalarField& y_; + Switch lowReCorrection_; + dimensionedScalar fwStar_; // Protected Member Functions - tmp chi() const; - - tmp fv1(const volScalarField& chi) const; - - tmp fv2 + //- Shielding function + virtual tmp psi ( const volScalarField& chi, const volScalarField& fv1 ) const; - tmp ft2(const volScalarField& chi) const; - - tmp Omega(const volTensorField& gradU) const; - - tmp Stilda - ( - const volScalarField& chi, - const volScalarField& fv1, - const volScalarField& Omega, - const volScalarField& dTilda - ) const; - - tmp r - ( - const volScalarField& nur, - const volScalarField& Omega, - const volScalarField& dTilda - ) const; - - tmp fw - ( - const volScalarField& Omega, - const volScalarField& dTilda - ) const; - - tmp psi + //- LES length scale + virtual tmp lengthScaleLES ( const volScalarField& chi, const volScalarField& fv1 @@ -180,9 +130,6 @@ protected: const volTensorField& gradU ) const; - void correctNut(const volScalarField& fv1); - virtual void correctNut(); - public: @@ -220,22 +167,8 @@ public: //- Re-read model coefficients if they have changed virtual bool read(); - //- Return the effective diffusivity for nuTilda - tmp DnuTildaEff() const; - - //- Return SGS kinetic energy - virtual tmp k() const; - - tmp nuTilda() const - { - return nuTilda_; - } - //- Return the LES field indicator virtual tmp LESRegion() const; - - //- Correct nuTilda and related properties - virtual void correct(); }; diff --git a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasIDDES/SpalartAllmarasIDDES.C b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasIDDES/SpalartAllmarasIDDES.C index a0a0535955..cd12c26e15 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasIDDES/SpalartAllmarasIDDES.C +++ b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasIDDES/SpalartAllmarasIDDES.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2015 OpenFOAM Foundation - Copyright (C) 2015-2020 OpenCFD Ltd. + Copyright (C) 2015-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -64,7 +64,7 @@ tmp SpalartAllmarasIDDES::ft const volScalarField& magGradU ) const { - return tanh(pow3(sqr(Ct_)*rd(this->nut_, magGradU))); + return tanh(pow3(sqr(Ct_)*this->r(this->nut_, magGradU, this->y_))); } @@ -74,36 +74,7 @@ tmp SpalartAllmarasIDDES::fl const volScalarField& magGradU ) const { - return tanh(pow(sqr(Cl_)*rd(this->nu(), magGradU), 10)); -} - - -template -tmp SpalartAllmarasIDDES::rd -( - const volScalarField& nur, - const volScalarField& magGradU -) const -{ - tmp tr - ( - min - ( - nur - /( - max - ( - magGradU, - dimensionedScalar("SMALL", magGradU.dimensions(), SMALL) - ) - *sqr(this->kappa_*this->y_) - ), - scalar(10) - ) - ); - tr.ref().boundaryFieldRef() == 0.0; - - return tr; + return tanh(pow(sqr(Cl_)*this->r(this->nu(), magGradU, this->y_), 10)); } @@ -113,7 +84,7 @@ tmp SpalartAllmarasIDDES::fdt const volScalarField& magGradU ) const { - return 1 - tanh(pow(Cdt1_*rd(this->nut_, magGradU), Cdt2_)); + return 1 - tanh(pow(Cdt1_*this->r(this->nut_, magGradU, this->y_), Cdt2_)); } @@ -130,7 +101,7 @@ tmp SpalartAllmarasIDDES::dTilda const volScalarField magGradU(mag(gradU)); const volScalarField psi(this->psi(chi, fv1)); - const volScalarField& lRAS(this->y_); + const volScalarField& lRAS = this->y_; const volScalarField lLES(psi*this->CDES_*this->delta()); const volScalarField alpha(this->alpha()); diff --git a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasIDDES/SpalartAllmarasIDDES.H b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasIDDES/SpalartAllmarasIDDES.H index 66199f664e..9f13f607b0 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasIDDES/SpalartAllmarasIDDES.H +++ b/src/TurbulenceModels/turbulenceModels/DES/SpalartAllmarasIDDES/SpalartAllmarasIDDES.H @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2016 OpenFOAM Foundation - Copyright (C) 2019 OpenCFD Ltd. + Copyright (C) 2019-2022 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -47,8 +47,8 @@ SourceFiles \*---------------------------------------------------------------------------*/ -#ifndef SpalartAllmarasIDDES_H -#define SpalartAllmarasIDDES_H +#ifndef Foam_SpalartAllmarasIDDES_H +#define Foam_SpalartAllmarasIDDES_H #include "SpalartAllmarasDES.H" #include "IDDESDelta.H" @@ -61,7 +61,7 @@ namespace LESModels { /*---------------------------------------------------------------------------*\ - Class SpalartAllmarasIDDES Declaration + Class SpalartAllmarasIDDES Declaration \*---------------------------------------------------------------------------*/ template @@ -75,14 +75,10 @@ class SpalartAllmarasIDDES const IDDESDelta& setDelta() const; tmp alpha() const; - tmp ft(const volScalarField& magGradU) const; - tmp fl(const volScalarField& magGradU) const; - tmp rd - ( - const volScalarField& nur, - const volScalarField& magGradU - ) const; + tmp ft(const volScalarField& magGradU) const; + + tmp fl(const volScalarField& magGradU) const; //- Delay function tmp fdt(const volScalarField& magGradU) const; @@ -105,9 +101,10 @@ protected: dimensionedScalar Cl_; dimensionedScalar Ct_; - // Fields - const IDDESDelta& IDDESDelta_; + //- IDDES delta + const IDDESDelta& IDDESDelta_; + // Protected Member Functions diff --git a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTIDDES/kOmegaSSTIDDES.C b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTIDDES/kOmegaSSTIDDES.C index 1909d02ebc..4714eed05c 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTIDDES/kOmegaSSTIDDES.C +++ b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTIDDES/kOmegaSSTIDDES.C @@ -107,8 +107,6 @@ tmp kOmegaSSTIDDES::rd } -// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // - template tmp kOmegaSSTIDDES::fdt ( @@ -119,6 +117,8 @@ tmp kOmegaSSTIDDES::fdt } +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + template tmp kOmegaSSTIDDES::dTilda ( diff --git a/src/TurbulenceModels/turbulenceModels/RAS/SpalartAllmaras/SpalartAllmaras.C b/src/TurbulenceModels/turbulenceModels/RAS/SpalartAllmaras/SpalartAllmaras.C index 3e38df486d..ec6dba6fb4 100644 --- a/src/TurbulenceModels/turbulenceModels/RAS/SpalartAllmaras/SpalartAllmaras.C +++ b/src/TurbulenceModels/turbulenceModels/RAS/SpalartAllmaras/SpalartAllmaras.C @@ -27,9 +27,6 @@ License \*---------------------------------------------------------------------------*/ #include "SpalartAllmaras.H" -#include "fvOptions.H" -#include "bound.H" -#include "wallDist.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -41,101 +38,14 @@ namespace RASModels // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // template -tmp SpalartAllmaras::chi() const -{ - return nuTilda_/this->nu(); -} - - -template -tmp SpalartAllmaras::fv1 +Foam::tmp SpalartAllmaras::dTilda ( - const volScalarField& chi + const volScalarField& chi, + const volScalarField& fv1, + const volTensorField& gradU ) const { - const volScalarField chi3(pow3(chi)); - - return chi3/(chi3 + pow3(Cv1_)); -} - - -template -tmp SpalartAllmaras::fv2 -( - const volScalarField::Internal& chi, - const volScalarField::Internal& fv1 -) const -{ - return scalar(1) - chi/(scalar(1) + chi*fv1); -} - - -template -tmp SpalartAllmaras::Stilda() -const -{ - const volScalarField chi(this->chi()); - - const volScalarField fv1(this->fv1(chi)); - - const volScalarField::Internal Omega - ( - ::sqrt(scalar(2))*mag(skew(fvc::grad(this->U_)().v())) - ); - - return - ( - max - ( - Omega + fv2(chi(), fv1())*nuTilda_()/sqr(kappa_*y_), - Cs_*Omega - ) - ); -} - - -template -tmp SpalartAllmaras::fw -( - const volScalarField::Internal& Stilda -) const -{ - const volScalarField::Internal r - ( - min - ( - nuTilda_ - /( - max - ( - Stilda, - dimensionedScalar(Stilda.dimensions(), SMALL) - ) - *sqr(kappa_*y_) - ), - scalar(10) - ) - ); - - const volScalarField::Internal g(r + Cw2_*(pow6(r) - r)); - - return - g*pow - ( - (scalar(1) + pow6(Cw3_))/(pow6(g) + pow6(Cw3_)), - scalar(1)/scalar(6) - ); -} - - -template -void SpalartAllmaras::correctNut() -{ - this->nut_ = nuTilda_*this->fv1(this->chi()); - this->nut_.correctBoundaryConditions(); - fv::options::New(this->mesh_).correct(this->nut_); - - BasicTurbulenceModel::correctNut(); + return this->y_; } @@ -154,7 +64,7 @@ SpalartAllmaras::SpalartAllmaras const word& type ) : - eddyViscosity> + SpalartAllmarasBase>> ( type, alpha, @@ -164,97 +74,7 @@ SpalartAllmaras::SpalartAllmaras phi, transport, propertiesName - ), - - sigmaNut_ - ( - dimensioned::getOrAddToDict - ( - "sigmaNut", - this->coeffDict_, - scalar(2)/scalar(3) - ) - ), - kappa_ - ( - dimensioned::getOrAddToDict - ( - "kappa", - this->coeffDict_, - 0.41 - ) - ), - - Cb1_ - ( - dimensioned::getOrAddToDict - ( - "Cb1", - this->coeffDict_, - 0.1355 - ) - ), - Cb2_ - ( - dimensioned::getOrAddToDict - ( - "Cb2", - this->coeffDict_, - 0.622 - ) - ), - Cw1_(Cb1_/sqr(kappa_) + (scalar(1) + Cb2_)/sigmaNut_), - Cw2_ - ( - dimensioned::getOrAddToDict - ( - "Cw2", - this->coeffDict_, - 0.3 - ) - ), - Cw3_ - ( - dimensioned::getOrAddToDict - ( - "Cw3", - this->coeffDict_, - 2.0 - ) - ), - Cv1_ - ( - dimensioned::getOrAddToDict - ( - "Cv1", - this->coeffDict_, - 7.1 - ) - ), - Cs_ - ( - dimensioned::getOrAddToDict - ( - "Cs", - this->coeffDict_, - 0.3 - ) - ), - - nuTilda_ - ( - IOobject - ( - "nuTilda", - this->runTime_.timeName(), - this->mesh_, - IOobject::MUST_READ, - IOobject::AUTO_WRITE - ), - this->mesh_ - ), - - y_(wallDist::New(this->mesh_).y()) + ) { if (type == typeName) { @@ -263,153 +83,6 @@ SpalartAllmaras::SpalartAllmaras } -// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // - -template -bool SpalartAllmaras::read() -{ - if (eddyViscosity>::read()) - { - sigmaNut_.readIfPresent(this->coeffDict()); - kappa_.readIfPresent(this->coeffDict()); - - Cb1_.readIfPresent(this->coeffDict()); - Cb2_.readIfPresent(this->coeffDict()); - Cw1_ = Cb1_/sqr(kappa_) + (scalar(1) + Cb2_)/sigmaNut_; - Cw2_.readIfPresent(this->coeffDict()); - Cw3_.readIfPresent(this->coeffDict()); - Cv1_.readIfPresent(this->coeffDict()); - Cs_.readIfPresent(this->coeffDict()); - - return true; - } - - return false; -} - - -template -tmp SpalartAllmaras::DnuTildaEff() const -{ - return tmp::New - ( - "DnuTildaEff", - (nuTilda_ + this->nu())/sigmaNut_ - ); -} - - -template -tmp SpalartAllmaras::k() const -{ - // (B:Eq. 4.50) - const scalar Cmu = 0.09; - - return tmp::New - ( - IOobject - ( - IOobject::groupName("k", this->alphaRhoPhi_.group()), - this->runTime_.timeName(), - this->mesh_ - ), - cbrt(this->fv1(this->chi())) - *nuTilda_ - *::sqrt(scalar(2)/Cmu) - *mag(symm(fvc::grad(this->U_))), - this->nut_.boundaryField().types() - ); -} - - -template -tmp SpalartAllmaras::epsilon() const -{ - // (B:Eq. 4.50) - const scalar Cmu = 0.09; - const dimensionedScalar nutSMALL(sqr(dimLength)/dimTime, SMALL); - - return tmp::New - ( - IOobject - ( - IOobject::groupName("epsilon", this->alphaRhoPhi_.group()), - this->runTime_.timeName(), - this->mesh_ - ), - pow(this->fv1(this->chi()), 0.5) - *pow(::sqrt(Cmu)*this->k(), 2) - /(nuTilda_ + this->nut_ + nutSMALL), - this->nut_.boundaryField().types() - ); -} - - -template -tmp SpalartAllmaras::omega() const -{ - // (P:p. 384) - const scalar betaStar = 0.09; - const dimensionedScalar k0(sqr(dimLength/dimTime), SMALL); - - return tmp::New - ( - IOobject - ( - IOobject::groupName("omega", this->alphaRhoPhi_.group()), - this->runTime_.timeName(), - this->mesh_ - ), - this->epsilon()/(betaStar*(this->k() + k0)), - this->nut_.boundaryField().types() - ); -} - - -template -void SpalartAllmaras::correct() -{ - if (!this->turbulence_) - { - return; - } - - { - // Construct local convenience references - const alphaField& alpha = this->alpha_; - const rhoField& rho = this->rho_; - const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_; - fv::options& fvOptions(fv::options::New(this->mesh_)); - - eddyViscosity>::correct(); - - const volScalarField::Internal Stilda(this->Stilda()); - - tmp nuTildaEqn - ( - fvm::ddt(alpha, rho, nuTilda_) - + fvm::div(alphaRhoPhi, nuTilda_) - - fvm::laplacian(alpha*rho*DnuTildaEff(), nuTilda_) - - Cb2_/sigmaNut_*alpha*rho*magSqr(fvc::grad(nuTilda_)) - == - Cb1_*alpha()*rho()*Stilda*nuTilda_() - - fvm::Sp(Cw1_*alpha()*rho()*fw(Stilda)*nuTilda_()/sqr(y_), nuTilda_) - + fvOptions(alpha, rho, nuTilda_) - ); - - nuTildaEqn.ref().relax(); - fvOptions.constrain(nuTildaEqn.ref()); - solve(nuTildaEqn); - fvOptions.correct(nuTilda_); - bound(nuTilda_, dimensionedScalar(nuTilda_.dimensions(), Zero)); - nuTilda_.correctBoundaryConditions(); - } - - // Update nut with latest available nuTilda - correctNut(); -} - - // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // } // End namespace RASModels diff --git a/src/TurbulenceModels/turbulenceModels/RAS/SpalartAllmaras/SpalartAllmaras.H b/src/TurbulenceModels/turbulenceModels/RAS/SpalartAllmaras/SpalartAllmaras.H index 5391c7d12c..3841682e7d 100644 --- a/src/TurbulenceModels/turbulenceModels/RAS/SpalartAllmaras/SpalartAllmaras.H +++ b/src/TurbulenceModels/turbulenceModels/RAS/SpalartAllmaras/SpalartAllmaras.H @@ -96,8 +96,6 @@ Note - The model is implemented without the trip-term since the model has almost always been used in fully turbulent applications rather than those where laminar-turbulent transition occurs. - - It has been argued that the \c ft2 term is not needed in the absence of the - trip-term, hence \c ft2 term is also not implementated. - The \c Stilda generation term should never be allowed to be zero or negative to avoid potential numerical issues and unphysical results for complex flows. To this end, a limiter proposed by Spalart (R:Note-1(b)) is applied @@ -112,11 +110,12 @@ SourceFiles \*---------------------------------------------------------------------------*/ -#ifndef SpalartAllmaras_H -#define SpalartAllmaras_H +#ifndef Foam_SpalartAllmaras_H +#define Foam_SpalartAllmaras_H #include "RASModel.H" #include "eddyViscosity.H" +#include "SpalartAllmarasBase.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -132,7 +131,7 @@ namespace RASModels template class SpalartAllmaras : - public eddyViscosity> + public SpalartAllmarasBase>> { // Private Member Functions @@ -145,55 +144,16 @@ class SpalartAllmaras protected: - // Protected Data - - // Model coefficients - - dimensionedScalar sigmaNut_; - dimensionedScalar kappa_; - dimensionedScalar Cb1_; - dimensionedScalar Cb2_; - dimensionedScalar Cw1_; - dimensionedScalar Cw2_; - dimensionedScalar Cw3_; - dimensionedScalar Cv1_; - dimensionedScalar Cs_; - - - // Fields - - //- Modified kinematic viscosity [m2/s] - volScalarField nuTilda_; - - //- Wall distance - // Note: different to wall distance in parent RASModel - // which is for near-wall cells only - const volScalarField::Internal& y_; - - // Protected Member Functions - tmp chi() const; - - tmp fv1(const volScalarField& chi) const; - - tmp fv2 + //- Length scale + virtual tmp dTilda ( - const volScalarField::Internal& chi, - const volScalarField::Internal& fv1 + const volScalarField& chi, + const volScalarField& fv1, + const volTensorField& gradU ) const; - tmp Stilda() const; - - tmp fw - ( - const volScalarField::Internal& Stilda - ) const; - - - //- Update nut with the latest available nuTilda - virtual void correctNut(); - public: @@ -224,27 +184,6 @@ public: //- Destructor virtual ~SpalartAllmaras() = default; - - - // Member Functions - - //- Re-read model coefficients if they have changed - virtual bool read(); - - //- Return the effective diffusivity for nuTilda - tmp DnuTildaEff() const; - - //- Return the (estimated) turbulent kinetic energy - virtual tmp k() const; - - //- Return the (estimated) turbulent kinetic energy dissipation rate - virtual tmp epsilon() const; - - //- Return the (estimated) specific dissipation rate - virtual tmp omega() const; - - //- Solve the turbulence equations and correct the turbulent viscosity - virtual void correct(); };