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-2023 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 "fvModels.H"
28 #include "fvConstraints.H"
29 #include "phaseSystem.H"
30 #include "dispersedDragModel.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace RASModels
37 {
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
41 template<class BasicMomentumTransportModel>
43 (
44  const alphaField& alpha,
45  const rhoField& rho,
46  const volVectorField& U,
47  const surfaceScalarField& alphaRhoPhi,
48  const surfaceScalarField& phi,
49  const viscosity& viscosity,
50  const word& type
51 )
52 :
53  kEpsilon<BasicMomentumTransportModel>
54  (
55  alpha,
56  rho,
57  U,
58  alphaRhoPhi,
59  phi,
60  viscosity,
61  type
62  ),
63 
64  gasTurbulencePtr_(nullptr),
65 
66  alphaInversion_
67  (
68  dimensioned<scalar>::lookupOrAddToDict
69  (
70  "alphaInversion",
71  this->coeffDict_,
72  0.3
73  )
74  ),
75 
76  Cp_
77  (
78  dimensioned<scalar>::lookupOrAddToDict
79  (
80  "Cp",
81  this->coeffDict_,
82  0.25
83  )
84  ),
85 
86  C4_
87  (
88  dimensioned<scalar>::lookupOrAddToDict
89  (
90  "C4",
91  this->coeffDict_,
92  this->C2_.value()
93  )
94  ),
95 
96  Cmub_
97  (
98  dimensioned<scalar>::lookupOrAddToDict
99  (
100  "Cmub",
101  this->coeffDict_,
102  0.6
103  )
104  )
105 {
106  if (type == typeName)
107  {
108  this->printCoeffs(type);
109  }
110 }
111 
112 
113 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
114 
115 template<class BasicMomentumTransportModel>
117 {
119  {
120  alphaInversion_.readIfPresent(this->coeffDict());
121  Cp_.readIfPresent(this->coeffDict());
122  C4_.readIfPresent(this->coeffDict());
123  Cmub_.readIfPresent(this->coeffDict());
124 
125  return true;
126  }
127  else
128  {
129  return false;
130  }
131 }
132 
133 
134 template<class BasicMomentumTransportModel>
137 {
138  if (!gasTurbulencePtr_)
139  {
140  const volVectorField& U = this->U_;
141 
142  const phaseModel& liquid =
143  refCast<const phaseModel>(this->properties());
144  const phaseSystem& fluid = liquid.fluid();
145  const phaseModel& gas = fluid.otherPhase(liquid);
146 
147  gasTurbulencePtr_ =
148  &U.db().lookupObject
149  <
151  >
152  (
154  (
155  momentumTransportModel::typeName,
156  gas.name()
157  )
158  );
159  }
160 
161  return *gasTurbulencePtr_;
162 }
163 
164 
165 template<class BasicMomentumTransportModel>
167 {
168  const phaseCompressibleMomentumTransportModel& gasTurbulence =
169  this->gasTurbulence();
170 
171  const phaseModel& liquid = refCast<const phaseModel>(this->properties());
172  const phaseSystem& fluid = liquid.fluid();
173  const phaseModel& gas = fluid.otherPhase(liquid);
174 
175  this->nut_ =
176  this->Cmu_*sqr(this->k_)/this->epsilon_
177  + Cmub_*gas.d()*gasTurbulence.alpha()
178  *(mag(this->U_ - gasTurbulence.U()));
179 
180  this->nut_.correctBoundaryConditions();
181  fvConstraints::New(this->mesh_).constrain(this->nut_);
182 }
183 
184 
185 template<class BasicMomentumTransportModel>
187 {
188  const phaseCompressibleMomentumTransportModel& gasTurbulence =
189  this->gasTurbulence();
190 
191  const phaseModel& liquid = refCast<const phaseModel>(this->properties());
192  const phaseSystem& fluid = liquid.fluid();
193  const phaseModel& gas = fluid.otherPhase(liquid);
194 
198 
199  volScalarField magUr(mag(this->U_ - gasTurbulence.U()));
200 
201  tmp<volScalarField> bubbleG
202  (
203  Cp_
204  *(
205  pow3(magUr)
206  + pow(drag.CdRe()*liquid.fluidThermo().nu()/gas.d(), 4.0/3.0)
207  *pow(magUr, 5.0/3.0)
208  )
209  *gas
210  /gas.d()
211  );
212 
213  return bubbleG;
214 }
215 
216 
217 template<class BasicMomentumTransportModel>
220 {
221  const volVectorField& U = this->U_;
222  const alphaField& alpha = this->alpha_;
223  const rhoField& rho = this->rho_;
224 
225  const momentumTransportModel& gasTurbulence = this->gasTurbulence();
226 
227  return
228  (
229  max(alphaInversion_ - alpha, scalar(0))
230  *rho
231  *min(gasTurbulence.epsilon()/gasTurbulence.k(), 1.0/U.time().deltaT())
232  );
233 }
234 
235 
236 template<class BasicMomentumTransportModel>
238 {
239  const alphaField& alpha = this->alpha_;
240  const rhoField& rho = this->rho_;
241 
242  const phaseCompressibleMomentumTransportModel& gasTurbulence =
243  this->gasTurbulence();
244 
245  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
246 
247  return
248  alpha*rho*bubbleG()
249  + phaseTransferCoeff*gasTurbulence.k()
250  - fvm::Sp(phaseTransferCoeff, this->k_);
251 }
252 
253 
254 template<class BasicMomentumTransportModel>
257 {
258  const alphaField& alpha = this->alpha_;
259  const rhoField& rho = this->rho_;
260 
261  const phaseCompressibleMomentumTransportModel& gasTurbulence =
262  this->gasTurbulence();
263 
264  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
265 
266  return
267  alpha*rho*this->C4_*this->epsilon_*bubbleG()/this->k_
268  + phaseTransferCoeff*gasTurbulence.epsilon()
269  - fvm::Sp(phaseTransferCoeff, this->epsilon_);
270 }
271 
272 
273 template<class BasicMomentumTransportModel>
275 {
277 }
278 
279 
280 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
281 
282 } // End namespace RASModels
283 } // End namespace Foam
284 
285 // ************************************************************************* //
Generic GeometricField class.
static word groupName(Name name, const word &group)
Continuous-phase k-epsilon model including bubble-generated turbulence.
Definition: LaheyKEpsilon.H:78
virtual tmp< fvScalarMatrix > epsilonSource() const
Source term for the epsilon equation.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
LaheyKEpsilon(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity, const word &type=typeName)
Construct from components.
Definition: LaheyKEpsilon.C:43
tmp< volScalarField > bubbleG() const
tmp< volScalarField > phaseTransferCoeff() const
virtual void correctNut()
Correct the eddy-viscosity nut.
virtual tmp< fvScalarMatrix > kSource() const
Source term for the k equation.
virtual bool read()
Read model coefficients if they have changed.
Standard k-epsilon turbulence model for incompressible and compressible flows including rapid distort...
Definition: kEpsilon.H:86
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: kEpsilon.C:227
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Definition: dictionaryIO.C:100
Generic dimensioned Type class.
Class to represent a interface between phases where one phase is considered dispersed within the othe...
BasicMomentumTransportModel::alphaField alphaField
BasicMomentumTransportModel::rhoField rhoField
Generic thermophysical properties class for a liquid in which the functions and coefficients for each...
Definition: liquid.H:53
Abstract base class for turbulence models (RAS, LES and laminar).
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
const volVectorField & U() const
Access function to velocity field.
virtual tmp< volScalarField > epsilon() const =0
Return the turbulence kinetic energy dissipation rate.
Templated abstract base class for multiphase compressible turbulence models.
const alphaField & alpha() const
Access function to phase fraction.
tmp< volScalarField > d() const
Return the Sauter-mean diameter.
Definition: phaseModel.C:181
const word & name() const
Return the name of this phase.
Definition: phaseModel.C:145
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:73
const phaseModel & otherPhase(const phaseModel &phase) const
Return the phase not given as an argument in a two-phase system.
Definition: phaseSystemI.H:172
const ModelType & lookupInterfacialModel(const phaseInterface &interface) const
Return a sub model for an interface.
A class for managing temporary objects.
Definition: tmp.H:55
Abstract base class for all fluid physical properties.
Definition: viscosity.H:50
A class for handling words, derived from string.
Definition: word.H:62
Info<< "Reading strained laminar flame speed field Su\n"<< endl;volScalarField Su(IOobject("Su", runTime.name(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field betav\n"<< endl;volScalarField betav(IOobject("betav", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Lobs\n"<< endl;volScalarField Lobs(IOobject("Lobs", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field CT\n"<< endl;volSymmTensorField CT(IOobject("CT", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Nv\n"<< endl;volScalarField Nv(IOobject("Nv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field nsv\n"<< 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
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const VolField< Type > &)
Namespace for OpenFOAM.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488