LaheyKEpsilon.C
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 \*---------------------------------------------------------------------------*/
25 
26 #include "LaheyKEpsilon.H"
28 #include "twoPhaseSystem.H"
29 #include "dragModel.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace RASModels
36 {
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
40 template<class BasicTurbulenceModel>
41 LaheyKEpsilon<BasicTurbulenceModel>::LaheyKEpsilon
42 (
43  const alphaField& alpha,
44  const rhoField& rho,
45  const volVectorField& U,
46  const surfaceScalarField& alphaRhoPhi,
47  const surfaceScalarField& phi,
48  const transportModel& transport,
49  const word& propertiesName,
50  const word& type
51 )
52 :
54  (
55  alpha,
56  rho,
57  U,
58  alphaRhoPhi,
59  phi,
60  transport,
61  propertiesName,
62  type
63  ),
64 
65  gasTurbulencePtr_(NULL),
66 
67  alphaInversion_
68  (
70  (
71  "alphaInversion",
72  this->coeffDict_,
73  0.3
74  )
75  ),
76 
77  Cp_
78  (
80  (
81  "Cp",
82  this->coeffDict_,
83  0.25
84  )
85  ),
86 
87  C3_
88  (
90  (
91  "C3",
92  this->coeffDict_,
93  this->C2_.value()
94  )
95  ),
96 
97  Cmub_
98  (
100  (
101  "Cmub",
102  this->coeffDict_,
103  0.6
104  )
105  )
106 {
107  if (type == typeName)
108  {
109  // Cannot correct nut yet: construction of the phases is not complete
110  // correctNut();
111 
112  this->printCoeffs(type);
113  }
114 }
115 
116 
117 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
118 
119 template<class BasicTurbulenceModel>
121 {
123  {
124  alphaInversion_.readIfPresent(this->coeffDict());
125  Cp_.readIfPresent(this->coeffDict());
126  C3_.readIfPresent(this->coeffDict());
127  Cmub_.readIfPresent(this->coeffDict());
128 
129  return true;
130  }
131  else
132  {
133  return false;
134  }
135 }
136 
137 
138 template<class BasicTurbulenceModel>
140 <
141  typename BasicTurbulenceModel::transportModel
142 >&
144 {
145  if (!gasTurbulencePtr_)
146  {
147  const volVectorField& U = this->U_;
148 
149  const transportModel& liquid = this->transport();
150  const twoPhaseSystem& fluid =
151  refCast<const twoPhaseSystem>(liquid.fluid());
152  const transportModel& gas = fluid.otherPhase(liquid);
153 
154  gasTurbulencePtr_ =
155  &U.db()
157  (
159  (
161  gas.name()
162  )
163  );
164  }
165 
166  return *gasTurbulencePtr_;
167 }
168 
169 
170 template<class BasicTurbulenceModel>
172 {
174  this->gasTurbulence();
175 
176  this->nut_ =
177  this->Cmu_*sqr(this->k_)/this->epsilon_
178  + Cmub_*gasTurbulence.transport().d()*gasTurbulence.alpha()
179  *(mag(this->U_ - gasTurbulence.U()));
180 
181  this->nut_.correctBoundaryConditions();
182 }
183 
184 
185 template<class BasicTurbulenceModel>
187 {
189  this->gasTurbulence();
190 
191  const transportModel& liquid = this->transport();
192  const twoPhaseSystem& fluid = refCast<const twoPhaseSystem>(liquid.fluid());
193  const transportModel& gas = fluid.otherPhase(liquid);
194 
195  volScalarField magUr(mag(this->U_ - gasTurbulence.U()));
196 
197  tmp<volScalarField> bubbleG
198  (
199  Cp_
200  *(
201  pow3(magUr)
202  + pow(fluid.drag(gas).CdRe()*liquid.nu()/gas.d(), 4.0/3.0)
203  *pow(magUr, 5.0/3.0)
204  )
205  *gas
206  /gas.d()
207  );
208 
209  return bubbleG;
210 }
211 
212 
213 template<class BasicTurbulenceModel>
216 {
217  const volVectorField& U = this->U_;
218  const alphaField& alpha = this->alpha_;
219  const rhoField& rho = this->rho_;
220 
221  const turbulenceModel& gasTurbulence = this->gasTurbulence();
222 
223  return
224  (
225  max(alphaInversion_ - alpha, scalar(0))
226  *rho
227  *min(gasTurbulence.epsilon()/gasTurbulence.k(), 1.0/U.time().deltaT())
228  );
229 }
230 
231 
232 template<class BasicTurbulenceModel>
234 {
235  const alphaField& alpha = this->alpha_;
236  const rhoField& rho = this->rho_;
237 
239  this->gasTurbulence();
240 
241  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
242 
243  return
244  alpha*rho*bubbleG()
245  + phaseTransferCoeff*gasTurbulence.k()
246  - fvm::Sp(phaseTransferCoeff, this->k_);
247 }
248 
249 
250 template<class BasicTurbulenceModel>
252 {
253  const alphaField& alpha = this->alpha_;
254  const rhoField& rho = this->rho_;
255 
257  this->gasTurbulence();
258 
259  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
260 
261  return
262  alpha*rho*this->C3_*this->epsilon_*bubbleG()/this->k_
263  + phaseTransferCoeff*gasTurbulence.epsilon()
264  - fvm::Sp(phaseTransferCoeff, this->epsilon_);
265 }
266 
267 
268 template<class BasicTurbulenceModel>
270 {
272 }
273 
274 
275 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
276 
277 } // End namespace RASModels
278 } // End namespace Foam
279 
280 // ************************************************************************* //
dimensionedScalar pow3(const dimensionedScalar &ds)
Continuous-phase k-epsilon model including bubble-generated turbulence.
Definition: LaheyKEpsilon.H:77
multiphaseSystem & fluid
Definition: createFields.H:10
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual tmp< fvScalarMatrix > kSource() const
virtual bool read()
Read model coefficients if they have changed.
A class for handling words, derived from string.
Definition: word.H:59
BasicTurbulenceModel::transportModel transportModel
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: kEpsilon.C:220
tmp< volScalarField > phaseTransferCoeff() const
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Namespace for OpenFOAM.
static word groupName(Name name, const word &group)
BasicTurbulenceModel::rhoField rhoField
Generic dimensioned Type class.
Standard k-epsilon turbulence model for incompressible and compressible flows including rapid distort...
Definition: kEpsilon.H:86
const phaseModel & otherPhase(const phaseModel &phase) const
Constant access the phase not given as an argument.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Class which solves the volume fraction equations for two phases.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual tmp< volScalarField > epsilon() const =0
Return the turbulence kinetic energy dissipation rate.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Templated abstract base class for multiphase compressible turbulence models.
Macros for easy insertion into run-time selection tables.
virtual tmp< fvScalarMatrix > epsilonSource() const
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
surfaceScalarField & phi
tmp< volScalarField > bubbleG() const
Abstract base class for turbulence models (RAS, LES and laminar).
const Time & time() const
Return time.
Definition: IOobject.C:251
static const word propertiesName
Default name of the turbulence properties dictionary.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
BasicTurbulenceModel::alphaField alphaField
dimensionedScalar deltaT() const
Return time step.
Definition: TimeState.C:79
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
virtual tmp< volScalarField > CdRe() const =0
Drag coefficient.
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:245
const dragModel & drag(const phaseModel &phase) const
Return the drag model for the given phase.
A class for managing temporary objects.
Definition: PtrList.H:118
U
Definition: pEqn.H:82