continuousGasKEpsilon.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 "continuousGasKEpsilon.H"
27 #include "fvOptions.H"
28 #include "twoPhaseSystem.H"
29 #include "dragModel.H"
30 #include "virtualMassModel.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace RASModels
37 {
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
41 template<class BasicTurbulenceModel>
42 continuousGasKEpsilon<BasicTurbulenceModel>::continuousGasKEpsilon
43 (
44  const alphaField& alpha,
45  const rhoField& rho,
46  const volVectorField& U,
47  const surfaceScalarField& alphaRhoPhi,
48  const surfaceScalarField& phi,
49  const transportModel& transport,
50  const word& propertiesName,
51  const word& type
52 )
53 :
55  (
56  alpha,
57  rho,
58  U,
59  alphaRhoPhi,
60  phi,
61  transport,
62  propertiesName,
63  type
64  ),
65 
66  liquidTurbulencePtr_(nullptr),
67 
68  nutEff_
69  (
70  IOobject
71  (
72  IOobject::groupName("nutEff", alphaRhoPhi.group()),
73  this->runTime_.timeName(),
74  this->mesh_,
77  ),
78  this->nut_
79  ),
80 
81  alphaInversion_
82  (
84  (
85  "alphaInversion",
86  this->coeffDict_,
87  0.7
88  )
89  )
90 {
91  if (type == typeName)
92  {
93  this->printCoeffs(type);
94  }
95 }
96 
97 
98 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
99 
100 template<class BasicTurbulenceModel>
102 {
104  {
105  alphaInversion_.readIfPresent(this->coeffDict());
106 
107  return true;
108  }
109  else
110  {
111  return false;
112  }
113 }
114 
115 
116 template<class BasicTurbulenceModel>
118 {
120 
121  const turbulenceModel& liquidTurbulence = this->liquidTurbulence();
122  const transportModel& gas = this->transport();
123  const twoPhaseSystem& fluid = refCast<const twoPhaseSystem>(gas.fluid());
124  const transportModel& liquid = fluid.otherPhase(gas);
125 
126  const virtualMassModel& virtualMass =
127  fluid.lookupSubModel<virtualMassModel>(gas, liquid);
128 
129  volScalarField thetal(liquidTurbulence.k()/liquidTurbulence.epsilon());
130  volScalarField rhodv(gas.rho() + virtualMass.Cvm()*liquid.rho());
131  volScalarField thetag((rhodv/(18*liquid.rho()*liquid.nu()))*sqr(gas.d()));
132  volScalarField expThetar
133  (
134  min
135  (
136  exp(min(thetal/thetag, scalar(50))),
137  scalar(1)
138  )
139  );
140  volScalarField omega((1 - expThetar)/(1 + expThetar));
141 
142  nutEff_ = omega*liquidTurbulence.nut();
143  fv::options::New(this->mesh_).correct(nutEff_);
144 }
145 
146 
147 template<class BasicTurbulenceModel>
148 const turbulenceModel&
150 {
151  if (!liquidTurbulencePtr_)
152  {
153  const volVectorField& U = this->U_;
154 
155  const transportModel& gas = this->transport();
156  const twoPhaseSystem& fluid =
157  refCast<const twoPhaseSystem>(gas.fluid());
158  const transportModel& liquid = fluid.otherPhase(gas);
159 
160  liquidTurbulencePtr_ =
162  (
164  (
166  liquid.name()
167  )
168  );
169  }
170 
171  return *liquidTurbulencePtr_;
172 }
173 
174 
175 template<class BasicTurbulenceModel>
178 {
179  volScalarField blend
180  (
181  max
182  (
183  min
184  (
185  (this->alpha_ - scalar(0.5))/(alphaInversion_ - 0.5),
186  scalar(1)
187  ),
188  scalar(0)
189  )
190  );
191 
192  return tmp<volScalarField>
193  (
194  new volScalarField
195  (
196  IOobject::groupName("nuEff", this->alphaRhoPhi_.group()),
197  blend*this->nut_
198  + (1.0 - blend)*rhoEff()*nutEff_/this->transport().rho()
199  + this->nu()
200  )
201  );
202 }
203 
204 
205 template<class BasicTurbulenceModel>
208 {
209  const transportModel& gas = this->transport();
210  const twoPhaseSystem& fluid = refCast<const twoPhaseSystem>(gas.fluid());
211  const transportModel& liquid = fluid.otherPhase(gas);
212 
213  const virtualMassModel& virtualMass =
214  fluid.lookupSubModel<virtualMassModel>(gas, liquid);
215 
216  return tmp<volScalarField>
217  (
218  new volScalarField
219  (
220  IOobject::groupName("rhoEff", this->alphaRhoPhi_.group()),
221  gas.rho() + (virtualMass.Cvm() + 3.0/20.0)*liquid.rho()
222  )
223  );
224 }
225 
226 
227 template<class BasicTurbulenceModel>
230 {
231  const volVectorField& U = this->U_;
232  const alphaField& alpha = this->alpha_;
233  const rhoField& rho = this->rho_;
234 
235  const turbulenceModel& liquidTurbulence = this->liquidTurbulence();
236 
237  return
238  (
239  max(alphaInversion_ - alpha, scalar(0))
240  *rho
241  *min
242  (
243  liquidTurbulence.epsilon()/liquidTurbulence.k(),
244  1.0/U.time().deltaT()
245  )
246  );
247 }
248 
249 
250 template<class BasicTurbulenceModel>
253 {
254  const turbulenceModel& liquidTurbulence = this->liquidTurbulence();
255  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
256 
257  return
258  phaseTransferCoeff*liquidTurbulence.k()
259  - fvm::Sp(phaseTransferCoeff, this->k_);
260 }
261 
262 
263 template<class BasicTurbulenceModel>
266 {
267  const turbulenceModel& liquidTurbulence = this->liquidTurbulence();
268  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
269 
270  return
271  phaseTransferCoeff*liquidTurbulence.epsilon()
272  - fvm::Sp(phaseTransferCoeff, this->epsilon_);
273 }
274 
275 
276 template<class BasicTurbulenceModel>
279 {
280  tmp<volScalarField> tk(this->k());
281 
283  (
285  (
286  IOobject
287  (
288  IOobject::groupName("R", this->alphaRhoPhi_.group()),
289  this->runTime_.timeName(),
290  this->mesh_,
293  ),
294  ((2.0/3.0)*I)*tk() - (nutEff_)*dev(twoSymm(fvc::grad(this->U_))),
295  tk().boundaryField().types()
296  )
297  );
298 }
299 
300 
301 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
302 
303 } // End namespace RASModels
304 } // End namespace Foam
305 
306 // ************************************************************************* //
BasicTurbulenceModel::transportModel transportModel
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:58
const modelType & lookupSubModel(const phasePair &key) const
Access a sub model between a phase pair.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
surfaceScalarField & phi
multiphaseSystem & fluid
Definition: createFields.H:11
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual tmp< fvScalarMatrix > epsilonSource() const
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual bool read()
Re-read model coefficients if they have changed.
label k
Boltzmann constant.
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.
Abstract base class for turbulence models (RAS, LES and laminar).
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
Class which solves the volume fraction equations for two phases.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensionedScalar exp(const dimensionedScalar &ds)
static const Identity< scalar > I
Definition: Identity.H:93
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)
virtual tmp< volScalarField > nuEff() const
Return the effective viscosity.
virtual tmp< fvScalarMatrix > kSource() const
const phaseModel & otherPhase(const phaseModel &phase) const
Return the phase not given as an argument.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
virtual tmp< volScalarField > Cvm() const =0
Return the virtual mass coefficient.
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:481
virtual tmp< volSymmTensorField > R() const
Return the Reynolds stress tensor.
BasicTurbulenceModel::alphaField alphaField
U
Definition: pEqn.H:72
tmp< volScalarField > phaseTransferCoeff() const
const Time & time() const
Return time.
Definition: IOobject.C:367
virtual tmp< volScalarField > epsilon() const =0
Return the turbulence kinetic energy dissipation rate.
virtual void correctNut()
Definition: kEpsilon.C:40
BasicTurbulenceModel::rhoField rhoField
virtual tmp< volScalarField > rhoEff() const
Return the effective density for the stress.
const turbulenceModel & liquidTurbulence() const
Return the turbulence model for the liquid phase.
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
A class for managing temporary objects.
Definition: PtrList.H:53
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:361
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
volScalarField & nu
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:53
Namespace for OpenFOAM.