PrinceBlanch.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) 2018-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 "PrinceBlanch.H"
28 #include "fvcGrad.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace diameterModels
36 {
37 namespace coalescenceModels
38 {
41  (
45  );
46 }
47 }
48 }
49 
51 
52 
53 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
54 
57 (
58  const populationBalanceModel& popBal,
59  const dictionary& dict
60 )
61 :
62  coalescenceModel(popBal, dict),
63  C1_(dimensionedScalar::lookupOrDefault("C1", dict, dimless, 0.356)),
64  h0_
65  (
66  dimensionedScalar::lookupOrDefault
67  (
68  "h0",
69  dict,
70  dimLength,
71  1e-4
72  )
73  ),
74  hf_
75  (
76  dimensionedScalar::lookupOrDefault
77  (
78  "hf",
79  dict,
80  dimLength,
81  1e-8
82  )
83  ),
84  turbulence_(dict.lookup("turbulence")),
85  buoyancy_(dict.lookup("buoyancy")),
86  laminarShear_(dict.lookup("laminarShear"))
87 {
88  if (laminarShear_)
89  {
90  shearStrainRate_.set
91  (
92  new volScalarField
93  (
94  IOobject
95  (
96  "shearStrainRate",
97  popBal_.time().name(),
98  popBal_.mesh()
99  ),
100  popBal_.mesh(),
102  (
103  "shearStrainRate",
105  Zero
106  )
107  )
108  );
109  }
110 }
111 
112 
113 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
114 
116 {
117  if (laminarShear_)
118  {
119  shearStrainRate_() =
120  sqrt(2.0)*mag(symm(fvc::grad(popBal_.continuousPhase().U())));
121  }
122 }
123 
124 
127 (
128  volScalarField& coalescenceRate,
129  const label i,
130  const label j
131 )
132 {
133  const phaseModel& continuousPhase = popBal_.continuousPhase();
134  const sizeGroup& fi = popBal_.sizeGroups()[i];
135  const sizeGroup& fj = popBal_.sizeGroups()[j];
137  popBal_.mesh().lookupObject<uniformDimensionedVectorField>("g");
138 
139  const dimensionedScalar rij(1/(1/fi.dSph() + 1/fj.dSph()));
140 
141  const volScalarField collisionEfficiency
142  (
143  exp
144  (
145  - sqrt
146  (
147  pow3(rij)*continuousPhase.rho()
148  /(16*popBal_.sigmaWithContinuousPhase(fi.phase()))
149  )
150  *log(h0_/hf_)
151  *cbrt(popBal_.continuousTurbulence().epsilon())/pow(rij, 2.0/3.0)
152  )
153  );
154 
155  if (turbulence_)
156  {
157  coalescenceRate +=
158  (
159  C1_*pi*sqr(fi.dSph() + fj.dSph())
160  *cbrt(popBal_.continuousTurbulence().epsilon())
161  *sqrt(pow(fi.dSph(), 2.0/3.0) + pow(fj.dSph(), 2.0/3.0))
162  )
163  *collisionEfficiency;
164  }
165 
166  if (buoyancy_)
167  {
168  const dimensionedScalar Sij(pi/4*sqr(fi.dSph() + fj.dSph()));
169 
170  coalescenceRate +=
171  (
172  Sij
173  *mag
174  (
175  sqrt
176  (
177  2.14*popBal_.sigmaWithContinuousPhase(fi.phase())
178  /(continuousPhase.rho()*fi.dSph())
179  + 0.505*mag(g)*fi.dSph()
180  )
181  - sqrt
182  (
183  2.14*popBal_.sigmaWithContinuousPhase(fi.phase())
184  /(continuousPhase.rho()*fj.dSph())
185  + 0.505*mag(g)*fj.dSph()
186  )
187  )
188  )
189  *collisionEfficiency;
190  }
191 
192  if (laminarShear_)
193  {
194  coalescenceRate +=
195  pow3(fi.dSph() + fj.dSph())/6
196  *shearStrainRate_()*collisionEfficiency;
197  }
198 }
199 
200 
201 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
const Time & time() const
Return time.
Definition: IOobject.C:318
Base class for coalescence models.
const populationBalanceModel & popBal_
Reference to the populationBalanceModel.
Model of Prince and Blanch (1990). The coalescence rate is calculated by.
Definition: PrinceBlanch.H:262
virtual void addToCoalescenceRate(volScalarField &coalescenceRate, const label i, const label j)
Add to coalescenceRate.
Definition: PrinceBlanch.C:127
virtual void precompute()
Precompute diameter independent expressions.
Definition: PrinceBlanch.C:115
PrinceBlanch(const populationBalanceModel &popBal, const dictionary &dict)
Definition: PrinceBlanch.C:57
Model for tracking the evolution of a dispersed phase size distribution due to coalescence (synonymou...
const fvMesh & mesh() const
Return reference to the mesh.
Single size class fraction field representing a fixed particle volume as defined by the user through ...
Definition: sizeGroup.H:102
const dimensionedScalar & dSph() const
Return representative spherical diameter of the sizeGroup.
Definition: sizeGroupI.H:51
const phaseModel & phase() const
Return const-reference to the phase.
Definition: sizeGroupI.H:37
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
const word & name() const
Return const reference to name.
virtual const volScalarField & rho() const =0
Return the density field.
Calculate the gradient of the given field.
const dimensionedScalar e
Elementary charge.
addToRunTimeSelectionTable(coalescenceModel, AdachiStuartFokkink, dictionary)
defineTypeNameAndDebug(AdachiStuartFokkink, 0)
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensionedScalar exp(const dimensionedScalar &ds)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimless
const dimensionSet dimLength
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
const dimensionSet dimVelocity
dimensionedScalar cbrt(const dimensionedScalar &ds)
dictionary dict