All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
nutWallFunctionFvPatchScalarField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
27 #include "momentumTransportModel.H"
28 #include "fvPatchFieldMapper.H"
29 #include "wallFvPatch.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  defineTypeNameAndDebug(nutWallFunctionFvPatchScalarField, 0);
37 }
38 
39 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
40 
42 {
43  if (!isA<wallFvPatch>(patch()))
44  {
46  << "Invalid wall function specification" << nl
47  << " Patch type for patch " << patch().name()
48  << " must be wall" << nl
49  << " Current patch type is " << patch().type() << nl << endl
50  << abort(FatalError);
51  }
52 }
53 
54 
56 (
57  Ostream& os
58 ) const
59 {
60  writeEntry(os, "Cmu", Cmu_);
61  writeEntry(os, "kappa", kappa_);
62  writeEntry(os, "E", E_);
63 }
64 
65 
66 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
67 
69 (
70  const fvPatch& p,
72 )
73 :
74  fixedValueFvPatchScalarField(p, iF),
75  Cmu_(0.09),
76  kappa_(0.41),
77  E_(9.8),
78  yPlusLam_(yPlusLam(kappa_, E_))
79 {
80  checkType();
81 }
82 
83 
85 (
87  const fvPatch& p,
89  const fvPatchFieldMapper& mapper
90 )
91 :
92  fixedValueFvPatchScalarField(ptf, p, iF, mapper),
93  Cmu_(ptf.Cmu_),
94  kappa_(ptf.kappa_),
95  E_(ptf.E_),
96  yPlusLam_(ptf.yPlusLam_)
97 {
98  checkType();
99 }
100 
101 
103 (
104  const fvPatch& p,
106  const dictionary& dict
107 )
108 :
109  fixedValueFvPatchScalarField(p, iF, dict),
110  Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
111  kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
112  E_(dict.lookupOrDefault<scalar>("E", 9.8)),
113  yPlusLam_(yPlusLam(kappa_, E_))
114 {
115  checkType();
116 }
117 
118 
120 (
122 )
123 :
124  fixedValueFvPatchScalarField(wfpsf),
125  Cmu_(wfpsf.Cmu_),
126  kappa_(wfpsf.kappa_),
127  E_(wfpsf.E_),
128  yPlusLam_(wfpsf.yPlusLam_)
129 {
130  checkType();
131 }
132 
133 
135 (
138 )
139 :
140  fixedValueFvPatchScalarField(wfpsf, iF),
141  Cmu_(wfpsf.Cmu_),
142  kappa_(wfpsf.kappa_),
143  E_(wfpsf.E_),
144  yPlusLam_(wfpsf.yPlusLam_)
145 {
146  checkType();
147 }
148 
149 
150 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
151 
154 (
155  const momentumTransportModel& turbModel,
156  const label patchi
157 )
158 {
159  return
160  refCast<const nutWallFunctionFvPatchScalarField>
161  (
162  turbModel.nut()().boundaryField()[patchi]
163  );
164 }
165 
166 
168 (
169  const scalar kappa,
170  const scalar E
171 )
172 {
173  scalar ypl = 11.0;
174 
175  for (int i=0; i<10; i++)
176  {
177  ypl = log(max(E*ypl, 1))/kappa;
178  }
179 
180  return ypl;
181 }
182 
183 
185 {
186  return yPlusLam_;
187 }
188 
189 
191 {
192  if (updated())
193  {
194  return;
195  }
196 
197  operator==(nut());
198 
199  fixedValueFvPatchScalarField::updateCoeffs();
200 }
201 
202 
204 {
206  writeLocalEntries(os);
207  writeEntry(os, "value", *this);
208 }
209 
210 
211 // ************************************************************************* //
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const dimensionedScalar & kappa
Coulomb constant: default SI units: [N.m2/C2].
dimensionedScalar log(const dimensionedScalar &ds)
scalar yPlusLam() const
Return the Y+ at the edge of the laminar sublayer.
virtual void writeLocalEntries(Ostream &) const
Write local wall function variables.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:280
Macros for easy insertion into run-time selection tables.
static const nutWallFunctionFvPatchScalarField & nutw(const momentumTransportModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
This boundary condition provides a turbulent kinematic viscosity condition when using wall functions...
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
Foam::fvPatchFieldMapper.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
static const char nl
Definition: Ostream.H:260
defineTypeNameAndDebug(combustionModel, 0)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
Abstract base class for turbulence models (RAS, LES and laminar).
virtual void checkType()
Check the type of the patch.
label patchi
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
nutWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
scalar yPlusLam_
Y+ at the edge of the laminar sublayer.
scalar nut
Namespace for OpenFOAM.