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