kLowReWallFunctionFvPatchScalarField.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) 2012-2024 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 
28 #include "momentumTransportModel.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 
36 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
37 
39 (
40  const fvPatch& p,
42  const dictionary& dict
43 )
44 :
45  fixedValueFvPatchField<scalar>(p, iF, dict),
46  Ceps2_(dict.lookupOrDefault<scalar>("Ceps2", dimless, 1.9))
47 {}
48 
49 
51 (
53  const fvPatch& p,
55  const fieldMapper& mapper
56 )
57 :
58  fixedValueFvPatchField<scalar>(ptf, p, iF, mapper),
59  Ceps2_(ptf.Ceps2_)
60 {}
61 
62 
64 (
67 )
68 :
69  fixedValueFvPatchField<scalar>(kwfpsf, iF),
70  Ceps2_(kwfpsf.Ceps2_)
71 {}
72 
73 
74 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
75 
77 {
78  if (updated())
79  {
80  return;
81  }
82 
83  const label patchi = patch().index();
84 
85  const momentumTransportModel& turbModel =
87 
90 
91  const scalarField& y = turbModel.yb()[patchi];
92 
93  const tmp<volScalarField> tk = turbModel.k();
94  const volScalarField& k = tk();
95 
96  const tmp<scalarField> tnuw = turbModel.nu(patchi);
97  const scalarField& nuw = tnuw();
98 
99  const scalar Cmu25 = pow025(nutw.Cmu());
100 
101  scalarField& kw = *this;
102 
103  // Set k wall values
104  forAll(kw, facei)
105  {
106  label celli = patch().faceCells()[facei];
107 
108  scalar uTau = Cmu25*sqrt(k[celli]);
109 
110  scalar yPlus = uTau*y[facei]/nuw[facei];
111 
112  if (yPlus > nutw.yPlusLam())
113  {
114  scalar Ck = -0.416;
115  scalar Bk = 8.366;
116  kw[facei] = Ck/nutw.kappa()*log(yPlus) + Bk;
117  }
118  else
119  {
120  scalar C = 11.0;
121  scalar Cf = (1.0/sqr(yPlus + C) + 2.0*yPlus/pow3(C) - 1.0/sqr(C));
122  kw[facei] = 2400.0/sqr(Ceps2_)*Cf;
123  }
124 
125  kw[facei] *= sqr(uTau);
126  }
127 
128  // Limit kw to avoid failure of the turbulence model due to division by kw
129  kw = max(kw, small);
130 
132 
133  // TODO: perform averaging for cells sharing more than one boundary face
134 }
135 
136 
138 {
139  writeEntry(os, "Ceps2", Ceps2_);
141 }
142 
143 
144 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
145 
147 (
150 );
151 
152 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
153 
154 } // End namespace Foam
155 
156 // ************************************************************************* //
scalar y
label k
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Graphite solid properties.
Definition: C.H:51
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Abstract base class for field mapping.
Definition: fieldMapper.H:48
This boundary condition supplies a fixed value constraint, and is the base class for a number of othe...
virtual void write(Ostream &) const
Write.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:88
const DimensionedField< Type, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:362
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:143
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:202
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:356
bool updated() const
Return true if the boundary condition has already been updated.
Definition: fvPatchField.H:374
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
label index() const
Return the index of this patch in the fvBoundaryMesh.
Definition: fvPatch.H:156
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:99
This boundary condition provides a turbulence kinetic energy wall function condition for low- and hig...
kLowReWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Abstract base class for turbulence models (RAS, LES and laminar).
const volScalarField::Boundary & yb() const
Return the near wall distance.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
virtual tmp< volScalarField > nu() const
Return the laminar viscosity.
This boundary condition provides a turbulent kinematic viscosity condition when using wall functions,...
static const nutWallFunctionFvPatchScalarField & nutw(const momentumTransportModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
const Type & lookupType(const word &group=word::null) const
Lookup and return the object of the given Type.
A class for managing temporary objects.
Definition: tmp.H:55
const scalar yPlus
const scalar uTau
label patchi
Namespace for OpenFOAM.
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
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimless
dimensionedScalar log(const dimensionedScalar &ds)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensionedScalar sqrt(const dimensionedScalar &ds)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
dimensionedScalar pow025(const dimensionedScalar &ds)
dictionary dict
volScalarField & p