atmEpsilonWallFunctionFvPatchScalarField.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2020 ENERCON GmbH
9 Copyright (C) 2020-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
31#include "turbulenceModel.H"
32#include "fvMatrix.H"
34
35// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
36
38(
39 const turbulenceModel& turbModel,
40 const List<scalar>& cornerWeights,
41 const fvPatch& patch,
42 scalarField& G0,
43 scalarField& epsilon0
44)
45{
46 const label patchi = patch.index();
47
48 const tmp<scalarField> tnutw = turbModel.nut(patchi);
49 const scalarField& nutw = tnutw();
50
51 const scalarField& y = turbModel.y()[patchi];
52
53 const tmp<scalarField> tnuw = turbModel.nu(patchi);
54 const scalarField& nuw = tnuw();
55
56 const tmp<volScalarField> tk = turbModel.k();
57 const volScalarField& k = tk();
58
59 const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchi];
60
61 const scalarField magGradUw(mag(Uw.snGrad()));
62
63 const scalar Cmu25 = pow025(wallCoeffs_.Cmu());
64 const scalar Cmu75 = pow(wallCoeffs_.Cmu(), 0.75);
65 const scalar kappa = wallCoeffs_.kappa();
66 const scalar yPlusLam = wallCoeffs_.yPlusLam();
67
68 const scalar t = db().time().timeOutputValue();
69 const scalarField z0(z0_->value(t));
70
71 #ifdef FULLDEBUG
72 for (const auto& z : z0)
73 {
74 if (z < VSMALL)
75 {
77 << "z0 field can only contain positive values. "
78 << "Please check input field z0."
79 << exit(FatalError);
80 }
81 }
82 #endif
83
85
86 // Set epsilon and G
87 forAll(nutw, facei)
88 {
89 const label celli = faceCells[facei];
90
91 const scalar yPlus = Cmu25*y[facei]*sqrt(k[celli])/nuw[facei];
92
93 const scalar w = cornerWeights[facei];
94
95 // (PGVB:Eq. 7, RH:Eq. 8)
96 scalar epsilonc =
97 w*Cmu75*pow(k[celli], 1.5)/(kappa*(y[facei] + z0[facei]));
98
99 scalar Gc =
100 w
101 *(nutw[facei] + nuw[facei])
102 *magGradUw[facei]
103 *Cmu25*sqrt(k[celli])
104 /(kappa*(y[facei] + z0[facei]));
105
106 if (lowReCorrection_ && yPlus < yPlusLam)
107 {
108 epsilonc = w*2.0*k[celli]*nuw[facei]/sqr(y[facei] + z0[facei]);
109 Gc = 0;
110 }
111
112 epsilon0[celli] += epsilonc;
113
114 G0[celli] += Gc;
115 }
116}
117
118
120(
121 Ostream& os
122) const
123{
124 os.writeEntryIfDifferent<bool>("lowReCorrection", false, lowReCorrection_);
125
126 if (z0_)
127 {
128 z0_->writeData(os);
129 }
130
131 wallCoeffs_.writeEntries(os);
132}
133
134
135// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
136
137Foam::atmEpsilonWallFunctionFvPatchScalarField::
138atmEpsilonWallFunctionFvPatchScalarField
139(
140 const fvPatch& p,
142)
143:
145 z0_(nullptr)
146{}
147
148
149Foam::atmEpsilonWallFunctionFvPatchScalarField::
150atmEpsilonWallFunctionFvPatchScalarField
151(
153 const fvPatch& p,
155 const fvPatchFieldMapper& mapper
156)
157:
159 z0_(ptf.z0_.clone(p.patch()))
160{}
161
162
163Foam::atmEpsilonWallFunctionFvPatchScalarField::
164atmEpsilonWallFunctionFvPatchScalarField
165(
166 const fvPatch& p,
168 const dictionary& dict
169)
170:
172 z0_(PatchFunction1<scalar>::New(p.patch(), "z0", dict))
173{}
174
175
176Foam::atmEpsilonWallFunctionFvPatchScalarField::
177atmEpsilonWallFunctionFvPatchScalarField
178(
180)
181:
183 z0_(ewfpsf.z0_.clone(this->patch().patch()))
184{}
185
186
187Foam::atmEpsilonWallFunctionFvPatchScalarField::
188atmEpsilonWallFunctionFvPatchScalarField
189(
192)
193:
195 z0_(ewfpsf.z0_.clone(this->patch().patch()))
196{}
197
198
199// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
200
202(
203 const fvPatchFieldMapper& m
204)
205{
207
208 if (z0_)
209 {
210 z0_->autoMap(m);
211 }
212}
213
214
216(
217 const fvPatchScalarField& ptf,
218 const labelList& addr
219)
220{
222
223 const auto& atmpsf =
224 refCast<const atmEpsilonWallFunctionFvPatchScalarField>(ptf);
225 if (z0_)
226 {
227 z0_->rmap(atmpsf.z0_(), addr);
228 }
229}
230
231
233(
234 Ostream& os
235) const
236{
238 writeLocalEntries(os);
239 writeEntry("value", os);
240}
241
242
243// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
244
245namespace Foam
246{
248 (
251 );
252}
253
254
255// ************************************************************************* //
scalar y
label k
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const Boundary & boundaryField() const
Return const-reference to the boundary field.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:77
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Definition: Ostream.H:251
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
scalar timeOutputValue() const
Return current time value.
Definition: TimeStateI.H:31
void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
void rmap(const atmBoundaryLayer &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
This boundary condition provides a wall constraint on the turbulent kinetic energy dissipation rate (...
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
autoPtr< PatchFunction1< scalar > > z0_
Surface roughness length field [m].
void writeLocalEntries(Ostream &) const
Write local wall function variables.
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void calculate(const turbulenceModel &turbulence, const List< scalar > &cornerWeights, const fvPatch &patch, scalarField &G, scalarField &epsilon)
Calculate the epsilon and G.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
This boundary condition provides wall functions for the turbulent kinetic energy dissipation rate (i....
const bool lowReCorrection_
Apply low-Re correction term (default = no)
wallFunctionCoefficients wallCoeffs_
Wall-function coefficients.
Smooth ATC in cells next to a set of patches supplied by type.
Definition: faceCells.H:59
virtual bool write()
Write the output fields.
A FieldMapper for finite-volume patch fields.
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:229
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:210
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:362
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
label index() const
Return the index of this patch in the fvBoundaryMesh.
Definition: fvPatch.H:203
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:113
const Time & time() const noexcept
Return time registry.
A class for managing temporary objects.
Definition: tmp.H:65
Abstract base class for turbulence models (RAS, LES and laminar).
const volVectorField & U() const
Access function to velocity field.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
const nearWallDist & y() const
Return the near wall distances.
scalar kappa() const noexcept
Return the object: kappa.
scalar yPlusLam() const noexcept
Return the object: yPlusLam.
scalar Cmu() const noexcept
Return the object: Cmu.
volScalarField & p
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
scalar yPlus
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
Definition: fvPatchField.H:676
Namespace for OpenFOAM.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
error FatalError
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dimensionedScalar pow025(const dimensionedScalar &ds)
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333