epsilonWallFunctionFvPatchScalarField.H
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-2018 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 Class
25  Foam::epsilonWallFunctionFvPatchScalarField
26 
27 Description
28  This boundary condition provides a turbulence dissipation wall constraint
29  for low- and high-Reynolds number turbulence models.
30 
31  The condition can be applied to wall boundaries for which it
32  - calculates \c epsilon and \c G
33  - specifies the near-wall epsilon value
34 
35  where
36 
37  \vartable
38  epsilon | turblence dissipation field
39  G | turblence generation field
40  \endvartable
41 
42  The model switches between laminar and turbulent functions based on the
43  laminar-to-turbulent y+ value derived from kappa and E.
44 
45 Usage
46  \table
47  Property | Description | Required | Default value
48  Cmu | model coefficient | no | 0.09
49  kappa | Von Karman constant | no | 0.41
50  E | model coefficient | no | 9.8
51  \endtable
52 
53  Example of the boundary condition specification:
54  \verbatim
55  <patchName>
56  {
57  type epsilonWallFunction;
58  }
59  \endverbatim
60 
61 See also
62  Foam::fixedInternalValueFvPatchField
63  Foam::omegaWallFunctionFvPatchScalarField
64 
65 SourceFiles
66  epsilonWallFunctionFvPatchScalarField.C
67 
68 \*---------------------------------------------------------------------------*/
69 
70 #ifndef epsilonWallFunctionFvPatchScalarField_H
71 #define epsilonWallFunctionFvPatchScalarField_H
72 
73 #include "fixedValueFvPatchField.H"
74 
75 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
76 
77 namespace Foam
78 {
79 
80 class turbulenceModel;
81 
82 /*---------------------------------------------------------------------------*\
83  Class epsilonWallFunctionFvPatchScalarField Declaration
84 \*---------------------------------------------------------------------------*/
85 
86 class epsilonWallFunctionFvPatchScalarField
87 :
88  public fixedValueFvPatchField<scalar>
89 {
90 protected:
91 
92  // Protected data
93 
94  //- Tolerance used in weighted calculations
95  static scalar tolerance_;
96 
97  //- Cmu coefficient
98  scalar Cmu_;
99 
100  //- Von Karman constant
101  scalar kappa_;
102 
103  //- E coefficient
104  scalar E_;
105 
106  //- y+ at the edge of the laminar sublayer
107  scalar yPlusLam_;
108 
109  //- Local copy of turbulence G field
110  scalarField G_;
111 
112  //- Local copy of turbulence epsilon field
114 
115  //- Initialised flag
116  bool initialised_;
117 
118  //- Master patch ID
119  label master_;
120 
121  //- List of averaging corner weights
123 
124 
125  // Protected Member Functions
126 
127  //- Check the type of the patch
128  virtual void checkType();
129 
130  //- Write local wall function variables
131  virtual void writeLocalEntries(Ostream&) const;
132 
133  //- Set the master patch - master is responsible for updating all
134  // wall function patches
135  virtual void setMaster();
136 
137  //- Create the averaging weights for cells which are bounded by
138  // multiple wall function faces
139  virtual void createAveragingWeights();
141  //- Helper function to return non-const access to an epsilon patch
143  (
144  const label patchi
145  );
147  //- Main driver to calculate the turbulence fields
148  virtual void calculateTurbulenceFields
149  (
151  scalarField& G0,
153  );
154 
155  //- Calculate the epsilon and G
156  virtual void calculate
157  (
158  const turbulenceModel& turbulence,
159  const List<scalar>& cornerWeights,
160  const fvPatch& patch,
161  scalarField& G,
163  );
164 
165  //- Return non-const access to the master patch ID
166  virtual label& master()
167  {
168  return master_;
169  }
170 
171 
172 public:
173 
174  //- Runtime type information
175  TypeName("epsilonWallFunction");
176 
177 
178  // Constructors
179 
180  //- Construct from patch and internal field
182  (
183  const fvPatch&,
185  );
186 
187  //- Construct from patch, internal field and dictionary
189  (
190  const fvPatch&,
192  const dictionary&
193  );
194 
195  //- Construct by mapping given
196  // epsilonWallFunctionFvPatchScalarField
197  // onto a new patch
199  (
201  const fvPatch&,
203  const fvPatchFieldMapper&
204  );
205 
206  //- Construct as copy
208  (
210  );
211 
212  //- Construct and return a clone
213  virtual tmp<fvPatchScalarField> clone() const
214  {
216  (
218  );
219  }
220 
221  //- Construct as copy setting internal field reference
223  (
226  );
227 
228  //- Construct and return a clone setting internal field reference
230  (
232  ) const
233  {
235  (
237  );
238  }
239 
240  //- Destructor
242  {}
243 
244 
245  // Member functions
246 
247  // Access
248 
249  //- Return non-const access to the master's G field
250  scalarField& G(bool init = false);
251 
252  //- Return non-const access to the master's epsilon field
253  scalarField& epsilon(bool init = false);
254 
255 
256  // Evaluation functions
258  //- Update the coefficients associated with the patch field
259  virtual void updateCoeffs();
260 
261  //- Update the coefficients associated with the patch field
262  virtual void updateWeightedCoeffs(const scalarField& weights);
263 
264  //- Manipulate matrix
265  virtual void manipulateMatrix(fvMatrix<scalar>& matrix);
266 
267  //- Manipulate matrix with given weights
268  virtual void manipulateMatrix
269  (
270  fvMatrix<scalar>& matrix,
271  const scalarField& weights
272  );
273 
274 
275  // I-O
276 
277  //- Write
278  virtual void write(Ostream&) const;
279 };
280 
281 
282 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
283 
284 } // End namespace Foam
285 
286 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
287 
288 #endif
289 
290 // ************************************************************************* //
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
virtual void writeLocalEntries(Ostream &) const
Write local wall function variables.
const dimensionedScalar G0
Conductance quantum: default SI units: [S].
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 any number of values (e...
Definition: dictionary.H:137
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:60
scalarField G_
Local copy of turbulence G field.
TypeName("epsilonWallFunction")
Runtime type information.
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
virtual void manipulateMatrix(fvMatrix< scalar > &matrix)
Manipulate matrix.
static scalar tolerance_
Tolerance used in weighted calculations.
const dimensionedScalar epsilon0
Electric constant: default SI units: [F/m].
Abstract base class for turbulence models (RAS, LES and laminar).
virtual label & master()
Return non-const access to the master patch ID.
virtual tmp< fvPatchScalarField > clone() const
Construct and return a clone.
virtual void calculateTurbulenceFields(const turbulenceModel &turbulence, scalarField &G0, scalarField &epsilon0)
Main driver to calculate the turbulence fields.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
scalar yPlusLam_
y+ at the edge of the laminar sublayer
Foam::fvPatchFieldMapper.
scalarField epsilon_
Local copy of turbulence epsilon field.
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
This boundary condition provides a turbulence dissipation wall constraint for low- and high-Reynolds ...
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:340
scalarField & G(bool init=false)
Return non-const access to the master&#39;s G field.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual epsilonWallFunctionFvPatchScalarField & epsilonPatch(const label patchi)
Helper function to return non-const access to an epsilon patch.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
Info<< "Reading field U\"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
Definition: createFields.H:94
label patchi
List< List< scalar > > cornerWeights_
List of averaging corner weights.
virtual void updateWeightedCoeffs(const scalarField &weights)
Update the coefficients associated with the patch field.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual void checkType()
Check the type of the patch.
virtual void setMaster()
Set the master patch - master is responsible for updating all.
A class for managing temporary objects.
Definition: PtrList.H:53
virtual void createAveragingWeights()
Create the averaging weights for cells which are bounded by.
epsilonWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Namespace for OpenFOAM.
scalarField & epsilon(bool init=false)
Return non-const access to the master&#39;s epsilon field.