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-2021 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 (
86  const fvPatch& p,
88  const dictionary& dict
89 )
90 :
91  fixedValueFvPatchScalarField(p, iF, dict),
92  Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
93  kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
94  E_(dict.lookupOrDefault<scalar>("E", 9.8)),
95  yPlusLam_(yPlusLam(kappa_, E_))
96 {
97  checkType();
98 }
99 
100 
102 (
104  const fvPatch& p,
106  const fvPatchFieldMapper& mapper
107 )
108 :
109  fixedValueFvPatchScalarField(ptf, p, iF, mapper),
110  Cmu_(ptf.Cmu_),
111  kappa_(ptf.kappa_),
112  E_(ptf.E_),
113  yPlusLam_(ptf.yPlusLam_)
114 {
115  checkType();
116 }
117 
118 
120 (
123 )
124 :
125  fixedValueFvPatchScalarField(wfpsf, iF),
126  Cmu_(wfpsf.Cmu_),
127  kappa_(wfpsf.kappa_),
128  E_(wfpsf.E_),
129  yPlusLam_(wfpsf.yPlusLam_)
130 {
131  checkType();
132 }
133 
134 
135 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
136 
139 (
140  const momentumTransportModel& turbModel,
141  const label patchi
142 )
143 {
144  return
145  refCast<const nutWallFunctionFvPatchScalarField>
146  (
147  turbModel.nut()().boundaryField()[patchi]
148  );
149 }
150 
151 
153 (
154  const scalar kappa,
155  const scalar E
156 )
157 {
158  scalar ypl = 11.0;
159 
160  for (int i=0; i<10; i++)
161  {
162  ypl = log(max(E*ypl, 1))/kappa;
163  }
164 
165  return ypl;
166 }
167 
168 
170 {
171  return yPlusLam_;
172 }
173 
174 
176 {
177  if (updated())
178  {
179  return;
180  }
181 
182  operator==(nut());
183 
184  fixedValueFvPatchScalarField::updateCoeffs();
185 }
186 
187 
189 {
191  writeLocalEntries(os);
192  writeEntry(os, "value", *this);
193 }
194 
195 
196 // ************************************************************************* //
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
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:156
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
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:62
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:260
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.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
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.