The k-epsilon-phit-f turbulence closure model for incompressible and compressible flows.
More...
|
| TypeName ("kEpsilonPhitF") |
|
| kEpsilonPhitF (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) |
|
virtual | ~kEpsilonPhitF ()=default |
|
virtual bool | read () |
|
tmp< volScalarField > | DkEff () const |
|
tmp< volScalarField > | DepsilonEff () const |
|
tmp< volScalarField > | DphitEff () const |
|
virtual tmp< volScalarField > | k () const |
|
virtual tmp< volScalarField > | epsilon () const |
|
virtual tmp< volScalarField > | phit () const |
|
virtual tmp< volScalarField > | f () const |
|
virtual void | correct () |
|
| eddyViscosity (const word &modelName, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName) |
|
virtual | ~eddyViscosity ()=default |
|
virtual tmp< volScalarField > | nut () const |
|
virtual tmp< scalarField > | nut (const label patchi) const |
|
virtual tmp< volScalarField > | k () const=0 |
|
virtual tmp< volSymmTensorField > | R () const |
|
virtual void | validate () |
|
| linearViscousStress (const word &modelName, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName) |
|
virtual | ~linearViscousStress ()=default |
|
virtual tmp< volSymmTensorField > | devRhoReff () const |
|
virtual tmp< volSymmTensorField > | devRhoReff (const volVectorField &U) const |
|
virtual tmp< fvVectorMatrix > | divDevRhoReff (volVectorField &U) const |
|
virtual tmp< fvVectorMatrix > | divDevRhoReff (const volScalarField &rho, volVectorField &U) const |
|
| TypeName ("RAS") |
|
| declareRunTimeSelectionTable (autoPtr, RASModel, dictionary,(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName),(alpha, rho, U, alphaRhoPhi, phi, transport, propertiesName)) |
|
| RASModel (const word &type, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName) |
|
virtual | ~RASModel ()=default |
|
const dimensionedScalar & | kMin () const |
|
const dimensionedScalar & | epsilonMin () const |
|
const dimensionedScalar & | omegaMin () const |
|
dimensionedScalar & | kMin () |
|
dimensionedScalar & | epsilonMin () |
|
dimensionedScalar & | omegaMin () |
|
virtual const dictionary & | coeffDict () const |
|
virtual tmp< volScalarField > | nuEff () const |
|
virtual tmp< scalarField > | nuEff (const label patchi) const |
|
virtual tmp< volScalarField > | omega () const |
|
template<class BasicTurbulenceModel>
class Foam::RASModels::kEpsilonPhitF< BasicTurbulenceModel >
The k-epsilon-phit-f turbulence closure model for incompressible and compressible flows.
The model is a three-transport-equation linear-eddy-viscosity turbulence closure model alongside an elliptic relaxation equation.
Input fields
k | : | Turbulent kinetic energy [m2/s2] |
epsilon | : | Turbulent kinetic energy dissipation rate [m2/s3] |
phit | : | Normalised wall-normal fluctuating velocity scale [-] |
f | : | Elliptic relaxation factor [1/s] |
Reference:
Standard model (Tag:LUU):
Laurence, D. R., Uribe, J. C., & Utyuzhnikov, S. V. (2005).
A robust formulation of the v2-f model.
Flow, Turbulence and Combustion, 73(3-4), 169–185.
DOI:10.1007/s10494-005-1974-8
The default model coefficients are (LUU:Eqs. 19-20):
kEpsilonPhitFCoeffs
{
includeNu true; // include nu in (LUU: Eq. 17), see Notes
Cmu 0.22; // Turbulent viscosity constant
Ceps1a 1.4; // Model constant for epsilon
Ceps1b 1.0; // Model constant for epsilon
Ceps1c 0.05; // Model constant for epsilon
Ceps2 1.9; // Model constant for epsilon
Cf1 1.4; // Model constant for f
Cf2 0.3; // Model constant for f
CL 0.25; // Model constant for L
Ceta 110.0; // Model constant for L
CT 6.0; // Model constant for T
sigmaK 1.0; // Turbulent Prandtl number for k
sigmaEps 1.3; // Turbulent Prandtl number for epsilon
sigmaPhit 1.0; // Turbulent Prandtl number for phit = sigmaK
}
- Note
- The name of the original variable replacing 'v2' is 'phi' (LUU:Eq. 14). However, the name 'phi' preexisted in OpenFOAM; therefore, this name was replaced by 'phit' herein.
Including nu
in DphitEff
even though it is not present in (LUU:Eq. 17) provided higher level of resemblance to benchmarks for the tests considered, particularly for the peak skin friction (yet, pressure-related predictions were unaffected). Users can switch off nu
in DphitEff
by using includeNu
entry in kEpsilonPhitFCoeffs
as shown above in order to follow the reference paper thereat. includeNu
is left true
by default. See GitLab issue #1560.
- Source files
-
- See also
- kEpsilon.C
Definition at line 127 of file kEpsilonPhitF.H.