All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
LaheyKEpsilon.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-2019 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::RASModels::LaheyKEpsilon
26 
27 Description
28  Continuous-phase k-epsilon model including bubble-generated turbulence.
29 
30  Reference:
31  \verbatim
32  Lahey Jr, R. T. (2005).
33  The simulation of multidimensional multiphase flows.
34  Nuclear Engineering and Design, 235(10), 1043-1060.
35  \endverbatim
36 
37  The default model coefficients are
38  \verbatim
39  LaheyKEpsilonCoeffs
40  {
41  Cmu 0.09;
42  C1 1.44;
43  C2 1.92;
44  C3 -0.33;
45  sigmak 1.0;
46  sigmaEps 1.3;
47  Cp 0.25;
48  Cmub 0.6;
49  alphaInversion 0.3;
50  }
51  \endverbatim
52 
53 SourceFiles
54  LaheyKEpsilon.C
55 
56 \*---------------------------------------------------------------------------*/
57 
58 #ifndef LaheyKEpsilon_H
59 #define LaheyKEpsilon_H
60 
61 #include "kEpsilon.H"
62 
63 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
64 
65 namespace Foam
66 {
67 namespace RASModels
68 {
69 
70 /*---------------------------------------------------------------------------*\
71  Class LaheyKEpsilon Declaration
72 \*---------------------------------------------------------------------------*/
73 
74 template<class BasicTurbulenceModel>
75 class LaheyKEpsilon
76 :
77  public kEpsilon<BasicTurbulenceModel>
78 {
79  // Private Data
80 
82  <
83  typename BasicTurbulenceModel::transportModel
84  > *gasTurbulencePtr_;
85 
86 
87  // Private Member Functions
88 
89  //- Return the turbulence model for the gas phase
91  <
92  typename BasicTurbulenceModel::transportModel
93  >&
94  gasTurbulence() const;
95 
96 
97 protected:
98 
99  // Protected data
100 
101  // Model coefficients
107 
108 
109  // Protected Member Functions
110 
111  virtual void correctNut();
114  virtual tmp<fvScalarMatrix> kSource() const;
115  virtual tmp<fvScalarMatrix> epsilonSource() const;
116 
117 
118 public:
120  typedef typename BasicTurbulenceModel::alphaField alphaField;
121  typedef typename BasicTurbulenceModel::rhoField rhoField;
122  typedef typename BasicTurbulenceModel::transportModel transportModel;
123 
124 
125  //- Runtime type information
126  TypeName("LaheyKEpsilon");
127 
128 
129  // Constructors
130 
131  //- Construct from components
133  (
134  const alphaField& alpha,
135  const rhoField& rho,
136  const volVectorField& U,
137  const surfaceScalarField& alphaRhoPhi,
138  const surfaceScalarField& phi,
139  const transportModel& transport,
140  const word& propertiesName = turbulenceModel::propertiesName,
141  const word& type = typeName
142  );
143 
144  //- Disallow default bitwise copy construction
145  LaheyKEpsilon(const LaheyKEpsilon&) = delete;
146 
147 
148  //- Destructor
149  virtual ~LaheyKEpsilon()
150  {}
151 
152 
153  // Member Functions
154 
155  //- Read model coefficients if they have changed
156  virtual bool read();
157 
158  //- Solve the turbulence equations and correct the turbulence viscosity
159  virtual void correct();
160 
161 
162  // Member Operators
163 
164  //- Disallow default bitwise assignment
165  void operator=(const LaheyKEpsilon&) = delete;
166 };
167 
168 
169 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
170 
171 } // End namespace RASModels
172 } // End namespace Foam
173 
174 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
175 
176 #ifdef NoRepository
177  #include "LaheyKEpsilon.C"
178 #endif
179 
180 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
181 
182 #endif
183 
184 // ************************************************************************* //
TypeName("LaheyKEpsilon")
Runtime type information.
BasicTurbulenceModel::rhoField rhoField
surfaceScalarField & phi
dimensionedScalar alphaInversion_
virtual bool read()
Read model coefficients if they have changed.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Continuous-phase k-epsilon model including bubble-generated turbulence.
Definition: LaheyKEpsilon.H:74
BasicTurbulenceModel::alphaField alphaField
Standard k-epsilon turbulence model for incompressible and compressible flows including rapid distort...
Definition: kEpsilon.H:83
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
BasicTurbulenceModel::transportModel transportModel
Templated abstract base class for multiphase compressible turbulence models.
tmp< volScalarField > phaseTransferCoeff() const
virtual ~LaheyKEpsilon()
Destructor.
virtual tmp< fvScalarMatrix > kSource() const
U
Definition: pEqn.H:72
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
tmp< volScalarField > bubbleG() const
LaheyKEpsilon(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName, const word &type=typeName)
Construct from components.
Definition: LaheyKEpsilon.C:42
void operator=(const LaheyKEpsilon &)=delete
Disallow default bitwise assignment.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
virtual tmp< fvScalarMatrix > epsilonSource() const
A class for managing temporary objects.
Definition: PtrList.H:53
Namespace for OpenFOAM.