diff --git a/src/TurbulenceModels/turbulenceModels/Base/kOmegaSST/kOmegaSSTBase.C b/src/TurbulenceModels/turbulenceModels/Base/kOmegaSST/kOmegaSSTBase.C index d750b52e1a..c15bbcf0cb 100644 --- a/src/TurbulenceModels/turbulenceModels/Base/kOmegaSST/kOmegaSSTBase.C +++ b/src/TurbulenceModels/turbulenceModels/Base/kOmegaSST/kOmegaSSTBase.C @@ -142,14 +142,31 @@ template tmp kOmegaSSTBase::epsilonByk ( - const volScalarField::Internal& F1, - const volScalarField::Internal& F2 + const volScalarField& F1, + const volTensorField& gradU ) const { return betaStar_*omega_(); } +template +tmp kOmegaSSTBase::GbyNu +( + const volScalarField::Internal& GbyNu0, + const volScalarField::Internal& F2, + const volScalarField::Internal& S2 +) const +{ + return min + ( + GbyNu0, + (c1_/a1_)*betaStar_*omega_() + *max(a1_*omega_(), b1_*F2*sqrt(S2)) + ); +} + + template tmp kOmegaSSTBase::kSource() const { @@ -181,9 +198,9 @@ tmp kOmegaSSTBase::omegaSource() const template tmp kOmegaSSTBase::Qsas ( - const volScalarField& S2, - const volScalarField& gamma, - const volScalarField& beta + const volScalarField::Internal& S2, + const volScalarField::Internal& gamma, + const volScalarField::Internal& beta ) const { return tmp @@ -422,13 +439,12 @@ void kOmegaSSTBase::correct() BasicEddyViscosityModel::correct(); - volScalarField divU(fvc::div(fvc::absolute(this->phi(), U))); + volScalarField::Internal divU(fvc::div(fvc::absolute(this->phi(), U))); tmp tgradU = fvc::grad(U); volScalarField S2(2*magSqr(symm(tgradU()))); - volScalarField GbyNu((tgradU() && dev(twoSymm(tgradU())))); - volScalarField G(this->GName(), nut*GbyNu); - tgradU.clear(); + volScalarField::Internal GbyNu0((tgradU() && dev(twoSymm(tgradU())))); + volScalarField::Internal G(this->GName(), nut*GbyNu0); // Update omega and G at the wall omega_.boundaryFieldRef().updateCoeffs(); @@ -439,10 +455,11 @@ void kOmegaSSTBase::correct() ); volScalarField F1(this->F1(CDkOmega)); + volScalarField F23(this->F23()); { - volScalarField gamma(this->gamma(F1)); - volScalarField beta(this->beta(F1)); + volScalarField::Internal gamma(this->gamma(F1)); + volScalarField::Internal beta(this->beta(F1)); // Turbulent frequency equation tmp omegaEqn @@ -451,20 +468,15 @@ void kOmegaSSTBase::correct() + fvm::div(alphaRhoPhi, omega_) - fvm::laplacian(alpha*rho*DomegaEff(F1), omega_) == - alpha*rho*gamma - *min - ( - GbyNu, - (c1_/a1_)*betaStar_*omega_*max(a1_*omega_, b1_*F23()*sqrt(S2)) - ) - - fvm::SuSp((2.0/3.0)*alpha*rho*gamma*divU, omega_) - - fvm::Sp(alpha*rho*beta*omega_, omega_) + alpha*rho*gamma*GbyNu(GbyNu0, F23(), S2()) + - fvm::SuSp((2.0/3.0)*alpha()*rho()*gamma*divU, omega_) + - fvm::Sp(alpha()*rho()*beta*omega_(), omega_) - fvm::SuSp ( - alpha*rho*(F1 - scalar(1))*CDkOmega/omega_, + alpha()*rho()*(F1() - scalar(1))*CDkOmega()/omega_(), omega_ ) - + Qsas(S2, gamma, beta) + + Qsas(S2(), gamma, beta) + omegaSource() + fvOptions(alpha, rho, omega_) ); @@ -484,13 +496,15 @@ void kOmegaSSTBase::correct() + fvm::div(alphaRhoPhi, k_) - fvm::laplacian(alpha*rho*DkEff(F1), k_) == - min(alpha*rho*G, (c1_*betaStar_)*alpha*rho*k_*omega_) - - fvm::SuSp((2.0/3.0)*alpha*rho*divU, k_) - - fvm::Sp(alpha*rho*betaStar_*omega_, k_) + alpha()*rho()*Pk(G) + - fvm::SuSp((2.0/3.0)*alpha()*rho()*divU, k_) + - fvm::Sp(alpha()*rho()*epsilonByk(F1, tgradU()), k_) + kSource() + fvOptions(alpha, rho, k_) ); + tgradU.clear(); + kEqn.ref().relax(); fvOptions.constrain(kEqn.ref()); solve(kEqn); diff --git a/src/TurbulenceModels/turbulenceModels/Base/kOmegaSST/kOmegaSSTBase.H b/src/TurbulenceModels/turbulenceModels/Base/kOmegaSST/kOmegaSSTBase.H index 8cf07622c8..6a38be2dc0 100644 --- a/src/TurbulenceModels/turbulenceModels/Base/kOmegaSST/kOmegaSSTBase.H +++ b/src/TurbulenceModels/turbulenceModels/Base/kOmegaSST/kOmegaSSTBase.H @@ -100,9 +100,6 @@ SourceFiles #ifndef kOmegaSSTBase_H #define kOmegaSSTBase_H -#include "RASModel.H" -#include "eddyViscosity.H" - // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // namespace Foam @@ -121,7 +118,7 @@ class kOmegaSSTBase // Disallow default bitwise copy construct and assignment kOmegaSSTBase(const kOmegaSSTBase&); - kOmegaSSTBase& operator=(const kOmegaSSTBase&); + void operator=(const kOmegaSSTBase&); protected: @@ -164,10 +161,10 @@ protected: // Protected Member Functions - tmp F1(const volScalarField& CDkOmega) const; - tmp F2() const; - tmp F3() const; - tmp F23() const; + virtual tmp F1(const volScalarField& CDkOmega) const; + virtual tmp F2() const; + virtual tmp F3() const; + virtual tmp F23() const; tmp blend ( @@ -179,6 +176,16 @@ protected: return F1*(psi1 - psi2) + psi2; } + tmp blend + ( + const volScalarField::Internal& F1, + const dimensionedScalar& psi1, + const dimensionedScalar& psi2 + ) const + { + return F1*(psi1 - psi2) + psi2; + } + tmp alphaK(const volScalarField& F1) const { return blend(F1, alphaK1_, alphaK2_); @@ -189,12 +196,18 @@ protected: return blend(F1, alphaOmega1_, alphaOmega2_); } - tmp beta(const volScalarField& F1) const + tmp beta + ( + const volScalarField::Internal& F1 + ) const { return blend(F1, beta1_, beta2_); } - tmp gamma(const volScalarField& F1) const + tmp gamma + ( + const volScalarField::Internal& F1 + ) const { return blend(F1, gamma1_, gamma2_); } @@ -212,8 +225,16 @@ protected: //- Return epsilon/k which for standard RAS is betaStar*omega virtual tmp epsilonByk ( - const volScalarField::Internal& F1, - const volScalarField::Internal& F2 + const volScalarField& F1, + const volTensorField& gradU + ) const; + + //- Return G/nu + virtual tmp GbyNu + ( + const volScalarField::Internal& GbyNu0, + const volScalarField::Internal& F2, + const volScalarField::Internal& S2 ) const; virtual tmp kSource() const; @@ -222,9 +243,9 @@ protected: virtual tmp Qsas ( - const volScalarField& S2, - const volScalarField& gamma, - const volScalarField& beta + const volScalarField::Internal& S2, + const volScalarField::Internal& gamma, + const volScalarField::Internal& beta ) const; @@ -326,7 +347,7 @@ public: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // #ifdef NoRepository -# include "kOmegaSSTBase.C" + #include "kOmegaSSTBase.C" #endif // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.H b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.H index cae207217a..80ee1b9311 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.H +++ b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDDES/kOmegaSSTDDES.H @@ -72,7 +72,7 @@ class kOmegaSSTDDES // Disallow default bitwise copy construct and assignment kOmegaSSTDDES(const kOmegaSSTDDES&); - kOmegaSSTDDES& operator=(const kOmegaSSTDDES&); + void operator=(const kOmegaSSTDDES&); protected: diff --git a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.C b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.C index a6e1b26b49..f0b35c3f0c 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.C +++ b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.C @@ -66,6 +66,30 @@ tmp kOmegaSSTDES::dTilda } +template +tmp kOmegaSSTDES::epsilonByk +( + const volScalarField& F1, + const volTensorField& gradU +) const +{ + volScalarField CDES(this->CDES(F1)); + return sqrt(this->k_())/dTilda(mag(gradU), CDES)()(); +} + + +template +tmp kOmegaSSTDES::GbyNu +( + const volScalarField::Internal& GbyNu0, + const volScalarField::Internal& F2, + const volScalarField::Internal& S2 +) const +{ + return GbyNu0; // Unlimited +} + + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // template @@ -150,96 +174,6 @@ bool kOmegaSSTDES::read() } -template -void kOmegaSSTDES::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_; - volScalarField& k = this->k_; - volScalarField& omega = this->omega_; - volScalarField& nut = this->nut_; - - DESModel::correct(); - - volScalarField divU(fvc::div(fvc::absolute(this->phi(), U))); - - tmp tgradU = fvc::grad(U); - volScalarField magGradU(mag(tgradU())); - volScalarField S2(2*magSqr(symm(tgradU()))); - volScalarField GbyNu((tgradU() && dev(twoSymm(tgradU())))); - volScalarField G(this->GName(), nut*GbyNu); - tgradU.clear(); - - // Update omega and G at the wall - omega.boundaryFieldRef().updateCoeffs(); - - volScalarField CDkOmega - ( - (2*this->alphaOmega2_)*(fvc::grad(k) & fvc::grad(omega))/omega - ); - - volScalarField F1(this->F1(CDkOmega)); - - { - volScalarField gamma(this->gamma(F1)); - volScalarField beta(this->beta(F1)); - - // Turbulent frequency equation - tmp omegaEqn - ( - fvm::ddt(alpha, rho, omega) - + fvm::div(alphaRhoPhi, omega) - - fvm::laplacian(alpha*rho*this->DomegaEff(F1), omega) - == - alpha*rho*gamma*GbyNu // Using unlimited GybNu - - fvm::SuSp((2.0/3.0)*alpha*rho*gamma*divU, omega) - - fvm::Sp(alpha*rho*beta*omega, omega) - - fvm::SuSp(alpha*rho*(F1 - scalar(1))*CDkOmega/omega, omega) - + this->omegaSource() - ); - - omegaEqn.ref().relax(); - - omegaEqn.ref().boundaryManipulate(omega.boundaryFieldRef()); - - solve(omegaEqn); - bound(omega, this->omegaMin_); - } - - { - volScalarField CDES(this->CDES(F1)); - volScalarField dTilda(this->dTilda(magGradU, CDES)); - - // Turbulent kinetic energy equation - tmp kEqn - ( - fvm::ddt(alpha, rho, k) - + fvm::div(alphaRhoPhi, k) - - fvm::laplacian(alpha*rho*this->DkEff(F1), k) - == - min(alpha*rho*G, (this->c1_*this->betaStar_)*alpha*rho*k*omega) - - fvm::SuSp((2.0/3.0)*alpha*rho*divU, k) - - fvm::Sp(alpha*rho*sqrt(k)/dTilda, k) // modified for DES - + this->kSource() - ); - - kEqn.ref().relax(); - solve(kEqn); - bound(k, this->kMin_); - } - - this->correctNut(S2); -} - - template tmp kOmegaSSTDES::LESRegion() const { diff --git a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.H b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.H index f9f95a11df..0f3cda9a9f 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.H +++ b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTDES/kOmegaSSTDES.H @@ -74,7 +74,7 @@ class kOmegaSSTDES // Disallow default bitwise copy construct and assignment kOmegaSSTDES(const kOmegaSSTDES&); - kOmegaSSTDES& operator=(const kOmegaSSTDES&); + void operator=(const kOmegaSSTDES&); protected: @@ -106,6 +106,21 @@ protected: const volScalarField& CDES ) const; + //- Return epsilon/k + virtual tmp epsilonByk + ( + const volScalarField& F1, + const volTensorField& gradU + ) const; + + //- Return G/nu + virtual tmp GbyNu + ( + const volScalarField::Internal& GbyNu0, + const volScalarField::Internal& F2, + const volScalarField::Internal& S2 + ) const; + public: @@ -144,9 +159,6 @@ public: //- Re-read model coefficients if they have changed virtual bool read(); - //- Solve the turbulence equations and correct the turbulence viscosity - virtual void correct(); - //- Return the LES field indicator virtual tmp LESRegion() const; }; @@ -159,7 +171,7 @@ public: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // #ifdef NoRepository -# include "kOmegaSSTDES.C" + #include "kOmegaSSTDES.C" #endif // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTIDDES/kOmegaSSTIDDES.H b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTIDDES/kOmegaSSTIDDES.H index 898776c835..e3ae10a2a3 100644 --- a/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTIDDES/kOmegaSSTIDDES.H +++ b/src/TurbulenceModels/turbulenceModels/DES/kOmegaSSTIDDES/kOmegaSSTIDDES.H @@ -85,7 +85,7 @@ class kOmegaSSTIDDES // Disallow default bitwise copy construct and assignment kOmegaSSTIDDES(const kOmegaSSTIDDES&); - kOmegaSSTIDDES& operator=(const kOmegaSSTIDDES&); + void operator=(const kOmegaSSTIDDES&); protected: @@ -160,7 +160,7 @@ public: // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // #ifdef NoRepository -# include "kOmegaSSTIDDES.C" + #include "kOmegaSSTIDDES.C" #endif // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSSTLM/kOmegaSSTLM.C b/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSSTLM/kOmegaSSTLM.C index 3bcdf7f648..5d14d97bf9 100644 --- a/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSSTLM/kOmegaSSTLM.C +++ b/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSSTLM/kOmegaSSTLM.C @@ -60,13 +60,13 @@ tmp kOmegaSSTLM::Pk template tmp kOmegaSSTLM::epsilonByk ( - const volScalarField::Internal& F1, - const volScalarField::Internal& F2 + const volScalarField& F1, + const volTensorField& gradU ) const { return min(max(gammaIntEff_, scalar(0.1)), scalar(1)) - *kOmegaSST::epsilonByk(F1, F2); + *kOmegaSST::epsilonByk(F1, gradU); } diff --git a/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSSTLM/kOmegaSSTLM.H b/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSSTLM/kOmegaSSTLM.H index 3d5acaef06..88626ff1ed 100644 --- a/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSSTLM/kOmegaSSTLM.H +++ b/src/TurbulenceModels/turbulenceModels/RAS/kOmegaSSTLM/kOmegaSSTLM.H @@ -164,8 +164,8 @@ protected: //- Modified form of the k-omega SST epsilon/k virtual tmp epsilonByk ( - const volScalarField::Internal& F1, - const volScalarField::Internal& F2 + const volScalarField& F1, + const volTensorField& gradU ) const; //- Freestream blending-function