From e38e5d67171471a1d6b345eee24803e5c2ace80e Mon Sep 17 00:00:00 2001
From: Kutalmis Bercin <kutalmis.bercin@esi-group.com>
Date: Tue, 22 Oct 2019 16:58:04 +0100
Subject: [PATCH] ENH: new RAS model: kEpsilonPhitF
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit
ENH: modify fWallFunction for kEpsilonPhitF model
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:
- Turbulent kinetic energy, \c k,
- Turbulent kinetic energy dissipation rate, \c epsilon,
- Normalised wall-normal fluctuating velocity scale, \c phit,
- Elliptic relaxation factor, \c f.
Reference:
\verbatim
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
\endverbatim
The default model coefficients are (LUU:Eqs. 19-20):
\verbatim
kEpsilonPhitFCoeffs
{
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
}
\endverbatim
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'
---
.../turbulentFluidThermoModels.C | 3 +
.../turbulentTransportModels.C | 3 +
.../turbulenceModels/Make/files | 1 +
.../RAS/kEpsilonPhitF/kEpsilonPhitF.C | 513 ++++++++++++++++++
.../RAS/kEpsilonPhitF/kEpsilonPhitF.H | 296 ++++++++++
.../RAS/kEpsilonPhitF/kEpsilonPhitFBase.C | 43 ++
.../RAS/kEpsilonPhitF/kEpsilonPhitFBase.H | 98 ++++
.../fWallFunctionFvPatchScalarField.C | 71 ++-
.../fWallFunctionFvPatchScalarField.H | 36 +-
9 files changed, 1022 insertions(+), 42 deletions(-)
create mode 100644 src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitF.C
create mode 100644 src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitF.H
create mode 100644 src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitFBase.C
create mode 100644 src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitFBase.H
diff --git a/src/TurbulenceModels/compressible/turbulentFluidThermoModels/turbulentFluidThermoModels.C b/src/TurbulenceModels/compressible/turbulentFluidThermoModels/turbulentFluidThermoModels.C
index 1887e54180a..86a19cd1336 100644
--- a/src/TurbulenceModels/compressible/turbulentFluidThermoModels/turbulentFluidThermoModels.C
+++ b/src/TurbulenceModels/compressible/turbulentFluidThermoModels/turbulentFluidThermoModels.C
@@ -94,6 +94,9 @@ makeRASModel(LRR);
#include "SSG.H"
makeRASModel(SSG);
+#include "kEpsilonPhitF.H"
+makeRASModel(kEpsilonPhitF);
+
// -------------------------------------------------------------------------- //
// LES models
diff --git a/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C b/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C
index 2d846cddc67..6004df24d32 100644
--- a/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C
+++ b/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C
@@ -90,6 +90,9 @@ makeRASModel(LRR);
#include "SSG.H"
makeRASModel(SSG);
+#include "kEpsilonPhitF.H"
+makeRASModel(kEpsilonPhitF);
+
// -------------------------------------------------------------------------- //
// LES models
diff --git a/src/TurbulenceModels/turbulenceModels/Make/files b/src/TurbulenceModels/turbulenceModels/Make/files
index 7a6ad662947..8a3c77f09e5 100644
--- a/src/TurbulenceModels/turbulenceModels/Make/files
+++ b/src/TurbulenceModels/turbulenceModels/Make/files
@@ -1,5 +1,6 @@
turbulenceModel.C
RAS/v2f/v2fBase.C
+RAS/kEpsilonPhitF/kEpsilonPhitFBase.C
LESdelta = LES/LESdeltas
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitF.C b/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitF.C
new file mode 100644
index 00000000000..2da4eee3edf
--- /dev/null
+++ b/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitF.C
@@ -0,0 +1,513 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2019 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "kEpsilonPhitF.H"
+#include "fvOptions.H"
+#include "bound.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace RASModels
+{
+
+// * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * * //
+
+template<class BasicTurbulenceModel>
+void kEpsilonPhitF<BasicTurbulenceModel>::correctNut()
+{
+ // (LUU:p. 173)
+ this->nut_ = Cmu_*phit_*k_*T_;
+ this->nut_.correctBoundaryConditions();
+ fv::options::New(this->mesh_).correct(this->nut_);
+
+ BasicTurbulenceModel::correctNut();
+}
+
+
+template<class BasicTurbulenceModel>
+volScalarField kEpsilonPhitF<BasicTurbulenceModel>::Ts() const
+{
+ // (LUU:Eq. 7)
+ return
+ max
+ (
+ k_/epsilon_,
+ CT_*sqrt
+ (
+ max
+ (
+ this->nu(),
+ dimensionedScalar(this->nu()().dimensions(), Zero)
+ )/epsilon_
+ )
+ );
+}
+
+
+template<class BasicTurbulenceModel>
+volScalarField kEpsilonPhitF<BasicTurbulenceModel>::Ls() const
+{
+ // (LUU:Eq. 7)
+ return
+ CL_*max
+ (
+ pow(k_, 1.5)/epsilon_,
+ Ceta_*pow025
+ (
+ pow3
+ (
+ max
+ (
+ this->nu(),
+ dimensionedScalar(this->nu()().dimensions(), Zero)
+ )
+ )/epsilon_
+ )
+ );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+template<class BasicTurbulenceModel>
+kEpsilonPhitF<BasicTurbulenceModel>::kEpsilonPhitF
+(
+ 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<RASModel<BasicTurbulenceModel>>
+ (
+ type,
+ alpha,
+ rho,
+ U,
+ alphaRhoPhi,
+ phi,
+ transport,
+ propertiesName
+ ),
+ kEpsilonPhitFBase(),
+
+ Cmu_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "Cmu",
+ this->coeffDict_,
+ 0.22
+ )
+ ),
+ Ceps1a_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "Ceps1a",
+ this->coeffDict_,
+ 1.4
+ )
+ ),
+ Ceps1b_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "Ceps1b",
+ this->coeffDict_,
+ 1.0
+ )
+ ),
+ Ceps1c_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "Ceps1c",
+ this->coeffDict_,
+ 0.05
+ )
+ ),
+ Ceps2_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "Ceps2",
+ this->coeffDict_,
+ 1.9
+ )
+ ),
+ Cf1_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "Cf1",
+ this->coeffDict_,
+ 1.4
+ )
+ ),
+ Cf2_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "Cf2",
+ this->coeffDict_,
+ 0.3
+ )
+ ),
+ CL_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "CL",
+ this->coeffDict_,
+ 0.25
+ )
+ ),
+ Ceta_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "Ceta",
+ this->coeffDict_,
+ 110.0
+ )
+ ),
+ CT_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "CT",
+ this->coeffDict_,
+ 6.0
+ )
+ ),
+ sigmaK_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "sigmaK",
+ this->coeffDict_,
+ 1.0
+ )
+ ),
+ sigmaEps_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "sigmaEps",
+ this->coeffDict_,
+ 1.3
+ )
+ ),
+ sigmaPhit_
+ (
+ dimensioned<scalar>::getOrAddToDict
+ (
+ "sigmaPhit",
+ this->coeffDict_,
+ 1.0
+ )
+ ),
+
+ k_
+ (
+ IOobject
+ (
+ IOobject::groupName("k", alphaRhoPhi.group()),
+ this->runTime_.timeName(),
+ this->mesh_,
+ IOobject::MUST_READ,
+ IOobject::AUTO_WRITE
+ ),
+ this->mesh_
+ ),
+ epsilon_
+ (
+ IOobject
+ (
+ IOobject::groupName("epsilon", alphaRhoPhi.group()),
+ this->runTime_.timeName(),
+ this->mesh_,
+ IOobject::MUST_READ,
+ IOobject::AUTO_WRITE
+ ),
+ this->mesh_
+ ),
+ phit_
+ (
+ IOobject
+ (
+ IOobject::groupName("phit", alphaRhoPhi.group()),
+ this->runTime_.timeName(),
+ this->mesh_,
+ IOobject::MUST_READ,
+ IOobject::AUTO_WRITE
+ ),
+ this->mesh_
+ ),
+ f_
+ (
+ IOobject
+ (
+ IOobject::groupName("f", alphaRhoPhi.group()),
+ this->runTime_.timeName(),
+ this->mesh_,
+ IOobject::MUST_READ,
+ IOobject::AUTO_WRITE
+ ),
+ this->mesh_
+ ),
+ T_
+ (
+ IOobject
+ (
+ "T",
+ this->runTime_.timeName(),
+ this->mesh_,
+ IOobject::NO_READ,
+ IOobject::NO_WRITE,
+ false
+ ),
+ this->mesh_,
+ dimensionedScalar(dimTime, Zero)
+ ),
+
+ phitMin_(dimensionedScalar("phitMin", phit_.dimensions(), SMALL)),
+ fMin_(dimensionedScalar("fMin", f_.dimensions(), SMALL)),
+ TMin_(dimensionedScalar("TMin", dimTime, SMALL)),
+ L2Min_(dimensionedScalar("L2Min", sqr(dimLength), SMALL))
+{
+ bound(k_, this->kMin_);
+ bound(epsilon_, this->epsilonMin_);
+ bound(phit_, phitMin_);
+ bound(f_, fMin_);
+
+ if (type == typeName)
+ {
+ this->printCoeffs(type);
+ }
+
+ if
+ (
+ mag(sigmaK_.value()) < VSMALL
+ || mag(sigmaEps_.value()) < VSMALL
+ || mag(sigmaPhit_.value()) < VSMALL
+ )
+ {
+ FatalErrorInFunction
+ << "Non-zero values are required for the model constants:" << nl
+ << "sigmaK = " << sigmaK_ << nl
+ << "sigmaEps = " << sigmaEps_ << nl
+ << "sigmaPhit = " << sigmaPhit_ << nl
+ << exit(FatalError);
+ }
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+template<class BasicTurbulenceModel>
+bool kEpsilonPhitF<BasicTurbulenceModel>::read()
+{
+ if (eddyViscosity<RASModel<BasicTurbulenceModel>>::read())
+ {
+ Cmu_.readIfPresent(this->coeffDict());
+ Ceps1a_.readIfPresent(this->coeffDict());
+ Ceps1b_.readIfPresent(this->coeffDict());
+ Ceps1c_.readIfPresent(this->coeffDict());
+ Ceps2_.readIfPresent(this->coeffDict());
+ Cf1_.readIfPresent(this->coeffDict());
+ Cf2_.readIfPresent(this->coeffDict());
+ CL_.readIfPresent(this->coeffDict());
+ Ceta_.readIfPresent(this->coeffDict());
+ CT_.readIfPresent(this->coeffDict());
+ sigmaK_.readIfPresent(this->coeffDict());
+ sigmaEps_.readIfPresent(this->coeffDict());
+ sigmaPhit_.readIfPresent(this->coeffDict());
+
+ return true;
+ }
+
+ return false;
+}
+
+
+template<class BasicTurbulenceModel>
+void kEpsilonPhitF<BasicTurbulenceModel>::correct()
+{
+ if (!this->turbulence_)
+ {
+ return;
+ }
+
+ // Construct local convenience 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_;
+ fv::options& fvOptions(fv::options::New(this->mesh_));
+
+ eddyViscosity<RASModel<BasicTurbulenceModel>>::correct();
+
+ const volScalarField::Internal divU
+ (
+ fvc::div(fvc::absolute(this->phi(), U))
+ );
+
+ tmp<volSymmTensorField> tS(symm(fvc::grad(U)));
+ volScalarField G(this->GName(), nut*(2.0*(dev(tS()) && tS())));
+ tS.clear();
+
+ T_ = Ts();
+ bound(T_, TMin_);
+
+ const volScalarField L2(type() + "L2", sqr(Ls()) + L2Min_);
+
+ const volScalarField::Internal Ceps1
+ (
+ "Ceps1",
+ Ceps1a_*(Ceps1b_ + Ceps1c_*sqrt(1.0/phit_()))
+ );
+
+
+ // Update epsilon and G at the wall
+ epsilon_.boundaryFieldRef().updateCoeffs();
+
+ // Turbulent kinetic energy dissipation rate equation (LUU:Eq. 4)
+ // k/T ~ epsilon
+ tmp<fvScalarMatrix> epsEqn
+ (
+ fvm::ddt(alpha, rho, epsilon_)
+ + fvm::div(alphaRhoPhi, epsilon_)
+ - fvm::laplacian(alpha*rho*DepsilonEff(), epsilon_)
+ ==
+ alpha()*rho()*Ceps1*G()/T_()
+ - fvm::SuSp
+ (
+ (2.0/3.0*Ceps1)*(alpha()*rho()*divU),
+ epsilon_
+ )
+ - fvm::Sp(alpha()*rho()*Ceps2_/T_(), epsilon_)
+ + fvOptions(alpha, rho, epsilon_)
+ );
+
+ epsEqn.ref().relax();
+ fvOptions.constrain(epsEqn.ref());
+ epsEqn.ref().boundaryManipulate(epsilon_.boundaryFieldRef());
+ solve(epsEqn);
+ fvOptions.correct(epsilon_);
+ bound(epsilon_, this->epsilonMin_);
+
+
+ // Turbulent kinetic energy equation (LUU:Eq. 3)
+ // epsilon/k ~ 1/Ts
+ tmp<fvScalarMatrix> kEqn
+ (
+ fvm::ddt(alpha, rho, k_)
+ + fvm::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()/T_(), k_)
+ + fvOptions(alpha, rho, k_)
+ );
+
+ kEqn.ref().relax();
+ fvOptions.constrain(kEqn.ref());
+ solve(kEqn);
+ fvOptions.correct(k_);
+ bound(k_, this->kMin_);
+
+
+ // Elliptic relaxation function equation (LUU:Eq. 18)
+ // All source terms are non-negative functions (LUU:p. 176)
+ tmp<fvScalarMatrix> fEqn
+ (
+ - fvm::laplacian(f_)
+ ==
+ - fvm::Sp(1.0/L2(), f_)
+ - (
+ (Cf1_ - 1.0)*(phit_() - 2.0/3.0)/T_()
+ -(Cf2_*G())/k_()
+ +(Cf2_*(2.0/3.0)*divU)
+ -(2.0*this->nu()*(fvc::grad(phit_) & fvc::grad(k_)))()/k_()
+ -(this->nu()*fvc::laplacian(phit_))()
+ )/L2()
+ );
+
+ fEqn.ref().relax();
+ solve(fEqn);
+ bound(f_, fMin_);
+
+
+ // Normalised wall-normal fluctuating velocity scale equation (LUU:Eq. 17)
+ // All source terms are non-negative functions (LUU:p. 176)
+ tmp<fvScalarMatrix> phitEqn
+ (
+ fvm::ddt(alpha, rho, phit_)
+ + fvm::div(alphaRhoPhi, phit_)
+ - fvm::laplacian(alpha*rho*DphitEff(), phit_)
+ ==
+ alpha()*rho()*f_()
+ - fvm::SuSp
+ (
+ alpha()*rho()*
+ (
+ G()/k_()
+ - (2.0/3.0)*divU
+ - (2.0*nut*(fvc::grad(phit_) & fvc::grad(k_)))()
+ /(k_()*sigmaPhit_*phit_())
+ )
+ , phit_
+ )
+ + fvOptions(alpha, rho, phit_)
+ );
+
+ phitEqn.ref().relax();
+ fvOptions.constrain(phitEqn.ref());
+ solve(phitEqn);
+ fvOptions.correct(phit_);
+ bound(phit_, phitMin_);
+
+ correctNut();
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace RASModels
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitF.H b/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitF.H
new file mode 100644
index 00000000000..18d07e82906
--- /dev/null
+++ b/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitF.H
@@ -0,0 +1,296 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2019 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 <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::RASModels::kEpsilonPhitF
+
+Group
+ grpRASTurbulence
+
+Description
+ 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:
+ - Turbulent kinetic energy, \c k,
+ - Turbulent kinetic energy dissipation rate, \c epsilon,
+ - Normalised wall-normal fluctuating velocity scale, \c phit,
+ - Elliptic relaxation factor, \c f.
+
+ Reference:
+ \verbatim
+ 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
+ \endverbatim
+
+ The default model coefficients are (LUU:Eqs. 19-20):
+ \verbatim
+ kEpsilonPhitFCoeffs
+ {
+ 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
+ }
+ \endverbatim
+
+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.
+
+SourceFiles
+ kEpsilonPhitF.C
+
+SeeAlso
+ v2f.H
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef kEpsilonPhitF_H
+#define kEpsilonPhitF_H
+
+#include "kEpsilonPhitFBase.H"
+#include "RASModel.H"
+#include "eddyViscosity.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace RASModels
+{
+
+/*---------------------------------------------------------------------------*\
+ Class kEpsilonPhitF Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class BasicTurbulenceModel>
+class kEpsilonPhitF
+:
+ public eddyViscosity<RASModel<BasicTurbulenceModel>>,
+ public kEpsilonPhitFBase
+{
+ // Private Member Functions
+
+ //- No copy construct
+ kEpsilonPhitF(const kEpsilonPhitF&) = delete;
+
+ //- No copy assignment
+ void operator=(const kEpsilonPhitF&) = delete;
+
+
+protected:
+
+ // Protected Data
+
+ // Model coefficients
+
+ dimensionedScalar Cmu_;
+ dimensionedScalar Ceps1a_;
+ dimensionedScalar Ceps1b_;
+ dimensionedScalar Ceps1c_;
+ dimensionedScalar Ceps2_;
+ dimensionedScalar Cf1_;
+ dimensionedScalar Cf2_;
+ dimensionedScalar CL_;
+ dimensionedScalar Ceta_;
+ dimensionedScalar CT_;
+ dimensionedScalar sigmaK_;
+ dimensionedScalar sigmaEps_;
+ dimensionedScalar sigmaPhit_;
+
+
+ // Fields
+
+ //- Turbulent kinetic energy [m2/s2]
+ volScalarField k_;
+
+ //- Turbulent kinetic energy dissipation rate [m2/s3]
+ volScalarField epsilon_;
+
+ //- Normalised wall-normal fluctuating velocity scale [-]
+ volScalarField phit_;
+
+ //- Elliptic relaxation factor [1/s]
+ volScalarField f_;
+
+ //- Turbulent time scale [s]
+ volScalarField T_;
+
+
+ // Bounding values
+
+ dimensionedScalar phitMin_;
+ dimensionedScalar fMin_;
+ dimensionedScalar TMin_;
+ dimensionedScalar L2Min_;
+
+
+ // Protected Member Functions
+
+ //- Update nut with the latest available k, phit, and T
+ virtual void correctNut();
+
+ //- Return the turbulent time scale, T
+ volScalarField Ts() const;
+
+ //- Return the turbulent length scale, L
+ volScalarField Ls() const;
+
+
+public:
+
+ typedef typename BasicTurbulenceModel::alphaField alphaField;
+ typedef typename BasicTurbulenceModel::rhoField rhoField;
+ typedef typename BasicTurbulenceModel::transportModel transportModel;
+
+
+ //- Runtime type information
+ TypeName("kEpsilonPhitF");
+
+
+ // Constructors
+
+ //- Construct from components
+ 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
+ );
+
+
+ //- Destructor
+ virtual ~kEpsilonPhitF() = default;
+
+
+ // Member Functions
+
+ //- Re-read model coefficients if they have changed
+ virtual bool read();
+
+ //- Return the effective diffusivity for k
+ tmp<volScalarField> DkEff() const
+ {
+ return tmp<volScalarField>
+ (
+ new volScalarField
+ (
+ "DkEff",
+ this->nut_/sigmaK_ + this->nu()
+ )
+ );
+ }
+
+ //- Return the effective diffusivity for epsilon
+ tmp<volScalarField> DepsilonEff() const
+ {
+ return tmp<volScalarField>
+ (
+ new volScalarField
+ (
+ "DepsilonEff",
+ this->nut_/sigmaEps_ + this->nu()
+ )
+ );
+ }
+
+ //- Return the effective diffusivity for phit
+ tmp<volScalarField> DphitEff() const
+ {
+ return tmp<volScalarField>
+ (
+ new volScalarField
+ (
+ "DphitEff",
+ this->nut_/sigmaPhit_ + this->nu()
+ )
+ );
+ }
+
+ //- Return the turbulent kinetic energy field
+ virtual tmp<volScalarField> k() const
+ {
+ return k_;
+ }
+
+ //- Return the turbulent kinetic energy dissipation rate field
+ virtual tmp<volScalarField> epsilon() const
+ {
+ return epsilon_;
+ }
+
+ //- Return the normalised wall-normal fluctuating velocity scale field
+ virtual tmp<volScalarField> phit() const
+ {
+ return phit_;
+ }
+
+ //- Return the elliptic relaxation factor field
+ virtual tmp<volScalarField> f() const
+ {
+ return f_;
+ }
+
+ //- Solve the transport equations and correct the turbulent viscosity
+ virtual void correct();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace RASModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+ #include "kEpsilonPhitF.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
\ No newline at end of file
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitFBase.C b/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitFBase.C
new file mode 100644
index 00000000000..e4a56e54bbe
--- /dev/null
+++ b/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitFBase.C
@@ -0,0 +1,43 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2019 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "kEpsilonPhitFBase.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace RASModels
+{
+
+ defineTypeNameAndDebug(kEpsilonPhitFBase, 0);
+
+} // End namespace RASModels
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitFBase.H b/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitFBase.H
new file mode 100644
index 00000000000..622bc3254cd
--- /dev/null
+++ b/src/TurbulenceModels/turbulenceModels/RAS/kEpsilonPhitF/kEpsilonPhitFBase.H
@@ -0,0 +1,98 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2019 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 <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::RASModels::kEpsilonPhitFBase
+
+Group
+ grpRASTurbulence
+
+Description
+ Abstract base-class for the \c k-epsilon-phit-f model to provide boundary
+ condition access to the \c phit and \c f fields.
+
+See also
+ Foam::RASModels::kEpsilonPhitF
+
+SourceFiles
+ kEpsilonPhitFBase.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef kEpsilonPhitFBase_H
+#define kEpsilonPhitFBase_H
+
+#include "RASModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace RASModels
+{
+
+/*---------------------------------------------------------------------------*\
+ Class kEpsilonPhitFBase Declaration
+\*---------------------------------------------------------------------------*/
+
+class kEpsilonPhitFBase
+{
+public:
+
+ //- Runtime type information
+ TypeName("kEpsilonPhitFBase");
+
+
+ // Constructors
+
+ kEpsilonPhitFBase()
+ {}
+
+
+ //- Destructor
+ virtual ~kEpsilonPhitFBase()
+ {}
+
+
+ // Member Functions
+
+ //- Return the normalised wall-normal fluctuating velocity scale field
+ virtual tmp<volScalarField> phit() const = 0;
+
+ //- Return the elliptic relaxation factor field
+ virtual tmp<volScalarField> f() const = 0;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace RASModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/fWallFunctions/fWallFunction/fWallFunctionFvPatchScalarField.C b/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/fWallFunctions/fWallFunction/fWallFunctionFvPatchScalarField.C
index da3175bf513..3573bcd3190 100644
--- a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/fWallFunctions/fWallFunction/fWallFunctionFvPatchScalarField.C
+++ b/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/fWallFunctions/fWallFunction/fWallFunctionFvPatchScalarField.C
@@ -29,6 +29,7 @@ License
#include "fWallFunctionFvPatchScalarField.H"
#include "nutWallFunctionFvPatchScalarField.H"
#include "v2f.H"
+#include "kEpsilonPhitF.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -111,7 +112,6 @@ void fWallFunctionFvPatchScalarField::updateCoeffs()
internalField().group()
)
);
- const v2fBase& v2fModel = refCast<const v2fBase>(turbModel);
const nutWallFunctionFvPatchScalarField& nutw =
nutWallFunctionFvPatchScalarField::nutw(turbModel, patchi);
@@ -121,43 +121,60 @@ void fWallFunctionFvPatchScalarField::updateCoeffs()
const tmp<scalarField> tnuw = turbModel.nu(patchi);
const scalarField& nuw = tnuw();
- const tmp<volScalarField> tk = turbModel.k();
- const volScalarField& k = tk();
+ scalarField& f = *this;
- const tmp<volScalarField> tepsilon = turbModel.epsilon();
- const volScalarField& epsilon = tepsilon();
+ if (isA<v2fBase>(turbModel))
+ {
+ const v2fBase& v2fModel = refCast<const v2fBase>(turbModel);
- const tmp<volScalarField> tv2 = v2fModel.v2();
- const volScalarField& v2 = tv2();
+ const tmp<volScalarField> tk = turbModel.k();
+ const volScalarField& k = tk();
- const scalar Cmu25 = pow025(nutw.Cmu());
- const scalar N = 6.0;
+ const tmp<volScalarField> tepsilon = turbModel.epsilon();
+ const volScalarField& epsilon = tepsilon();
- scalarField& f = *this;
+ const tmp<volScalarField> tv2 = v2fModel.v2();
+ const volScalarField& v2 = tv2();
- // Set f wall values
- forAll(f, facei)
- {
- const label celli = patch().faceCells()[facei];
+ const scalar Cmu25 = pow025(nutw.Cmu());
+ const scalar N = 6.0;
- const scalar uTau = Cmu25*sqrt(k[celli]);
+ // Set f wall values
+ forAll(f, facei)
+ {
+ const label celli = patch().faceCells()[facei];
- const scalar yPlus = uTau*y[facei]/nuw[facei];
+ const scalar uTau = Cmu25*sqrt(k[celli]);
- if (nutw.yPlusLam() < yPlus)
- {
- const scalar v2c = v2[celli];
- const scalar epsc = epsilon[celli];
- const scalar kc = k[celli];
+ const scalar yPlus = uTau*y[facei]/nuw[facei];
- f[facei] = N*v2c*epsc/(sqr(kc) + ROOTVSMALL);
- f[facei] /= sqr(uTau) + ROOTVSMALL;
- }
- else
- {
- f[facei] = 0.0;
+ if (nutw.yPlusLam() < yPlus)
+ {
+ const scalar v2c = v2[celli];
+ const scalar epsc = epsilon[celli];
+ const scalar kc = k[celli];
+
+ f[facei] = N*v2c*epsc/(sqr(kc) + ROOTVSMALL);
+ f[facei] /= sqr(uTau) + ROOTVSMALL;
+ }
+ else
+ {
+ f[facei] = 0.0;
+ }
}
}
+ else if (isA<kEpsilonPhitFBase>(turbModel))
+ {
+ // (LUU:p. 176)
+ f = 0.0;
+ }
+ else
+ {
+ FatalErrorInFunction
+ << "The RAS model is neither the v2f nor kEpsilonPhitF model. "
+ << "Therefore, fWallFunction is not usable." << nl
+ << exit(FatalError);
+ }
fixedValueFvPatchField<scalar>::updateCoeffs();
diff --git a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/fWallFunctions/fWallFunction/fWallFunctionFvPatchScalarField.H b/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/fWallFunctions/fWallFunction/fWallFunctionFvPatchScalarField.H
index 2462d8226fd..88fca25ff1e 100644
--- a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/fWallFunctions/fWallFunction/fWallFunctionFvPatchScalarField.H
+++ b/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/fWallFunctions/fWallFunction/fWallFunctionFvPatchScalarField.H
@@ -32,24 +32,24 @@ Group
Description
This boundary condition provides a wall constraint on the elliptic
- relaxation factor, \c f, which is executed in the \c v2-f eddy viscosity
- turbulence model. The condition is applicable for low- and high-Reynolds
- number turbulent flow cases.
-
- For \c f, the viscous sublayer and log-law region blending approaches are
- claimed to be inviable (Popovac and Hanjalić (2007), p. 194). Therefore,
- the only boundary condition blending mode is the stepwise mode
- where the viscous sublayer and log-law region contributions switch over
- \c y+ value derived from the \c kappa and \c E.
+ relaxation factor, \c f, which is executed in the \c v2-f based
+ turbulence closure models for low- and high-Reynolds number
+ turbulent flow cases.
Reference:
\verbatim
- Remark on the blending approach:
- Popovac, M., & Hanjalić, K. (2007).
- Compound wall treatment for RANS computation of complex
- turbulent flows and heat transfer.
- Flow, turbulence and combustion, 78(2), 177-202.
- doi:10.1007/s10494-006-9067-x
+ Remark on the blending approach for f (tag:PH):
+ Popovac, M., & Hanjalić, K. (2007).
+ Compound wall treatment for RANS computation of complex
+ turbulent flows and heat transfer.
+ Flow, turbulence and combustion, 78(2), 177-202.
+ DOI:10.1007/s10494-006-9067-x
+
+ Wall-boundary expression for f in kEpsilonPhitF 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
\endverbatim
Usage
@@ -69,6 +69,12 @@ Note
the specified \c nutWallFunction in order to ensure that each patch
possesses the same set of values for these coefficients.
+ For \c f, the viscous and inertial sublayer blending approaches were
+ claimed to be inviable (PH:p. 194). Therefore, the only blending mode
+ for the v2-f model is the stepwise mode where the viscous and inertial
+ sublayer contributions switch over a sublayer-intersection value of
+ \c y+ estimated from the \c kappa and \c E.
+
See also
Foam::fixedValueFvPatchField
--
GitLab