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-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 "continuousGasKEpsilon.H"
27 #include "fvModels.H"
28 #include "fvConstraints.H"
29 #include "phaseSystem.H"
30 #include "dragModel.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace RASModels
38 {
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
42 template<class BasicMomentumTransportModel>
44 (
45  const alphaField& alpha,
46  const rhoField& rho,
47  const volVectorField& U,
48  const surfaceScalarField& alphaRhoPhi,
49  const surfaceScalarField& phi,
50  const viscosity& viscosity,
51  const word& type
52 )
53 :
54  kEpsilon<BasicMomentumTransportModel>
55  (
56  alpha,
57  rho,
58  U,
59  alphaRhoPhi,
60  phi,
61  viscosity,
62  type
63  ),
64 
65  liquidTurbulencePtr_(nullptr),
66 
67  nutEff_
68  (
69  IOobject
70  (
71  this->groupName("nutEff"),
72  this->runTime_.name(),
73  this->mesh_,
74  IOobject::READ_IF_PRESENT,
75  IOobject::AUTO_WRITE
76  ),
77  this->nut_
78  ),
79 
80  alphaInversion_
81  (
82  dimensioned<scalar>::lookupOrAddToDict
83  (
84  "alphaInversion",
85  this->coeffDict_,
86  0.7
87  )
88  )
89 {
90  if (type == typeName)
91  {
92  this->printCoeffs(type);
93  }
94 }
95 
96 
97 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
98 
99 template<class BasicMomentumTransportModel>
101 {
103  {
104  alphaInversion_.readIfPresent(this->coeffDict());
105 
106  return true;
107  }
108  else
109  {
110  return false;
111  }
112 }
113 
114 
115 template<class BasicMomentumTransportModel>
117 {
119 
120  const momentumTransportModel& liquidTurbulence = this->liquidTurbulence();
121 
122  const phaseModel& gas = refCast<const phaseModel>(this->properties());
123  const phaseSystem& fluid = gas.fluid();
124  const phaseModel& liquid = fluid.otherPhase(gas);
125 
130 
131  volScalarField thetal(liquidTurbulence.k()/liquidTurbulence.epsilon());
132  volScalarField rhodv(gas.rho() + virtualMass.Cvm()*liquid.rho());
133  volScalarField thetag
134  (
135  (rhodv/(18*liquid.rho()*liquid.thermo().nu()))*sqr(gas.d())
136  );
137  volScalarField expThetar
138  (
139  min
140  (
141  exp(min(thetal/thetag, scalar(50))),
142  scalar(1)
143  )
144  );
145  volScalarField omega((1 - expThetar)/(1 + expThetar));
146 
147  nutEff_ = omega*liquidTurbulence.nut();
148  fvConstraints::New(this->mesh_).constrain(nutEff_);
149 }
150 
151 
152 template<class BasicMomentumTransportModel>
155 {
156  if (!liquidTurbulencePtr_)
157  {
158  const volVectorField& U = this->U_;
159 
160  const phaseModel& gas = refCast<const phaseModel>(this->properties());
161  const phaseSystem& fluid = gas.fluid();
162  const phaseModel& liquid = fluid.otherPhase(gas);
163 
164  liquidTurbulencePtr_ =
165  &U.db().lookupType<momentumTransportModel>(liquid.name());
166  }
167 
168  return *liquidTurbulencePtr_;
169 }
170 
171 
172 template<class BasicMomentumTransportModel>
175 {
176  volScalarField blend
177  (
178  max
179  (
180  min
181  (
182  (this->alpha_ - scalar(0.5))/(alphaInversion_ - 0.5),
183  scalar(1)
184  ),
185  scalar(0)
186  )
187  );
188 
189  return volScalarField::New
190  (
191  this->groupName("nuEff"),
192  blend*this->nut_
193  + (1.0 - blend)*rhoEff()*nutEff_/this->rho_
194  + this->nu()
195  );
196 }
197 
198 
199 template<class BasicMomentumTransportModel>
202 {
203  const phaseModel& gas = refCast<const phaseModel>(this->properties());
204  const phaseSystem& fluid = gas.fluid();
205  const phaseModel& liquid = fluid.otherPhase(gas);
206 
211 
212  return volScalarField::New
213  (
214  this->groupName("rhoEff"),
215  gas.rho() + (virtualMass.Cvm() + 3.0/20.0)*liquid.rho()
216  );
217 }
218 
219 
220 template<class BasicMomentumTransportModel>
223 {
224  const volVectorField& U = this->U_;
225  const alphaField& alpha = this->alpha_;
226  const rhoField& rho = this->rho_;
227 
228  const momentumTransportModel& liquidTurbulence = this->liquidTurbulence();
229 
230  return
231  (
232  max(alphaInversion_ - alpha, scalar(0))
233  *rho
234  *min
235  (
236  liquidTurbulence.epsilon()/liquidTurbulence.k(),
237  1.0/U.time().deltaT()
238  )
239  );
240 }
241 
242 
243 template<class BasicMomentumTransportModel>
246 {
247  const momentumTransportModel& liquidTurbulence = this->liquidTurbulence();
248  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
249 
250  return
251  phaseTransferCoeff*liquidTurbulence.k()
252  - fvm::Sp(phaseTransferCoeff, this->k_);
253 }
254 
255 
256 template<class BasicMomentumTransportModel>
259 {
260  const momentumTransportModel& liquidTurbulence = this->liquidTurbulence();
261  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
262 
263  return
264  phaseTransferCoeff*liquidTurbulence.epsilon()
265  - fvm::Sp(phaseTransferCoeff, this->epsilon_);
266 }
267 
268 
269 template<class BasicMomentumTransportModel>
272 {
273  tmp<volScalarField> tk(this->k());
274 
276  (
277  this->groupName("R"),
278  ((2.0/3.0)*I)*tk() - (nutEff_)*dev(twoSymm(fvc::grad(this->U_)))
279  );
280 }
281 
282 
283 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
284 
285 } // End namespace RASModels
286 } // End namespace Foam
287 
288 // ************************************************************************* //
label k
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &)
Return a temporary field constructed from name,.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
virtual tmp< fvScalarMatrix > epsilonSource() const
virtual tmp< volSymmTensorField > sigma() const
Return the stress tensor [m^2/s^2].
virtual tmp< volScalarField > nuEff() const
Return the effective viscosity.
tmp< volScalarField > phaseTransferCoeff() const
virtual tmp< volScalarField > rhoEff() const
Return the effective density for the stress.
virtual tmp< fvScalarMatrix > kSource() const
continuousGasKEpsilon(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.
const momentumTransportModel & liquidTurbulence() const
Return the turbulence model for the liquid phase.
virtual bool read()
Re-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 correctNut()
Definition: kEpsilon.C:41
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
virtual const word & name() const
Return the name of the liquid.
Generic thermophysical properties class for a liquid in which the functions and coefficients for each...
Definition: liquid.H:53
scalar rho(scalar p, scalar T) const
Liquid density [kg/m^3].
Definition: liquidI.H:26
Abstract base class for turbulence models (RAS, LES and laminar).
virtual tmp< volScalarField > nut() const =0
Return the turbulence viscosity.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
virtual tmp< volScalarField > epsilon() const =0
Return the turbulence kinetic energy dissipation rate.
const phaseSystem & fluid() const
Return the system to which this phase belongs.
Definition: phaseModel.C:127
tmp< volScalarField > d() const
Return the Sauter-mean diameter.
Definition: phaseModel.C:145
virtual const volScalarField & rho() const =0
Return the density field.
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:111
const ModelType & lookupInterfacialModel(const phaseInterface &interface) const
Return a sub model for an interface.
A class for managing temporary objects.
Definition: tmp.H:55
virtual tmp< volScalarField > Cvm() const =0
Return the virtual mass coefficient.
Abstract base class for all fluid physical properties.
Definition: viscosity.H:50
A class for handling words, derived from string.
Definition: word.H:62
const scalar omega
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< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const VolField< Type > &)
Namespace for OpenFOAM.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
static const Identity< scalar > I
Definition: Identity.H:93
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488