LaheyKEpsilon.C
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 \*---------------------------------------------------------------------------*/
25 
26 #include "LaheyKEpsilon.H"
27 #include "fvOptions.H"
28 #include "phaseSystem.H"
29 #include "dragModel.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace RASModels
36 {
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
40 template<class BasicMomentumTransportModel>
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& type
50 )
51 :
53  (
54  alpha,
55  rho,
56  U,
57  alphaRhoPhi,
58  phi,
59  transport,
60  type
61  ),
62 
63  gasTurbulencePtr_(nullptr),
64 
65  alphaInversion_
66  (
68  (
69  "alphaInversion",
70  this->coeffDict_,
71  0.3
72  )
73  ),
74 
75  Cp_
76  (
78  (
79  "Cp",
80  this->coeffDict_,
81  0.25
82  )
83  ),
84 
85  C3_
86  (
88  (
89  "C3",
90  this->coeffDict_,
91  this->C2_.value()
92  )
93  ),
94 
95  Cmub_
96  (
98  (
99  "Cmub",
100  this->coeffDict_,
101  0.6
102  )
103  )
104 {
105  if (type == typeName)
106  {
107  this->printCoeffs(type);
108  }
109 }
110 
111 
112 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
113 
114 template<class BasicMomentumTransportModel>
116 {
118  {
119  alphaInversion_.readIfPresent(this->coeffDict());
120  Cp_.readIfPresent(this->coeffDict());
121  C3_.readIfPresent(this->coeffDict());
122  Cmub_.readIfPresent(this->coeffDict());
123 
124  return true;
125  }
126  else
127  {
128  return false;
129  }
130 }
131 
132 
133 template<class BasicMomentumTransportModel>
135 <
136  typename BasicMomentumTransportModel::transportModel
137 >&
139 {
140  if (!gasTurbulencePtr_)
141  {
142  const volVectorField& U = this->U_;
143 
144  const transportModel& liquid = this->transport();
145  const phaseSystem& fluid = liquid.fluid();
146  const transportModel& gas = fluid.otherPhase(liquid);
147 
148  gasTurbulencePtr_ =
149  &U.db().lookupObject
150  <
152  >
153  (
155  (
156  momentumTransportModel::typeName,
157  gas.name()
158  )
159  );
160  }
161 
162  return *gasTurbulencePtr_;
163 }
164 
165 
166 template<class BasicMomentumTransportModel>
168 {
170  gasTurbulence = this->gasTurbulence();
171 
172  this->nut_ =
173  this->Cmu_*sqr(this->k_)/this->epsilon_
174  + Cmub_*gasTurbulence.transport().d()*gasTurbulence.alpha()
175  *(mag(this->U_ - gasTurbulence.U()));
176 
177  this->nut_.correctBoundaryConditions();
178  fv::options::New(this->mesh_).correct(this->nut_);
179 }
180 
181 
182 template<class BasicMomentumTransportModel>
184 {
186  gasTurbulence = this->gasTurbulence();
187 
188  const transportModel& liquid = this->transport();
189  const phaseSystem& fluid = liquid.fluid();
190  const transportModel& gas = fluid.otherPhase(liquid);
191 
192  const dragModel& drag = fluid.lookupSubModel<dragModel>(gas, liquid);
193 
194  volScalarField magUr(mag(this->U_ - gasTurbulence.U()));
195 
196  tmp<volScalarField> bubbleG
197  (
198  Cp_
199  *(
200  pow3(magUr)
201  + pow(drag.CdRe()*liquid.thermo().nu()/gas.d(), 4.0/3.0)
202  *pow(magUr, 5.0/3.0)
203  )
204  *gas
205  /gas.d()
206  );
207 
208  return bubbleG;
209 }
210 
211 
212 template<class BasicMomentumTransportModel>
215 {
216  const volVectorField& U = this->U_;
217  const alphaField& alpha = this->alpha_;
218  const rhoField& rho = this->rho_;
219 
220  const momentumTransportModel& gasTurbulence = this->gasTurbulence();
221 
222  return
223  (
224  max(alphaInversion_ - alpha, scalar(0))
225  *rho
226  *min(gasTurbulence.epsilon()/gasTurbulence.k(), 1.0/U.time().deltaT())
227  );
228 }
229 
230 
231 template<class BasicMomentumTransportModel>
233 {
234  const alphaField& alpha = this->alpha_;
235  const rhoField& rho = this->rho_;
236 
238  gasTurbulence = this->gasTurbulence();
239 
240  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
241 
242  return
243  alpha*rho*bubbleG()
244  + phaseTransferCoeff*gasTurbulence.k()
245  - fvm::Sp(phaseTransferCoeff, this->k_);
246 }
247 
248 
249 template<class BasicMomentumTransportModel>
252 {
253  const alphaField& alpha = this->alpha_;
254  const rhoField& rho = this->rho_;
255 
257  gasTurbulence = 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 BasicMomentumTransportModel>
270 {
272 }
273 
274 
275 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
276 
277 } // End namespace RASModels
278 } // End namespace Foam
279 
280 // ************************************************************************* //
Generic thermophysical properties class for a liquid in which the functions and coefficients for each...
Definition: liquid.H:50
Templated abstract base class for multiphase compressible turbulence models.
const transportModel & transport() const
Access function to incompressible transport model.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual bool read()
Read model coefficients if they have changed.
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const modelType & lookupSubModel(const phasePair &key) const
Return a sub model between a phase pair.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Generic dimensioned Type class.
LaheyKEpsilon(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.
Definition: LaheyKEpsilon.C:42
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Continuous-phase k-epsilon model including bubble-generated turbulence.
Definition: LaheyKEpsilon.H:74
Info<< "Reading strained laminar flame speed field Su\"<< endl;volScalarField Su(IOobject("Su", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field betav\"<< endl;volScalarField betav(IOobject("betav", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Lobs\"<< endl;volScalarField Lobs(IOobject("Lobs", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field CT\"<< endl;volSymmTensorField CT(IOobject("CT", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Nv\"<< endl;volScalarField Nv(IOobject("Nv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field nsv\"<< endl;volSymmTensorField nsv(IOobject("nsv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);IOdictionary PDRProperties(IOobject("PDRProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE));autoPtr< PDRDragModel > drag
Definition: createFields.H:192
phi
Definition: pEqn.H:104
Standard k-epsilon turbulence model for incompressible and compressible flows including rapid distort...
Definition: kEpsilon.H:83
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:67
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
const phaseModel & otherPhase(const phaseModel &phase) const
Return the phase not given as an argument in a two-phase system.
Definition: phaseSystemI.H:112
BasicMomentumTransportModel::alphaField alphaField
Definition: kEpsilon.H:115
phaseSystem & fluid
Definition: createFields.H:11
tmp< volScalarField > phaseTransferCoeff() const
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
Abstract base class for turbulence models (RAS, LES and laminar).
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
virtual tmp< fvScalarMatrix > kSource() const
dimensionedScalar pow3(const dimensionedScalar &ds)
const volVectorField & U() const
Access function to velocity field.
U
Definition: pEqn.H:72
const Time & time() const
Return time.
Definition: IOobject.C:328
BasicMomentumTransportModel::rhoField rhoField
Definition: kEpsilon.H:116
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
BasicMomentumTransportModel::transportModel transportModel
Definition: kEpsilon.H:117
tmp< volScalarField > bubbleG() const
dimensioned< scalar > mag(const dimensioned< Type > &)
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
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: kEpsilon.C:217
virtual tmp< volScalarField > CdRe() const =0
Drag coefficient.
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
Definition: fvOptions.C:101
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:322
virtual tmp< volScalarField > epsilon() const =0
Return the turbulence kinetic energy dissipation rate.
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:53
Namespace for OpenFOAM.