mixtureKEpsilon.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2013-2016 OpenFOAM Foundation
9 Copyright (C) 2019-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Class
28 Foam::RASModels::mixtureKEpsilon
29
30Group
31 grpRASTurbulence
32
33Description
34 Mixture k-epsilon turbulence model for two-phase gas-liquid systems
35
36 The basic structure of the model is based on:
37 \verbatim
38 Behzadi, A., Issa, R. I., & Rusche, H. (2004).
39 Modelling of dispersed bubble and droplet flow at high phase fractions.
40 Chemical Engineering Science, 59(4), 759-770.
41 \endverbatim
42
43 but an effective density for the gas is used in the averaging and an
44 alternative model for bubble-generated turbulence from:
45 \verbatim
46 Lahey Jr, R. T. (2005).
47 The simulation of multidimensional multiphase flows.
48 Nuclear Engineering and Design, 235(10), 1043-1060.
49 \endverbatim
50
51 The default model coefficients are
52 \verbatim
53 mixtureKEpsilonCoeffs
54 {
55 Cmu 0.09;
56 C1 1.44;
57 C2 1.92;
58 C3 C2;
59 sigmak 1.0;
60 sigmaEps 1.3;
61 }
62 \endverbatim
63
64SourceFiles
65 mixtureKEpsilon.C
66
67\*---------------------------------------------------------------------------*/
68
69#ifndef mixtureKEpsilon_H
70#define mixtureKEpsilon_H
71
72#include "RASModel.H"
73#include "eddyViscosity.H"
74
75// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
76
77namespace Foam
78{
79namespace RASModels
80{
81
82/*---------------------------------------------------------------------------*\
83 Class mixtureKEpsilon Declaration
84\*---------------------------------------------------------------------------*/
85
86template<class BasicTurbulenceModel>
88:
89 public eddyViscosity<RASModel<BasicTurbulenceModel>>
90{
91 // Private data
92
93 mutable mixtureKEpsilon<BasicTurbulenceModel> *liquidTurbulencePtr_;
94
95
96 // Private Member Functions
97
98 //- No copy construct
99 mixtureKEpsilon(const mixtureKEpsilon&) = delete;
100
101 //- No copy assignment
102 void operator=(const mixtureKEpsilon&) = delete;
103
104 //- Return the turbulence model for the other phase
105 mixtureKEpsilon<BasicTurbulenceModel>& liquidTurbulence() const;
106
107
108protected:
109
110 // Protected data
111
112 // Model coefficients
121
122 // Fields
126
127 // Mixture fields
133
134
135 // Protected Member Functions
136
138
140 (
141 volScalarField& vsf,
142 const volScalarField& refVsf
143 ) const;
144
145 void initMixtureFields();
146
147 virtual void correctNut();
148
149 tmp<volScalarField> Ct2() const;
150
154
156 (
157 const volScalarField& fc,
158 const volScalarField& fd
159 ) const;
160
162 (
163 const volScalarField& fc,
164 const volScalarField& fd
165 ) const;
166
168 (
169 const surfaceScalarField& fc,
170 const surfaceScalarField& fd
171 ) const;
172
174 virtual tmp<fvScalarMatrix> kSource() const;
175 virtual tmp<fvScalarMatrix> epsilonSource() const;
176
177 //- Return the effective diffusivity for k
178 tmp<volScalarField> DkEff(const volScalarField& nutm) const
179 {
181 (
183 (
184 "DkEff",
185 nutm/sigmak_
186 )
187 );
188 }
189
190 //- Return the effective diffusivity for epsilon
192 {
194 (
196 (
197 "DepsilonEff",
198 nutm/sigmaEps_
199 )
200 );
201 }
202
203public:
205 typedef typename BasicTurbulenceModel::alphaField alphaField;
206 typedef typename BasicTurbulenceModel::rhoField rhoField;
207 typedef typename BasicTurbulenceModel::transportModel transportModel;
208
209
210 //- Runtime type information
211 TypeName("mixtureKEpsilon");
212
213
214 // Constructors
215
216 //- Construct from components
218 (
219 const alphaField& alpha,
220 const rhoField& rho,
221 const volVectorField& U,
222 const surfaceScalarField& alphaRhoPhi,
223 const surfaceScalarField& phi,
224 const transportModel& transport,
225 const word& propertiesName = turbulenceModel::propertiesName,
226 const word& type = typeName
227 );
228
229
230 //- Destructor
231 virtual ~mixtureKEpsilon() = default;
232
233
234 // Member Functions
235
236 //- Re-read model coefficients if they have changed
237 virtual bool read();
238
239 //- Return the turbulence kinetic energy
240 virtual tmp<volScalarField> k() const
241 {
242 return k_;
243 }
244
245 //- Return the turbulence kinetic energy dissipation rate
246 virtual tmp<volScalarField> epsilon() const
247 {
248 return epsilon_;
249 }
250
251 //- Solve the turbulence equations and correct the turbulence viscosity
252 virtual void correct();
253};
254
255
256// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
257
258} // End namespace RASModels
259} // End namespace Foam
260
261// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
262
263#ifdef NoRepository
264 #include "mixtureKEpsilon.C"
265#endif
266
267// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
268
269#endif
270
271// ************************************************************************* //
surfaceScalarField & phi
Mixture k-epsilon turbulence model for two-phase gas-liquid systems.
BasicTurbulenceModel::alphaField alphaField
virtual tmp< fvScalarMatrix > epsilonSource() const
tmp< volScalarField > DepsilonEff(const volScalarField &nutm) const
Return the effective diffusivity for epsilon.
BasicTurbulenceModel::rhoField rhoField
tmp< volScalarField > Ct2() const
tmp< volScalarField > rhogEff() const
tmp< volScalarField > mixU(const volScalarField &fc, const volScalarField &fd) const
virtual ~mixtureKEpsilon()=default
Destructor.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
tmp< volScalarField > rhom() const
tmp< volScalarField > bubbleG() const
autoPtr< volScalarField > rhom_
tmp< volScalarField > DkEff(const volScalarField &nutm) const
Return the effective diffusivity for k.
void correctInletOutlet(volScalarField &vsf, const volScalarField &refVsf) const
autoPtr< volScalarField > km_
autoPtr< volScalarField > Ct2_
virtual tmp< fvScalarMatrix > kSource() const
BasicTurbulenceModel::transportModel transportModel
tmp< surfaceScalarField > mixFlux(const surfaceScalarField &fc, const surfaceScalarField &fd) const
wordList epsilonBoundaryTypes(const volScalarField &epsilon) const
tmp< volScalarField > rholEff() const
autoPtr< volScalarField > epsilonm_
tmp< volScalarField > mix(const volScalarField &fc, const volScalarField &fd) const
TypeName("mixtureKEpsilon")
Runtime type information.
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
virtual bool read()
Re-read model coefficients if they have changed.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:58
A class for managing temporary objects.
Definition: tmp.H:65
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from Foam::string.
Definition: word.H:68
U
Definition: pEqn.H:72
Namespace for OpenFOAM.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
volScalarField & alpha
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:73