LuoSvendsen.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) 2019-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 "LuoSvendsen.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace diameterModels
36 {
37 namespace binaryBreakupModels
38 {
41  (
45  );
46 }
47 }
48 }
49 
50 
51 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
52 
54 (
55  const populationBalanceModel& popBal,
56  const dictionary& dict
57 )
58 :
59  binaryBreakupModel(popBal, dict),
60  gammaUpperReg2by11_(),
61  gammaUpperReg5by11_(),
62  gammaUpperReg8by11_(),
63  C4_(dimensionedScalar::lookupOrDefault("C4", dict, dimless, 0.923)),
64  beta_(dimensionedScalar::lookupOrDefault("beta", dict, dimless, 2.05)),
65  minEddyRatio_
66  (
67  dimensionedScalar::lookupOrDefault("minEddyRatio", dict, dimless, 11.4)
68  ),
69  kolmogorovLengthScale_
70  (
71  IOobject
72  (
73  "kolmogorovLengthScale",
74  popBal_.time().name(),
75  popBal_.mesh()
76  ),
77  popBal_.mesh(),
79  (
80  "kolmogorovLengthScale",
81  dimLength,
82  Zero
83  )
84  )
85 {
86  List<Tuple2<scalar, scalar>> gammaUpperReg2by11Table;
87  List<Tuple2<scalar, scalar>> gammaUpperReg5by11Table;
88  List<Tuple2<scalar, scalar>> gammaUpperReg8by11Table;
89 
90  gammaUpperReg2by11Table.append(Tuple2<scalar, scalar>(0.0, 1.0));
91  gammaUpperReg5by11Table.append(Tuple2<scalar, scalar>(0.0, 1.0));
92  gammaUpperReg8by11Table.append(Tuple2<scalar, scalar>(0.0, 1.0));
93 
94  for (scalar z = 1e-2; z <= 10.0; z = z + 1e-2)
95  {
96  Tuple2<scalar, scalar> gamma2by11
97  (
98  z,
99  incGammaRatio_Q(2.0/11.0, z)
100  );
101 
102  Tuple2<scalar, scalar> gamma5by11
103  (
104  z,
105  incGammaRatio_Q(5.0/11.0, z)
106  );
107 
108  Tuple2<scalar, scalar> gamma8by11
109  (
110  z,
111  incGammaRatio_Q(8.0/11.0, z)
112  );
113 
114  gammaUpperReg2by11Table.append(gamma2by11);
115  gammaUpperReg5by11Table.append(gamma5by11);
116  gammaUpperReg8by11Table.append(gamma8by11);
117  }
118 
119  gammaUpperReg2by11_ =
121  (
122  "gamma2by11",
124  linearInterpolationWeights::typeName,
125  autoPtr<TableReader<scalar>>(nullptr),
126  gammaUpperReg2by11Table
127  );
128 
129  gammaUpperReg5by11_ =
131  (
132  "gamma5by11",
134  linearInterpolationWeights::typeName,
135  autoPtr<TableReader<scalar>>(nullptr),
136  gammaUpperReg5by11Table
137  );
138 
139  gammaUpperReg8by11_ =
141  (
142  "gamma8by11",
144  linearInterpolationWeights::typeName,
145  autoPtr<TableReader<scalar>>(nullptr),
146  gammaUpperReg8by11Table
147  );
148 }
149 
150 
151 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
152 
154 {
155  kolmogorovLengthScale_ =
156  pow025
157  (
158  pow3
159  (
160  popBal_.continuousPhase().fluidThermo().nu()
161  )
162  /popBal_.continuousTurbulence().epsilon()
163  );
164 }
165 
166 
167 void
169 (
170  volScalarField& binaryBreakupRate,
171  const label i,
172  const label j
173 )
174 {
175  const phaseModel& continuousPhase = popBal_.continuousPhase();
176  const sizeGroup& fi = popBal_.sizeGroups()[i];
177  const sizeGroup& fj = popBal_.sizeGroups()[j];
178 
179  const dimensionedScalar cf
180  (
181  pow(fi.x()/fj.x(), 2.0/3.0) + pow((1 - fi.x()/fj.x()), 2.0/3.0) - 1
182  );
183 
184  const volScalarField b
185  (
186  12*cf*popBal_.sigmaWithContinuousPhase(fi.phase())
187  /(
188  beta_*continuousPhase.rho()*pow(fj.dSph(), 5.0/3.0)
189  *pow(popBal_.continuousTurbulence().epsilon(), 2.0/3.0)
190  )
191  );
192 
193  const volScalarField xiMin(minEddyRatio_*kolmogorovLengthScale_/fj.dSph());
194 
195  const volScalarField tMin(b/pow(xiMin, 11.0/3.0));
196 
197  volScalarField integral(3/(11*pow(b, 8.0/11.0)));
198 
199  forAll(integral, celli)
200  {
201  integral[celli] *=
202  2*pow(b[celli], 3.0/11.0)*tgamma(5.0/11.0)
203  *(
204  gammaUpperReg5by11_->value(b[celli])
205  - gammaUpperReg5by11_->value(tMin[celli])
206  );
207  }
208 
209  binaryBreakupRate +=
210  C4_*(1 - popBal_.alphas())/fj.x()
211  *cbrt
212  (
213  popBal_.continuousTurbulence().epsilon()
214  /sqr(fj.dSph())
215  )
216  *integral;
217 }
218 
219 
220 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:178
Macros for creating standard TableReader-s.
Definition: TableReader.H:51
A 2-tuple for storing two objects of different types.
Definition: Tuple2.H:66
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Base class for binary breakup models that provide a breakup rate between a size class pair directly,...
Model of Luo and Svendsen (1996). The breakup rate is calculated by.
Definition: LuoSvendsen.H:231
virtual void precompute()
Precompute diameter independent expressions.
Definition: LuoSvendsen.C:153
virtual void addToBinaryBreakupRate(volScalarField &binaryBreakupRate, const label i, const label j)
Add to binary breakupRate.
Definition: LuoSvendsen.C:169
LuoSvendsen(const populationBalanceModel &popBal, const dictionary &dict)
Definition: LuoSvendsen.C:54
Model for tracking the evolution of a dispersed phase size distribution due to coalescence (synonymou...
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 dimensionedScalar & x() const
Return representative volume of the sizeGroup.
Definition: sizeGroupI.H:58
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
virtual const volScalarField & rho() const =0
Return the density field.
volScalarField & b
Definition: createFields.H:25
addToRunTimeSelectionTable(binaryBreakupModel, LehrMilliesMewes, dictionary)
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
const doubleScalar e
Definition: doubleScalar.H:106
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)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
const dimensionSet dimless
const dimensionSet dimLength
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar cbrt(const dimensionedScalar &ds)
scalar incGammaRatio_Q(const scalar a, const scalar x)
Normalised upper incomplete gamma function.
Definition: incGamma.C:221
dimensionedScalar pow025(const dimensionedScalar &ds)
dictionary dict