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-2018 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 "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>
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_(nullptr),
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  this->printCoeffs(type);
110  }
111 }
112 
113 
114 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
115 
116 template<class BasicTurbulenceModel>
118 {
120  {
121  alphaInversion_.readIfPresent(this->coeffDict());
122  Cp_.readIfPresent(this->coeffDict());
123  C3_.readIfPresent(this->coeffDict());
124  Cmub_.readIfPresent(this->coeffDict());
125 
126  return true;
127  }
128  else
129  {
130  return false;
131  }
132 }
133 
134 
135 template<class BasicTurbulenceModel>
137 <
138  typename BasicTurbulenceModel::transportModel
139 >&
141 {
142  if (!gasTurbulencePtr_)
143  {
144  const volVectorField& U = this->U_;
145 
146  const transportModel& liquid = this->transport();
147  const twoPhaseSystem& fluid =
148  refCast<const twoPhaseSystem>(liquid.fluid());
149  const transportModel& gas = fluid.otherPhase(liquid);
150 
151  gasTurbulencePtr_ =
152  &U.db()
154  (
156  (
158  gas.name()
159  )
160  );
161  }
162 
163  return *gasTurbulencePtr_;
164 }
165 
166 
167 template<class BasicTurbulenceModel>
169 {
171  this->gasTurbulence();
172 
173  this->nut_ =
174  this->Cmu_*sqr(this->k_)/this->epsilon_
175  + Cmub_*gasTurbulence.transport().d()*gasTurbulence.alpha()
176  *(mag(this->U_ - gasTurbulence.U()));
177 
178  this->nut_.correctBoundaryConditions();
179  fv::options::New(this->mesh_).correct(this->nut_);
180 
181  BasicTurbulenceModel::correctNut();
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  const dragModel& drag = fluid.lookupSubModel<dragModel>(gas, liquid);
196 
197  volScalarField magUr(mag(this->U_ - gasTurbulence.U()));
198 
199  tmp<volScalarField> bubbleG
200  (
201  Cp_
202  *(
203  pow3(magUr)
204  + pow(drag.CdRe()*liquid.nu()/gas.d(), 4.0/3.0)
205  *pow(magUr, 5.0/3.0)
206  )
207  *gas
208  /gas.d()
209  );
210 
211  return bubbleG;
212 }
213 
214 
215 template<class BasicTurbulenceModel>
218 {
219  const volVectorField& U = this->U_;
220  const alphaField& alpha = this->alpha_;
221  const rhoField& rho = this->rho_;
222 
223  const turbulenceModel& gasTurbulence = this->gasTurbulence();
224 
225  return
226  (
227  max(alphaInversion_ - alpha, scalar(0))
228  *rho
229  *min(gasTurbulence.epsilon()/gasTurbulence.k(), 1.0/U.time().deltaT())
230  );
231 }
232 
233 
234 template<class BasicTurbulenceModel>
236 {
237  const alphaField& alpha = this->alpha_;
238  const rhoField& rho = this->rho_;
239 
241  this->gasTurbulence();
242 
243  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
244 
245  return
246  alpha*rho*bubbleG()
247  + phaseTransferCoeff*gasTurbulence.k()
248  - fvm::Sp(phaseTransferCoeff, this->k_);
249 }
250 
251 
252 template<class BasicTurbulenceModel>
254 {
255  const alphaField& alpha = this->alpha_;
256  const rhoField& rho = this->rho_;
257 
259  this->gasTurbulence();
260 
261  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
262 
263  return
264  alpha*rho*this->C3_*this->epsilon_*bubbleG()/this->k_
265  + phaseTransferCoeff*gasTurbulence.epsilon()
266  - fvm::Sp(phaseTransferCoeff, this->epsilon_);
267 }
268 
269 
270 template<class BasicTurbulenceModel>
272 {
274 }
275 
276 
277 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
278 
279 } // End namespace RASModels
280 } // End namespace Foam
281 
282 // ************************************************************************* //
const modelType & lookupSubModel(const phasePair &key) const
Access a sub model between a phase pair.
BasicTurbulenceModel::rhoField rhoField
surfaceScalarField & phi
multiphaseSystem & fluid
Definition: createFields.H:11
const volVectorField & U() const
Access function to velocity field.
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)
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Generic dimensioned Type class.
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
Abstract base class for turbulence models (RAS, LES and laminar).
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:182
Class which solves the volume fraction equations for two phases.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
const transportModel & transport() const
Access function to incompressible transport model.
BasicTurbulenceModel::alphaField alphaField
Standard k-epsilon turbulence model for incompressible and compressible flows including rapid distort...
Definition: kEpsilon.H:83
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
BasicTurbulenceModel::transportModel transportModel
Templated abstract base class for multiphase compressible turbulence models.
const phaseModel & otherPhase(const phaseModel &phase) const
Return the phase not given as an argument.
tmp< volScalarField > phaseTransferCoeff() const
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
virtual tmp< fvScalarMatrix > kSource() const
dimensionedScalar pow3(const dimensionedScalar &ds)
U
Definition: pEqn.H:72
const Time & time() const
Return time.
Definition: IOobject.C:360
virtual tmp< volScalarField > epsilon() const =0
Return the turbulence kinetic energy dissipation rate.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
tmp< volScalarField > bubbleG() const
LaheyKEpsilon(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName, const word &type=typeName)
Construct from components.
Definition: LaheyKEpsilon.C:42
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:221
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:354
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:53
Namespace for OpenFOAM.