omegaWallFunctionFvPatchScalarField.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2016 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::omegaWallFunctionFvPatchScalarField
26 
27 Group
28  grpWallFunctions
29 
30 Description
31  This boundary condition provides a wall constraint on turbulnce specific
32  dissipation, omega for both low and high Reynolds number turbulence models.
33 
34  The near-wall omega may be either blended between the viscous region and
35  logarithmic region values using:
36 
37  \f[
38  \omega = sqrt(\omega_{vis}^2 + \omega_{log}^2)
39  \f]
40 
41  where
42 
43  \vartable
44  \omega_{vis} | omega in viscous region
45  \omega_{log} | omega in logarithmic region
46  \endvartable
47 
48  see eq.(15) of:
49  \verbatim
50  Menter, F., Esch, T.
51  "Elements of Industrial Heat Transfer Prediction"
52  16th Brazilian Congress of Mechanical Engineering (COBEM),
53  Nov. 2001
54  \endverbatim
55 
56  or switched between these values based on the laminar-to-turbulent y+ value
57  derived from kappa and E. Recent tests have shown that the standard
58  switching method provides more accurate results for 10 < y+ < 30 when used
59  with high Reynolds number wall-functions and both methods provide accurate
60  results when used with continuous wall-functions. Based on this the
61  standard switching method is used by default.
62 
63 Usage
64  \table
65  Property | Description | Required | Default value
66  Cmu | Model coefficient | no | 0.09
67  kappa | von Karman constant | no | 0.41
68  E | Model coefficient | no | 9.8
69  beta1 | Model coefficient | no | 0.075
70  blended | Blending switch | no | false
71  \endtable
72 
73  Example of the boundary condition specification:
74  \verbatim
75  <patchName>
76  {
77  type omegaWallFunction;
78  }
79  \endverbatim
80 
81 See also
82  Foam::fixedInternalValueFvPatchField
83  Foam::epsilonWallFunctionFvPatchScalarField
84 
85 SourceFiles
86  omegaWallFunctionFvPatchScalarField.C
87 
88 \*---------------------------------------------------------------------------*/
89 
90 #ifndef omegaWallFunctionFvPatchScalarField_H
91 #define omegaWallFunctionFvPatchScalarField_H
92 
93 #include "fixedValueFvPatchField.H"
94 
95 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
96 
97 namespace Foam
98 {
99 
100 class turbulenceModel;
101 
102 /*---------------------------------------------------------------------------*\
103  Class omegaWallFunctionFvPatchScalarField Declaration
104 \*---------------------------------------------------------------------------*/
105 
106 class omegaWallFunctionFvPatchScalarField
107 :
108  public fixedValueFvPatchField<scalar>
109 {
110 protected:
111 
112  // Protected data
113 
114  //- Tolerance used in weighted calculations
115  static scalar tolerance_;
116 
117  //- Cmu coefficient
118  scalar Cmu_;
119 
120  //- Von Karman constant
121  scalar kappa_;
122 
123  //- E coefficient
124  scalar E_;
125 
126  //- beta1 coefficient
127  scalar beta1_;
128 
129  //- beta1 coefficient
130  Switch blended_;
131 
132  //- y+ at the edge of the laminar sublayer
133  scalar yPlusLam_;
134 
135  //- Local copy of turbulence G field
136  scalarField G_;
137 
138  //- Local copy of turbulence omega field
140 
141  //- Initialised flag
142  bool initialised_;
144  //- Master patch ID
145  label master_;
146 
147  //- List of averaging corner weights
149 
150 
151  // Protected Member Functions
153  //- Check the type of the patch
154  virtual void checkType();
156  //- Write local wall function variables
157  virtual void writeLocalEntries(Ostream&) const;
159  //- Set the master patch - master is responsible for updating all
160  // wall function patches
161  virtual void setMaster();
162 
163  //- Create the averaging weights for cells which are bounded by
164  // multiple wall function faces
165  virtual void createAveragingWeights();
166 
167  //- Helper function to return non-const access to an omega patch
169  (
170  const label patchi
171  );
172 
173  //- Main driver to calculate the turbulence fields
174  virtual void calculateTurbulenceFields
175  (
177  scalarField& G0,
178  scalarField& omega0
179  );
180 
181  //- Calculate the omega and G
182  virtual void calculate
183  (
185  const List<scalar>& cornerWeights,
186  const fvPatch& patch,
187  scalarField& G,
189  );
190 
191  //- Return non-const access to the master patch ID
192  virtual label& master()
193  {
194  return master_;
195  }
196 
197 
198 public:
199 
200  //- Runtime type information
201  TypeName("omegaWallFunction");
202 
203 
204  // Constructors
205 
206  //- Construct from patch and internal field
208  (
209  const fvPatch&,
211  );
212 
213  //- Construct from patch, internal field and dictionary
215  (
216  const fvPatch&,
218  const dictionary&
219  );
220 
221  //- Construct by mapping given
222  // omegaWallFunctionFvPatchScalarField
223  // onto a new patch
225  (
227  const fvPatch&,
230  );
231 
232  //- Construct as copy
234  (
236  );
237 
238  //- Construct and return a clone
239  virtual tmp<fvPatchScalarField> clone() const
240  {
242  (
244  );
245  }
246 
247  //- Construct as copy setting internal field reference
249  (
252  );
253 
254  //- Construct and return a clone setting internal field reference
256  (
258  ) const
259  {
261  (
263  );
264  }
265 
266 
267  // Member functions
268 
269  // Access
270 
271  //- Return non-const access to the master's G field
272  scalarField& G(bool init = false);
273 
274  //- Return non-const access to the master's omega field
275  scalarField& omega(bool init = false);
277 
278  // Evaluation functions
279 
280  //- Update the coefficients associated with the patch field
281  virtual void updateCoeffs();
282 
283  //- Update the coefficients associated with the patch field
284  virtual void updateWeightedCoeffs(const scalarField& weights);
285 
286  //- Manipulate matrix
287  virtual void manipulateMatrix(fvMatrix<scalar>& matrix);
288 
289  //- Manipulate matrix with given weights
290  virtual void manipulateMatrix
291  (
292  fvMatrix<scalar>& matrix,
293  const scalarField& weights
294  );
295 
296 
297  // I-O
298 
299  //- Write
300  virtual void write(Ostream&) const;
301 };
302 
303 
304 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
305 
306 } // End namespace Foam
307 
308 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
309 
310 #endif
311 
312 // ************************************************************************* //
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 updateCoeffs()
Update the coefficients associated with the patch field.
const dimensionedScalar G0
Conductance quantum: default SI units: [S].
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
omegaWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
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
virtual tmp< fvPatchScalarField > clone() const
Construct and return a clone.
virtual void updateWeightedCoeffs(const scalarField &weights)
Update the coefficients associated with the patch field.
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
virtual void createAveragingWeights()
Create the averaging weights for cells which are bounded by.
virtual void checkType()
Check the type of the patch.
Abstract base class for turbulence models (RAS, LES and laminar).
scalarField & G(bool init=false)
Return non-const access to the master&#39;s G field.
virtual omegaWallFunctionFvPatchScalarField & omegaPatch(const label patchi)
Helper function to return non-const access to an omega patch.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Foam::fvPatchFieldMapper.
static scalar tolerance_
Tolerance used in weighted calculations.
virtual void calculateTurbulenceFields(const turbulenceModel &turbulence, scalarField &G0, scalarField &omega0)
Main driver to calculate the turbulence fields.
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
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:340
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
List< List< scalar > > cornerWeights_
List of averaging corner weights.
scalarField G_
Local copy of turbulence G field.
scalarField & omega(bool init=false)
Return non-const access to the master&#39;s omega field.
label patchi
virtual void writeLocalEntries(Ostream &) const
Write local wall function variables.
This boundary condition provides a wall constraint on turbulnce specific dissipation, omega for both low and high Reynolds number turbulence models.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual void setMaster()
Set the master patch - master is responsible for updating all.
A class for managing temporary objects.
Definition: PtrList.H:53
TypeName("omegaWallFunction")
Runtime type information.
virtual label & master()
Return non-const access to the master patch ID.
virtual void manipulateMatrix(fvMatrix< scalar > &matrix)
Manipulate matrix.
scalar yPlusLam_
y+ at the edge of the laminar sublayer
scalarField omega_
Local copy of turbulence omega field.
Namespace for OpenFOAM.
virtual void calculate(const turbulenceModel &turbulence, const List< scalar > &cornerWeights, const fvPatch &patch, scalarField &G, scalarField &omega)
Calculate the omega and G.