zhyang-dev/openfoam
1
0
Fork 0
mirror of https://develop.openfoam.com/Development/openfoam.git synced 2025-11-28 03:28:01 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

ENH: enforce consistent boundness on turbulence models.

Browse Source 
- remove epsilonSmall, omegaSmall
- k0/epsilon0/omega0 become kMin/epsilonMin/omegaMin
- add qMin/zetaMin for consistency

These files still need some attention:
    dynOneEqEddy.C
    NonlinearKEShih.C
    settlingFoam

BUG: incompressible::LESModels:dynOneEqEddy::correct()
- avoid tmp field destruction for consistency with the compressible
  version

Possible TODO:
   - set kMin to zero (instead of SMALL) and introduce kSmall
     to avoid division by zero
This commit is contained in:
Mark Olesen
2010-03-24 09:07:53 +01:00
parent 240587029a
commit 0c8fb634f0
45 changed files with 374 additions and 338 deletions
Download Patch File Download Diff File Expand all files Collapse all files

21
applications/solvers/combustion/PDRFoam/PDRModels/turbulence/PDRkEpsilon/PDRkEpsilon.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -121,9 +121,8 @@ PDRkEpsilon::PDRkEpsilon
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh_
autoCreateK("k", mesh_)
),
epsilon_
(
IOobject
@ -134,9 +133,8 @@ PDRkEpsilon::PDRkEpsilon
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh_
autoCreateEpsilon("epsilon", mesh_)
),
mut_
(
IOobject
@ -147,9 +145,8 @@ PDRkEpsilon::PDRkEpsilon
IOobject::NO_READ,
IOobject::NO_WRITE
),
Cmu_*rho_*sqr(k_)/(epsilon_ + epsilonSmall_)
autoCreateMut("mut", mesh_)
),
alphat_
(
IOobject
@ -163,7 +160,9 @@ PDRkEpsilon::PDRkEpsilon
autoCreateAlphat("alphat", mesh_)
)
{
mut_ = Cmu_*rho_*sqr(k_)/(epsilon_ + epsilonSmall_);
bound(epsilon_, epsilonMin_);
mut_ = Cmu_*rho_*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
alphat_ = mut_/Prt_;
@ -250,7 +249,7 @@ void PDRkEpsilon::correct()
if (!turbulence_)
{
// Re-calculate viscosity
mut_ = rho_*Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
mut_ = rho_*Cmu_*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
// Re-calculate thermal diffusivity
@ -302,7 +301,7 @@ void PDRkEpsilon::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -320,7 +319,7 @@ void PDRkEpsilon::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity
mut_ = rho_*Cmu_*sqr(k_)/epsilon_;

14
applications/solvers/multiphase/settlingFoam/kEpsilon.H
View File

@ -5,8 +5,8 @@ if (turbulence)
y.correct();
}
dimensionedScalar k0("k0", k.dimensions(), SMALL);
dimensionedScalar epsilon0("epsilon0", epsilon.dimensions(), SMALL);
dimensionedScalar kMin("kMin", k.dimensions(), SMALL);
dimensionedScalar epsilonMin("epsilonMin", epsilon.dimensions(), SMALL);
volScalarField divU = fvc::div(phi/fvc::interpolate(rho));
@ -15,7 +15,7 @@ if (turbulence)
tgradU.clear();
volScalarField Gcoef =
Cmu*k/sigmak*(g & fvc::grad(rho))/(epsilon + epsilon0);
Cmu*k/sigmak*(g & fvc::grad(rho))/(epsilon + epsilonMin);
#include "wallFunctions.H"
@ -40,7 +40,7 @@ if (turbulence)
epsEqn.relax();
epsEqn.solve();
bound(epsilon, epsilon0);
bound(epsilon, epsilonMin);
// Turbulent kinetic energy equation
@ -60,11 +60,13 @@ if (turbulence)
- fvm::Sp(rho*epsilon/k, k)
);
bound(k, k0);
//FIXME: why no kEqn.relax()?
bound(k, kMin);
//- Re-calculate viscosity
mut = rho*Cmu*sqr(k)/(epsilon + epsilon0);
mut = rho*Cmu*sqr(k)/(epsilon + epsilonMin);
#include "wallViscosity.H"
}

16
src/finiteVolume/cfdTools/general/bound/bound.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -30,11 +30,11 @@ License
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
void Foam::bound(volScalarField& vsf, const dimensionedScalar& vsf0)
void Foam::bound(volScalarField& vsf, const dimensionedScalar& lowerBound)
{
scalar minVsf = min(vsf).value();
const scalar minVsf = min(vsf).value();
if (minVsf < vsf0.value())
if (minVsf < lowerBound.value())
{
Info<< "bounding " << vsf.name()
<< ", min: " << minVsf
@ -47,13 +47,13 @@ void Foam::bound(volScalarField& vsf, const dimensionedScalar& vsf0)
max
(
vsf.internalField(),
fvc::average(max(vsf, vsf0))().internalField()
*pos(-vsf.internalField())
fvc::average(max(vsf, lowerBound))().internalField()
* pos(-vsf.internalField())
),
vsf0.value()
lowerBound.value()
);
vsf.boundaryField() = max(vsf.boundaryField(), vsf0.value());
vsf.boundaryField() = max(vsf.boundaryField(), lowerBound.value());
}
}

12
src/finiteVolume/cfdTools/general/bound/bound.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -23,11 +23,13 @@ License
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
InNamespace
Foam::bound
Foam
Description
Bound the given scalar field if it has gone unbounded. Used extensively
in RAS and LES turbulence models.
Bound the given scalar field if it has gone unbounded.
Used extensively in RAS and LES turbulence models, but also of use
within solvers.
SourceFiles
bound.C
@ -49,7 +51,7 @@ namespace Foam
//- Bound the given scalar field if it has gone unbounded.
// Used extensively in RAS and LES turbulence models.
void bound(volScalarField& vsf, const dimensionedScalar& vsf0);
void bound(volScalarField&, const dimensionedScalar& lowerBound);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

4
src/turbulenceModels/compressible/LES/DeardorffDiffStress/DeardorffDiffStress.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -132,7 +132,7 @@ void DeardorffDiffStress::correct(const tmp<volTensorField>& tgradU)
}
K = 0.5*tr(B_);
bound(K, k0());
bound(K, kMin_);
updateSubGridScaleFields(K);
}

8
src/turbulenceModels/compressible/LES/LESModel/LESModel.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -79,11 +79,11 @@ LESModel::LESModel
printCoeffs_(lookupOrDefault<Switch>("printCoeffs", false)),
coeffDict_(subOrEmptyDict(type + "Coeffs")),
k0_("k0", dimVelocity*dimVelocity, SMALL),
kMin_("kMin", dimVelocity*dimVelocity, SMALL),
delta_(LESdelta::New("delta", U.mesh(), *this))
{
readIfPresent("k0", k0_);
readIfPresent("kMin", kMin_);
// Force the construction of the mesh deltaCoeffs which may be needed
// for the construction of the derived models and BCs
@ -167,7 +167,7 @@ bool LESModel::read()
coeffDict_ <<= *dictPtr;
}
readIfPresent("k0", k0_);
readIfPresent("kMin", kMin_);
delta_().read(*this);

16
src/turbulenceModels/compressible/LES/LESModel/LESModel.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -83,7 +83,7 @@ protected:
Switch printCoeffs_;
dictionary coeffDict_;
dimensionedScalar k0_;
dimensionedScalar kMin_;
autoPtr<LESdelta> delta_;
@ -170,16 +170,16 @@ public:
return coeffDict_;
}
//- Return the value of k0 which k is not allowed to be less than
const dimensionedScalar& k0() const
//- Return the lower allowable limit for k (default: SMALL)
const dimensionedScalar& kMin() const
{
return k0_;
return kMin_;
}
//- Allow k0 to be changed
dimensionedScalar& k0()
//- Allow kMin to be changed
dimensionedScalar& kMin()
{
return k0_;
return kMin_;
}
//- Access function to filter width

14
src/turbulenceModels/compressible/LES/dynOneEqEddy/dynOneEqEddy.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -74,10 +74,10 @@ dimensionedScalar dynOneEqEddy::ce_(const volSymmTensorField& D) const
pow(KK + filter_(k_), 1.5)/(2*delta()) - filter_(pow(k_, 1.5))/delta();
volScalarField ee =
2*delta()*ck_(D)*
(
2*delta()*ck_(D)
*(
filter_(sqrt(k_)*magSqr(D))
- 2*sqrt(KK + filter_(k_))*magSqr(filter_(D))
- 2*sqrt(KK + filter_(k_))*magSqr(filter_(D))
);
return average(ee*mm)/average(mm*mm);
@ -129,7 +129,11 @@ void dynOneEqEddy::correct(const tmp<volTensorField>& tgradU)
- fvm::Sp(ce_(D)*rho()*sqrt(k_)/delta(), k_)
);
bound(k_, dimensionedScalar("0", k_.dimensions(), 1.0e-10));
//FIXME: why not this?
// kEqn.relax();
// kEqn.solve();
bound(k_, kMin_);
updateSubGridScaleFields(D);
}

4
src/turbulenceModels/compressible/LES/lowReOneEqEddy/lowReOneEqEddy.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -117,7 +117,7 @@ void lowReOneEqEddy::correct(const tmp<volTensorField>& tgradU)
- fvm::Sp(ce_*rho()*sqrt(k_)/delta(), k_)
);
bound(k_, k0());
bound(k_, kMin_);
updateSubGridScaleFields();
}

4
src/turbulenceModels/compressible/LES/oneEqEddy/oneEqEddy.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -107,7 +107,7 @@ void oneEqEddy::correct(const tmp<volTensorField>& tgradU)
kEqn.relax();
kEqn.solve();
bound(k_, k0());
bound(k_, kMin_);
updateSubGridScaleFields();
}

18
src/turbulenceModels/compressible/RAS/LRR/LRR.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -229,7 +229,9 @@ LRR::LRR
<< exit(FatalError);
}
mut_ = Cmu_*rho_*sqr(k_)/(epsilon_ + epsilonSmall_);
bound(epsilon_, epsilonMin_);
mut_ = Cmu_*rho_*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
alphat_ = mut_/Prt_;
@ -324,7 +326,7 @@ void LRR::correct()
if (!turbulence_)
{
// Re-calculate viscosity
mut_ = rho_*Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
mut_ = rho_*Cmu_*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
// Re-calculate thermal diffusivity
@ -359,7 +361,7 @@ void LRR::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Reynolds stress equation
@ -406,15 +408,15 @@ void LRR::correct()
R_.dimensions(),
symmTensor
(
k0_.value(), -GREAT, -GREAT,
k0_.value(), -GREAT,
k0_.value()
kMin_.value(), -GREAT, -GREAT,
kMin_.value(), -GREAT,
kMin_.value()
)
)
);
k_ = 0.5*tr(R_);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

18
src/turbulenceModels/compressible/RAS/LaunderGibsonRSTM/LaunderGibsonRSTM.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -258,7 +258,9 @@ LaunderGibsonRSTM::LaunderGibsonRSTM
<< exit(FatalError);
}
mut_ = Cmu_*rho_*sqr(k_)/(epsilon_ + epsilonSmall_);
bound(epsilon_, epsilonMin_);
mut_ = Cmu_*rho_*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
alphat_ = mut_/Prt_;
@ -357,7 +359,7 @@ void LaunderGibsonRSTM::correct()
if (!turbulence_)
{
// Re-calculate viscosity
mut_ = rho_*Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
mut_ = rho_*Cmu_*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
// Re-calculate thermal diffusivity
@ -397,7 +399,7 @@ void LaunderGibsonRSTM::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Reynolds stress equation
@ -453,15 +455,15 @@ void LaunderGibsonRSTM::correct()
R_.dimensions(),
symmTensor
(
k0_.value(), -GREAT, -GREAT,
k0_.value(), -GREAT,
k0_.value()
kMin_.value(), -GREAT, -GREAT,
kMin_.value(), -GREAT,
kMin_.value()
)
)
);
k_ == 0.5*tr(R_);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate turbulent viscosity

14
src/turbulenceModels/compressible/RAS/LaunderSharmaKE/LaunderSharmaKE.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -187,7 +187,9 @@ LaunderSharmaKE::LaunderSharmaKE
autoCreateAlphat("alphat", mesh_)
)
{
mut_ = rho_*Cmu_*fMu()*sqr(k_)/(epsilon_ + epsilonSmall_);
bound(epsilon_, epsilonMin_);
mut_ = rho_*Cmu_*fMu()*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
alphat_ = mut_/Prt_;
@ -275,7 +277,7 @@ void LaunderSharmaKE::correct()
if (!turbulence_)
{
// Re-calculate viscosity
mut_ == rho_*Cmu_*fMu()*sqr(k_)/(epsilon_ + epsilonSmall_);
mut_ == rho_*Cmu_*fMu()*sqr(k_)/epsilon_;
// Re-calculate thermal diffusivity
alphat_ = mut_/Prt_;
@ -320,7 +322,7 @@ void LaunderSharmaKE::correct()
epsEqn().relax();
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -338,11 +340,11 @@ void LaunderSharmaKE::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity
mut_ == Cmu_*fMu()*rho_*sqr(k_)/(epsilon_ + epsilonSmall_);
mut_ == Cmu_*fMu()*rho_*sqr(k_)/epsilon_;
// Re-calculate thermal diffusivity

16
src/turbulenceModels/compressible/RAS/RASModel/RASModel.C
View File

@ -81,11 +81,9 @@ RASModel::RASModel
printCoeffs_(lookupOrDefault<Switch>("printCoeffs", false)),
coeffDict_(subOrEmptyDict(type + "Coeffs")),
k0_("k0", dimVelocity*dimVelocity, SMALL),
epsilon0_("epsilon0", k0_.dimensions()/dimTime, SMALL),
epsilonSmall_("epsilonSmall", epsilon0_.dimensions(), SMALL),
omega0_("omega0", dimless/dimTime, SMALL),
omegaSmall_("omegaSmall", omega0_.dimensions(), SMALL),
kMin_("kMin", dimVelocity*dimVelocity, SMALL),
epsilonMin_("epsilonMin", kMin_.dimensions()/dimTime, SMALL),
omegaMin_("omegaMin", dimless/dimTime, SMALL),
y_(mesh_)
{
@ -219,11 +217,9 @@ bool RASModel::read()
coeffDict_ <<= *dictPtr;
}
k0_.readIfPresent(*this);
epsilon0_.readIfPresent(*this);
epsilonSmall_.readIfPresent(*this);
omega0_.readIfPresent(*this);
omegaSmall_.readIfPresent(*this);
kMin_.readIfPresent(*this);
epsilonMin_.readIfPresent(*this);
omegaMin_.readIfPresent(*this);
return true;
}

74
src/turbulenceModels/compressible/RAS/RASModel/RASModel.H
View File

@ -93,19 +93,13 @@ protected:
scalar yPlusLam_;
//- Lower limit of k
dimensionedScalar k0_;
dimensionedScalar kMin_;
//- Lower limit of epsilon
dimensionedScalar epsilon0_;
//- Small epsilon value used to avoid divide by zero
dimensionedScalar epsilonSmall_;
dimensionedScalar epsilonMin_;
//- Lower limit for omega
dimensionedScalar omega0_;
//- Small omega value used to avoid divide by zero
dimensionedScalar omegaSmall_;
dimensionedScalar omegaMin_;
//- Near wall distance boundary field
nearWallDist y_;
@ -185,66 +179,40 @@ public:
// Access
//- Return the lower allowable limit for k
const dimensionedScalar& k0() const
//- Return the lower allowable limit for k (default: SMALL)
const dimensionedScalar& kMin() const
{
return k0_;
return kMin_;
}
//- Return the lower allowable limit for epsilon
const dimensionedScalar& epsilon0() const
//- Return the lower allowable limit for epsilon (default: SMALL)
const dimensionedScalar& epsilonMin() const
{
return epsilon0_;
return epsilonMin_;
}
//- Return the value of epsilonSmall which is added to epsilon when
// calculating mut
const dimensionedScalar& epsilonSmall() const
//- Return the lower allowable limit for omega (default: SMALL)
const dimensionedScalar& omegaMin() const
{
return epsilonSmall_;
return omegaMin_;
}
//- Return the lower allowable limit for omega
const dimensionedScalar& omega0() const
//- Allow kMin to be changed
dimensionedScalar& kMin()
{
return omega0_;
return kMin_;
}
//- Return the value of omegaSmall which is added to omega when
// calculating mut
const dimensionedScalar& omegaSmall() const
//- Allow epsilonMin to be changed
dimensionedScalar& epsilonMin()
{
return omegaSmall_;
return epsilonMin_;
}
//- Allow k0 to be changed
dimensionedScalar& k0()
//- Allow omegaMin to be changed
dimensionedScalar& omegaMin()
{
return k0_;
}
//- Allow epsilon0 to be changed
dimensionedScalar& epsilon0()
{
return epsilon0_;
}
//- Allow epsilonSmall to be changed
dimensionedScalar& epsilonSmall()
{
return epsilonSmall_;
}
//- Allow omega0 to be changed
dimensionedScalar& omega0()
{
return omega0_;
}
//- Allow omegaSmall to be changed
dimensionedScalar& omegaSmall()
{
return omegaSmall_;
return omegaMin_;
}
//- Return the near wall distances

14
src/turbulenceModels/compressible/RAS/RNGkEpsilon/RNGkEpsilon.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -186,7 +186,9 @@ RNGkEpsilon::RNGkEpsilon
autoCreateAlphat("alphat", mesh_)
)
{
mut_ = Cmu_*rho_*sqr(k_)/(epsilon_ + epsilonSmall_);
bound(epsilon_, epsilonMin_);
mut_ = Cmu_*rho_*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
alphat_ = mut_/Prt_;
@ -276,7 +278,7 @@ void RNGkEpsilon::correct()
if (!turbulence_)
{
// Re-calculate viscosity
mut_ = rho_*Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
mut_ = rho_*Cmu_*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
// Re-calculate thermal diffusivity
@ -327,7 +329,7 @@ void RNGkEpsilon::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -339,12 +341,12 @@ void RNGkEpsilon::correct()
- fvm::laplacian(DkEff(), k_)
==
G - fvm::SuSp(2.0/3.0*rho_*divU, k_)
- fvm::Sp(rho_*(epsilon_)/k_, k_)
- fvm::Sp(rho_*epsilon_/k_, k_)
);
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

12
src/turbulenceModels/compressible/RAS/kEpsilon/kEpsilon.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -168,7 +168,9 @@ kEpsilon::kEpsilon
autoCreateAlphat("alphat", mesh_)
)
{
mut_ = Cmu_*rho_*sqr(k_)/(epsilon_ + epsilonSmall_);
bound(epsilon_, epsilonMin_);
mut_ = Cmu_*rho_*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
alphat_ = mut_/Prt_;
@ -257,7 +259,7 @@ void kEpsilon::correct()
if (!turbulence_)
{
// Re-calculate viscosity
mut_ = rho_*Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
mut_ = rho_*Cmu_*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
// Re-calculate thermal diffusivity
@ -300,7 +302,7 @@ void kEpsilon::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -318,7 +320,7 @@ void kEpsilon::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

8
src/turbulenceModels/compressible/RAS/kOmegaSST/kOmegaSST.C
View File

@ -258,12 +258,14 @@ kOmegaSST::kOmegaSST
autoCreateAlphat("alphat", mesh_)
)
{
bound(omega_, omegaMin_);
mut_ =
(
a1_*rho_*k_
/ max
(
a1_*(omega_ + omegaSmall_),
a1_*omega_,
F2()*sqrt(magSqr(symm(fvc::grad(U_))))
)
);
@ -422,7 +424,7 @@ void kOmegaSST::correct()
omegaEqn().boundaryManipulate(omega_.boundaryField());
solve(omegaEqn);
bound(omega_, omega0_);
bound(omega_, omegaMin_);
// Turbulent kinetic energy equation
tmp<fvScalarMatrix> kEqn
@ -438,7 +440,7 @@ void kOmegaSST::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

16
src/turbulenceModels/compressible/RAS/realizableKE/realizableKE.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -69,7 +69,7 @@ tmp<volScalarField> realizableKE::rCmu
volScalarField As = sqrt(6.0)*cos(phis);
volScalarField Us = sqrt(S2/2.0 + magSqr(skew(gradU)));
return 1.0/(A0_ + As*Us*k_/(epsilon_ + epsilonSmall_));
return 1.0/(A0_ + As*Us*k_/epsilon_);
}
@ -200,10 +200,10 @@ realizableKE::realizableKE
autoCreateAlphat("alphat", mesh_)
)
{
bound(k_, k0_);
bound(epsilon_, epsilon0_);
bound(k_, kMin_);
bound(epsilon_, epsilonMin_);
mut_ = rCmu(fvc::grad(U_))*rho_*sqr(k_)/(epsilon_ + epsilonSmall_);
mut_ = rCmu(fvc::grad(U_))*rho_*sqr(k_)/epsilon_;
mut_.correctBoundaryConditions();
alphat_ = mut_/Prt_;
@ -341,7 +341,7 @@ void realizableKE::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -353,12 +353,12 @@ void realizableKE::correct()
- fvm::laplacian(DkEff(), k_)
==
G - fvm::SuSp(2.0/3.0*rho_*divU, k_)
- fvm::Sp(rho_*(epsilon_)/k_, k_)
- fvm::Sp(rho_*epsilon_/k_, k_)
);
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity
mut_ = rCmu(gradU, S2, magS)*rho_*sqr(k_)/epsilon_;

10
src/turbulenceModels/incompressible/LES/DeardorffDiffStress/DeardorffDiffStress.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -122,15 +122,15 @@ void DeardorffDiffStress::correct(const tmp<volTensorField>& tgradU)
forAll(B_, celli)
{
B_[celli].component(symmTensor::XX) =
max(B_[celli].component(symmTensor::XX), k0().value());
max(B_[celli].component(symmTensor::XX), kMin_.value());
B_[celli].component(symmTensor::YY) =
max(B_[celli].component(symmTensor::YY), k0().value());
max(B_[celli].component(symmTensor::YY), kMin_.value());
B_[celli].component(symmTensor::ZZ) =
max(B_[celli].component(symmTensor::ZZ), k0().value());
max(B_[celli].component(symmTensor::ZZ), kMin_.value());
}
K = 0.5*tr(B_);
bound(K, k0());
bound(K, kMin_);
updateSubGridScaleFields(K);
}

8
src/turbulenceModels/incompressible/LES/LESModel/LESModel.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -78,10 +78,10 @@ LESModel::LESModel
printCoeffs_(lookupOrDefault<Switch>("printCoeffs", false)),
coeffDict_(subOrEmptyDict(type + "Coeffs")),
k0_("k0", dimVelocity*dimVelocity, SMALL),
kMin_("kMin", dimVelocity*dimVelocity, SMALL),
delta_(LESdelta::New("delta", U.mesh(), *this))
{
readIfPresent("k0", k0_);
readIfPresent("kMin", kMin_);
// Force the construction of the mesh deltaCoeffs which may be needed
// for the construction of the derived models and BCs
@ -167,7 +167,7 @@ bool LESModel::read()
delta_().read(*this);
readIfPresent("k0", k0_);
readIfPresent("kMin", kMin_);
return true;
}

16
src/turbulenceModels/incompressible/LES/LESModel/LESModel.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -83,7 +83,7 @@ protected:
Switch printCoeffs_;
dictionary coeffDict_;
dimensionedScalar k0_;
dimensionedScalar kMin_;
autoPtr<LESdelta> delta_;
@ -171,16 +171,16 @@ public:
return delta_();
}
//- Return the value of k0 which k is not allowed to be less than
const dimensionedScalar& k0() const
//- Return the lower allowable limit for k (default: SMALL)
const dimensionedScalar& kMin() const
{
return k0_;
return kMin_;
}
//- Allow k0 to be changed
dimensionedScalar& k0()
//- Allow kMin to be changed
dimensionedScalar& kMin()
{
return k0_;
return kMin_;
}

10
src/turbulenceModels/incompressible/LES/LRRDiffStress/LRRDiffStress.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -132,15 +132,15 @@ void LRRDiffStress::correct(const tmp<volTensorField>& tgradU)
forAll(B_, celli)
{
B_[celli].component(symmTensor::XX) =
max(B_[celli].component(symmTensor::XX), k0().value());
max(B_[celli].component(symmTensor::XX), kMin_.value());
B_[celli].component(symmTensor::YY) =
max(B_[celli].component(symmTensor::YY), k0().value());
max(B_[celli].component(symmTensor::YY), kMin_.value());
B_[celli].component(symmTensor::ZZ) =
max(B_[celli].component(symmTensor::ZZ), k0().value());
max(B_[celli].component(symmTensor::ZZ), kMin_.value());
}
K = 0.5*tr(B_);
bound(K, k0());
bound(K, kMin_);
updateSubGridScaleFields(K);
}

14
src/turbulenceModels/incompressible/LES/dynOneEqEddy/dynOneEqEddy.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -80,10 +80,10 @@ dimensionedScalar dynOneEqEddy::ce(const volSymmTensorField& D) const
pow(KK + filter_(k_), 1.5)/(2*delta()) - filter_(pow(k_, 1.5))/delta();
volScalarField ee =
2*delta()*ck(D)
2*delta()*ck(D)
*(
filter_(sqrt(k_)*magSqr(D))
- 2*sqrt(KK + filter_(k_))*magSqr(filter_(D))
filter_(sqrt(k_)*magSqr(D))
- 2*sqrt(KK + filter_(k_))*magSqr(filter_(D))
);
dimensionedScalar mmmm = average(magSqr(mm));
@ -135,8 +135,10 @@ dynOneEqEddy::dynOneEqEddy
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void dynOneEqEddy::correct(const tmp<volTensorField>& gradU)
void dynOneEqEddy::correct(const tmp<volTensorField>& tgradU)
{
const volTensorField& gradU = tgradU();
GenEddyVisc::correct(gradU);
volSymmTensorField D = symm(gradU);
@ -156,7 +158,7 @@ void dynOneEqEddy::correct(const tmp<volTensorField>& gradU)
kEqn.relax();
kEqn.solve();
bound(k_, k0());
bound(k_, kMin_);
updateSubGridScaleFields(D);
}

24
src/turbulenceModels/incompressible/LES/kOmegaSSTSAS/kOmegaSSTSAS.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2008-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2008-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -263,8 +263,7 @@ kOmegaSSTSAS::kOmegaSSTSAS
)
),
omega0_("omega0", dimless/dimTime, SMALL),
omegaSmall_("omegaSmall", dimless/dimTime, SMALL),
omegaMin_("omegaMin", dimless/dimTime, SMALL),
y_(mesh_),
Cmu_
(
@ -324,6 +323,8 @@ kOmegaSSTSAS::kOmegaSSTSAS
mesh_
)
{
bound(omega_, omegaMin_);
updateSubGridScaleFields(magSqr(symm(fvc::grad(U))));
printCoeffs();
@ -346,9 +347,8 @@ void kOmegaSSTSAS::correct(const tmp<volTensorField>& gradU)
volVectorField gradK = fvc::grad(k_);
volVectorField gradOmega = fvc::grad(omega_);
volScalarField L = sqrt(k_)/(pow025(Cmu_)*(omega_ + omegaSmall_));
volScalarField CDkOmega =
(2.0*alphaOmega2_)*(gradK & gradOmega)/(omega_ + omegaSmall_);
volScalarField L = sqrt(k_)/(pow025(Cmu_)*omega_);
volScalarField CDkOmega = (2.0*alphaOmega2_)*(gradK & gradOmega)/omega_;
volScalarField F1 = this->F1(CDkOmega);
volScalarField G = nuSgs_*2.0*S2;
@ -368,14 +368,12 @@ void kOmegaSSTSAS::correct(const tmp<volTensorField>& gradU)
kEqn.relax();
kEqn.solve();
}
bound(k_, k0());
bound(k_, kMin_);
volScalarField grad_omega_k = max
(
magSqr(gradOmega)/
sqr(omega_ + omegaSmall_),
magSqr(gradK)/
sqr(k_ + k0())
magSqr(gradOmega)/sqr(omega_),
magSqr(gradK)/sqr(k_)
);
// Turbulent frequency equation
@ -397,7 +395,7 @@ void kOmegaSSTSAS::correct(const tmp<volTensorField>& gradU)
+ FSAS_
*max
(
dimensionedScalar("zero",dimensionSet(0, 0 , -2, 0, 0),0. ),
dimensionedScalar("zero",dimensionSet(0, 0, -2, 0, 0), 0.0),
zetaTilda2_*kappa_*S2*(L/Lvk2(S2))
- 2.0/alphaPhi_*k_*grad_omega_k
)
@ -406,7 +404,7 @@ void kOmegaSSTSAS::correct(const tmp<volTensorField>& gradU)
omegaEqn.relax();
omegaEqn.solve();
}
bound(omega_, omega0_);
bound(omega_, omegaMin_);
updateSubGridScaleFields(S2);
}

5
src/turbulenceModels/incompressible/LES/kOmegaSSTSAS/kOmegaSSTSAS.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2008-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2008-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -103,8 +103,7 @@ protected:
dimensionedScalar zetaTilda2_;
dimensionedScalar FSAS_;
dimensionedScalar omega0_;
dimensionedScalar omegaSmall_;
dimensionedScalar omegaMin_;
wallDist y_;
dimensionedScalar Cmu_;

4
src/turbulenceModels/incompressible/LES/locDynOneEqEddy/locDynOneEqEddy.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -152,7 +152,7 @@ void locDynOneEqEddy::correct(const tmp<volTensorField>& gradU)
kEqn.relax();
kEqn.solve();
bound(k_, k0());
bound(k_, kMin_);
updateSubGridScaleFields(D, KK);
}

4
src/turbulenceModels/incompressible/LES/oneEqEddy/oneEqEddy.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -113,7 +113,7 @@ void oneEqEddy::correct(const tmp<volTensorField>& gradU)
kEqn.relax();
kEqn.solve();
bound(k_, k0());
bound(k_, kMin_);
updateSubGridScaleFields();
}

20
src/turbulenceModels/incompressible/RAS/LRR/LRR.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -186,9 +186,6 @@ LRR::LRR
autoCreateNut("nut", mesh_)
)
{
nut_ = Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
nut_.correctBoundaryConditions();
if (couplingFactor_.value() < 0.0 || couplingFactor_.value() > 1.0)
{
FatalErrorIn
@ -201,6 +198,11 @@ LRR::LRR
<< exit(FatalError);
}
bound(epsilon_, epsilonMin_);
nut_ = Cmu_*sqr(k_)/epsilon_;
nut_.correctBoundaryConditions();
printCoeffs();
}
@ -320,7 +322,7 @@ void LRR::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Reynolds stress equation
@ -368,15 +370,15 @@ void LRR::correct()
R_.dimensions(),
symmTensor
(
k0_.value(), -GREAT, -GREAT,
k0_.value(), -GREAT,
k0_.value()
kMin_.value(), -GREAT, -GREAT,
kMin_.value(), -GREAT,
kMin_.value()
)
)
);
k_ = 0.5*tr(R_);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

10
src/turbulenceModels/incompressible/RAS/LamBremhorstKE/LamBremhorstKE.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -140,7 +140,9 @@ LamBremhorstKE::LamBremhorstKE
autoCreateLowReNut("nut", mesh_)
)
{
nut_ = Cmu_*fMu_*sqr(k_)/(epsilon_ + epsilonSmall_);
bound(epsilon_, epsilonMin_);
nut_ = Cmu_*fMu_*sqr(k_)/epsilon_;
nut_.correctBoundaryConditions();
printCoeffs();
@ -260,7 +262,7 @@ void LamBremhorstKE::correct()
epsEqn().relax();
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -276,7 +278,7 @@ void LamBremhorstKE::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

16
src/turbulenceModels/incompressible/RAS/LaunderGibsonRSTM/LaunderGibsonRSTM.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -224,7 +224,9 @@ LaunderGibsonRSTM::LaunderGibsonRSTM
autoCreateNut("nut", mesh_)
)
{
nut_ = Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
bound(epsilon_, epsilonMin_);
nut_ = Cmu_*sqr(k_)/epsilon_;
nut_.correctBoundaryConditions();
if (couplingFactor_.value() < 0.0 || couplingFactor_.value() > 1.0)
@ -362,7 +364,7 @@ void LaunderGibsonRSTM::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Reynolds stress equation
@ -419,15 +421,15 @@ void LaunderGibsonRSTM::correct()
R_.dimensions(),
symmTensor
(
k0_.value(), -GREAT, -GREAT,
k0_.value(), -GREAT,
k0_.value()
kMin_.value(), -GREAT, -GREAT,
kMin_.value(), -GREAT,
kMin_.value()
)
)
);
k_ == 0.5*tr(R_);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate turbulent viscosity

10
src/turbulenceModels/incompressible/RAS/LaunderSharmaKE/LaunderSharmaKE.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -146,7 +146,9 @@ LaunderSharmaKE::LaunderSharmaKE
autoCreateLowReNut("nut", mesh_)
)
{
nut_ = Cmu_*fMu()*sqr(k_)/(epsilonTilda_ + epsilonSmall_);
bound(epsilonTilda_, epsilonMin_);
nut_ = Cmu_*fMu()*sqr(k_)/epsilonTilda_;
nut_.correctBoundaryConditions();
printCoeffs();
@ -256,7 +258,7 @@ void LaunderSharmaKE::correct()
epsEqn().relax();
solve(epsEqn);
bound(epsilonTilda_, epsilon0_);
bound(epsilonTilda_, epsilonMin_);
// Turbulent kinetic energy equation
@ -272,7 +274,7 @@ void LaunderSharmaKE::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

12
src/turbulenceModels/incompressible/RAS/LienCubicKE/LienCubicKE.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -170,6 +170,8 @@ LienCubicKE::LienCubicKE
autoCreateEpsilon("epsilon", mesh_)
),
//FIXME - epsilon is not bounded
gradU_(fvc::grad(U)),
eta_(k_/epsilon_*sqrt(2.0*magSqr(0.5*(gradU_ + gradU_.T())))),
ksi_(k_/epsilon_*sqrt(2.0*magSqr(0.5*(gradU_ - gradU_.T())))),
@ -226,7 +228,9 @@ LienCubicKE::LienCubicKE
)
)
{
nut_ = Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_) + C5viscosity_;
bound(epsilon_, epsilonMin_);
nut_ = Cmu_*sqr(k_)/epsilon_ + C5viscosity_;
nut_.correctBoundaryConditions();
printCoeffs();
@ -350,7 +354,7 @@ void LienCubicKE::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -367,7 +371,7 @@ void LienCubicKE::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

8
src/turbulenceModels/incompressible/RAS/LienCubicKELowRe/LienCubicKELowRe.C
View File

@ -280,10 +280,12 @@ LienCubicKELowRe::LienCubicKELowRe
)
)
{
bound(epsilon_, epsilonMin_);
nut_ = Cmu_
* (scalar(1) - exp(-Am_*yStar_))
/ (scalar(1) - exp(-Aepsilon_*yStar_) + SMALL)
* sqr(k_)/(epsilon_ + epsilonSmall_)
* sqr(k_)/epsilon_
// cubic term C5, implicit part
+ max
(
@ -432,7 +434,7 @@ void LienCubicKELowRe::correct()
# include "wallDissipationI.H"
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -449,7 +451,7 @@ void LienCubicKELowRe::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

10
src/turbulenceModels/incompressible/RAS/LienLeschzinerLowRe/LienLeschzinerLowRe.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -180,9 +180,11 @@ LienLeschzinerLowRe::LienLeschzinerLowRe
autoCreateLowReNut("nut", mesh_)
)
{
bound(epsilon_, epsilonMin_);
nut_ = Cmu_*(scalar(1) - exp(-Am_*yStar_))
/(scalar(1) - exp(-Aepsilon_*yStar_) + SMALL)*sqr(k_)
/(epsilon_ + epsilonSmall_);
/(epsilon_);
nut_.correctBoundaryConditions();
printCoeffs();
@ -322,7 +324,7 @@ void LienLeschzinerLowRe::correct()
#include "wallDissipationI.H"
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -339,7 +341,7 @@ void LienLeschzinerLowRe::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

32
src/turbulenceModels/incompressible/RAS/NonlinearKEShih/NonlinearKEShih.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -144,6 +144,7 @@ NonlinearKEShih::NonlinearKEShih
)
),
//FIXME: should be named 'kappa_' or 'kappa'?
kappa_
(
dimensioned<scalar>::lookupOrAddToDict
@ -174,6 +175,7 @@ NonlinearKEShih::NonlinearKEShih
IOobject::AUTO_WRITE
),
mesh_
//FIXME: what about autoCreateK("k", mesh_)
),
epsilon_
@ -187,15 +189,32 @@ NonlinearKEShih::NonlinearKEShih
IOobject::AUTO_WRITE
),
mesh_
//FIXME: what about autoCreateK("epsilon", mesh_)
),
//FIXME: epsilon is not bounded
gradU_(fvc::grad(U)),
eta_(k_/epsilon_*sqrt(2.0*magSqr(0.5*(gradU_ + gradU_.T())))),
ksi_(k_/epsilon_*sqrt(2.0*magSqr(0.5*(gradU_ - gradU_.T())))),
Cmu_(2.0/(3.0*(A1_ + eta_ + alphaKsi_*ksi_))),
fEta_(A2_ + pow(eta_, 3.0)),
nut_("nut", Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_)),
// FIXME: epsilon is not bounded
nut_("nut", Cmu_*sqr(k_)/(epsilon_ + epsilonMin_)),
// FIXME: why not use the following?
// nut_
// (
// IOobject
// (
// "nut",
// runTime_.timeName(),
// mesh_,
// IOobject::NO_READ,
// IOobject::AUTO_WRITE
// ),
// autoCreateNut("nut", mesh_)
// ),
nonlinearStress_
(
@ -215,6 +234,11 @@ NonlinearKEShih::NonlinearKEShih
)
)
{
bound(epsilon_, epsilonMin_);
//FIXME: could use this
// nut_ = Cmu_*sqr(k_)/epsilon_;
#include "wallNonlinearViscosityI.H"
printCoeffs();
@ -343,7 +367,7 @@ void NonlinearKEShih::correct()
#include "wallDissipationI.H"
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -360,7 +384,7 @@ void NonlinearKEShih::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

16
src/turbulenceModels/incompressible/RAS/RASModel/RASModel.C
View File

@ -80,11 +80,9 @@ RASModel::RASModel
printCoeffs_(lookupOrDefault<Switch>("printCoeffs", false)),
coeffDict_(subOrEmptyDict(type + "Coeffs")),
k0_("k0", dimVelocity*dimVelocity, SMALL),
epsilon0_("epsilon0", k0_.dimensions()/dimTime, SMALL),
epsilonSmall_("epsilonSmall", epsilon0_.dimensions(), SMALL),
omega0_("omega0", dimless/dimTime, SMALL),
omegaSmall_("omegaSmall", omega0_.dimensions(), SMALL),
kMin_("kMin", dimVelocity*dimVelocity, SMALL),
epsilonMin_("epsilonMin", kMin_.dimensions()/dimTime, SMALL),
omegaMin_("omegaMin", dimless/dimTime, SMALL),
y_(mesh_)
{
@ -213,11 +211,9 @@ bool RASModel::read()
coeffDict_ <<= *dictPtr;
}
k0_.readIfPresent(*this);
epsilon0_.readIfPresent(*this);
epsilonSmall_.readIfPresent(*this);
omega0_.readIfPresent(*this);
omegaSmall_.readIfPresent(*this);
kMin_.readIfPresent(*this);
epsilonMin_.readIfPresent(*this);
omegaMin_.readIfPresent(*this);
return true;
}

74
src/turbulenceModels/incompressible/RAS/RASModel/RASModel.H
View File

@ -87,19 +87,13 @@ protected:
dictionary coeffDict_;
//- Lower limit of k
dimensionedScalar k0_;
dimensionedScalar kMin_;
//- Lower limit of epsilon
dimensionedScalar epsilon0_;
//- Small epsilon value used to avoid divide by zero
dimensionedScalar epsilonSmall_;
dimensionedScalar epsilonMin_;
//- Lower limit for omega
dimensionedScalar omega0_;
//- Small omega value used to avoid divide by zero
dimensionedScalar omegaSmall_;
dimensionedScalar omegaMin_;
//- Near wall distance boundary field
nearWallDist y_;
@ -176,66 +170,40 @@ public:
// Access
//- Return lower allowable limit for k
const dimensionedScalar& k0() const
//- Return lower allowable limit for k (default: SMALL)
const dimensionedScalar& kMin() const
{
return k0_;
return kMin_;
}
//- Return the lower allowable limit for epsilon
const dimensionedScalar& epsilon0() const
//- Return the lower allowable limit for epsilon (default: SMALL)
const dimensionedScalar& epsilonMin() const
{
return epsilon0_;
return epsilonMin_;
}
//- Return the value of epsilonSmall which is added to epsilon when
// calculating nut
const dimensionedScalar& epsilonSmall() const
//- Return the lower allowable limit for omega (default: SMALL)
const dimensionedScalar& omegaMin() const
{
return epsilonSmall_;
return omegaMin_;
}
//- Return the lower allowable limit for omega0
const dimensionedScalar& omega0() const
//- Allow kMin to be changed
dimensionedScalar& kMin()
{
return omega0_;
return kMin_;
}
//- Return the value of omegaSmall which is added to omega when
// calculating nut
const dimensionedScalar& omegaSmall() const
//- Allow epsilonMin to be changed
dimensionedScalar& epsilonMin()
{
return omegaSmall_;
return epsilonMin_;
}
//- Allow k0 to be changed
dimensionedScalar& k0()
//- Allow omegaMin to be changed
dimensionedScalar& omegaMin()
{
return k0_;
}
//- Allow epsilon0 to be changed
dimensionedScalar& epsilon0()
{
return epsilon0_;
}
//- Allow epsilonSmall to be changed
dimensionedScalar& epsilonSmall()
{
return epsilonSmall_;
}
//- Allow omega0 to be changed
dimensionedScalar& omega0()
{
return omega0_;
}
//- Allow omegaSmall to be changed
dimensionedScalar& omegaSmall()
{
return omegaSmall_;
return omegaMin_;
}
//- Return the near wall distances

10
src/turbulenceModels/incompressible/RAS/RNGkEpsilon/RNGkEpsilon.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -155,7 +155,9 @@ RNGkEpsilon::RNGkEpsilon
autoCreateNut("nut", mesh_)
)
{
nut_ = Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
bound(epsilon_, epsilonMin_);
nut_ = Cmu_*sqr(k_)/epsilon_;
nut_.correctBoundaryConditions();
printCoeffs();
@ -274,7 +276,7 @@ void RNGkEpsilon::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -290,7 +292,7 @@ void RNGkEpsilon::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

10
src/turbulenceModels/incompressible/RAS/kEpsilon/kEpsilon.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -128,7 +128,9 @@ kEpsilon::kEpsilon
autoCreateNut("nut", mesh_)
)
{
nut_ = Cmu_*sqr(k_)/(epsilon_ + epsilonSmall_);
bound(epsilon_, epsilonMin_);
nut_ = Cmu_*sqr(k_)/epsilon_;
nut_.correctBoundaryConditions();
printCoeffs();
@ -237,7 +239,7 @@ void kEpsilon::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -254,7 +256,7 @@ void kEpsilon::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

12
src/turbulenceModels/incompressible/RAS/kOmega/kOmega.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -137,7 +137,9 @@ kOmega::kOmega
autoCreateNut("nut", mesh_)
)
{
nut_ = k_/(omega_ + omegaSmall_);
bound(omega_, omegaMin_);
nut_ = k_/omega_;
nut_.correctBoundaryConditions();
printCoeffs();
@ -246,7 +248,7 @@ void kOmega::correct()
omegaEqn().boundaryManipulate(omega_.boundaryField());
solve(omegaEqn);
bound(omega_, omega0_);
bound(omega_, omegaMin_);
// Turbulent kinetic energy equation
@ -263,11 +265,11 @@ void kOmega::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity
nut_ = k_/(omega_ + omegaSmall_);
nut_ = k_/omega_;
nut_.correctBoundaryConditions();
}

8
src/turbulenceModels/incompressible/RAS/kOmegaSST/kOmegaSST.C
View File

@ -236,12 +236,14 @@ kOmegaSST::kOmegaSST
autoCreateNut("nut", mesh_)
)
{
bound(omega_, omegaMin_);
nut_ =
(
a1_*k_
/ max
(
a1_*(omega_ + omegaSmall_),
a1_*omega_,
F2()*mag(symm(fvc::grad(U_)))
)
);
@ -376,7 +378,7 @@ void kOmegaSST::correct()
omegaEqn().boundaryManipulate(omega_.boundaryField());
solve(omegaEqn);
bound(omega_, omega0_);
bound(omega_, omegaMin_);
// Turbulent kinetic energy equation
tmp<fvScalarMatrix> kEqn
@ -392,7 +394,7 @@ void kOmegaSST::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity

18
src/turbulenceModels/incompressible/RAS/qZeta/qZeta.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -126,6 +126,9 @@ qZeta::qZeta
)
),
qMin_("qMin", dimVelocity, SMALL),
zetaMin_("zetaMin", dimVelocity/dimTime, SMALL),
k_
(
IOobject
@ -193,7 +196,11 @@ qZeta::qZeta
autoCreateNut("nut", mesh_)
)
{
nut_ = Cmu_*fMu()*sqr(k_)/(epsilon_ + epsilonSmall_);
bound(epsilon_, epsilonMin_);
bound(q_, qMin_);
bound(zeta_, zetaMin_);
nut_ = Cmu_*fMu()*sqr(k_)/epsilon_;
nut_.correctBoundaryConditions();
printCoeffs();
@ -263,6 +270,9 @@ bool qZeta::read()
sigmaZeta_.readIfPresent(coeffDict());
anisotropic_.readIfPresent("anisotropic", coeffDict());
qMin_.readIfPresent(*this);
zetaMin_.readIfPresent(*this);
return true;
}
else
@ -302,7 +312,7 @@ void qZeta::correct()
zetaEqn().relax();
solve(zetaEqn);
bound(zeta_, epsilon0_/(2*sqrt(k0_)));
bound(zeta_, zetaMin_);
// q equation
@ -318,7 +328,7 @@ void qZeta::correct()
qEqn().relax();
solve(qEqn);
bound(q_, sqrt(k0_));
bound(q_, qMin_);
// Re-calculate k and epsilon

34
src/turbulenceModels/incompressible/RAS/qZeta/qZeta.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -66,6 +66,11 @@ class qZeta
dimensionedScalar sigmaZeta_;
Switch anisotropic_;
//- Lower limit of q
dimensionedScalar qMin_;
//- Lower limit of zeta
dimensionedScalar zetaMin_;
// Fields
@ -107,6 +112,33 @@ public:
// Member Functions
// Access
//- Return lower allowable limit for q (default: SMALL)
const dimensionedScalar& qMin() const
{
return qMin_;
}
//- Return the lower allowable limit for zeta (default: SMALL)
const dimensionedScalar& zetaMin() const
{
return zetaMin_;
}
//- Allow qMin to be changed
dimensionedScalar& qMin()
{
return qMin_;
}
//- Allow zetaMin to be changed
dimensionedScalar& zetaMin()
{
return zetaMin_;
}
//- Return the turbulence viscosity
virtual tmp<volScalarField> nut() const
{

14
src/turbulenceModels/incompressible/RAS/realizableKE/realizableKE.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -69,7 +69,7 @@ tmp<volScalarField> realizableKE::rCmu
volScalarField As = sqrt(6.0)*cos(phis);
volScalarField Us = sqrt(S2/2.0 + magSqr(skew(gradU)));
return 1.0/(A0_ + As*Us*k_/(epsilon_ + epsilonSmall_));
return 1.0/(A0_ + As*Us*k_/epsilon_);
}
@ -178,10 +178,10 @@ realizableKE::realizableKE
autoCreateNut("nut", mesh_)
)
{
bound(k_, k0_);
bound(epsilon_, epsilon0_);
bound(k_, kMin_);
bound(epsilon_, epsilonMin_);
nut_ = rCmu(fvc::grad(U_))*sqr(k_)/(epsilon_ + epsilonSmall_);
nut_ = rCmu(fvc::grad(U_))*sqr(k_)/epsilon_;
nut_.correctBoundaryConditions();
printCoeffs();
@ -303,7 +303,7 @@ void realizableKE::correct()
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, epsilon0_);
bound(epsilon_, epsilonMin_);
// Turbulent kinetic energy equation
@ -319,7 +319,7 @@ void realizableKE::correct()
kEqn().relax();
solve(kEqn);
bound(k_, k0_);
bound(k_, kMin_);
// Re-calculate viscosity