turbulenceModel.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) 2013-2015 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::turbulenceModel
26 
27 Group
28  grpTurbulence
29 
30 Description
31  Abstract base class for turbulence models (RAS, LES and laminar).
32 
33 SourceFiles
34  turbulenceModel.C
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #ifndef turbulenceModel_H
39 #define turbulenceModel_H
40 
41 #include "IOdictionary.H"
42 #include "primitiveFieldsFwd.H"
43 #include "volFieldsFwd.H"
44 #include "surfaceFieldsFwd.H"
45 #include "fvMatricesFwd.H"
46 #include "nearWallDist.H"
47 #include "geometricOneField.H"
48 
49 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50 
51 namespace Foam
52 {
53 
54 // Forward declarations
55 class fvMesh;
56 
57 /*---------------------------------------------------------------------------*\
58  Class turbulenceModel Declaration
59 \*---------------------------------------------------------------------------*/
60 
61 class turbulenceModel
62 :
63  public IOdictionary
64 {
65 
66 protected:
67 
68  // Protected data
69 
70  const Time& runTime_;
71  const fvMesh& mesh_;
72 
75  const surfaceScalarField& phi_;
76 
77  //- Near wall distance boundary field
79 
80 
81 private:
82 
83  // Private Member Functions
84 
85  //- Disallow default bitwise copy construct
87 
88  //- Disallow default bitwise assignment
89  void operator=(const turbulenceModel&);
90 
91 
92 public:
93 
94  //- Runtime type information
95  TypeName("turbulenceModel");
96 
97  //- Default name of the turbulence properties dictionary
98  static const word propertiesName;
99 
100 
101  // Constructors
102 
103  //- Construct from components
105  (
106  const volVectorField& U,
108  const surfaceScalarField& phi,
109  const word& propertiesName
110  );
111 
112 
113  //- Destructor
114  virtual ~turbulenceModel()
115  {}
116 
117 
118  // Member Functions
119 
120  //- Read model coefficients if they have changed
121  virtual bool read() = 0;
123  const Time& time() const
124  {
125  return runTime_;
126  }
128  const fvMesh& mesh() const
129  {
130  return mesh_;
131  }
132 
133  //- Const access to the coefficients dictionary
134  virtual const dictionary& coeffDict() const = 0;
135 
136  //- Helper function to return the name of the turbulence G field
137  inline word GName() const
138  {
139  return word(type() + ":G");
140  }
141 
142  //- Access function to velocity field
143  inline const volVectorField& U() const
144  {
145  return U_;
146  }
147 
148  //- Access function to phase flux field
149  inline const surfaceScalarField& alphaRhoPhi() const
150  {
151  return alphaRhoPhi_;
152  }
153 
154  //- Return the volumetric flux field
155  virtual tmp<surfaceScalarField> phi() const;
156 
157  //- Return the near wall distances
158  const nearWallDist& y() const
159  {
160  return y_;
161  }
162 
163  //- Return the laminar viscosity
164  virtual tmp<volScalarField> nu() const = 0;
165 
166  //- Return the laminar viscosity on patch
167  virtual tmp<scalarField> nu(const label patchi) const = 0;
168 
169  //- Return the turbulence viscosity
170  virtual tmp<volScalarField> nut() const = 0;
171 
172  //- Return the turbulence viscosity on patch
173  virtual tmp<scalarField> nut(const label patchi) const = 0;
174 
175  //- Return the effective viscosity
176  virtual tmp<volScalarField> nuEff() const = 0;
177 
178  //- Return the effective viscosity on patch
179  virtual tmp<scalarField> nuEff(const label patchi) const = 0;
180 
181  //- Return the laminar dynamic viscosity
182  virtual tmp<volScalarField> mu() const = 0;
183 
184  //- Return the laminar dynamic viscosity on patch
185  virtual tmp<scalarField> mu(const label patchi) const = 0;
186 
187  //- Return the turbulence dynamic viscosity
188  virtual tmp<volScalarField> mut() const = 0;
189 
190  //- Return the turbulence dynamic viscosity on patch
191  virtual tmp<scalarField> mut(const label patchi) const = 0;
192 
193  //- Return the effective dynamic viscosity
194  virtual tmp<volScalarField> muEff() const = 0;
195 
196  //- Return the effective dynamic viscosity on patch
197  virtual tmp<scalarField> muEff(const label patchi) const = 0;
198 
199  //- Return the turbulence kinetic energy
200  virtual tmp<volScalarField> k() const = 0;
201 
202  //- Return the turbulence kinetic energy dissipation rate
203  virtual tmp<volScalarField> epsilon() const = 0;
204 
205  //- Return the Reynolds stress tensor
206  virtual tmp<volSymmTensorField> R() const = 0;
207 
208  //- Validate the turbulence fields after construction
209  // Update derived fields as required
210  virtual void validate();
211 
212  //- Solve the turbulence equations and correct the turbulence viscosity
213  virtual void correct() = 0;
214 };
215 
216 
217 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
218 
219 } // End namespace Foam
220 
221 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
222 
223 #endif
224 
225 // ************************************************************************* //
virtual void validate()
Validate the turbulence fields after construction.
word GName() const
Helper function to return the name of the turbulence G field.
const surfaceScalarField & phi_
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
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual ~turbulenceModel()
Destructor.
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
const volVectorField & U() const
Access function to velocity field.
virtual tmp< volScalarField > muEff() const =0
Return the effective dynamic viscosity.
virtual tmp< volSymmTensorField > R() const =0
Return the Reynolds stress tensor.
virtual tmp< volScalarField > mut() const =0
Return the turbulence dynamic viscosity.
Generic GeometricField class.
TypeName("turbulenceModel")
Runtime type information.
const fvMesh & mesh_
Abstract base class for turbulence models (RAS, LES and laminar).
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
const surfaceScalarField & alphaRhoPhi_
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
virtual tmp< surfaceScalarField > phi() const
Return the volumetric flux field.
Forward declarations of the specialisations of Field<T> for scalar, vector and tensor.
const surfaceScalarField & alphaRhoPhi() const
Access function to phase flux field.
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
const Time & time() const
Distance calculation for cells with face on a wall. Searches pointNeighbours to find closest...
Definition: nearWallDist.H:51
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
virtual tmp< volScalarField > nuEff() const =0
Return the effective viscosity.
virtual const dictionary & coeffDict() const =0
Const access to the coefficients dictionary.
Forward declarations of fvMatrix specializations.
label patchi
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
virtual tmp< volScalarField > epsilon() const =0
Return the turbulence kinetic energy dissipation rate.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
nearWallDist y_
Near wall distance boundary field.
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
A class for managing temporary objects.
Definition: PtrList.H:54
virtual bool read()=0
Read model coefficients if they have changed.
virtual tmp< volScalarField > mu() const =0
Return the laminar dynamic viscosity.
const nearWallDist & y() const
Return the near wall distances.
Namespace for OpenFOAM.
const fvMesh & mesh() const
const volVectorField & U_