Go to the documentation of this file.
43 template<
class BasicTurbulenceModel>
44 tmp<volScalarField> kL<BasicTurbulenceModel>::Cmu()
const
47 return (0.556 + 0.108*Rt_)/(1.0 + 0.308*Rt_ + 0.00837*
sqr(Rt_));
51 template<
class BasicTurbulenceModel>
52 tmp<volScalarField> kL<BasicTurbulenceModel>::CmuPrime()
const
55 return 0.556/(1.0 + 0.277*Rt_);
59 template<
class BasicTurbulenceModel>
60 tmp<volScalarField> kL<BasicTurbulenceModel>::nutPrime()
const
63 return CmuPrime()*
sqrt(k_)*L_;
67 template<
class BasicTurbulenceModel>
68 tmp<volScalarField> kL<BasicTurbulenceModel>::epsilonCanopy()
const
71 this->mesh_.template findObject<volScalarField>(
"plantCd");
73 this->mesh_.template findObject<volScalarField>(
"leafAreaDensity");
78 const auto& Cd = *CdPtr;
79 const auto& LAD = *LADPtr;
82 return Cd*LAD*
mag(
U)*k_;
90 this->runTime_.timeName(),
102 template<
class BasicTurbulenceModel>
103 tmp<volScalarField> kL<BasicTurbulenceModel>::epsilon()
const
106 tmp<volScalarField> tepsilonCanopy = epsilonCanopy();
109 tmp<volScalarField> tepsilonPlain =
pow3(Cmu0_)*
pow(k_, 1.5)/L_;
112 tmp<volScalarField> tepsilon =
max(tepsilonPlain, tepsilonCanopy);
120 template<
class BasicTurbulenceModel>
121 tmp<volScalarField> kL<BasicTurbulenceModel>::canopyHeight()
const
123 const auto* canopyHeightPtr =
124 this->mesh_.template findObject<volScalarField>(
"canopyHeight");
128 const auto& canopyHeight = *canopyHeightPtr;
137 this->runTime_.timeName(),
149 template<
class BasicTurbulenceModel>
150 tmp<volScalarField> kL<BasicTurbulenceModel>::L()
const
156 tmp<volScalarField> tLcanopy = kappa_*canopyHeight();
160 return max(Lcanopy, Lplain);
164 template<
class BasicTurbulenceModel>
165 void kL<BasicTurbulenceModel>::stratification(
const volScalarField& fVB)
167 tmp<volScalarField> tLg =
L();
170 tmp<volScalarField> tcanopyHeight = canopyHeight();
173 tmp<volScalarField> tLcanopy = kappa_*canopyHeight;
176 const scalar Cmu0 = Cmu0_.value();
177 const scalar CbStable = CbStable_.value();
178 const scalar CbUnstable = CbUnstable_.value();
182 if (y_[i] > canopyHeight[i])
187 const scalar Lb = CbStable*
sqrt(k_[i])/
sqrt(fVB[i]);
202 Rt_[i] -
sqr(Rt_[i] + 1.0)/(Rt_[i] - 1.0)
208 *
sqrt(1.0 -
pow(Cmu0, 6.0)*
pow(CbUnstable, -2.0)*Rt_[i]);
224 template<
class BasicTurbulenceModel>
227 this->nut_ = Cmu()*
sqrt(k_)*L_;
228 this->nut_.correctBoundaryConditions();
231 BasicTurbulenceModel::correctNut();
235 template<
class BasicTurbulenceModel>
248 template<
class BasicTurbulenceModel>
257 const word& propertiesName,
344 this->runTime_.timeName(),
356 this->runTime_.timeName(),
369 this->runTime_.timeName(),
380 bound(k_, this->kMin_);
382 if (
type == typeName)
384 this->printCoeffs(
type);
391 template<
class BasicTurbulenceModel>
396 kappa_.readIfPresent(this->coeffDict());
397 sigmak_.readIfPresent(this->coeffDict());
398 beta_.readIfPresent(this->coeffDict());
399 Cmu0_.readIfPresent(this->coeffDict());
400 Lmax_.readIfPresent(this->coeffDict());
401 CbStable_.readIfPresent(this->coeffDict());
402 CbUnstable_.readIfPresent(this->coeffDict());
411 template<
class BasicTurbulenceModel>
414 if (!this->turbulence_)
477 bound(k_, this->kMin_);
480 Rt_ = fBV*
sqr(k_/tepsilon);
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const vector L(dict.get< vector >("L"))
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
A class for handling words, derived from Foam::string.
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
void clear() const noexcept
BasicTurbulenceModel::transportModel transportModel
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
const dimensionedScalar G
Newtonian constant of gravitation.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
BasicTurbulenceModel::rhoField rhoField
Templated abstract base class for RAS turbulence models.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to database if not present.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
bool read(const char *buf, int32_t &val)
Same as readInt32.
static const wallDist & New(const fvMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
virtual tmp< fvScalarMatrix > kSource() const
Add explicit source for turbulent kinetic energy.
virtual void correctNut()
Correct the turbulence viscosity.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAll(list, i)
Loop across all elements in list.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
virtual bool read()
Re-read model coefficients if they have changed.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
dimensionedScalar pow3(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
A one-equation (turbulent kinetic energy k) turbulence closure model for incompressible and compressi...
Bound the given scalar field if it has gone unbounded.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
BasicTurbulenceModel::alphaField alphaField
Generic dimensioned Type class.
GeometricField< vector, fvPatchField, volMesh > volVectorField
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
const dimensionedScalar e
Elementary charge.
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Eddy viscosity turbulence model base class.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
const dimensionSet dimVolume(pow3(dimLength))
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
static const gravity & New(const Time &runTime)
Construct on Time.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)