diff --git a/src/TurbulenceModels/compressible/turbulentFluidThermoModels/turbulentFluidThermoModels.C b/src/TurbulenceModels/compressible/turbulentFluidThermoModels/turbulentFluidThermoModels.C index 4d94dedf1..24f6f157d 100644 --- a/src/TurbulenceModels/compressible/turbulentFluidThermoModels/turbulentFluidThermoModels.C +++ b/src/TurbulenceModels/compressible/turbulentFluidThermoModels/turbulentFluidThermoModels.C @@ -57,6 +57,9 @@ makeRASModel(SpalartAllmaras); #include "kEpsilon.H" makeRASModel(kEpsilon); +#include "RNGkEpsilon.H" +makeRASModel(RNGkEpsilon); + #include "realizableKE.H" makeRASModel(realizableKE); diff --git a/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C b/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C index 506fbbdd0..791ea15db 100644 --- a/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C +++ b/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C @@ -55,6 +55,9 @@ makeRASModel(SpalartAllmaras); #include "kEpsilon.H" makeRASModel(kEpsilon); +#include "RNGkEpsilon.H" +makeRASModel(RNGkEpsilon); + #include "realizableKE.H" makeRASModel(realizableKE); diff --git a/src/TurbulenceModels/turbulenceModels/RAS/RNGkEpsilon/RNGkEpsilon.C b/src/TurbulenceModels/turbulenceModels/RAS/RNGkEpsilon/RNGkEpsilon.C new file mode 100644 index 000000000..563d33ebf --- /dev/null +++ b/src/TurbulenceModels/turbulenceModels/RAS/RNGkEpsilon/RNGkEpsilon.C @@ -0,0 +1,322 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +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 "RNGkEpsilon.H" +#include "bound.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace RASModels +{ + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +template +void RNGkEpsilon::correctNut() +{ + this->nut_ = Cmu_*sqr(k_)/epsilon_; + this->nut_.correctBoundaryConditions(); +} + + +template +tmp RNGkEpsilon::kSource() const +{ + return tmp + ( + new fvScalarMatrix + ( + k_, + dimVolume*this->rho_.dimensions()*k_.dimensions() + /dimTime + ) + ); +} + + +template +tmp RNGkEpsilon::epsilonSource() const +{ + return tmp + ( + new fvScalarMatrix + ( + epsilon_, + dimVolume*this->rho_.dimensions()*epsilon_.dimensions() + /dimTime + ) + ); +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template +RNGkEpsilon::RNGkEpsilon +( + const alphaField& alpha, + const rhoField& rho, + const volVectorField& U, + const surfaceScalarField& alphaRhoPhi, + const surfaceScalarField& phi, + const transportModel& transport, + const word& propertiesName, + const word& type +) +: + eddyViscosity > + ( + type, + alpha, + rho, + U, + alphaRhoPhi, + phi, + transport, + propertiesName + ), + + Cmu_ + ( + dimensioned::lookupOrAddToDict + ( + "Cmu", + this->coeffDict_, + 0.0845 + ) + ), + C1_ + ( + dimensioned::lookupOrAddToDict + ( + "C1", + this->coeffDict_, + 1.42 + ) + ), + C2_ + ( + dimensioned::lookupOrAddToDict + ( + "C2", + this->coeffDict_, + 1.68 + ) + ), + C3_ + ( + dimensioned::lookupOrAddToDict + ( + "C3", + this->coeffDict_, + -0.33 + ) + ), + sigmak_ + ( + dimensioned::lookupOrAddToDict + ( + "sigmak", + this->coeffDict_, + 0.71942 + ) + ), + sigmaEps_ + ( + dimensioned::lookupOrAddToDict + ( + "sigmaEps", + this->coeffDict_, + 0.71942 + ) + ), + eta0_ + ( + dimensioned::lookupOrAddToDict + ( + "eta0", + this->coeffDict_, + 4.38 + ) + ), + beta_ + ( + dimensioned::lookupOrAddToDict + ( + "beta", + this->coeffDict_, + 0.012 + ) + ), + + k_ + ( + IOobject + ( + IOobject::groupName("k", U.group()), + this->runTime_.timeName(), + this->mesh_, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + this->mesh_ + ), + epsilon_ + ( + IOobject + ( + IOobject::groupName("epsilon", U.group()), + this->runTime_.timeName(), + this->mesh_, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + this->mesh_ + ) +{ + bound(k_, this->kMin_); + bound(epsilon_, this->epsilonMin_); + + if (type == typeName) + { + correctNut(); + this->printCoeffs(type); + } +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template +bool RNGkEpsilon::read() +{ + if (eddyViscosity >::read()) + { + Cmu_.readIfPresent(this->coeffDict()); + C1_.readIfPresent(this->coeffDict()); + C2_.readIfPresent(this->coeffDict()); + C3_.readIfPresent(this->coeffDict()); + sigmak_.readIfPresent(this->coeffDict()); + sigmaEps_.readIfPresent(this->coeffDict()); + eta0_.readIfPresent(this->coeffDict()); + beta_.readIfPresent(this->coeffDict()); + + return true; + } + else + { + return false; + } +} + + +template +void RNGkEpsilon::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& nut = this->nut_; + + eddyViscosity >::correct(); + + volScalarField divU(fvc::div(fvc::absolute(this->phi(), U))); + + tmp tgradU = fvc::grad(U); + volScalarField S2((tgradU() && dev(twoSymm(tgradU())))); + tgradU.clear(); + + volScalarField G(this->GName(), nut*S2); + + volScalarField eta(sqrt(mag(S2))*k_/epsilon_); + volScalarField eta3(eta*sqr(eta)); + + volScalarField R + ( + ((eta*(-eta/eta0_ + scalar(1)))/(beta_*eta3 + scalar(1))) + ); + + // Update epsilon and G at the wall + epsilon_.boundaryField().updateCoeffs(); + + // Dissipation equation + tmp epsEqn + ( + fvm::ddt(alpha, rho, epsilon_) + + fvm::div(alphaRhoPhi, epsilon_) + - fvm::laplacian(alpha*rho*DepsilonEff(), epsilon_) + == + (C1_ - R)*alpha*rho*G*epsilon_/k_ + - fvm::SuSp(((2.0/3.0)*C1_ + C3_)*alpha*rho*divU, epsilon_) + - fvm::Sp(C2_*alpha*rho*epsilon_/k_, epsilon_) + + epsilonSource() + ); + + epsEqn().relax(); + + epsEqn().boundaryManipulate(epsilon_.boundaryField()); + + solve(epsEqn); + bound(epsilon_, this->epsilonMin_); + + + // Turbulent kinetic energy equation + + tmp kEqn + ( + fvm::ddt(alpha, rho, k_) + + fvm::div(alphaRhoPhi, k_) + - fvm::Sp(fvc::ddt(alpha, rho) + fvc::div(alphaRhoPhi), k_) + - fvm::laplacian(alpha*rho*DkEff(), k_) + == + alpha*rho*G + - fvm::SuSp((2.0/3.0)*alpha*rho*divU, k_) + - fvm::Sp(alpha*rho*epsilon_/k_, k_) + + kSource() + ); + + kEqn().relax(); + solve(kEqn); + bound(k_, this->kMin_); + + correctNut(); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/TurbulenceModels/turbulenceModels/RAS/RNGkEpsilon/RNGkEpsilon.H b/src/TurbulenceModels/turbulenceModels/RAS/RNGkEpsilon/RNGkEpsilon.H new file mode 100644 index 000000000..315563713 --- /dev/null +++ b/src/TurbulenceModels/turbulenceModels/RAS/RNGkEpsilon/RNGkEpsilon.H @@ -0,0 +1,221 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +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::RASModels::RNGkEpsilon + +Group + grpRASTurbulence + +Description + Renormalization group k-epsilon turbulence model for incompressible and + compressible flows. + + Reference: + \verbatim + Yakhot, V. S. A. S. T. B. C. G., Orszag, S. A., Thangam, S., + Gatski, T. B., & Speziale, C. G. (1992). + Development of turbulence models for shear flows + by a double expansion technique. + Physics of Fluids A: Fluid Dynamics (1989-1993), 4(7), 1510-1520. + + For the RDT-based compression term: + El Tahry, S. H. (1983). + k-epsilon equation for compressible reciprocating engine flows. + Journal of Energy, 7(4), 345-353. + \endverbatim + + The default model coefficients correspond to the following: + \verbatim + RNGkEpsilonCoeffs + { + Cmu 0.0845; + C1 1.42; + C2 1.68; + C3 -0.33; + sigmak 0.71942; + sigmaEps 0.71942; + eta0 4.38; + beta 0.012; + } + \endverbatim + +SourceFiles + RNGkEpsilon.C + +\*---------------------------------------------------------------------------*/ + +#ifndef RNGkEpsilon_H +#define RNGkEpsilon_H + +#include "RASModel.H" +#include "eddyViscosity.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace RASModels +{ + +/*---------------------------------------------------------------------------*\ + Class RNGkEpsilon Declaration +\*---------------------------------------------------------------------------*/ + +template +class RNGkEpsilon +: + public eddyViscosity > +{ + // Private Member Functions + + // Disallow default bitwise copy construct and assignment + RNGkEpsilon(const RNGkEpsilon&); + RNGkEpsilon& operator=(const RNGkEpsilon&); + + +protected: + + // Protected data + + // Model coefficients + + dimensionedScalar Cmu_; + dimensionedScalar C1_; + dimensionedScalar C2_; + dimensionedScalar C3_; + dimensionedScalar sigmak_; + dimensionedScalar sigmaEps_; + dimensionedScalar eta0_; + dimensionedScalar beta_; + + + // Fields + + volScalarField k_; + volScalarField epsilon_; + + + // Protected Member Functions + + virtual void correctNut(); + virtual tmp kSource() const; + virtual tmp epsilonSource() const; + + +public: + + typedef typename BasicTurbulenceModel::alphaField alphaField; + typedef typename BasicTurbulenceModel::rhoField rhoField; + typedef typename BasicTurbulenceModel::transportModel transportModel; + + + //- Runtime type information + TypeName("RNGkEpsilon"); + + + // Constructors + + //- Construct from components + RNGkEpsilon + ( + const alphaField& alpha, + const rhoField& rho, + const volVectorField& U, + const surfaceScalarField& alphaRhoPhi, + const surfaceScalarField& phi, + const transportModel& transport, + const word& propertiesName = turbulenceModel::propertiesName, + const word& type = typeName + ); + + + //- Destructor + virtual ~RNGkEpsilon() + {} + + + // Member Functions + + //- Re-read model coefficients if they have changed + virtual bool read(); + + //- Return the effective diffusivity for k + tmp DkEff() const + { + return tmp + ( + new volScalarField + ( + "DkEff", + (this->nut_/sigmak_ + this->nu()) + ) + ); + } + + //- Return the effective diffusivity for epsilon + tmp DepsilonEff() const + { + return tmp + ( + new volScalarField + ( + "DepsilonEff", + (this->nut_/sigmaEps_ + this->nu()) + ) + ); + } + + //- Return the turbulence kinetic energy + virtual tmp k() const + { + return k_; + } + + //- Return the turbulence kinetic energy dissipation rate + virtual tmp epsilon() const + { + return epsilon_; + } + + //- Solve the turbulence equations and correct the turbulence viscosity + virtual void correct(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace RASModels +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository +# include "RNGkEpsilon.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/TurbulenceModels/turbulenceModels/RAS/kEpsilon/kEpsilon.C b/src/TurbulenceModels/turbulenceModels/RAS/kEpsilon/kEpsilon.C index f37b906a2..1543385ff 100644 --- a/src/TurbulenceModels/turbulenceModels/RAS/kEpsilon/kEpsilon.C +++ b/src/TurbulenceModels/turbulenceModels/RAS/kEpsilon/kEpsilon.C @@ -245,7 +245,6 @@ void kEpsilon::correct() ( fvm::ddt(alpha, rho, epsilon_) + fvm::div(alphaRhoPhi, epsilon_) - - fvm::Sp(fvc::ddt(alpha, rho) + fvc::div(alphaRhoPhi), epsilon_) - fvm::laplacian(alpha*rho*DepsilonEff(), epsilon_) == C1_*alpha*rho*G*epsilon_/k_ @@ -268,7 +267,6 @@ void kEpsilon::correct() ( fvm::ddt(alpha, rho, k_) + fvm::div(alphaRhoPhi, k_) - - fvm::Sp(fvc::ddt(alpha, rho) + fvc::div(alphaRhoPhi), k_) - fvm::laplacian(alpha*rho*DkEff(), k_) == alpha*rho*G