mixtureKEpsilon.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::RASModels::mixtureKEpsilon
26 
27 Description
28  Mixture k-epsilon turbulence model for two-phase gas-liquid systems
29 
30  The basic structure of the model is based on:
31  \verbatim
32  Behzadi, A., Issa, R. I., & Rusche, H. (2004).
33  Modelling of dispersed bubble and droplet flow at high phase fractions.
34  Chemical Engineering Science, 59(4), 759-770.
35  \endverbatim
36 
37  but an effective density for the gas is used in the averaging and an
38  alternative model for bubble-generated turbulence from:
39  \verbatim
40  Lahey Jr, R. T. (2005).
41  The simulation of multidimensional multiphase flows.
42  Nuclear Engineering and Design, 235(10), 1043-1060.
43  \endverbatim
44 
45  The default model coefficients are
46  \verbatim
47  mixtureKEpsilonCoeffs
48  {
49  Cmu 0.09;
50  C1 1.44;
51  C2 1.92;
52  C3 C2;
53  Cp 0.25; // Bubble generated turbulence
54  alphap 0.3; // Gas phase fraction below which
55  // bubble generated turbulence is included
56  sigmak 1.0;
57  sigmaEps 1.3;
58  }
59  \endverbatim
60 
61 SourceFiles
62  mixtureKEpsilon.C
63 
64 \*---------------------------------------------------------------------------*/
65 
66 #ifndef mixtureKEpsilon_H
67 #define mixtureKEpsilon_H
68 
69 #include "RASModel.H"
70 #include "eddyViscosity.H"
71 
72 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
73 
74 namespace Foam
75 {
76 namespace RASModels
77 {
78 
79 /*---------------------------------------------------------------------------*\
80  Class mixtureKEpsilon Declaration
81 \*---------------------------------------------------------------------------*/
82 
83 template<class BasicMomentumTransportModel>
84 class mixtureKEpsilon
85 :
86  public eddyViscosity<RASModel<BasicMomentumTransportModel>>
87 {
88  // Private Data
89 
91  *liquidTurbulencePtr_;
92 
93 
94  // Private Member Functions
95 
96  //- Return the turbulence model for the other phase
97  mixtureKEpsilon<BasicMomentumTransportModel>& liquidTurbulence() const;
98 
99 
100 protected:
101 
102  // Protected data
103 
104  // Model coefficients
114 
115  // Fields
119 
120  // Mixture fields
126 
127 
128  // Protected Member Functions
129 
131 
132  void correctInletOutlet
133  (
134  volScalarField& vsf,
135  const volScalarField& refVsf
136  ) const;
137 
138  void initMixtureFields();
139 
140  virtual void correctNut();
141 
142  tmp<volScalarField> Ct2() const;
143 
146  tmp<volScalarField> rhom() const;
147 
149  (
150  const volScalarField& fc,
151  const volScalarField& fd
152  ) const;
153 
155  (
156  const volScalarField& fc,
157  const volScalarField& fd
158  ) const;
159 
161  (
162  const surfaceScalarField& fc,
163  const surfaceScalarField& fd
164  ) const;
165 
167  virtual tmp<fvScalarMatrix> kSource() const;
168  virtual tmp<fvScalarMatrix> epsilonSource() const;
169 
170  //- Return the effective diffusivity for k
171  tmp<volScalarField> DkEff(const volScalarField& nutm) const
172  {
173  return volScalarField::New
174  (
175  "DkEff",
176  nutm/sigmak_
177  );
178  }
179 
180  //- Return the effective diffusivity for epsilon
182  {
183  return volScalarField::New
184  (
185  "DepsilonEff",
186  nutm/sigmaEps_
187  );
188  }
189 
190 public:
192  typedef typename BasicMomentumTransportModel::alphaField alphaField;
193  typedef typename BasicMomentumTransportModel::rhoField rhoField;
194  typedef typename BasicMomentumTransportModel::transportModel transportModel;
195 
196 
197  //- Runtime type information
198  TypeName("mixtureKEpsilon");
199 
200 
201  // Constructors
202 
203  //- Construct from components
205  (
206  const alphaField& alpha,
207  const rhoField& rho,
208  const volVectorField& U,
209  const surfaceScalarField& alphaRhoPhi,
210  const surfaceScalarField& phi,
211  const transportModel& transport,
212  const word& type = typeName
213  );
214 
215  //- Disallow default bitwise copy construction
216  mixtureKEpsilon(const mixtureKEpsilon&) = delete;
217 
218 
219  //- Destructor
220  virtual ~mixtureKEpsilon()
221  {}
222 
223 
224  // Member Functions
225 
226  //- Re-read model coefficients if they have changed
227  virtual bool read();
228 
229  //- Return the turbulence kinetic energy
230  virtual tmp<volScalarField> k() const
231  {
232  return k_;
233  }
234 
235  //- Return the turbulence kinetic energy dissipation rate
236  virtual tmp<volScalarField> epsilon() const
237  {
238  return epsilon_;
239  }
240 
241  //- Solve the turbulence equations and correct the turbulence viscosity
242  virtual void correct();
243 
244 
245  // Member Operators
246 
247  //- Disallow default bitwise assignment
248  void operator=(const mixtureKEpsilon&) = delete;
249 };
250 
251 
252 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
253 
254 } // End namespace RASModels
255 } // End namespace Foam
256 
257 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
258 
259 #ifdef NoRepository
260  #include "mixtureKEpsilon.C"
261 #endif
262 
263 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
264 
265 #endif
266 
267 // ************************************************************************* //
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
tmp< surfaceScalarField > mixFlux(const surfaceScalarField &fc, const surfaceScalarField &fd) const
tmp< volScalarField > mixU(const volScalarField &fc, const volScalarField &fd) const
autoPtr< volScalarField > epsilonm_
TypeName("mixtureKEpsilon")
Runtime type information.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
tmp< volScalarField > DepsilonEff(const volScalarField &nutm) const
Return the effective diffusivity for epsilon.
autoPtr< volScalarField > Ct2_
tmp< volScalarField > mix(const volScalarField &fc, const volScalarField &fd) const
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
autoPtr< volScalarField > rhom_
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
wordList epsilonBoundaryTypes(const volScalarField &epsilon) const
autoPtr< volScalarField > km_
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:49
tmp< volScalarField > rhom() const
tmp< volScalarField > bubbleG() const
tmp< volScalarField > rhogEff() const
A class for handling words, derived from string.
Definition: word.H:59
BasicMomentumTransportModel::alphaField alphaField
virtual bool read()
Re-read model coefficients if they have changed.
virtual tmp< fvScalarMatrix > epsilonSource() const
phi
Definition: correctPhi.H:3
BasicMomentumTransportModel::transportModel transportModel
virtual tmp< fvScalarMatrix > kSource() const
void operator=(const mixtureKEpsilon &)=delete
Disallow default bitwise assignment.
tmp< volScalarField > Ct2() const
BasicMomentumTransportModel::rhoField rhoField
U
Definition: pEqn.H:72
Mixture k-epsilon turbulence model for two-phase gas-liquid systems.
tmp< volScalarField > DkEff(const volScalarField &nutm) const
Return the effective diffusivity for k.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
mixtureKEpsilon(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &type=typeName)
Construct from components.
A class for managing temporary objects.
Definition: PtrList.H:53
void correctInletOutlet(volScalarField &vsf, const volScalarField &refVsf) const
tmp< volScalarField > rholEff() const
Namespace for OpenFOAM.
virtual ~mixtureKEpsilon()
Destructor.