Simulation of droplets impact on curved surfaces with lattice Boltzmann method

Simulation of droplets impact on curved surfaces with lattice Boltzmann method

International Journal of Heat and Mass Transfer 55 (2012) 6938–6943 Contents lists available at SciVerse ScienceDirect International Journal of Heat...

863KB Sizes 1 Downloads 61 Views

International Journal of Heat and Mass Transfer 55 (2012) 6938–6943

Contents lists available at SciVerse ScienceDirect

International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt

Simulation of droplets impact on curved surfaces with lattice Boltzmann method Shengqiang Shen, Feifei Bi, Yali Guo ⇑ Key Lab. of Ocean Energy Utilization and Energy Conservation of Ministry of Education, School of Energy and Power Engineering, Dalian University of Technology, Dalian 116024, China

a r t i c l e

i n f o

Article history: Received 5 April 2012 Received in revised form 3 July 2012 Accepted 4 July 2012 Available online 8 August 2012 Keywords: Lattice Boltzmann method Droplet impacting Curved surface

a b s t r a c t Several dimensionless parameters are studied to describe their effects on the deformation of a droplet after impact on a 2D round surface by using lattice Boltzmann implementation of pseudo-potential model. Four typical deformation process can be found: moving, spreading, nucleating and falling. In addition, in some special cases, part splashing is involved. It is observed that impact velocity of droplet has a significant influence on the droplet impacting dynamics. With the increasing of the impact velocity, different states have been found during the process. Moreover, when the surface is hydrophobic, splash occurs. Ó 2012 Elsevier Ltd. All rights reserved.

1. Introduction The droplets impacting on solid surface exist widely in the nature and modern industry, such as in metallurgy, chemical engineering, aerospace, materials processing and so on. To study this process can contribute to industrial processes, such as the liquid fuel combustion which spurted from nozzle spay on the wall of combustor, spray evaporation in cooling or heat pump system, and droplets dripping process of horizontal tube falling film evaporation in desalination plant. Therefore, research on droplets impact on solid surface attracts great interests of the researchers. Chandra et al. [1] studied the droplet collision with a heated solid surface, and found that the droplet would rebound if the temperature of solid surface was higher than Leidenfronst’s temperature. The droplet impact on different plat and inclined solid surfaces and the effects of various factors during the impacting process were studied by Pasandideh-Fard and Kang [2–4]. It was found that contact angle was an important factor in the impacting process. Sikalo et al. [5–7] obtained the dynamic contact angle of various droplets collision with plat and inclined solid surfaces from the experimental data of droplets impacting on different surfaces. In recent years, numerical simulation was draw more and more attention in solving complex problems because of the lack of experimental technique to deal with the complex measurement. Prashant et al. [8] carried out the experiments that a droplet impacted on a flat surface and simulated it using VOF model in 2D. Agreements were obtained from the computational results and

⇑ Corresponding author. Tel.: +86 411 84708774. E-mail address: [email protected] (Y. Guo). 0017-9310/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ijheatmasstransfer.2012.07.007

experimental results, which proved the feasibility to study the droplet impinging process with VOF. Considering the relationship between the dynamic angles and dynamic lines, Bussmann et al. [9] studied the flow of a droplet impact on plat and inclined surface with a 3D VOF method. Based on the VOF model, Li and coworkers [10] simulated a droplet collision with a heated solid surface, and compared the result with the experimental results. As a new method, lattice Boltzmann method (LBM) becomes a promising numerical algorithm for its parallelism, simple structure, simplicity in coding and the straightforward incorporation of microscopic interactions. Dupis [11] used a three-dimensional lattice Boltzmann model to investigate the spreading of mesoscopic droplets on homogeneous and heterogeneous surfaces. Yan [12] reported a new numerical scheme of lattice Boltzmann method which combined the existing model of Inamuro [13] and Briant [14–16] for calculating liquid droplet behavior on particle wetting surface typically for the liquid–gas system of a large density ratio. Elham [17] presented a new algorithm to simulate 2D dynamic wetting phenomena using the single phase lattice Boltzmann method. The simulation results were compared with analytical solutions and experiments and a good agreement was achieved. A LBM based on single-phase free surface tracking algorithm had been developed by Xing [18], which was used to simulate the droplet deformation under surface tension effect, the droplet spreading on solid surface and the capillary rise under various wetting properties. Li and coworkers [19] used the lattice Boltzmann pseudopotential method to simulate the droplets impacting on solid surfaces, and the simulation results were consistent with the previous theory prediction and experimental results. Until now, most of these studies always focused on plat and inclined solid surfaces. However, the curved solid surfaces are more common and

S. Shen et al. / International Journal of Heat and Mass Transfer 55 (2012) 6938–6943

complex, this paper will focus on the simulation of a droplet impact on a curved solid surface.

m0k ¼

s X q k¼1

m sk

!,

k k

s X q

k

k¼1

sk

Hence, given fik ; qk and mk can be calculated, and then so can kðeqÞ and fi .

2. Numerical method

6939

! ð7Þ

m0k ; meq k

2.1. Pseudo-potential model (S-C model) 2.2. Interaction potential The lattice Boltzmann approach solves the Navier–Stokes equation by the following evolution of partial distribution function fi on a regular, d-dimensional lattice formed of sites x. The D2Q9 lattice topology we use here has the following velocity vectors, see Fig. 1.

½e0 ; e1 ; e2 ; e3 ; e4 ; e5 ; e6 ; e7 ; e8    0 1 0 1 0 1 1 1 1 ¼ 0 0 1 0 1 1 1 1 1

F kr ðxÞ ¼ uk ðxÞ

ð2Þ

fik ðx; tÞ

where is the particle distribution function with particle velocity ei at the point x and time t, Dt is the time step of the simulation, s is the single relaxation time that controls the rate of approach to equilibrium. The index k is the total number of components, running from 1 to S. The density of the kth fluid is deP fined as qk ðx; tÞ ¼ 9i¼0 fik ðx; tÞ; the fluid velocity of the kth mk is deP9 fined as qk mk ¼ i¼0 ei fik ðx; tÞ. On the right-hand side of Eq. (2), we adopt the BGK [20] collision term and choose the equilibrium distribution function as follows:

1

Xki ðx; tÞ ¼  ½fik ðx; tÞ  fikðeqÞ ðx; tÞ

ð3Þ

s

" kðeqÞ

fi

ðx; tÞ ¼ qk W i 1 þ

# 2

2 eq er meq ðer meq k k Þ mk þ c2s 2c4s 2c2s

8 > < 4=9; i ¼ 0 W i ¼ 1=9; i ¼ 1; 2; 3; 4 > : 1=36; i ¼ 5; 6; 7; 8 0 qk meq k ¼ qk mk þ sk F k :

ð4Þ

ð6Þ

ð8Þ

where Gkk0 ðx; x0 Þ is a Green’s function and satisfies Gkk0 ðx; x0 Þ ¼ Gkk0 ðx0 ; xÞ:uk ðxÞ is a function of x through its dependence of qk, uk ðxÞ ¼ 1  expðqk Þ. If only homogeneous isotropic interaction between the nearest neighbors are considered Gkk0 ðx; x0 Þ can be reduced to the following symmetric matrix with constant element,

8 0 > < G; jx  x j ¼ 1 pffiffiffi 0 Gkk0 ðx; x Þ ¼ G=4; jx  x0 j ¼ 2 > : 0; otherwise

ð9Þ

Here G is the strength of the interaction potential between different components. By varying G we can control the surface tension between the fluids. At large enough G two mixed fluids separate and a well defined interface forms between the fluids of order or two or more lattice spacing depending on the strength of G. 2.3. Fluid–solid interaction The model, developed by Martys and Chen [22] to describe the interaction between a fluid and a wall, is implemented as follows:

F t ðxÞ ¼ uðxÞ ð5Þ

s XX Gkk0 ðx; x0 Þuðx0 Þðx0  xÞ x0 k¼1

ð1Þ

The lattice Boltzmann dynamic equation is described as the following equation:

fik ðx þ ei Dt; t þ DtÞ ¼ fik ðx; tÞ ¼ Xki ðx; tÞ ði ¼ 0; 1; 2 . . . 8Þ

To model surface tension in multicomponent fluids or fluids that are characterized by a nonideal gas equation of state, the interaction potential Fr can be incorporated as [21]:

X

Gti sðx þ ei Þei

8 > < Gt ; ci ¼ c pffiffiffi Gti ¼ Gt =4; ci ¼ 2c > : 0; ci ¼ 0

ð10Þ

ð11Þ

where m0k is a common velocity on top of which, an extra component-specific velocity due to interparticle interaction is added for each component. Fk is the total interparticle force acting on kth component, including surface tension force Fr, interaction force between fluid and solid surface Ft, and gravitational force Fg. Thus in this model, interaction of particles is given through the equilibrium distribution. To conserve momentum at each collision in the absence of interaction forces m0k has to satisfy the equation:

where s = 0 or 1 for fluid or solid, respectively. Gti is the fluid–solid interaction potential parameter. By adjusting Gti (positive for nonweting fluid and negative for wetting fluid) for each fluid we control which fluid preferentially wets a surface. A contact angle is defined as the angle at which a liquid/vapor interface meets a solid surface. Typically, a fluid is regarded as hydrophilic for small contact angle and hydrophobic if the contact angle is large. We find that reasonably well defined contact angles may be obtained by adjusting Gti for each fluid such that one of the fluids wets a surface.

Fig. 1. Discrete velocity set of two-dimensional 9-velocity (D2Q9) model.

Fig. 2. The initial and boundary conditions in domain.

6940

S. Shen et al. / International Journal of Heat and Mass Transfer 55 (2012) 6938–6943

Fig. 3. Time evolution of droplet shape for We = 14 and Gti1 = 0.175 (t⁄ = tv/D).

3. Results and discussion

q0 = q1 = 1.0; G01 = G10 = G = 0.33, Gti0 = 0.018, Gti1 = 0.175, s0 = s1 = 1.0, g1 = 0.00252 and g0 = 0.00003. The Reynolds num-

3.1. Initial and boundary conditions

ber, Weber number and Bond number are 41.4, 14, and 72.16. In the above equations, subscripts 0 and 1 respectively indicate gas and liquid phase, and g is the gravitational constant. At the initial condition, the Bond number is a constant for the fixed gravitational force and surface tension, whose value tells that the gravitational force is dominant over the surface tension. From the simulation results we can find that the deforming process after the droplet impacting on the pipe surface with initial velocity can be divided into four stages: moving, spreading, nucleating at the bottom of the pipe and dripping, as shown in Fig. 3. The first moving stage is presented in (a) and (b), which just takes only a very short period after impacting. The upper portion of the droplet can hold circular under the action of inertia force and surface tension, and only a very small part of the bottom contacts with the surface. Then from the (c) and (d), it can be seen that the rims form and corona flows appear around the rims. At t⁄ = 6.52, where the dimensionless time t⁄ is defined as tv/D, the droplet attains the maximum radial extension. After that, the droplet spreads along the solid surface as the rims grow, and the corona flows turn into film flows and move with the rims. As the contact area increases, the droplet finally converges at the bottom at t⁄ = 26.08. More and more liquid converge down and gradually gather into the nucleus again at t⁄ = 45.65. This type of droplet is called nucleating droplet. The kinetic energy of the nucleating droplet overcomes the effect of surface tension and wall friction, and therefore, drips down as shown in (h) and (i).

The test physical problem selected in a 2D, as shown in Fig. 2, a square domain with Lx  Ly = 300  300 is considered. The center of a pipe with radius R0 = 100 is at the center of the computational domain. A liquid droplet with diameter R1 = 46 is placed in contact with the top midline of pipe surface (y = 250), and collides with the surface with an impact velocity v at initial time (t = 0). A non-slip boundary condition is imposed to implement bounce back reflection on the particles arriving at a node adjacent to the surface. A dry pipe surface is used as the solid substrates, which is assumed to be homogeneous and smooth. In the S-C model, Gti can be used to relate with h, that is the contact angle. By varying Gti we can control the wettabliable of the pipe surface, and the periodic boundary condition is used on all the sides of the domain. Two dimensionless parameters that describe the dynamic impacting process, Weber number and Bond number are defined: We = qm2D/r, Bo = qgD2/r, where q, D, r and m are the density, diameter, coefficient of surface tension, and velocity, respectively. 3.2. Results and discussion 3.2.1. The process of a droplet impact on the surface The dynamic behavior of a droplet impinging on a pipe surface is investigated. The dimensional parameters of the case are:

S. Shen et al. / International Journal of Heat and Mass Transfer 55 (2012) 6938–6943

6941

Fig. 4. Time evolution of droplet shape for We = 3.5 and Gti1 = 0.175 (t⁄ = tv/D).

Fig. 5. Time evolution of droplet shape for We = 33 and Gti1 = 0.175 (t⁄ = tv/D).

3.2.2. The effects of impact velocity In order to investigate the effects of the impact velocity on the dynamic process, the droplet collisions on the pipe surface with various We numbers are simulated. Figs. 4 and 5 show the time evolution of droplet shape for We = 3.5 and We = 33, respectively, and in both the processes Bo = 72.16. From Fig. 4, it is seen that at t⁄ = 4.89(t = 3000), the rims are thicker than that at t⁄ = 9.78(t = 3000) for We = 14 in Fig. 3, which due to the comparable velocity of the peripheral rims and film flows close to the surface. The vertical diameter of the nucleating droplet is elongated downward and reaches its maximum at t⁄ = 29.35. Meanwhile, the velocity of the nucleating droplet decreases to zero because of the depletion of kinetic energy. The initial impact energy of

droplet and potential energy are dissipated in overcoming viscous flow and in producing new surface area. Then the nucleating droplet recoils upwards, and eventually, reaches an equilibrium state after oscillation for a quite time at t⁄ = 40.76. In the case of We = 33(Fig. 5), the thickness of the rims at t⁄ = 15.02 is smaller than that that at t⁄ = 9.78(t = 3000) for We = 14 in Fig. 3, that is, the amount of converging droplet increase as an increase of impact velocity, so the spreading velocity increases as an increase of impact velocity. Also, splashes do occur around the contact line of both side at t⁄ = 16.5. We have calculated the droplet impacts on a pipe surface for various We numbers and Gti numbers, and the droplet evolutions at every condition are obtained. Based on the evolution, we classify

6942

S. Shen et al. / International Journal of Heat and Mass Transfer 55 (2012) 6938–6943

Fig. 6. A droplet impact on different surfaces with different velocities.

Fig. 7. Time evolution of droplet shape for hydrophilic surface when We = 14 and Gti1 = 0.3 (t⁄ = tv/D).

Fig. 8. Time evolution of droplet shape for hydrophobic surface when We = 14 and Gti1 = 0.01 (t⁄ = tv/D).

S. Shen et al. / International Journal of Heat and Mass Transfer 55 (2012) 6938–6943

the results into three types of states, as shown in Fig. 6. It is seen that, at Gti1 = 0.175 the equilibrium states appear in the region of We number lower than Wemin, and the splashing occurs at high We number over Wemax. From the results, we know the two critical points are Wemin = 13 and Wemax = 23 when Gti1 = 0.175, and the dripping occurs between the two points. When Gti1equals other values, there are the corresponding Wemin and Wemax numbers as shown in Fig. 6. 3.2.3. The effects of characteristics of surface In order to investigate the effect of the wettability on the droplet evolution after impact, the droplet impacts on the hydrophilic wall (Gti1 = 0.3, 0.25) and on the hydrophobic wall (Gti1 = 0.1, 0.01) are simulated for various We numbers. Figs. 7 and 8 show the time evolution of droplet impact on a hydrophilic and a hydrophobic pipe surfaces for We = 14, Gti1 = 0.3 and We = 14, Gti1 = 0.01 respectively. It can be seen that the rims are thicker than for Gti1 = 0.175 in Fig. 3, and the oscillation behavior is not as obvious as that in Fig. 4, because the energy is dissipated greatly for hydrophilicity on the surface, and the left is too small to recoil. On the other hand, in the case of the hydrophobic wall Gti1 = 0.01 (Fig. 8), the droplet splashes partly before attaining the maximum radial extension. From Figs. 3, 7 and 8, the dynamic processes vary greatly due to the wettability of the surface. Fig. 6 represents the droplet evolutions at different Gti1 and We numbers combined conditions. For Gti1 = 0.3 and 0.25, the droplet only converges at the bottom at different Weber numbers, which due to stronger interaction between fluid and solid wall. At Gti1 = 0.2, the droplet drips only at We = 20, and converges or splashes at the value more or smaller than 20, respectively. When the value of Gti1 decreases to 0.175, 0.15 and 0.1, the interaction force becomes weak. The droplet converges at low Weber numbers, drips fall down when Weber number rises over a certain value, and splashes at high Weber numbers. The condition for Gti1 = 0.05 is similar with Gti1 = 0.2, the droplet drips only at We = 7.5. For hydrophobic wall Gti1 = 0.1 and 0.01, the interaction force is so weak that the droplet splashes even at smaller velocity. 4. Conclusions The lattice Boltzmann algorithm is applied to simulate the droplet impacting on pipe surface and the influences of physical factors on the droplet evolution are analyzed. From the results, the following conclusions can be obtained: (1) The process of droplets impacting on a pipe surface can be described into four stages, moving, spreading, nucleating, and dripping. And under some particular circumstances, the part-splashing phenomenon will occur. (2) Initial impacting velocity plays a significant role in impact dynamic. When the droplet impacts on the same surface, there are two critical velocities, Wemin and Wemax. The droplet drips down at the value between them; converge at the value lower than Wemin and splashes partly at the value larger than Wemax.

6943

(3) For the hydrophilic surfaces, the droplet only converge at the bottom at different velocities, with the decrease of wettability of the surfaces, the droplet will drip down and splash at larger velocity, and for the hydrophobic surfaces, the droplet will splash easily even at small velocity.

Acknowledgement This work is supported by Chinese National Natural Foundation Projects (Nos. 50976016 and 50906005). References [1] S. Chandra, C.T. Averdisian, On the collision of a droplet with a solid surface, Proc. Roy. Soc. Lond. A 432 (1991) 13–14. [2] M. Pasandideh-Fard, Y.M. Qiao, S. Chandra, J. Mostaghimi, Capillary effects during droplet impact on a solid surface, Phys. Fluids 8 (1996) 3. [3] T. Mao, D.C.S. Kuhn, H. Tran, Spread and rebound of liquid droplets upon impact on flat surfaces, AIChE J. 43 (1997) 2169–2179. [4] B.S. Kang, S.H. Lee, On the dynamic behavior of a liquid droplet impacting upon an inclined surface, Exp. Fluids 29 (2000) 380–387. [5] Š. Šikalo, C. Tropea, E.N. Ganic, Dynamic wetting angle of a spreading droplet, Exp. Therm. Fluid Sci. 29 (2005) 795–802. [6] Š. Šikalo, C. Tropea, E.N. Ganic, Impact of droplets onto inclined surfaces, J. Colloid Interf. Sci. 286 (2005) 661–669. [7] Š. Šikalo, M. Marengo, C. Tropea, E.N. Ganic, Analysis of impact of droplets on horizontal surfaces, Exp. Therm. Fluid Sci. 25 (2005) 503–510. [8] R. Prashant, Gunjal, V. Vivek, Ranade, V. Raghunath, Chaudhari, Dynamics of drop impact on solid surface: experiments and VOF simulations, AIChE J. 51 (2005) 59–78. [9] M. Bussmann, S. Afkhami, Drop impact simulation with a velocity-dependent contact angle, Chem. Eng. Sci. 62 (2007) 7214–7224. [10] S.Q. Shen, Y. Li, Numerical simulation of droplet impacting on isothermal flat solid surface, J. Eng. Thermophys. 30 (2009) 2116–2118. [11] A. Dupuis, J.M. Yeomans, Lattice Boltzmann modeling of droplets on chemically heterogeneous surfaces, Future Gener. Comp. Sy. 20 (2004) 993– 1001. [12] Y.Y. Yan, Y.Q. Zu, A lattice Boltzmann method for incompressible two-phase flows on partial wetting surface with large density ratio, J. Comput. Phys. 227 (2007) 763–775. [13] T. Inamuro, T. Ogata, S. Tajima, N. Konishi, A lattice Boltzmann method for incompressible two-phase flows with large density differences, J. Comput. Phys. 198 (2004) 628–644. [14] A.J. Briant, P. Papatzacos, J.M. Yeomans, Lattice Boltzmann simulations of contact line motion in a liquid–gas system, Philos. Trans. Roy. Soc. Lond. A 360 (2002) 485–495. [15] A.J. Briant, A.J. Wagner, J.M. Yeomans, Lattice Boltzmann simulations of contact line motion: I. Liquid–gas systems, Phys. Rev. E69 (2004) 031602. [16] A.J. Briant, J.M. Yeomans, Lattice Boltzmann simulations of contact line motion: II. Binary fluids, Phys. Rev. E 69 (2004) 031603. [17] A. Elham, K. Carolin, Lattice Boltzmann method for dynamic wetting problems, J. Colloid Interf. Sci. 335 (2009) 84–93. [18] Xiuqing Xing, David Lee Butter, Chun Yang, A lattice Boltzmann based singlephase method for modelling surface tension and wetting, Comput. Mater. Sci. 39 (2007) 282–290. [19] S.L. Quan, S. Li, A simulation of impact of droplets on solid surface by using the lattice Boltzmann method, Chinese J. Comput. Mech. 26 (2009) 627–632. [20] Y.H. Qian, D. D’Humieres, P. Lallemand, Lattice BGK model for Navier–Stokes equation, Europhys. Lett. 17 (1992) 479–484. [21] X. Shan, H. Chen, Simulation of non-ideal gases and liquid–gas phase transitions by the lattice Boltzmann equation, Phys. Rev. E 49 (1994) 2941– 2948. [22] N.S. Matrys, H. Chen, Simulation of multicomponent fluids in complex threedimensional geometries by the lattice Boltzmann method, Phys. Rev. E 53 (1996) 743–750.