Modeling droplet coalescence kinetics in microfluidic devices using population balances

Chemical Engineering Science 201 (2019) 475–483

Contents lists available at ScienceDirect

Chemical Engineering Science journal homepage: www.elsevier.com/locate/ces

Modeling droplet coalescence kinetics in microfluidic devices using population balances Yhan O’Neil Williams ⇑, Nelmary Roas-Escalona, Gieberth Rodríguez-Lopez, Andrea Villa-Torrealba, Jhoan Toro-Mendoza ⇑ Laboratorio de Dispersiones e Interfases, Centro de Estudios Interdisciplinarios de la Física (CEIF), Instituto Venezolano de Investigaciones Científicas (IVIC), Caracas 1020A, Venezuela

h i g h l i g h t s

g r a p h i c a l a b s t r a c t


Accurate coalescence kinetics

modeling of microfluidic experiments.
Testing of a wide range of physicochemical conditions using population balances.
Determination of the effects of interfacial mobility and hydrodynamic interactions.
Potential application of the methodology as a rapid diagnostic tool.

a r t i c l e

i n f o

Article history: Received 20 October 2018 Received in revised form 6 January 2019 Accepted 16 February 2019 Available online 7 March 2019 Keywords: O/W emulsion Coalescence kinetics Coalescence time Population balance model

a b s t r a c t The coalescence kinetics of oil-in-water emulsions in a wide range of properties and flow microfluidic conditions is quantified. The conditions were chosen in order to mimic situations found in industrial processes involving liquid-liquid dispersions. With that aim, a numerical scheme based on population balance equations is proposed, applied, and validated by comparison to microfluidic experiments reported in the literature. A coalescence efficiency model accounting for colloidal and hydrodynamic interactions, and interface mobility is incorporated using the Smoluchowski collision kernel. The latter assures the accurate estimation of the droplet size evolution which governs the interfacial area and rate of mass transfer. Besides, the combined effect of interfacial tension and oil viscosity on the coalescence kinetics is properly quantified with one single fitting parameter. From the kinetics, the estimated coalescence time increases as the shear rate and volume fraction of the dispersed phase diminish. The close agreement of our results with the experimental findings substantiates the accuracy and wider application of the methodology here described as a diagnostic tool beneficial to industrial process design and control. Ó 2019 Published by Elsevier Ltd.

1. Introduction At the large scale industrial applications, liquid-liquid dispersions called emulsions may be subjected to flows, for example, during mixing or through pipes. In such a flow, droplets are ⇑ Corresponding authors. E-mail addresses: [email protected] (Y.O’Neil Williams), jhoantoro@ gmail.com (J. Toro-Mendoza). https://doi.org/10.1016/j.ces.2019.02.040 0009-2509/Ó 2019 Published by Elsevier Ltd.

sheared and collide, causing coalescence (Tsouris and Tavlarides, 1994). The efficiency of the coalescence process is determined by the droplet collision rate, fluid flow type/regime, and colloidal and hydrodynamic interactions (Chesters, 1991). In fact, both the collision rate and the inter-droplet interactions are explicit functions of the droplet size which accordingly changes as a result of coalescence (Klink et al., 2011). Since transport properties, such as the effective viscosity of an emulsion, and mass transfer rates between the fluids depend on the droplet size distribution, it is

476

Y.O’Neil Williams et al. / Chemical Engineering Science 201 (2019) 475–483

useful to predict how flow conditions change droplet sizes (Shardt et al., 2013). Droplet interactions are also critical in microfluidic devices where individual droplets can be formed and manipulated at will (Teh et al., 2008). In particular, monodisperse droplets are generated by focusing one liquid current into another using the adequate channel geometry, for example a T-junction or other flow-focusing devices (Gu et al., 2011; Fu et al., 2012). Further downstream, the use of the appropriate channel geometries drives droplets together leading to coalescence (Xu et al., 2011). The study of this process is advantageous compared to conventional bulk measurement methods since precise control of the factors affecting the process as well as the direct visualization of the coalescence events is possible (Bremond and Bibette, 2012). These events in microfluidics are usually symmetric between two isolated droplets which undergo head on collisions under a constant interaction force. Despite all the advantages that microfluidics possess as a diagnostic screening platform, only a few reports of coalescence dynamics in emulsions with a high dispersed phase fraction exist (Krebs et al., 2013, 2012a,b; Muijlwijk et al., 2017). Typically, industrial coalescers process concentrated emulsions with dispersed phase as high as 60–80% (Abdurahman et al., 2013). Consequently, the effect of multiple interactions is strong when compared to more dilute systems. The reduced use of microfluidics as diagnostic tool for some industrial process conditions may be due to the fact that these devices are frequently used to study the specific behavior of single droplets or a few of them. Nevertheless, microfluidics are suitable to observe and quantify the time evolution of liquid dispersions under different physicochemical conditions as well. Therefore, any prediction tool reproducing those conditions is, in principle, feasibly useful too. For the mentioned challenges, we propose the use of a population balance framework to describe the coalescence kinetics as a complement of microfluidic experiments. For that, we use microfluidic experiments as a benchmark to compare the results obtained from our modeling since they are widely accepted as an effective diagnostic tool. In fact, both methods combined are quite useful as shown by Lazzari et al. where a microfluidic platform combined with population balance equation PBE modeling revealed relevant insights into the nanocrystal formation process accounting for surface ligands (Lazzari et al., 2018). Our proposed methodology offers the screening of process conditions which are crucial for the optimum design of chemical processing systems without the purchase of pumps, optical microscopy and other costly materials which form a part of the integrated microfluidic setup. In this study, we present the results of population balance modeling to predict the collective droplet coalescence kinetics of dense oil in water emulsions without surfactant. We incorporate a drainage based coalescence efficiency model which accounts for both colloidal forces and interface mobility. Besides, an auxiliary scheme for the hydrodynamic interactions due to collisions between different sized droplets is proposed to rationalize the experimental observations. Then, we apply the methodology to study the influence of oil type and viscosity on the coalescence kinetics. Ten systems were studied: hexadecane, three mineral oils and five silicone oils represented the dispersed phase, covering a range of 3.5–106 mPa s, and, pure water was the continuous phase and, finally, mineral oil-inbrine emulsions. The droplets were subjected to shear rates that are large enough to neglect Brownian motion but sufficiently small so that the droplet coalescence is the dominant mechanism. In this sense, all the selected flow conditions assured that the capillary number Ca of the system was  1. From the collective kinetic prediction, we estimate coalescence times for the experimental flow conditions.

2. Theoretical description 2.1. Population balance The PBE is an integro-differential equation used to model the evolution of the particle size distribution through the conservation of the number of dispersed particles. The particle size distribution function is characterized by a number density function NDF that is a function of variables in the particle state space, which consists of external and internal coordinates (Ramkrishna, 2000)

@nðn; X; t Þ þ r  ðhUjninÞ ¼ Sn ðn; X; tÞ: @t

ð1Þ

where (1) nðn; X; t Þ is the dispersed phase number density function, with X (space coordinates) and n(dispersed phase size) being the external and internal coordinates, respectively. hUjni is the dispersed phase velocity field conditional to the particle size, and Sn includes all the source terms for the PBE. Sn includes the birth B and death functions D due to droplet breakup and coalescence and is given by

Sn ¼ BB þ BC  DB  DC :

ð2Þ

In order to determine the evolution of droplet sizes, instead of solving for the entire size distribution the corresponding moments of the NDF can be solved lowering the computational cost (Marchisio et al., 2003). Homogenizers have distinct zones where local shear forces can change dramatically. On the contrary, in microfluidics, velocity fluctuations are minor due to very low shear rates. Thus, we treated the flow conditions of the matching experiments as a homogeneous system to avoid the complexities associated with including spatial variations. For all that, in this work the moments of the NDF are solved using the length based expressions from the analytical solution for a spatially-homogeneous, transient PBE (McCoy and Madras, 2003). The length-based formulation used here follow the corresponding expressions by Bhutani et al. (2016). The steady-state solution of the dimensionless particle number density according to McCoy and Madras is  12 ~ n ð1Þ ¼ 2e n Sm3 =b ð1=m0 ð0ÞÞ where m0 ð0Þ is the initial value of the moment 0 of the NDF, e S breakage frequency and b represents ~ (1) can assume arbitrary the aggregation frequency. Although n values, conditions for this work represented a system where coa~ (1)< 1. n ~ (1)> 1 would replescence is dominant, in which case n resent conditions dominated by the breakup of droplets. To ~ (1) values for the experimental conestablish the corresponding n ditions modeled in this article e S was not taken as 0, but rather, for simplicity we used e S ¼ 1 following the work form Silva et al. (2008). In their work a breakage dominant case was to be studied. For that they established the aggregation constant as b ¼ 1 and ~ (1) was given arbitrary values to represent the dominant cases n of droplet breakup and aggregation, respectively. The equations used to calculate the aggregation frequency are presented in the following section. In order to follow the droplet size evolution the moments are evaluated over a given time lapse by a solution to the first few moments wherein only the dispersed phase diameter n is needed as an input from the population balance equation. As few as four moments is considered sufficient for an accurate prediction th

(Marchisio et al., 2003; Marchisio and Fox, 2005). The k moment of the number density function is defined as

Z mk ¼

nk nðnÞdn:

ð3Þ

To calculate n, the Sauter Mean Diameter SMD is used and can be defined in terms of the moments as

477

Y.O’Neil Williams et al. / Chemical Engineering Science 201 (2019) 475–483

m3 : m2

d32 ¼

ð4Þ

The volume fraction of the dispersed phase /d is related to the third moment of the number density distribution function as

  6 : m3 ¼ /

ð5Þ

p

The corresponding equations to estimate the first four moments of the NDF is presented in the Appendix from the formulation of the analytical problem of McCoy and Madras (2003), Bhutani et al. (2016). 2.2. Collision frequency We quantify the rate at which the number of droplets changes due to droplet-droplet collision leading to coalescence via Smoluchowski scheme. For that, droplet centers are assumed to move along streamlines and coalescence occurs when the distance between the droplets is less than the sum of their radii (Klink et al., 2011). Neglecting hydrodynamic interactions and other inter-droplet forces and assuming that all droplet collisions result in coalescence, the collision kernel between two droplets with volumes v i and v j in simple shear flow is

Asmol



v i; v j



¼

3 c_  1=3 v i þ v 1=3 : j p

ð6Þ

where c_ is the shear rate. While the widely known Smoluchowski kernel assures conservation of mass, it weights evenly the different collisions between disparate sized droplets. It implies that the information of the individual coalescence events is lost but facilitates the accounting of collisions leading to coalescence. To accurately model the droplet population dynamics of an emulsion subjected to a shear flow field only a fraction of the collisions (denoted as collision efficiency e in the present study) given by Eq. (6) result in successful collisions leading to coalescence. e should include contributions from hydrodynamic effects, colloidal forces and flow properties (Chesters, 1991). In particular, from hydrodynamic analysis van de Ven concluded that a pair of spheres would not touch in the absence of any attractive forces (Van de Ven, 1989). Furthermore, some studies have revealed that the collision efficiency is dependent on the size ratio of the two droplets or aggregates that collide. These studies concluded that the efficiency of equal sized droplet collisions can be significantly greater than the efficiency for unequal sized droplet collisions (Huang and Hellums, 1993a,b; Tandon and Diamond, 1997). For all that, we present a methodology to predict the collision efficiency e taking into account the interplay of shear flow hydrodynamics, colloidal forces and the droplet size ratio. To include effects due to these factors on the population dynamics, a modified collision kernel is defined as

  b V i ; V j ¼ eAsmol ;

ð7Þ

  where b V i ; V j is the coalescence kernel. e is made of a contribution due to hydrodynamic interactions and coalescence probability due to colloidal and flow properties and the droplet size ratio, a1 ¼ e. In this study a is an estimate of the probability of coalescence during droplet collisions under the action of given local flow conditions based on the film drainage model where a is established as a function of van der Waals attraction, the shear rate, interfacial tension and the viscosity ratio. The starting point is to relate a to the ratio of the mean coalescence time tcoal to the mean contact time t cont by a ¼ exp ðt coal =t cont Þ (Liao and Lucas, 2010; Sajjadi et al., 2013). The average contact time between two droplets during a collision is given by tcont ¼ p=2c_ . Instead, the coalescence time is determined by the interface mobility and the shape of the

droplets in the contact area. For pure liquid-liquid systems the film flow is controlled by interface motion and by disperse phase viscosity (Podgórska, 2005). The deformation of droplets depends on their size: we can expect that small droplets remain spherical, due to high capillary pressure, while larger ones flatten when approaching each other (Davis et al., 1989). In that regard, we propose the use of a drainage model that elucidates deformable, partially-mobile drops in viscous simple shear. Hence, we use

2

!1 3



 8pcD2eq 3 3 l 5; a ¼ exp 4C 3 d Ca2 lc A

ð8Þ

where A is the Hamaker constant, ld the dispersed phase viscosity and C 3 a constant taking values 6 1(Chesters, 1991; Ivanov et al., 1999). The capillary number is expressed as Ca ¼ lc c_ Deq =c, with lc the continuous phase viscosity, c_ the shear rate, Deq the equivalent droplet diameter and c the interfacial tension (Christopher   et al., 2009). Deq ¼ 2di dj = di þ dj , where di and dj represent the diameters of the colliding droplet pair. In this work, C 3 was taken as

5:5  102

and

A ¼ 4:40  1020

J,

6:29  1020

J

and

21

4:90  10 J for silicone oils, mineral oils and hexadecane, respectively (Israelachvili, 2011). Due to the relatively high dispersed phase volume and the reduced dimensions of the channel, we assumed that the droplets would follow an Arrhenius-type dependence on the droplet disparity. The latter is established through the size ratio . ¼ di =dj . The size dependence model 1 is expressed as





1 ¼ 237.3  588.2 þ 487.  135 exp 

dj di

ð9Þ

A similar functional form of hydrodynamic interactions was established by Tandon and Diamond (1997) where platelets and their aggregates were modeled as unequal-sized spheres. In their work the hydrodynamic efficiency function had the form

h ¼ X 1 ð.ÞðC A ÞY

1 ð.Þ

, where X 1 ð.Þ=0:2258.3  0:1579.2 þ 0:7356.þ 0:1077 and Y 1 ð.Þ=1:1866.3 þ 2:7509.2  2:0556. þ 0:6411, where C A is a non-dimensional attractive force parameter. 3. Results and discussion

The experimental results of hexadecane/water dispersions by Krebs et al. (2012a) were used for comparison and model validation. These experiments were performed in a microfluidic device made of glass with a rectangular cross-section. All channels have a uniform depth of 45 lm. The width of the channels is 100 lm, except for the chamber in the center of the chip with a width of 500 lm and a length of 5 mm. The physical properties of the oils used in these experiments are presented in Table 1. The values are taken from (Krebs et al., 2013, 2012a).

Table 1 Oil viscosities

l and oil/water interfacial tensions c at 293 K.

Oil Hexadecane Mineral 1 Mineral 2 Mineral 3 Silicone 1 PMX-200 Silicone 2 PMX-200 Silicone 3 PMX-200 Silicone 4 PMX-200 Silicone 5 PMX-200

5 10 20 50 100

l (mPa s)

c (mN/m)

3:5  0:02 7:7  0:03 19:6  0:04 70:6  0:03 5:9  0:04 10:7  0:02 21:1  0:05 52:1  0:03 106:4  0:02

50:5  0:3 22:9  0:8 23:8  0:5 25:0  1:1 45:1  1:2 47:9  0:6 49:9  0:6 49:2  0:4 47:4  0:7

478

Y.O’Neil Williams et al. / Chemical Engineering Science 201 (2019) 475–483

3.1. Model prediction for evolution of droplet size The mean droplet sizes in the experiments are all established as hydraulic diameters dh , however, we established them as a mean Sauter diameter SMD. In order to establish the solution to the number density function the mean droplet sizes and the dispersed phase fraction / are used. The mean hydraulic diameter of the droplets is given by dh ¼ 4A=P, where A and P are the cross sectional area and perimeter of the droplet, respectively. In microfluidics it is a common practice when the droplets are more like disks than spheres to express the average droplet sizes as hydraulic diameters. This occurs when the droplets are of sizes comparable to the channel height. Average shear rates were determined from the mean flow velocities used in the experiments and the uniform depth of the channel. Flow rates ranged from 0.26 to 8 cm s1 in the experiments presented, and the uniform depth of the channel was 45 lm. Hence, using the relation velocity/channel height the shear rate used in the predictions spanned 57.78–1778.78 s1. In order to predict the collision frequency (b) Eqs. (6)–(9) are used for a given shear rate and dispersed phase fraction. We calculate droplet volumes from the evolution of the SMD over the experimental time of 2 s. The coalescence kinetics were calculated from the events between colliding droplets as follows. We express the mean relative volume as the ratio between the evolving mean volume V f of droplets to the mean initial droplet volume V 1 as u ¼ V f =V 1 . u is related to the mean number of events N c ¼ u  1. N c indicates the number of coalescence events that led to the formation of the larger droplet. The coalescence rate is defined from a combination of two droplets V i and V j as rij ¼ N c =Ac t exp , where Ac is the monitored area of the collision chamber, and t exp the experimental time. Droplet volumes V i and V j refer to droplets that have not undergone coalescence ði; j ¼ 1Þ or are formed through coalescence i; j > 1 from i  1 coalescence events. The initial SMD of the droplets entering the chamber was d1 ¼ 65  5:5 lm and /  0; 35. The SMD of a droplet with volume V 4 , which was formed through coalescence of droplets with d1 ¼ 65 lm is only 76 lm. Moreover, 10 combinations of initial droplet pair sizes are possible when considering the error (65  5:5 lm). Therefore, colliding droplets in the range V 1  V 4 spanned a range of 60–76 lm. Since the moments are placed in function of the SMD, all relative volumes are calculated using D1 ¼ 65 lm. The relatively small range of SMD of the droplets is probably due to the fact that colliding droplet pairs are expressed as an effective diameter. From the experimental results it is observed that the droplet coalescence rate decreases as droplets formed through coalescence and those who have not undergone coalescence events collide. A larger scatter is obtained with the increase in disparity. Therefore, to model the kinetics observed we test the developed droplet size dependence model presented in Eq. (9). r11 rates are between equal sized droplets. Nonetheless, the average droplet size contains a standard deviation. For instance in Fig. 1 the initial mean droplet size according to Krebs et al. was 65  5:5 lm (Krebs et al., 2012a). Therefore the d1 spanned a range of 60–71 lm. Taking the mean droplet size 65 lm we calculate the volume of a droplet and use Eq. (6) for a given shear rate. Using the droplet and flow properties we predict the coalescence efficiency as predicted by Eq. (8). For that, an equivalent diameter is to be established, but for same sized droplets the Deq ¼ d1 . The same procedure is repeated for the limits represented by the standard deviation. Now to determine the moments of the distribution, only the average is used.Therefore the initial droplet volume is always taken as the average regardless of the droplet sizes being evaluated. For Eq. (9) di = dj .

r 14 rates are determined similarly to r11 . From the experiments Krebs et al. found that for the system presented in Fig. 1, a droplet with an initial size of 65 lm achieves a size of 76 lm undergoing 4 coalescence events the d4 (Krebs et al., 2012a). With that in mind we established equivalent diameters from that value and the range established for the d1 . The corresponding volumes of the colliding droplets d1 and d4 again were used in Eq. (6). Droplet and flow characteristics were again evaluated in Eq. (8). In this case the equivalent diameter is established using the range of sizes for d1 and the size of d4 as established in the experiment. The hydrodynamic disparity function Eq. (9) is used in all the cases where di ¼ d1 and dj ¼ d4 . Three collision frequencies are calculated, when d1 = 65 lm, d1 = 60 lm and d1 = 71 lm, respectively. From the ranges of droplet sizes presented for the initial droplet size d1 and the final droplet size represented by d4 , we established through a mean that the droplet sizes d2 ranged from 65–73 lm and d3 71–74 lm. The small range of droplet sizes are due to the fact that the newly formed droplets due to coalescence were expressed as an equivalent diameter. The collision frequencies of the rates r12 and r 13 are determined using the range of d1 and the corresponding range of d2 and d3 , respectively. Again droplet volumes are calculated and used to calculate the collision frequency using Eq. (6). Coalescence efficiency is estimated using Eq. (8) and the hydrodynamic interaction via Eq. (9). di ¼ d1 in all cases. Once the corresponding collision frequencies are found for each experimental shear rate and dispersed phase volume fraction, ~ (1) is determined using A.11. Python scripts are used to evaluate n the corresponding equations for the moments of the NDF A.10– A.15. The relation between the moments is used to find the evolution of the average droplet size over the 2 s lapse as expressed in Eq. (4). In Fig. 1 the span of possible values of the coalescence kinetics predicted in this work are represented by the lines. The solid lines are the kinetics when the initial droplet size is taken as the average, while the dashed lines represent the kinetics from the limits of the standard deviation of the initial droplet size. The symbols represent the experimental results. It can be seen that the coalescence kinetics are adequately predicted by the solution to the number density function. Also, it is observed that r ij increases with v y at a constant / owing to the increase in the momentum of the colliding droplets. Furthermore, initially at low velocities the coalescence rates increase almost linearly, suggesting that almost all collisions lead to coalescence. However, at greater velocities a far less steep increase is observed owing to a decrease in drainage efficiency a. At larger flow velocities droplets have less time to interact, thus, leading to a less efficient film drainage, and rupture. Besides, coalescence efficiency of larger droplets can also be reduced for an increased probability to deform in the collision contact area. As is with the experimental results, rij decreases with the disparity between droplet sizes (d1  d4 Þ possibly due to droplet concentrations. That is, a droplet born due to binary collision further collides with its parent drops, which are larger in number. The observed scatter points are well described by the droplet ranges established as d2 and d3 , but the limit described by the r 14 deviates slightly from the experimental observations. This scatter can be explained by the occurrence of multiple simultaneous collisions involving more than two droplets not predicted by the binary collision model. 3.2. Effect of volume fraction on coalescence kinetics From the experimental data Bins described in Ref. (Krebs et al., 2012a), we can estimate the effect of / for a pair of equally sized colliding droplets. The mean coalescence time was taken as the time between each event in Fig. 2 as a function of v y . Since the

479

Y.O’Neil Williams et al. / Chemical Engineering Science 201 (2019) 475–483

 Fig. 1. Coalescence rates r ij of the hexadecane/water system for i ¼ ð1 6 j 6 4Þ as a function of the mean dispersed phase flux h/i v y . droplet velocity in the collision chamber.

vy

is the mean y component of the

Fig. 2. (a) Coalescence times for a 66:5 lm colliding droplet pair as a function of flow velocity for different mean dispersed phase fraction. (b) Coalescence times for colliding droplet pairs as a function of flow velocity at different mean dispersed phase fractions.

solution to the PBE used here requires a constant aggregation rate, for any given v y and /, the duration of the events are equal. It can be seen from Fig. 2a-b that the average coalescence time decreases with both flow velocity and /. Fig. 2a. presents the results for a fixed size and variable / . Instead, in Fig. 2b. both parameters were changed. It is evident that a minimal difference is observed in both cases. The latter indicates that the number of neighboring droplets governs the mean time between each event. In other words, the presence of more neighboring droplets may change the timescale of coalescing droplets as their presence may accelerate/decelerate the droplet approach velocity, a behavior similar to the wake effect caused on trailing droplets due to the drag force (Leal, 2004; Yang et al., 2001). In hindsight, we observe that at low flow velocities the

mean time between each event is within the experimental error predicted by the trajectory analysis in Ref. (Krebs et al., 2012a). However, as flow velocities increase, the difference between the experimental values and those presented here becomes important. Nevertheless, the timescale of the events found are entirely possible and fall in the same order of magnitude of the experiments, but are better than more sophisticated theoretic models, for example those presented by Prince and Blanch (1990a). Interestingly, from the experimental results reported from (Krebs et al., 2012a), the coalescence time did not present a dependence on /, but rather, slightly on the size of the colliding droplets. In retrospect, the coalescence time predicted here estimated as the mean time between each event is based on the overall collective kinetics, while the

480

Y.O’Neil Williams et al. / Chemical Engineering Science 201 (2019) 475–483

experimental results were obtained using the trajectory analysis of thousands of colliding droplet pairs. 3.3. Short time behavior Both long and short time behavior in the evolution of the droplet size distribution have been analyzed for several initial systems undergoing simple shear. At short times, for the case when the initial drop size distribution is log-normal, Eq. (1). reduces to

 ln

Nt N0



  

4/c_ 1 þ 3 exp ð2 lnð2ÞrÞ t; ¼ 4 p

ð10Þ

where r is the standard deviation of the initial average droplet size, and N 0 and N 1 represent the initial and total number concentration, respectively. To evaluate the coalescence kinetics in the short time regime, we substitute the total number concentration of droplets for the moment 0 by Nt ¼ m0 ðt Þ. In Fig. 3 the evolution of the total number of droplets over 2 s is plotted. Excluding the brief period at initial times when t < 0:25 s, we notice that the number concentration of all coalescence rates evolves exponentially. We also observe that the initial change in N t is reduced as the kinetics are governed by more disparate collisions at t < 0:25 s. When the population of droplets moves towards larger diameters, the size of the colliding droplet deviates from the initial average leading to a reduction in the rate of coalescence at a given flow velocity owing to a reduced probability of collision among large droplets in comparison to small droplets. The steep initial linear slope observed can be estimated from ð4/c_ =pÞð1 þ 3 expð2 lnð2ÞrÞ=4, taking the standard deviation from the average initial size 65  5 lm. The gradient of the slope represents a first order approximation of coalescence kinetics over a small time period.

3.4. Predicting the effects of viscosity ratio on the coalescence kinetics An additional validation to our methodology was made against mineral oil-in-brine experiments reported by Krebs et al. (2012c) and the results are shown in Fig. 4. In this case, the initial growth observed when droplets undergo gentle collisions is followed by a shift after achieving a maximum. This transition indicates a less favorable condition for coalescence, despite the fact that the increase in shear rate should promote an increase in the number of collisions. Our model predicts a continued reduction in coalescence rate even after a flux rate of 1.6. Interestingly, at that point the experimental results convey an opposite tendency, which was not widely discussed by Krebs et al. (2012c). The current comparison highlights this deviation from the expected behavior around 1:6 cm s1. It can be explained in terms of an important modification of the interfacial mobility due to a combined effect of an increase in the flow velocity and the presence of salt. On one hand, both lc and c values were accordingly included in the calculations, although these were not enough explain the observed tendency. On the other hand, the increase in the flow velocity is insufficient to generate the required local shear stresses to enhance the coalescence rate. To adequately predict this transition it is necessary to derive a function for interfacial mobility which would depend on specific salt interface interactions and flow conditions. Moreover, in surfactant-free systems, the addition of salt is known to decrease the rate of film drainage due to the formation of interfacial tension gradients (Prince and Blanch, 1990b). In Fig. 5 the kinetics as a function of Ca and k ¼ ld =lc exhibit a reduction in the coalescence rate as k increases for both mineral and silicone oil. However, this difference is enhanced when Ca > 1  103 . At very low Ca the hydrodynamic forces responsible for the coliisions are small, leading to a more efficient coalescence process. However, for larger k, successful coalescence conditions

Fig. 3. Short time evolution of the colliding droplet concentration as a function of flow velocity at h/i ¼ 0:35.

481

Y.O’Neil Williams et al. / Chemical Engineering Science 201 (2019) 475–483

Based on the trends observed in Fig. 5 we extended the shear rates used to predict the capillary number at which no coalescence events would be observed, critical capillary number Cac . Besides, we aimed to validate the value of the constant C 3 and the hydrodynamic disparity function used as to ensure the underlying physics of the results. As can be observed in Fig. 6, the trends shown in Fig. 5 are maintained, wherein the Cac is smaller for the more viscous oils. The previous, again suggests that smaller k favors droplet coalescence. Interestingly, Hu et al. (2000) observed that coalescence will occur for Ca below a critical value Cac , and found that Cac for coalescence in the absence of additives scales as

r11[mm-2s-1]

102

Exp PBE Major Error Region

101

0

1

2

<φ>[cm s-1] Fig. 4. Coalescence rates r ij of mineral oil in brine system for i ¼ ð1Þ as a function of  the mean dispersed phase flux h/i v y . v y is the mean y component of the droplet velocity in the collision chamber. The arrow highlight the deviation from the expected behavior.

requires that Ca to be still smaller. In accordance, we observe that depending on the value of /; r ij may increase monotonously or pass through a maximum in a given range of c_ . The presence of the maximum can be considered as an indication of a transition between two limiting behaviors: from partial to zero interfacial mobility. Since the interface mobility acts as an intermediary between the external and internal flow, consequently, coalescence is directly affected too. Hence, for small k, the droplet interface is more mobile, leading to lower coalescence times. Inversely, for large k, droplets behave like nearly rigid spheres, and for an increase in c_ , interface mobility would hinder drainage of the separating film. The occurrence of the maximum at lower Ca as k increases, confirms that coalescence is favored by a low ld . In addition, being Ca a measure of droplet deformability, the results of Fig. 5 suggest that deformation will slow down film drainage, thus, reducing the coalescence rate. However, an increase in ld would imply a less deformable droplet which counteracts the effect of the increase in Ca.

Cac / ðkÞ0:410:06 ð2aÞ0:820:03 , where a is the droplet radius. From Fig. 6 we observe that the results presented here are in good qualitative agreement with the values predicted from the theoretic scaling relations even when the range of the properties studied in this work is greater than those conducted by their experiments. The latter, substantiates the values of the constant C 3 used in Eq. (8) and the size dependence model presented in Eq. (9). In summary, the coalescence kinetics were measured under laminar flow conditions, as encountered in some relevant industrial applications, for example oil-water separators. Moreover, we advocate that both the proposed model and methodology is capable of studying the phase separation process and coalescence times of oils with different physical properties. 4. Conclusions
This study presents the modeling of demulsification kinetics of oil in water emulsions in a microfluidic device by a droplet number density solution to the population balance equation. The effect of droplet size ratio, dispersed phase volume fraction, shear rate, and viscosity ratio on the coalescence kinetics of oilwater emulsions were analyzed via 10 emulsion systems with different physical properties.
The results prove that the kinetics for disparate collisions can be modeled using an Arrhenius-type size dependence model. We observe that the coalescence rate increases with shear rate, but, after an initial steep aggregation rate, the coalescence rate is reduced at large shear rates. Likewise, we observe that as interface mobility of the droplets decreases the coalescence rate is reduced as predicted by theory.
It is shown that this analysis can also be applied to estimate the average coalescence times of colliding droplets in an emulsion. These estimated coalescence times become smaller than the

r11[mm-2 s-1]

(a)

(b)

102

102

101

101

λ 5.92 10.7 21.1 52.1 106.4

λ 7.71 19.6 70.6

100

0

1

100

2

3

4

5

0

Ca x10

1

2

3

3

Fig. 5. (a) Coalescence rates of colliding droplet pairs as a function of the capillary number for mineral oils for / ¼ 0:26 and d1 ¼ 64:8 lm. (b) Coalescence rates of colliding droplet pairs as a function of the capillary number for silicone oils for / ¼ 0:26 and d1 ¼ 66:9 lm.

482

Y.O’Neil Williams et al. / Chemical Engineering Science 201 (2019) 475–483

(a)

300

λ

r11[mm-2 s-1]

250

1.42 0.97 0.57

7.71 19.6 70.6

200

(b)

500

Cac x 102

λ Cac x 102 5.92 10.7 21.1 52.1 106.4

400 300

150

1.54 1.21 0.91 0.63 0.47

200

100 50

100

0

0 0

1

2

3

4

5

6

7

0

Ca x10

1

2

3

2

Fig. 6. (a) Coalescence rates of colliding droplet pairs as a function of the capillary number for mineral oils for / ¼ 0:26 and d1 ¼ 64:8 lm. (b) Coalescence rates of colliding droplet pairs as a function of the capillary number for silicone oils for / ¼ 0:26 and d1 ¼ 66:9 lm. Marked values in the graph are from the theoretic scaling relation.

experimentally measured as shear rate increases. In hindsight, it is observed that according to the model, droplet size has less of an effect on the coalescence time compared to the dispersed phase fraction, a behavior opposite to the reported from the experiments.
The sum of all the aspects discussed in this article make us uphold the proposed model and methodology as a viable tool for screening the coalescence kinetics in process conditions found in relevant industrial applications. Conflict of Interest

where m is the average number formed after breakup, aðnÞ the breakage frequency, e S the constant coefficient in power law breakage frequency, b the breakage daughter distribution function, and b the aggregation kernel. An analytical solution to the population balance Eq. (1) exists when conditions are spatially-homogeneous, which is given as

  nn ðn; tÞ ¼ k1 n2 exp k2 n3 :

ð:7Þ

k1 and k2 are given by

~ 2n ðtÞ k1 ¼ 3n

m0 ð0Þ2 ; m3

ð:8Þ

The authors declared that there is no conflict of interest.

~ n ðt Þ k2 ¼ n

Acknowledgments

m0 ð0Þ ; m3

ð:9Þ

where

Thanks are given to Ph.D. Yeni Sánchez for helpful discussions and criticisms. This research was financially supported by project IVIC-1013.

 3 ~n ð1Þtanh n ~ n ð1Þ m0 ð20Þbt 1þn ~ n ð1Þ4 ~ n ðt Þ ¼ n   5: n ~ n ð1Þ m0 ð20Þbt ~ n ð1Þ þ tanh n n

Appendix A

~ n ð1Þ is given in terms of the breakage and aggreThe constant n gation frequency as

The length based version of the problem developed by McCoy and Madras (2003) can be written as follows for an initial number density function given by

 m0 ð0Þ m0 ð0Þ ; nn ðn; 0Þ ¼ 3n2 exp n3 m3 m3 2



ð:1Þ

and breakage and aggregation kernels given by

m ¼ 2;

ð:2Þ

aðnÞ ¼ e Sn3 ;

ð:3Þ

3n

n31

~ n ð1Þ ¼ n

; n < n1 ;

ð:5Þ

b ¼ aggregation frequency;

ð:6Þ

!12 

 1 : m0 ð0Þ

~n ðtÞ; m0 ðtÞ ¼ m0 ð0Þn

ð:4Þ

bðnjn1 Þ ¼ 0; nin1 ;

2e Sm3 b

ð:10Þ

ð:11Þ

The moments of the NDF given in this section as Eq. (7) were needed for the calculation of the Sauter Mean Diameter, d32 for verification. The first four length-based moments of this function can be written as

m1 ðtÞ ¼

2

bðnjn1 Þ ¼

2

m2 ðtÞ ¼ and

ð:12Þ

C

  4 ; 3

ð:13Þ

C 5=3

  5 ; 3

ð:14Þ

k1 4=3 3k2

k1 3k2

Y.O’Neil Williams et al. / Chemical Engineering Science 201 (2019) 475–483

m3 ðt Þ ¼ m3 ;

ð:15Þ

where C is the standard gamma function. References Abdurahman, N.H., Azhari, N.H., Yunus, Y.M., 2013. Formulation and evaluation of water-continuous emulsion of heavy crude oil prepared for pipeline transportation. IJESIT 2 (5). Bhutani, G., Brito-Parada, P.R., Cilliers, J.J., 2016. Polydispersed flow modelling using population balances in an adaptive mesh finite element framework. Comput. Chem. Eng. 87, 208–225. Bremond, N., Bibette, J., 2012. Exploring emulsion science with microfluidics. Soft Matter 8 (41), 10549–10559. Chesters, A.K., 1991. Modelling of coalescence processes in fluid-liquid dispersions: a review of current understanding. Chem. Eng. Res. Des. 69 (A4), 259–270. Christopher, G.F., Bergstein, J., End, N.B., Poon, M., Nguyen, C., Anna, S.L., 2009. Coalescence and splitting of confined droplets at microfluidic junctions. Lab Chip 9 (8), 1102–1109. Davis, R.H., Schonberg, J.A., Rallison, J.M., 1989. The lubrication force between two viscous drops. Phys. Fluids 1 (1), 77–81. Fu, T., Wu, Y., Ma, Y., Li, H.Z., 2012. Droplet formation and breakup dynamics in microfluidic flow-focusing devices: from dripping to jetting. Chem. Eng. Sci. 84, 207–217. Gu, H., Duits, M.H.G., Mugele, F., 2011. Droplets formation and merging in twophase flow microfluidics. Int. J. Mol. Sci 12 (4), 2572–2597. Huang, P.Y., Hellums, J.D., 1993a. Aggregation and disaggregation kinetics of human blood platelets: Part i. development and validation of a population balance method. Biophys. J. 65 (1), 334–343. Huang, P.Y., Hellums, J.D., 1993b. Aggregation and disaggregation kinetics of human blood platelets: Part ii. shear-induced platelet aggregation. Biophys. J. 65 (1), 344. Hu, Y.T., Pine, D.J., Leal, G., 2000. Drop deformation, breakup, and coalescence with compatibilizer. Phys. Fluids 12 (3), 484–489. Israelachvili, J.N., 2011. Intermolecular and Surface Forces. Academic Press. Ivanov, I.B., Danov, K.D., Kralchevsky, P.A., 1999. Flocculation and coalescence of micron-size emulsion droplets. Colloids Surf., A 152 (1-2), 161–182. Klink, I.M., Phillips, R.J., Dungan, S.R., 2011. Effect of emulsion drop-size distribution upon coalescence in simple shear flow: A population balance study. J. Colloid Interface Sci. 353 (2), 467–475. Krebs, T., Schroen, K., Boom, R., 2012a. A microfluidic method to study demulsification kinetics. Lab Chip 12 (6), 1060–1070. Krebs, T., Schroën, K., Boom, R., 2012b. Coalescence dynamics of surfactantstabilized emulsions studied with microfluidics. Soft Matter 8 (41), 10650– 10657. Krebs, T., Schroen, C.P.G.H., Boom, R.M., 2012c. A microfluidic study of oil-water separation kinetics. Adv. Fluid Mech. IX 74, 427–438.

483

Krebs, T., Schroën, C., Boom, R.M., 2013. Coalescence kinetics of oil-in-water emulsions studied with microfluidics. Fuel 106, 327–334. Lazzari, S., Theiler, P.M., Shen, Y., Coley, C.W., Stemmer, A., Jensen, K.F., 2018. Ligand-mediated nanocrystal growth. Langmuir 34 (10), 3307–3315. Leal, L.G., 2004. Flow induced coalescence of drops in a viscous fluid. Phys. Fluids 16 (6), 1833–1851. Liao, Y., Lucas, D., 2010. A literature review on mechanisms and models for the coalescence process of fluid particles. Chem. Eng. Sci. 65 (10), 2851–2864. Marchisio, D.L., Fox, R.O., 2005. Solution of population balance equations using the direct quadrature method of moments. J. Aerosol Sci 36 (1), 43–73. Marchisio, D.L., Pikturna, J.T., Fox, R.O., Vigil, R.D., Barresi, A.A., 2003. Quadrature method of moments for population-balance equations. AlChE J. 49 (5), 1266– 1276. McCoy, B.J., Madras, G., 2003. Analytical solution for a population balance equation with aggregation and fragmentation. Chem. Eng. Sci. 58 (13), 3049–3051. Muijlwijk, K., Colijn, I., Harsono, H., Krebs, T., Berton-Carabin, C., Schroën, K., 2017. Coalescence of protein-stabilised emulsions studied with microfluidics. Food Hydrocoll. 70, 96–104. Podgórska, W., 2005. Scale-up effects in coalescing dispersions-comparison of liquid-liquid systems differing in interface mobility. Chem. Eng. Sci. 60 (8), 2115–2125. Prince, M.J., Blanch, H.W., 1990a. Bubble coalescence and break-up in air-sparged bubble columns. AIChE J. 36 (10), 1485–1499. Prince, M.J., Blanch, H.W., 1990b. Transition electrolyte concentrations for bubble coalescence. AIChE J. 36 (9), 1425–1429. Ramkrishna, D., 2000. Population Balances: Theory and Applications to Particulate Systems in Engineering. Elsevier. Sajjadi, B., Abdul Raman, A.A., Ehsan Shah, R.S.S-R., Ibrahim, S., 2013. Review on applicable breakup/coalescence models in turbulent liquid-liquid flows. Rev. Chem. Eng. 29 (3). Shardt, O., Derksen, J.J., Mitra, S.K., 2013. Simulations of droplet coalescence in simple shear flow. Langmuir 29 (21), 6201–6212. Silva, LF., Damian, R.B., Lage, P., 2008. Implementation and analysis of numerical solution of the population balance equation in cfd packages. Comput. Chem. Eng. 32 (12), 2933–2945. Tandon, P., Diamond, S.L., 1997. Hydrodynamic effects and receptor interactions of platelets and their aggregates in linear shear flow. Biophys. J. 73 (5), 2819– 2835. Teh, S-Y., Lin, R., Hung, L-H., Lee, A.P., 2008. Droplet microfluidics. Lab Chip 8 (2), 198–220. Tsouris, C., Tavlarides, L.L., 1994. Breakage and coalescence models for drops in turbulent dispersions. AlChE J. 40 (3), 395–406. Van de Ven, T.G.M., 1989. Colloidal Hydrodynamics. Academic Press. Xu, B., Nguyen, N-T., Neng Wong, T., 2011. Droplet coalescence in microfluidic systems. Micro Nanosyst. 3 (2), 131–136. Yang, H., Park, C.C., Hu, Y.T., Leal, L.G., 2001. The coalescence of two equal-sized drops in a two-dimensional linear flow. Phys. Fluids 13 (5), 1087–1106.