v2WallFunctionFvPatchScalarField.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-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 
28 #include "momentumTransportModel.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace RASModels
36 {
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
41 (
42  const fvPatch& p,
44 )
45 :
47 {}
48 
49 
51 (
53  const fvPatch& p,
55  const fvPatchFieldMapper& mapper
56 )
57 :
58  fixedValueFvPatchField<scalar>(ptf, p, iF, mapper)
59 {}
60 
61 
63 (
64  const fvPatch& p,
66  const dictionary& dict
67 )
68 :
70 {}
71 
72 
74 (
76 )
77 :
79 {}
80 
81 
83 (
84  const v2WallFunctionFvPatchScalarField& v2wfpsf,
86 )
87 :
89 {}
90 
91 
92 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
93 
95 {
96  if (updated())
97  {
98  return;
99  }
100 
101  const label patchi = patch().index();
102 
103  const momentumTransportModel& turbModel =
105  (
107  (
108  momentumTransportModel::typeName,
109  internalField().group()
110  )
111  );
112 
114  nutWallFunctionFvPatchScalarField::nutw(turbModel, patchi);
115 
116  const scalarField& y = turbModel.y()[patchi];
117 
118  const tmp<volScalarField> tk = turbModel.k();
119  const volScalarField& k = tk();
120 
121  const tmp<scalarField> tnuw = turbModel.nu(patchi);
122  const scalarField& nuw = tnuw();
123 
124  const scalar Cmu25 = pow025(nutw.Cmu());
125 
126  scalarField& v2 = *this;
127 
128  // Set v2 wall values
129  forAll(v2, facei)
130  {
131  label celli = patch().faceCells()[facei];
132 
133  scalar uTau = Cmu25*sqrt(k[celli]);
134 
135  scalar yPlus = uTau*y[facei]/nuw[facei];
136 
137  if (yPlus > nutw.yPlusLam())
138  {
139  scalar Cv2 = 0.193;
140  scalar Bv2 = -0.94;
141  v2[facei] = Cv2/nutw.kappa()*log(yPlus) + Bv2;
142  }
143  else
144  {
145  scalar Cv2 = 0.193;
146  v2[facei] = Cv2*pow4(yPlus);
147  }
148 
149  v2[facei] *= sqr(uTau);
150  }
151 
153 
154  // TODO: perform averaging for cells sharing more than one boundary face
155 }
156 
157 
158 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
159 
161 (
164 );
165 
166 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
167 
168 } // End namespace RASModels
169 } // End namespace Foam
170 
171 // ************************************************************************* //
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
const char *const group
Group name for atomic constants.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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)
label index() const
Return the index of this patch in the fvBoundaryMesh.
Definition: fvPatch.H:173
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
bool updated() const
Return true if the boundary condition has already been updated.
Definition: fvPatchField.H:356
This boundary condition provides a turbulence stress normal to streamlines wall function condition fo...
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
dimensionedScalar pow025(const dimensionedScalar &ds)
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
const nearWallDist & y() const
Return the near wall distances.
label k
Boltzmann constant.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Macros for easy insertion into run-time selection tables.
scalar uTau
static const nutWallFunctionFvPatchScalarField & nutw(const momentumTransportModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
scalar y
This boundary condition provides a turbulent kinematic viscosity condition when using wall functions...
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:93
static word groupName(Name name, const word &group)
Foam::fvPatchFieldMapper.
v2WallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:325
makePatchTypeField(fvPatchScalarField, fWallFunctionFvPatchScalarField)
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
Abstract base class for turbulence models (RAS, LES and laminar).
label patchi
scalar yPlus
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:229
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensionedScalar pow4(const dimensionedScalar &ds)
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:167
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
Namespace for OpenFOAM.
const DimensionedField< Type, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:332