89,99 €
Theory and Application of Multiphase Lattice Boltzmann Methods presents a comprehensive review of all popular multiphase Lattice Boltzmann Methods developed thus far and is aimed at researchers and practitioners within relevant Earth Science disciplines as well as Petroleum, Chemical, Mechanical and Geological Engineering. Clearly structured throughout, this book will be an invaluable reference on the current state of all popular multiphase Lattice Boltzmann Methods (LBMs). The advantages and disadvantages of each model are presented in an accessible manner to enable the reader to choose the model most suitable for the problems they are interested in. The book is targeted at graduate students and researchers who plan to investigate multiphase flows using LBMs.
Throughout the text most of the popular multiphase LBMs are analyzed both theoretically and through numerical simulation. The authors present many of the mathematical derivations of the models in greater detail than is currently found in the existing literature. The approach to understanding and classifying the various models is principally based on simulation compared against analytical and observational results and discovery of undesirable terms in the derived macroscopic equations and sometimes their correction. A repository of FORTRAN codes for multiphase LBM models is also provided.
Sie lesen das E-Book in den Legimi-Apps auf:
Seitenzahl: 518
Veröffentlichungsjahr: 2015
Cover
Title Page
Copyright
Preface
About the companion website
Chapter 1: Introduction
1.1 History of the Lattice Boltzmann method
1.2 The Lattice Boltzmann method
1.3 Multiphase LBM
1.4 Comparison of models
1.5 Units in this book and parameter conversion
1.6 Appendix: Einstein summation convention
1.7 Use of the Fortran code in the book
Chapter 2: Single-component multiphase Shan–Chen-type model
2.1 Introduction
2.2 Typical equations of state
2.3 Thermodynamic consistency
2.4 Analytical surface tension
2.5 Contact angle
2.6 Capillary rise
2.7 Parallel flow and relative permeabilities
2.8 Forcing term in the SC model
2.9 Multirange pseudopotential (Inter-particle Force Model B)
2.10 Conclusions
2.11 Appendix A: Analytical solution for layered multiphase flow in a channel
2.12 Appendix B: FORTRAN code to simulate single component multiphase droplet contacting a wall, as shown in Figure 2.7(c)
Chapter 3: Shan and Chen-type multi-component multiphase models
3.1 Multi-component multiphase SC LBM
3.2 Derivation of the pressure
3.3 Determining and the surface tension
3.4 Contact angle
3.5 Flow through capillary tubes
3.6 Layered two-phase flow in a 2D channel
3.7 Pressure or velocity boundary conditions
3.8 Displacement in a 3D porous medium
Chapter 4: Rothman–Keller multiphase Lattice Boltzmann model
4.1 Introduction
4.2 RK color-gradient model
4.3 Theoretical analysis (Chapman–Enskog expansion)
4.4 Layered two-phase flow in a 2D channel
4.5 Interfacial tension and isotropy of the RK model
4.6 Drainage and capillary filling
4.7 MRT RK model
4.8 Contact angle
4.9 Tests of inlet/outlet boundary conditions
4.10 Immiscible displacements in porous media
4.11 Appendix A
4.12 Appendix B
Chapter 5: Free-energy-based multiphase Lattice Boltzmann model
5.1 Swift free-energy based single-component multiphase LBM
5.2 Chapman–Enskog expansion
5.3 Issue of Galilean invariance
5.4 Phase separation
5.5 Contact angle
5.6 Swift free-energy-based multi-component multiphase LBM
5.7 Appendix
Chapter 6: Inamuro's multiphase Lattice Boltzmann model
6.1 Introduction
6.2 Droplet collision
6.3 Appendix
Chapter 7: He–Chen–Zhang multiphase Lattice Boltzmann model
7.1 Introduction
7.2 HCZ model
7.3 Chapman–Enskog analysis
7.4 Surface tension and phase separation
7.5 Layered two-phase flow in a channel
7.6 Rayleigh–Taylor instability
7.7 Contact angle
7.8 Capillary rise
7.9 Geometric scheme to specify the contact angle and its hysteresis
7.10 Oscillation of an initially ellipsoidal droplet
7.11 Appendix A
7.12 Appendix B: 2D code
7.13 Appendix C: 3D code
Chapter 8: Axisymmetric multiphase HCZ model
8.1 Introduction
8.2 Methods
8.3 The Laplace law
8.4 Oscillation of an initially ellipsoidal droplet
8.5 Cylindrical liquid column break
8.6 Droplet collision
8.7 A revised axisymmetric HCZ model (Huang et al. 2014a)
8.8 Bubble rise
8.9 Conclusion
8.10 Appendix A: Chapman–Enskog analysis
Chapter 9: Extensions of the HCZ model for high-density ratio two-phase flows
9.1 Introduction
9.2 Model I (Lee and Lin 2005)
9.3 Model II (Amaya-Bower and Lee 2010)
9.4 Model III (Lee and Liu 2010)
9.5 Model IV
9.6 Numerical tests for different models
9.7 Conclusions
9.8 Appendix A: Analytical solutions for layered two-phase flow in a channel
9.9 Appendix B: 2D code based on Amaya-Bower and Lee (2010)
Chapter 10: Axisymmetric high-density ratio two-phase LBMs (extension of the HCZ model)
10.1 Introduction
10.2 The model based on Lee and Lin (2005)
10.3 Axisymmetric model based on Lee and Liu (2010)
References
Index
End User License Agreement
xi
xii
xiii
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
Cover
Table of Contents
Preface
Begin Reading
Chapter 1: Introduction
Figure 1.1 Discrete velocity models (a) D2Q7, (b) D2Q9, (c) D3Q15, and (d) D3Q19.
Figure 1.2 Citations in each year for the SC model (Shan and Chen 1993) (data from
Web of Science
).
Figure 1.3 Content of the book.
Chapter 2: Single-component multiphase Shan–Chen-type model
Figure 2.1 plot of van der Waals equations of state for . van der Waals constants atm and (Atkins 1978), and atm L .
l
and
g
represent the coexistence liquid and gas states at K.
Figure 2.2 diagram of C–S EOS (Eq. (2.30)). The constants are lu/(mu ), /mu, and lu/(ts tu). Curves of tu, tu, and tu are plotted.
Figure 2.3 plot of the SC equations of state. with and .
Figure 2.4 , and the integration as functions of . The C–S EOS (Eq. (2.30)) with is used. At the equilibrium state, and coexist.
Figure 2.5 Contact angle.
Figure 2.6 Different contact angles obtained through adjusting the parameter . The R–K EOS (Eq. (2.29)) with , and was used in the simulations.
Figure 2.7 Different equilibrium contact angles obtained by adjusting the parameter . C–S EOS with was used in the simulations. (a) , (b) , (c) , (d) . The contour of is shown in all frames. The unit of the density is mu/lu. The coexisting densities are approximately and .
Figure 2.8 Capillary rise. (a) Capillary rise at ts, (b) equilibrium rise height at ts, and (c) the height of liquid inside and outside (on the left boundary) the capillary tube as a function of time. The R–K EOS with was used in the simulations. The coexisting densities are approximately , and mu/lu.
Figure 2.9 Co-current immiscible two-phase flow in a 2D channel. The wetting phase (Fluid 2) flows along upper and lower plate while non-wetting phase (Fluid 1) flows in the center region.
Figure 2.10 Velocity profiles in the 2D channel with . The solid lines are analytical solutions and the squares are the LBM results. In (a) and (b) the wetting phase is less viscous (), while in (c) and (d) the wetting phase is more viscous (). lu/ts was applied only on (a) and (c), the central non-wetting Fluid 1 in Figure 2.9, and (b) and (d), the external wetting Fluid 2 in Figure 2.9, respectively. Velocity
u
has units lu/ts.
Figure 2.11 Relative permeabilities () as function of wetting saturation for two-phase flow in a 2D channel, . Lines show analytical solution and matched points from simulations.
Figure 2.12 Coexistence curve for the C–S EOS. The solid line is the analytical solution obtained from the Maxwell construction (corresponding to thermodynamic consistency). The numerical results were obtained from Scheme V coupled with Inter-particle Force Model A and the C–S EOS ().
Figure 2.13 Results obtained from the SC model with C-S-EOS and forcing Scheme III for different values of . (a) Reduced density of liquid, (b) reduced density of gas, (c) surface tension, and (d) liquid-to-gas density ratio as a function of . Reprinted from Huang et al. (2011a), copyright (2011), with permission from APS.
Figure 2.14 Results obtained from the original SC model (Scheme I) for different values of . (a) Reduced density of liquid, (b) reduced density of gas, (c) surface tension, and (d) liquid-to-gas density ratio as a function of for C–S EOS. Reprinted from Huang et al. (2011a), copyright (2011), with permission from APS.
Figure 2.15 Results obtained by the SC model in combination with Scheme V for different values of . (a) Reduced density of liquid, (b) reduced density of gas, (c) surface tension, and (d) density ratio as a function of for C–S EOS. Reprinted from Huang et al. (2011a), copyright (2011), with permission from APS.
Figure 2.16 Stability of forcing Schemes I, III, and V as a function of . The region to the upper right of each line is the corresponding stable region of each scheme. The C–S EOS is used with , and . Reprinted from Huang et al. (2011a), copyright (2011), with permission from APS.
Figure 2.17 Surface tension comparison between the SC model with forcing Scheme III and the analytical solution for different and . In the simulations, the inter-particle force including next-nearest neighbors for the C–S EOS is used and . Reprinted from Huang et al. (2011a), copyright (2011), with permission from APS.
Figure 2.18 Coexistence curve for the C–S EOS. The solid line is the analytical solution obtained from the Maxwell construction (corresponding to thermodynamic consistency). The numerical results were obtained from the original SC model with . Reprinted from Huang et al. (2011a), copyright (2011), with permission from APS.
Chapter 3: Shan and Chen-type multi-component multiphase models
Figure 3.1 Scaled component densities (left
y
axis) and scaled lattice surface tension (right-hand
y
axis) as a function of the scaled cohesion parameter . The simulations were carried out for two sets of initial densities. The vertical line depicts the value used in the rest of the study. depicts the critical value above which stable phase separation is possible (see text). Reprinted from Huang et al. (2007), copyright (2007), with permission from APS.
Figure 3.2 Simulations of different contact angles for multi-component fluids interacting with a surface (parameters , , , and measured contact angles for each case are listed in Table 3.1). The density of the first fluid is shown in dark grey. Simulation domain is . Images represent 24,000 time steps from an initial condition of a rectangle of the first fluid surrounded by the second fluid. In each region the fluids have density 2 and dissolved density 0.06. for both fluids. Reprinted from Huang et al. 2007, copyright (2007), with permission from APS.
Figure 3.3 3D simulations of different contact angles for multi-component fluids interacting with a surface (, , value is labeled in each frame). The density of the first fluid is shown in dark grey. Simulation domain is . Images represent 50,000 time steps from an initial condition of a cube of the first fluid surrounded by the second fluid. In each region, the fluids have density 2 and dissolved density 0.06. for both fluids.
Figure 3.4 Contact angles for two fluids interacting with a surface, versus the value of the adhesion parameter, (). Simulation domain is . Results represent 24,000 time steps from an initial condition of a (“measured 1”) and a (“measured 2”) rectangle of the first fluid surrounded by the second fluid. (a) in each region, the fluids have density 2 and dissolved density 0.06. (b) in each region, the fluids have density 8 and dissolved density 0.24. for both fluids. Reprinted from Huang et al. (2007), copyright (2007), with permission from APS.
Figure 3.5 Difference (in degrees) between the contact angle computed using Eq. (3.33) and the contact angle measured using the algorithm described in the text as a function of and (). Black dots show the parameter values used in simulations. Reprinted from Huang et al. (2007), copyright (2007), with permission from APS.
Figure 3.7 Primary drainage by fluid displacement simulations in a capillary tube with a square-shaped cross-section. (a) , , no displacement; (b) , , the hole is penetrated and the wetting fluid is displaced.
Figure 3.6 Primary drainage curve obtained by fluid displacement simulations in capillary tube with square-shaped cross-section. Reproduced by Jianlin Zhao, China University of Petroleum, Eastern China.
Figure 3.8 Velocity profiles comparison for a case of , , , , and when the body force is applied only to either Fluid 1 (a) or Fluid 2 (b). The numerical solution is obtained with the SC model. Source: Huang et al. 2011b. Reproduced with permission of Elsevier.
Figure 3.9 Dynamic distributions of the non-wetting phase during primary drainage for the GB1b porous medium system (Pan et al. 2004). The time interval between each frame is 3000 time steps. The transparent dark grey isosurface indicates the interface between the pore space and the solid phase. Dark grey and the pore space that is not dark grey represent the non-wetting and wetting fluid, respectively (produced by Jianlin Zhao, China University of Petroleum, Eastern China).
Chapter 4: Rothman–Keller multiphase Lattice Boltzmann model
Figure 4.1 Velocity profiles for cases with identical densities and applied to both fluids. (a) , , and ; (b) , , and .
Figure 4.2 Velocity profiles for cases with identical densities and applied to both fluids. The original RK approach is (a) , , and ; (b) , , and .
Figure 4.3 Velocity profiles for cases with identical densities and applied to both fluids. (a) , , and ; (b) , , and .
Figure 4.4 Comparison of the velocity profiles obtained from Cases A to H and the corresponding analytical solutions. The parameters of Cases A to H are illustrated in Table 4.2. ‘LBM’ and ‘LBM (revised)’ denote the original RK model and the RK model eliminating unwanted terms, respectively.
Figure 4.5 Comparison of the velocity profiles obtained from Cases I and J and the corresponding analytical solutions. (a) Case I: , ; (b) Case J: , . In these simulations and the unwanted terms are eliminated.
Figure 4.6 (a) Isotropy and the magnitude of maximum spurious current as functions of with . (b) Isotropy and the magnitude of maximum spurious current as functions of
A
with .
Figure 4.7 (Color online) Snapshots of the injection of a non-wetting fluid into two parallel capillary tubes with different , which is the pressure difference between the inlet and outlet. The widths of the left and right tubes are 24 lu and 32 lu, respectively. (a) , ts, (b) , i.e., , ts, and (c) , ts. The black, dark grey and light grey represent the solid, non-wetting fluid, and wetting fluid respectively, .
Figure 4.8 (a) Schematic diagram for liquid filling a capillary tube that initially contains gas.
H
is the width of the capillary tube and
l
is the length of the liquid column that has penetrated the capillary. The light grey and black are the non-wetting component (gas, Fluid 1) and the solid, respectively. (b) The length of the liquid column (lu) that has penetrated the capillary as a function of time. The initial penetration is about 15 lu at .
Figure 4.9 Equilibrium contact angles (measured from Fluid 1, dark grey one) for the present MRT RK model (, ). (a) , , , (b) , , (c) , , , and (d) , , .
Figure 4.10 Contact angles obtained from the present MRT RK model (measured from Fluid 1) compared with the analytical solution. Viscosity ratios , , , and are simulated.
Figure 4.11 Equilibrium of the coflow. Fluid 1 (dark grey, non-wetting) and Fluid 2 (light grey, wetting) are injected from the upper boundary (two fluids have an identical kinematic viscosity because ). (a) Inlet velocities for Fluids 1 and 2 are and , respectively. (b) Inlet velocities for Fluids 1 and 2 are and , respectively.
Figure 4.12 From stable displacement to capillary fingering (). (a) Case 7 (log ), (b) a case of , , , , , and (c) Case 13 (log). The grey and black represent the non-wetting fluid (Fluid 1) and solid, respectively. The white area is occupied by wetting fluid (Fluid 2).
Figure 4.13 Viscous fingering obtained from typical combination of different
M
and
Ca
. (a) Case 2 (log , log ), (b) Case 1 (log , log ), (c) Case 6 (log , log). The grey and black represent the non-wetting fluid (Fluid 1) and solid, respectively. The white area is occupied by wetting fluid (Fluid 2). Source: Huang et al. 2014b. Reproduced with permission of Elsevier.
Figure 4.14 Plot of the constant saturation domains for the simulations of the porous medium with different viscosity ratios and capillary numbers. The thick gray dashed lines represent the boundaries approximately drawn to separate the three domains. Source: Huang et al. 2014b. Reproduced with permission of Elsevier.
Chapter 5: Free-energy-based multiphase Lattice Boltzmann model
Figure 5.2 Coexistence curve of the van der Waals fluid. In the EOS (i.e., Eq. (5.7)), and . The corresponding critical density and temperature are and , respectively.
Figure 5.1 (a) Equilibrium state of a droplet immersed in gas. In the simulation, and . The equilibrium densities of liquid and gas are and , respectively. The maximum magnitude of the spurious current is lu/ts. (b) Equilibrium snapshot of the coexisting liquid and gas with flat interfaces for . Densities above and below the middle horizontal line are initialized as and , respectively. The equilibrium densities are and .
Figure 5.3 EOS with , , and . The values at points 1 and 2 are 4.08 and 2.92, respectively, and are at the limits of the non-physical portion of the EOS. In order to simulate spontaneous phase separation, the initial densities should be selected between these values. The equilibrium and are about 4.50 and 2.55, respectively.
Figure 5.4 Phase separation (a) ts, (b) ts, (c) ts, and (d) ts. In this case initially there is a small-amplitude perturbation in the density field (Eq. (5.82)), , . Density contours are shown in (b)–(d). The densities correspond to the greyscales shown in the key and are identical in all images.
Figure 5.5 Different contact angles obtained through adjusting the parameter . The computational domain is , , and . The , measured , and calculated are (a) , , , (b) , , , (c) , , , and (d) , , .
Figure 5.6 The equilibrium contact angle as a function of the density gradient at the wall, i.e., . The circles are angles measured from LBM results. The solid lines are calculated from Eq. (5.87): (a) and (b) . In the LBM simulation, the computational domain is and .
Chapter 6: Inamuro's multiphase Lattice Boltzmann model
Figure 6.1
p
–
v
diagram of the EOS (Eq. (6.10)). Here and are analogous to density and molar volume
v
, respectively.
Figure 6.2 Two-droplet collision (coalescence mode): lu/ts, . (a) Dimensionless time , (b) , (c) , and (d) .
Chapter 7: He–Chen–Zhang multiphase Lattice Boltzmann model
Figure 7.1 Equilibrium state of a droplet immersed in gas. In the simulation , and . The maximum velocity magnitude measured is lu/ts.
Figure 7.2 Layered two-phase flow simulation with . Fluid in the center of the channel is liquid with . There is gas flow near the channel walls with . (a) and (b) .
Figure 7.3 Layered two-phase flow simulation with . Fluid in the center of the channel is liquid with . There is gas flow near the channel walls with . (a) , and (b) .
Figure 7.4 Interface evolution for Rayleigh–Taylor instability simulation with lu/ts, , and . The non-dimensional time, which is normalized by , is labeled in each Figure Fifteen equal-interval contours at are drawn in the figures.
Figure 7.5 Interface evolution for Rayleigh–Taylor instability simulation with lu/ts, , and . The non-dimensional time, which is normalized by , is labeled in the Figure Fifteen equal-interval contours are drawn in each figure.
Figure 7.6 (a) The normalized positions and (b) velocities of the spike tip and bubble front as a function of time (both the above simulated cases and are shown, ). Time and length are normalized by and
W
, respectively.
Figure 7.7 Evolution of a droplet contacting with a wall (because , the expected equilibrium contact angle is about ). (a) , initially a semi-circle droplet is initialized, (b) , (c) , and (d) . In the simulation , and .
Figure 7.8 Equilibrium contact angles (initial condition is a semi-circular droplet as shown in Figure 7.7(a)). (a) , measured , (b) , measured , (c) , measured , and (d) , measured . In the simulation , and .
Figure 7.9 The equilibrium contact angle as a function of . The circles show results measured from LBM simulations, including the simulations in Figure 7.8, labeled as
a
,
b
,
c
, and
d
. The solid line connects points and in the coordinate system. In the simulation , and .
Figure 7.10 Capillary rise. (a) Capillary rise at ts, (b) equilibrium rise height at ts, and (c) the height of liquid inside and outside (on the left boundary) the capillary tube as a function of time. The coexisting densities are , and .
Figure 7.11 Schematic of the contact angle in 3D situations. Reprinted from Wang et al. (2013), copyright (2013), with permission from APS.
Figure 7.12 Different contact angles are achieved by putting different specified into Eq. (7.59). Reprinted from (Wang et al. (2013), copyright (2013), with permission from APS.
Figure 7.13 Four typical motion modes of contact points obtained by specifying the contact angle hysteresis effect in the HCZ model. . Reprinted from Wang et al. (2013), copyright (2013), with permission from APS.
Figure 7.14 Movement of a 3D hemispherical droplet, which is driven by a shear flow at (a) ts, (b) ts, and (c) ts. The hysteresis window is . The shapes of contact lines at different times are shown in the right-hand column. Reprinted from Wang et al. (2013), copyright (2013), with permission from APS.
Figure 7.15 Radii of the oscillatory spheroid as a function of time for Case I. The dashed and solid lines represent the radius measured from the center of the droplet in the
z
and
x
directions, respectively. Quantities are in lattice units.
Chapter 8: Axisymmetric multiphase HCZ model
Figure 8.1 Geometry of computational domain in the coordinates. is the axisymmetric axis.
Figure 8.2 The pattern of spurious currents (case of ). The maximum velocity magnitude is about 0.0075 lu/ts.
Figure 8.3 The surface tension as a function of in the axisymmetric HCZ model.
Figure 8.4 Drop shape at different times for Case III (parameters in Table 8.1). (a) ts, (b) ts, and (c) ts. The horizontal and vertical radii for the times shown in (a), (b), and (c) are also labeled in Figure 8.5. The left-hand column shows the shape of the droplet in an axisymmetric plane; the right-hand column is the corresponding 3D view of the droplet.
Figure 8.5 Radii of the oscillatory spheroid as a function of time for Cases I to IV. In each frame the dashed and solid lines represent the radius measured from the center of the droplet in the horizontal and vertical directions, respectively. Quantities are in lattice units.
Figure 8.6 The column of cylindrical liquid breaks up into two drops under an axisymmetric disturbance (case of lu). (a) , (b) ts, (c) ts, (d) ts, (e) ts, and (f) ts.
Figure 8.7 The final equilibrium states of liquid break up under an axisymmetric disturbance for four different cases. (a) lu, (b) lu, (c) lu, and (d) lu.
Figure 8.8 Theoretical and simulated (plotted as points from Table 8.2) main and satellite drop sizes as function of wave number .
Figure 8.9 Evolution of coalescence of two droplets (case of and ). (a) , (b) , (c) , (d) , (e) , (f) , (g) , and (h) . . Axes are in lattice units.
Figure 8.10 Evolution of droplet collision (case of and , in separation regime). The MRT version is used in the simulation. (a) , (b) , (c) , (d) , (e) , (f) , (g) , and (h) .
Figure 8.11 3D view of Figure 8.10 (case of and ). (a) , (b) , (c) , (d) , (e) , (f) , (g) , and (h) .
Figure 8.12 Evolution of droplet collision (case of and ). In the MRT simulation the gradient is calculated using a second-order finite difference scheme. (a) , (b) , (c) , (d) , (e) , (f) , (g) , and (h) .
Figure 8.13 Evolution of droplet collision (case of and , BGK simulation). (a) , (b) , (c) , (d) , (e) , (f) , (g) , and (h) .
Figure 8.14 Comparison of LBM results (Amaya-Bower and Lee (2010), left-hand column) and experimental data (Bhaga and Weber 1981, right-hand column). The Figure is copied from Table 2 of Amaya-Bower and Lee (2010). The Eötvös number and Morton number of cases A4, A5, and A6 (upper, middle, and lower rows) are listed in Table 8.3. Source: Amaya-Bower and Lee (2010). Reproduced with permission of Elsevier.
Figure 8.16 Comparison of terminal bubble shapes observed in experiments and predicted by LBM (A1–A9, for various
Eo
and
Mo
, Table 8.3, density ratios are 15.5). Left, middle, and right columns represent experiment results (Bhaga and Weber 1981), LBM results, and the descriptions of the bubble shape (Bhaga and Weber 1981), respectively. *Density ratio 3.
Figure 8.15 The top left and top right photographs are the experiment results (Bhaga and Weber 1981) for Cases A4 and A9, respectively. (a1) Bubble shape evolution of Case A4 with
revised
surface tension calculation, (a2) Case A4 with
original
surface tension calculation, (b1) Case A9 with
revised
surface tension calculation, and (b2) Case A9 with
original
surface tension calculation. The parameters in Cases A4 and A9 are listed in Table 8.3.
Figure 8.17 Comparison of terminal bubble shapes observed in Cases A5 and A6 with different density ratios. The left column represents experiment results (Bhaga and Weber 1981); the middle and right columns are LBM results with density ratio 3 and density ratio 15, respectively.
Figure 8.18 Terminal bubble wakes in (a) Case B1 (upper row), (b) Case B2 (middle row), and (c) Case B3 (bottom row). The left, middle, and right columns are results from experiment (Bhaga and Weber 1981), numerical work (Hua and Lou 2007), and our LBM, respectively. Source: Huang et al. 2014. Reproduced with permission of Elsevier.
Chapter 9: Extensions of the HCZ model for high-density ratio two-phase flows
Figure 9.1 The pressure as a function of (Eq. (8.9) with , and ).
Figure 9.2 Schematic diagram of a droplet splashing on a thin liquid film.
Figure 9.3 Snapshots of a droplet splashing on a thin liquid film with , and using Model I. The labeled non-dimensional time .
Figure 9.4 Velocity fields after 2,000,000 time steps. (a) Case 2 (Model I), (b) Case 4 (Model II), (c) Case 6 (Model III), and (d) Case 8 (Model ). Values of are magnified by , and times in (a), (b), (c), and (d), respectively.
Figure 9.5 Pressure fields for (a) Case 2 (Model I), (b) Case 4 (Model II), (c) Case 6 (Model III), and (d) Case 8 (Model ).
Figure 9.6 Immiscible parallel two-phase flow in a 2D channel.
Figure 9.7 Velocity profile comparison between analytical solution and LBM results for immiscible parallel two-phase flow when using Model I. The body force is applied to both liquid and gas. (a) Density ratio 50. (b) Density ratio 200. In Analytical Solution I the interface thickness is assumed to be zero. In Analytical Solution II a finite interface thickness is considered. The thicknesses are measured from the density profile (not shown), and the interface thicknesses in Cases (a) and (b) are 6 lu and 8 lu, respectively.
Figure 9.8 Velocity profile comparison between Analytical Solution I and LBM results using Model I for fine mesh. In LBM simulations 300 lu is used to represent the width of the channel and is applied. (a) Density ratio 50 (). (b) Density ratio 200 (). In the analytical solutions the interface thickness is assumed to be zero.
Figure 9.9 Couette flow velocity profiles for Model III. Case A: lower boundary moving with . Case B: upper boundary moving with . The density ratio is 50.
Figure 9.10 Density contours and velocity field for Model III. Case (a): density ratio 500 with velocities of the two plates lu/ts.
Figure 9.11 Density contours and velocity field for Model III. Case (b): density ratio 50 with velocities of the two plates lu/ts.
Figure 9.12 Droplet centriod velocity as function of time for Model III. (a) Density ratio 500. (b) Density ratio 50. The non-dimensional velocity is .
Figure 9.13 Schematic diagram for Analytical Solution II. Dark grey and light grey denote two different phases. (a) A layered two-phase flow with channel width
b
and interface thickness
W
. (b) Analytical Solution I for a channel width with interface at . (c) The velocity inside the interface region is assumed to be .
Chapter 10: Axisymmetric high-density ratio two-phase LBMs (extension of the HCZ model)
Figure 10.1 LBM simulation result for and .
Figure 10.2 3D view of Figure 10.1 (, ).
Figure 10.3 The measured spreading radius as a function of when and . The straight line is the power law .
Figure 10.4 Snapshots of two-droplet head-on collision for Case A (, ) (a) compared with the experimental result and (b) in Qian and Law (1997).
Figure 10.5 Snapshots of two-droplet head-on collision for Case B (, ) (a) compared with the experimental result and (b) in Qian and Law (1997).
Figure 10.6 Snapshots of two-droplet head-on collision for Case C (, ) (a) compared with the experimental result and (b) in Qian and Law (1997).
Figure 10.7 Gas entrapment in LBM simulation (Case A: , , ). Dimensionless times indicated.
Figure 10.8 Entrapped gas dissolution in LBM simulation (Case A: , , ). Dimensionless times indicated.
Figure 10.9 Bubble rise (Case I: , , ). The top right is the experimental snapshot for terminal bubble shape (Figure3(d) in Bhaga and Weber (1981)). The normalized time () is shown beside each bubble. The vertical black lines represent the axisymmetric axes.
Figure 10.10 (a) Mass of the bubble as a function time and (b) bubble rise velocity as a function of time for Case I with , , and .
Figure 10.11 Bubble rise (Case II: , , ). The normalized time is shown beside each bubble. The vertical black lines represent the axisymmetric axes.
Figure 10.12 (a) Mass of the bubble as a function time. (b) Bubble rise velocity as a function of time in Case II with , , and .
Chapter 1: Introduction
Table 1.1 Overview of the weighting coefficients and sound speeds.
Table 1.2 Citations of main articles about Lattice Boltzmann multiphase models.
Table 1.3 Comparison of Lattice Boltzmann multiphase models.
Table 1.4 Units in Lattice Boltzmann methods.
Table 1.5 Properties of an experimental multiphase system.
Chapter 2: Single-component multiphase Shan–Chen-type model
Table 2.1 Parameters in EOS (following Yuan and Schaefer (2006)).
Table 2.2 Critical properties in our study for vdW, R–K, and C–S EOS.
Table 2.3 Overview of the forcing schemes.
Table 2.4 Density ratios and surface tension obtained from the SC model with Schemes I and V for (C–S EOS is used).
Table 2.5 Density ratio for different surface tensions (SC model with forcing Scheme III, , C–S EOS).
Chapter 3: Shan and Chen-type multi-component multiphase models
Table 3.1 Adhesion parameters and contact angles (, degrees) for fluid 1 ( and ).
Table 3.2 Performance of the SC model for simulation of a layered two-phase flow with different viscosity ratios.
Chapter 4: Rothman–Keller multiphase Lattice Boltzmann model
Table 4.1 Units of parameters.
Table 4.2 Parameters for the RK model simulations of parallel flows with different densities (, , is applied to both fluids).
Table 4.3 Parameters for the RK model simulations of parallel flows with different densities (, is applied to both fluids).
Table 4.4 Interfacial tension as a function of for different viscosity ratios.
Table 4.5 Spurious current in the MRT and BGK RK simulations (cases of contact angle with different viscosity ratios ().
Table 4.6 Displacement in porous medium cases with different viscosity ratios
M
and capillary numbers
Ca
.
Chapter 7: He–Chen–Zhang multiphase Lattice Boltzmann model
Table 7.1 Theoretical and simulated oscillation period
T
of ellipsoidal drops ( lu).
Chapter 8: Axisymmetric multiphase HCZ model
Table 8.1 Theoretical and simulated oscillation period
T
of ellipsoidal drops ( lu and lu).
Table 8.2 Measured equilibrium drop sizes for cases of cylindrical liquid break up with lu.
Table 8.3 Cases for bubble rising (, , lu, tube diameter ).
Chapter 9: Extensions of the HCZ model for high-density ratio two-phase flows
Table 9.1 Cases of an initial stationary drop inside a box with theoretical and .
Chapter 10: Axisymmetric high-density ratio two-phase LBMs (extension of the HCZ model)
Table 10.2 Parameters in LBM simulations (in all these simulations mu/, lu).
Table 10.3 Parameters in bubble rise simulations (cases with , , , , expected , interface thickness 5 lu).
Haibo Huang
Department of Modern Mechanics,
University of Science and Technology of China
Michael C. Sukop
Department of Earth and Environment
Florida International University
Xi-Yun Lu
Department of Modern Mechanics
University of Science and Technology of China
This edition first published 2015 © 2015 by John Wiley & Sons, Ltd
Registered office: John Wiley & Sons, Ltd, The Atrium, Southern Gate, Chichester, West Sussex, PO19 8SQ, UK
Editorial offices: 9600 Garsington Road, Oxford, OX4 2DQ, UK
The Atrium, Southern Gate, Chichester, West Sussex, PO19 8SQ, UK
111 River Street, Hoboken, NJ 07030-5774, USA
For details of our global editorial offices, for customer services and for information about how to apply for permission to reuse the copyright material in this book please see our website at www.wiley.com/wiley-blackwell.
The right of the author to be identified as the author of this work has been asserted in accordance with the UK Copyright, Designs and Patents Act 1988.
All rights reserved. No part of this publication may be reproduced, stored in a retrieval system, or transmitted, in any form or by any means, electronic, mechanical, photocopying, recording or otherwise, except as permitted by the UK Copyright, Designs and Patents Act 1988, without the prior permission of the publisher.
Designations used by companies to distinguish their products are often claimed as trademarks. All brand names and product names used in this book are trade names, service marks, trademarks or registered trademarks of their respective owners. The publisher is not associated with any product or vendor mentioned in this book.
Limit of Liability/Disclaimer of Warranty: While the publisher and author(s) have used their best efforts in preparing this book, they make no representations or warranties with respect to the accuracy or completeness of the contents of this book and specifically disclaim any implied warranties of merchantability or fitness for a particular purpose. It is sold on the understanding that the publisher is not engaged in rendering professional services and neither the publisher nor the author shall be liable for damages arising herefrom. If professional advice or other expert assistance is required, the services of a competent professional should be sought.
Library of Congress Cataloging-in-Publication Data
Huang, Haibo (Engineering professor)
Multiphase lattice Boltzmann methods : theory and application / Haibo Huang, Michael C. Sukop, and Xi-Yun Lu.
pages cm
Includes bibliographical references and index.
ISBN 978-1-118-97133-8 (cloth)
1. Lattice Boltzmann methods. 2. Multiphase flow. 3. Fluid dynamics. 4. Fluid dynamics–Mathematical models. I. Sukop, Michael C. II. Lu, Xi-Yun. III. Title.
TA357.5.M84H83 2015
530.13′2–dc23
2015004408
A catalogue record for this book is available from the British Library.
Wiley also publishes its books in a variety of electronic formats. Some content that appears in print may not be available in electronic books.
Cover image: Water-oil distribution inside porous media: © Haibo Huang
When we began working with the multiphase Lattice Boltzmann method (LBM), we always asked which of the available models was the best or most suitable for the particular multiphase conditions we investigated and what the advantages and disadvantages of each model were. Over the years we have gained experience with these issues and in this book we share these experiences to hopefully make it easier for researchers to choose an appropriate model.
In this book most of the popular multiphase LBMs are analyzed both theoretically and through numerical simulation. We present many of the mathematical derivations of the models in greater detail than they can be found in the existing literature. Our approach to understanding and classifying the various models is principally based on simulation compared against analytical and observational results and discovery of undesirable terms in the derived macroscopic equations and sometimes their correction. However, we are not numerical analysts, and this limits our ability to fully characterize the numerics of the models. Our hope is that readers can understand some of the advantages and disadvantages of each model quickly and choose the model most suitable for the problems they are interested in.
More LBM books appear with each passing year. As far as we know, this book is the first one addressing all of the popular multiphase models. The book is intended as a reference book on the current state of multiphase LBMs and is targeted at graduate students and researchers who plan to investigate multiphase flows using LBMs.
This book grew out of roughly a decade of teaching and research at the University of Science and Technology of China and Florida International University, USA. We are grateful for the support we received from the National Natural Science Foundation of China (Grant Nos. 11172297 and 11472269), the Program for New Century Excellent Talents in University, Ministry of Education, China (NCET-12-0506), and Alexander von Humboldt Foundation, Germany. We also acknowledge funding from the United States National Science Foundation. This material is based upon work supported by the National Science Foundation under Grant No. 0440253. Any opinions, findings, and conclusions or recommendations are those of the authors and do not necessarily reflect the views of the National Science Foundation. In addition we thank our universities, departments, colleagues, and students for enabling us conduct our research and pursue the creation of this book. Finally, Jianlin Zhao, China University of Petroleum, Eastern China, produced two of the book's figures and Matlab code for the Maxwell construction, while Alejandro Garcia of Florida International University assisted by producing gfortran versions of the book's code, which are available on the internet.
Haibo Huang, Michael C. Sukop, Xi-Yun Lu
University of Science and Technology of China, Hefei, China
Florida International University, Miami, Florida
November 2014
This book is accompanied by a companion website: www.wiley.com/go/huang/ boltzmann
The website includes:
Compaq Fortran and GNU Fortran (GFortran) codes for downloading
Multiphase fluid phenomena and flows occur when two or more fluids that do not readily mix (such as air and water) share an interface. Multiphase fluid interactions are nearly ubiquitous in natural and industrial processes. Multiphase phenomena and flows can involve single component multiphase fluids, e.g., water and its own vapor, and multi-component multiphase fluids, e.g., oil/water. Some practical examples of multiphase fluid problems are the recovery and enhanced recovery of petroleum resources from reservoirs, non-aqueous phase liquid contamination of groundwater, soil water behavior, surface wetting phenomena, fuel cell operation, and the movement and evolution of clouds.
Computational fluid dynamics (CFD) has become very important in fluid flow studies. The Lattice Boltzmann method (LBM) has developed very quickly in the last two decades and has become a novel and powerful CFD tool – particularly for multiphase flows. The LBM has some major advantages compared to traditional CFD methods. First, it originates from Boltzmann's kinetic molecular dynamics – a more foundational level than normal continuum approaches. The LBM is able to recover the traditional macroscopic scale continuity and Navier–Stokes (N–S) equations, which are discretized and solved numerically in the common CFD methods. In the LBM, the more fundamental Boltzmann equation is directly discretized. Alternatively, the LBM can be viewed from its discrete-particle, more molecular-dynamics-like lattice gas origins. Second, in the LBM the pressure is usually related to the density through an ideal gas equation of state (for single-phase flow) or through a non-ideal van der Waals-like equation of state for some types of complex multiphase fluids. The pressure fields can be obtained directly once the density field is known. Hence, the Poisson equation – which can be computationally expensive – does not have to be solved in the LBM. The third advantage of the LBM is that the method is easy to parallelize due to the locality of much of the computation. Finally, no-slip boundary condition can be easily handled by simple bounce-back scheme.
The LBM has had great success in studies of single-phase flows, with commercial software known as POWERFLOW (Exa Corporation, https://www.exa.com/), based on the LBM, appearing about ten years ago. In contrast, multiphase LBMs are still undergoing development and there are many multiphase Lattice Boltzmann models available.
LBMs trace their roots to cellular automata, which were originally conceived by Stanislaw Ulam and John von Neumann in the 1940s. Cellular automata consist of a discretization of space on which individual cells exist in a particular state (say 0 or 1), and update their state at each time step according to a rule that takes as input the states of some set of the cell's neighbors. Sukop and Thorne (2006) provide an introduction to cellular automata. Wolfram (1983, 2002) studied simple cellular automata systematically and inspired some of the earliest application to fluids, leading to the first paper to propose a lattice gas cellular automaton (LGCA) for the N–S equations (Frisch et al. 1986). The use of a triangular grid restored some of the symmetry required to properly simulate fluids. Rothman and Zaleski (1997), Wolf-Gladrow (2000), Succi (2001), and Sukop and Thorne (2006) all provide instructive information on this model and the extensions that appeared. All of the LGCA models suffer from inherent defects, however, in particular the lack of Galilean invariance for fast flows and statistical noise (Qian et al. 1992, Wolf-Gladrow 2000). These are explicit particle-based Boolean models that include the random fluctuations that one would expect at a molecular level of gas simulation and hence required extensive averaging to recover the smooth behavior expected at macroscopic scales.
A second major step towards the modern LBM was taken by McNamara and Zanetti (1988), who dispensed with the individual particles of the LGCAs and replaced them with an averaged but still directionally discrete distribution function. This completely eliminated the statistical noise of the LGCA. A major simplification was introduced by Qian et al. (1992): the collision matrix of Higuera et al. (1989) is replaced by a single relaxation time, leading to the Bhatnagar, Gross, and Krook (BGK) model. After that, the LBM developed very quickly. Sukop and Thorne (2006) showed that there were fewer than 20 papers on the topic in 1992; more than 600 were published in 2013.
Later Lallemand and Luo (2000) and Luo (1998) showed that the LBM can be derived from the continuous Boltzmann equation (Boltzmann 1964/1995). Hence, it can be considered as a special discretized form of the Boltzmann equation (Nourgaliev et al. 2003). From the Chapman–Enskog expansion (Wolf-Gladrow 2000), the governing continuity and N–S equations can be recovered from the LBM. Without solving Poisson's equation, the pressure field can be obtained directly from the density distributions.
Today, the use of LBM spans a broad variety of disciplines. For example, an overview of the LBM for material science and engineering can be found in Raabe (2004). Application of the LBM to biophysics can be found in Boyd et al. (2005) and Sun et al. (2003).
The LBM can be derived from the BGK approximation of the Boltzmann equation (He and Luo 1997),
where is the single-particle distribution function in the phase space , and is the Maxwell–Boltzmann distribution function. x is the position vector, is the microscopic velocity, F(x,t) is a body force, and is the relaxation time, which determines the kinematic viscosity.
In the lattice BGK method, a discrete distribution function is introduced to represent the fluid. This distribution function satisfies the following Lattice Boltzmann equation (He and Luo 1997):
where is the density distribution function related to the discrete velocity direction i and is a relaxation time, which is related to the kinematic viscosity by , where is the sound speed. is the source term added into the standard Lattice Boltzmann equation. The equilibrium distribution function can be calculated as (Luo 1998)
In Eqs (1.2)) and (1.3)) the are the discrete velocities, as defined below, and s are weights, as given in Table 1.1. is the macroscopic density and u is the macroscopic velocity vector. Discrete velocity models are usually specified as DnQm, where n is the space dimension and m is the number of velocities. The popular 2D and 3D discrete velocity models are D2Q7, D2Q9, D3Q15, and D3Q19, which are shown in Figure 1.1.
Table 1.1 Overview of the weighting coefficients and sound speeds.
Model
D2Q7
(
),
(
),
D2Q9
(
),
(
),
(
)
D3Q15
(
),
(
),
(
)
D3Q19
(
),
(
),
(
)
Figure 1.1 Discrete velocity models (a) D2Q7, (b) D2Q9, (c) D3Q15, and (d) D3Q19.
For the D2Q7 model (Frisch et al. 1986), the discrete velocities are
For the D2Q9 model, the discrete velocities are given by (Qian et al. 1992)
For the D3Q15 model (Wolf-Gladrow 2000), the velocities are
For the D3Q19 model (Wolf-Gladrow 2000), they are
In the above equations, c is the lattice speed and is defined as . Here, we define 1 lattice unit () as 1 lu, 1 time step () as 1 ts, and 1 mass unit as 1 mu. There are other velocity models available, for example the D3Q27 model (He and Luo 1997), but we do not use them in simulations in this book.
In Eq. (1.3) s are weighting coefficients that can be derived theoretically (He and Luo 1997). can be derived from
where when , otherwise and we use the Einstein summation convention as detailed in the appendix to this chapter. Hence, or . As a detailed example, the computation of for the D2Q9 model is given in the following (calculation of each term from to is shown):
while the contribution from () is
In Eq. (1.3) is the density of the fluid, which can be obtained from:
This is simply the sum of the , revealing them as portions of the overall density associated with one of the discrete velocity directions. For , the macroscopic fluid velocity is given by
or in terms of the vector components of u as
which means the discrete velocities weighted by the directional densities. Application examples for viscous single-phase flow can be found in Yu et al. (2003), Dünweg and Ladd (2009), Aidun and Clausen (2010), and many others. A discussion on the H theorem in the context of the LBM can be found in Succi et al. (2002).
Numerous macroscopic numerical methods have been developed for solving the two-phase N–S equations (Scardovelli and Zaleski 1999), such as the front-tracking method, the volume-of-fluid (VOF) method, the level set method, and so on. The first three methods are the most popular ones. However, the front-tracking method is usually not able to simulate interface coalescence or break-up (Liu et al. 2012; Scardovelli and Zaleski 1999). In the VOF and level set methods an interface reconstruction step or interface reinitialization is usually required, which may be non-physical or complex to implement (Liu et al. 2012). In addition, numerical instability may appear when the VOF and level set methods are applied to simulate surface-tension-dominated flows in complex geometries (Scardovelli and Zaleski 1999).
Compared to common CFD methods, the LBM has many advantages (Chen and Doolen 1998). First, it is based on the molecular kinetic theory (Luo 1998). At the macroscopic scale it is able to recover N–S equations. Second, for single-phase flow simulations it usually involves an ideal-gas equation of state. Hence, it is not necessary to solve a Poisson equation for the pressure in the LBM. This saves significant computer central processing unit (CPU) time compared with common CFD methods. Third, it is easy to program and parallelize with much of the computational burden local to a node.
In the last decade LBM has become a numerically robust and efficient technique for simulating both single-phase and multiphase fluids (Chen and Doolen 1998; Guo and Shu 2013; He et al. 1999; Lee and Lin 2005; Rothman and Keller 1988; Shan and Chen 1993; Swift et al. 1995). Compared with conventional methods for multiphase flows, LBM usually automatically maintains sharp interfaces, and explicit interface tracking is not needed (Házi et al. 2002; Inamuro et al. 2004; Sankaranarayanan et al. 2003).
There
