# Deformation and break-up of viscoelastic droplets Using Lattice Boltzmann Models

###### Abstract

We investigate the break-up of Newtonian/viscoelastic droplets in a viscoelastic/Newtonian matrix under the hydrodynamic conditions of a confined shear flow. Our numerical approach is based on a combination of Lattice-Boltzmann models (LBM) and Finite Difference (FD) schemes. LBM are used to model two immiscible fluids with variable viscosity ratio (i.e. the ratio of the droplet to matrix viscosity); FD schemes are used to model viscoelasticity, and the kinetics of the polymers is introduced using constitutive equations for viscoelastic fluids with finitely extensible non-linear elastic dumbbells with Peterlin’s closure (FENE-P). We study both strongly and weakly confined cases to highlight the role of matrix and droplet viscoelasticity in changing the droplet dynamics after the startup of a shear flow. Simulations provide easy access to quantities such as droplet deformation and orientation and will be used to quantitatively predict the critical Capillary number at which the droplet breaks, the latter being strongly correlated to the formation of multiple neckings at break-up. This study complements our previous investigation on the role of droplet viscoelasticity (A. Gupta & M. Sbragaglia, Phys. Rev. E 90, 023305 (2014)), and is here further extended to the case of matrix viscoelasticity.

###### keywords:

Droplet Microfluidics, Deformation and Break-up, Viscoelasticity, Lattice Boltzmann Models;^{†}

^{†}volume:

00 \journalnameProcedia IUTAM \runauthAnupam Gupta \jidpiutam

IUTAM Symposium on Multiphase flows with phase change: challenges and opportunities, Hyderabad, India (December 08 – December 11, 2014)

agupta@roma2.infn.it University of “Tor Vergata”

\deftermDroplet/matrix velocity \defitem\deftermDroplet/matrix pressure \defitem\deftermDynamic shear viscosity (matrix (M), fluid solvent (f), droplet (d)) \defitem\deftermViscosity ratio between the dispersed (droplet) and continuum (matrix) phase \defitem\deftermGap spacing, wall velocity, shear rate \defitem\deftermDroplet radius \defitem, , \deftermPolymer conformation tensor, polymer relaxation time, polymer viscosity \defitem\deftermFENE-P potential, maximum elongation of the polymers \defitemCa, De\deftermCapillary number, Deborah number \defitem\deftermDroplet emulsion time

## 1 Introduction

Emulsion properties are largely determined by their microstructure which can be tuned and designed for a huge variety of applications Christophher07 . In particular, deformation and break-up of dispersed droplets determine the emulsion rheology Larson . Droplet deformation and break-up in Newtonian fluids have been extensively studied in the literature Taylor34 ; Grace ; Stone . The effect of an unconfined shear flow on droplets of one fluid suspended in another immiscible fluid was first considered long time ago by Taylor Taylor34 : he estimated the largest stable droplet radius by balancing the surface stresses due to interfacial tension and viscous stresses due to shear. A dimensionless measure of this balance is provided by the Capillary number , where is the dynamic viscosity of the fluid matrix, the shear rate, the droplet radius at rest and the surface tension. Break-up occurs at a critical Capillary number when the viscous forces overcome the surface forces. The problem of droplet deformation and break-up under confined shear flow between two parallel plates has also been addressed in a series of theoretical and experimental papers (see Shapira ; Sibillo06 ; Vananroye07 ; Janssen10 and references therein). It was suggested that under confined conditions, a uniform shear flow can be exploited to generate quasi monodisperse emulsions by controlled break-up at near-critical conditions Sibillo06 ; Renardy07 . The properties of confined droplets that contain viscoelastic components are less studied Cardinaels09 ; Minale10 ; Cardinaels11 and the critical conditions for break-up have been rarely explored so far. Recent experiments suggest that viscoelasticity changes profoundly the critical Capillary numbers in confined conditions Cardinaels11 . Complementing these kind of results with the help of numerical simulations would be of extreme interest. Simulations provide easier access to quantities such as droplet deformation and orientation as well as the velocity flow field and pressure field inside and outside the droplet. The goal of this paper is to use numerical simulations to characterize the idealized problem of a Newtonian/viscoelastic droplet subject to simple shear in a confined viscoelastic/Newtonian matrix.

## 2 Theoretical Model

Our numerical approach is based on a hybrid combination of Lattice-Boltzmann models (LBM) and finite difference (FD) schemes, the former used to model two immiscible fluids with variable viscosity ratio, and the latter used to model viscoelasticity using the FENE-P constitutive equations. LBM have already been used to model droplet deformation problems Xi99 ; VanDerSman08 ; Komrakovaa13 ; Liuetal12 and also viscoelastic flows Onishi2 ; Malaspinas10 . The approach we use has already been studied and validated in a dedicated work SbragagliaGuptaScagliarini , where we have provided evidence that the model is able to capture quantitatively rheological properties of dilute suspensions as well as deformation and orientation of single droplets in confined shear flows. We just recall here the relevant continuum equations which are integrated in both the droplet (d) and the matrix (M) phases. In the droplet phase we integrate both the NS (Navier-Stokes) and FENE-P reference equations:

(1) | |||||

(2) |

Here, is the dynamic viscosity of the fluid, the viscosity parameter for the FENE-P solute, the polymer relaxation time, the solvent density, the solvent pressure, the transpose of , the polymer-conformation tensor, the identity tensor, the FENE-P potential that ensures finite extensibility, and is the maximum possible extension of the polymers bird . In the outer matrix (M) phase (indicated with a prime), the corresponding equations are

(3) | |||||

(4) |

with the solvent matrix shear viscosity. In all the cases, immiscibility between the droplet phase and the matrix phase is introduced using the so-called “Shan-Chen” model SbragagliaGuptaScagliarini . In all the numerical simulations presented in this paper, we work with unitary viscosity ratio, defined in terms of the total (solvent+polymer) shear viscosity. In particular, when studying matrix viscoelasticity (MV), we will choose a case with in equation (1) and and polymer concentration ; for the simulations with droplet viscoelasticity (DV), we will choose a case with in equation (3) with and polymer concentration . The degree of viscoelasticity is computed from the Deborah number

(5) |

where Ca is always computed in the matrix phase while the Deborah number is computed in either the matrix or the droplet phase, dependently on the case studied. In Eq. (5), is the first normal stress difference which develops in homogeneous steady shear. Solving the constitutive equations for such a hydrodynamic problem, , , , both the polymer feedback stress and the first normal stress difference for the FENE-P model bird ; Lindner03 follow (primed variables replace non-primed variables for matrix phases)

(6) |

(7) |

In the Oldroyd-B limit () we can use the asymptotic expansion of the hyperbolic functions and we get

(8) |

Equation (8) shows that De is clearly dependent on the ratio between the polymer relaxation time and the droplet emulsion time

(9) |

As evident from Eq. (6), the model supports a thinning effect at large shear, although such effect will not be important in our calculations, all the numerical simulations being performed with fluid pairs with nearly constant shear viscosities. In the following sections, we report the Deborah number based on the definition (8), as we estimated the difference between (8) and (5) to be at maximum of a few percent for the values of considered. Also, we focus mainly on the droplet deformation and break-up problems, being the quantitative benchmarks against known analytical results for the rheology of dilute suspensions bird ; Herrchen present in another dedicated methodological publication SbragagliaGuptaScagliarini . In a previous study SbragagliaGuptaPRE we investigated the role of droplet viscoelasticity: a non trivial interplay between confinement and viscoelasticity has emerged. With the use of numerical simulations we had the opportunity to change separately the viscosity ratio of the Newtonian phases, the maximum extension of the polymers, and the degree of viscoelasticity, thus allowing for a systematic analysis of the viscoelastic effects while keeping the shear viscosity of the droplet fixed to the reference Newtonian case. In particular, by increasing the finite extensibility of the polymers, it was observed that the resistance against elongation may be enough to prevent both droplet elongation and subsequent triple break-up, thus altering significantly the critical Capillary number for viscoelastic droplets under confinement. In this paper, we push the analysis a bit further and we propose a comparative study between matrix and droplet viscoelasticity. To simplify matters, we will also keep the maximum elongation of the polymers fixed to , since we have exhaustively treated the importance of in our previous study SbragagliaGuptaPRE .

## 3 Results

In all the cases discussed in this section, a spherical droplet with radius is initially placed halfway between the walls. The critical Capillary number is computed by identifying the pre-critical () and the post-critical wall velocity (), i.e. the largest (smallest) wall velocity for which the droplet is stable (breaks). All the simulations described refer to cases with polymer relaxation times ranging in the interval lbu and finite extensibility parameter , corresponding to Deborah numbers ranging in the interval . The numerical simulations have been carried out in three dimensional domains . The droplet radius and the vertical gap have been changed in the ranges lattice cells and lattice cells to achieve different confinement ratios . The stream-flow (x) direction is resolved with lattice cells, depending on the droplet elongation properties, while the transverse-flow (y) direction is resolved with lattice cells. Periodic conditions are applied in the stream-flow and in the transverse-flow directions. The droplet is subjected to a linear shear flow , , with the shear introduced with two opposite velocities in the stream-flow direction () at the upper () and lower wall ().

In Figs. 2 and 3 we show the time history for droplets in post-critical conditions at changing confinement and viscoelasticity. Fig. 2 refers to cases with lower confinement ratio and including matrix and droplet viscoelasticity. For each case we consider three representative snapshots showing (i) the initial droplet deformation (ii) the droplet deformation prior to break-up (iii) and the droplet in post break-up conditions. We use the droplet emulsion time (9) as a unit of time. For the Newtonian case we find , which is different from the usual unconfined result Grace ; Janssen10 . This can be attributed to the finite Reynolds number () of our simulations RenardyCristini01 . This fact said, we observe that both droplet and matrix viscoelasticity do not have an important influence on the critical Capillary numbers for break-up. In the case of droplet viscoelasticity we find a small stabilization that increases the critical Capillary number by some percent; wheres matrix viscoelasticity is not producing any visible effect on the critical Capillary number. It must be emphasized that we tuned the polymeric viscosity in such a way to reproduce always a unitary viscosity ratio between the droplet phase and the matrix phase. Fig. 3 is essentially the counterpart of Fig. 2 for an increased confinement ratio. A series of hints are given by the visual inspection of the droplet shapes and the associated critical Capillary numbers. First, the confinement ratio is already large enough to stabilize long droplet shapes in the Newtonian case, thus triggering the emergence of triple break-up Janssen10 . Such droplet shapes would be unstable in unconfined flows: confinement makes them stable and the droplet can sustain larger Capillary numbers before break-up. Break-up mechanism itself changes, as the droplet can reach a minimum length at which a Rayleigh-Plateau instability Janssen10 develops at the interface and breaks the droplet in equally sized daughter droplets (Panel (c)). This fact is known from the literature Janssen10 ; Cardinaels11 and lends further support SbragagliaGuptaScagliarini to the validity of the numerical approach. Second, and more interestingly, the role of matrix and droplet viscoelasticity seems opposite. Droplet viscoelasticity reduces droplet elongation and higher Capillary numbers are needed to break the droplet. The break-up process still leads to the formation of multiple neckings but the degree of monodispersity of the resulting daughter droplets gets affected at the Deborah number studied SbragagliaGuptaPRE . On the other hand, we observe that matrix viscoelasticity completely suppresses the formation of multiple neckings and the break-up process looks much more similar to the unbounded case. We remark that the finite extensibility parameter has been kept fixed to . In another study SbragagliaGuptaPRE , we investigated systematically the importance of the finite extensibility parameter for the case of droplet viscoelasticity. It has been found that increasing could lead to a situation where elongated droplet shapes cannot not be stable anymore due to the net increase of the polymer elongational viscosity, which actually increases at increasing . In such a case, also with droplet viscoelasticity the critical Capillary number decreases with respect to the Newtonian case.

In Fig. 4 we show the dimensionless droplet elongation as a function of time for several values of De and Ca, at fixed confinement ratio , comparing both matrix and droplet viscoelasticity with the Newtonian case. Since the shape of highly deformed and confined droplets deviates from an ellipsoid, we estimated the droplet elongation from the projection of the droplet length () in the velocity direction. In Panel (a) of Fig. 4 we report results for a given Capillary number . If compared with the Newtonian case, the maximum elongation of the droplet is suppressed in the case of droplet viscoelasticity, while is enhanced in the case of matrix viscoelasticity, which is a signature that the critical Capillary number for matrix viscoelasticity is smaller than the Newtonian one. This happens for Deborah number just above unity, whereas results with small Deborah numbers are actually very close to the Newtonian case. In Panel (b) we report the pre-critical and post-critical time history for both Newtonian and viscoelastic droplets. We notice that the maximum dimensionless elongation achieved before break-up, , is larger for the Newtonian case compared to the case of matrix viscoelasticity, which indeed breaks at smaller Capillary number. On the other hand, the case with droplet viscoelasticity achieves a maximum elongation before break-up that is roughly doubled with respect to the Newtonian case. Also the critical Capillary number is roughly doubled (see Fig. 5).

Panel (a) of Fig. 5 summarizes and extends the findings of Figs. 2-4 to other confinement ratios and degrees of viscoelasticity in both the matrix and droplet phases. We report data for weakly viscoelastic systems () and also data with Deborah number just above unity (). As already noticed elsewhere SbragagliaGuptaPRE , for Newtonian droplets the role of confinement is almost insignificant up to , whereas for larger confinement ratios a monotonous increase of is observed. The emergence of the up-turn in is intimately connected to the change of the break-up mechanism. For small Deborah numbers the curve vs. does not change much for droplet viscoelasticity, whereas some decrease in can be readily seen for the case of matrix viscoelasticity. The black open circles indicate situations where ternary break-up is observed. Actually, the tendency of viscoelasticity to alter and change the stability properties of confined droplets is perceptibly more pronounced in the case of matrix viscoelasticity than droplet viscoelasticity. For Deborah number just above unity () both matrix and droplet viscoelasticity alter significantly the critical Capillary number at break-up and the changes are more pronounced and amplified at larger confinement ratios. In Panel (b) of Fig. 5, we report the maximum dimensionless elongation of the droplet, , as a function of the confinement ratio. It is clear that the trends for and are quite similar.

## 4 Conclusions

The effects of viscoelasticity and geometrical confinement on droplet break-up in a shear flow have been studied. Our analysis strongly benefited of numerical simulations, where we could model immiscible fluids in the presence of viscoelastic behaviour in either the droplet or the matrix phase. Numerical simulations offer great flexibility and easy access to hydrodynamic quantities and droplet interface dynamics, thus revealing particularly suited for the study at hand. We have found that the effect of viscoelasticity is rather insignificant in unbounded cases, whereas it gets amplified with confinement. In particular, viscoelasticity reduces the capability of micro-confined shear flows to generate monodisperse emulsions. This is a feature that we found in a previous study on droplet viscoelasticity SbragagliaGuptaPRE and we confirm in presence of matrix viscoelasticity. At small Deborah number, the tendency of viscoelasticity to alter and change the stability properties of confined droplets is more pronounced in the case of matrix viscoelasticity, if compared to the case of droplet viscoelasticity.

For future research it would be extremely interesting to repeat some of the numerical simulations at changing the polymer concentration and/or at changing the geometry of the system (i.e. T-shaped channels, flow-focusing devices) to reveal other interesting features on the dynamics of viscoelastic fluids in complex geometries.

## 5 Acknowledgments

We are particularly grateful to F. Bonacccorso for technical support. We kindly acknowledge funding from the European Research Council under the European Community’s Seventh Framework Programme (FP7/2007-2013) / ERC Grant Agreement N. 279004. A. Gupta acknowledges R. Pandit and S. S. Ray for fruitful discussions during his visits to IISc Bangalore (February 2014, August 2014) and ICTS-TIFR Bangalore (October 2013).

## References

- (1) G.F. Christopher & S.L. Anna J Phys D Appl Phys 2007 40:R319–R336
- (2) R.G. Larson, The Structure and Rheology of Complex Fluids, Oxford University Press, New York (1999)
- (3) G.I. Taylor, Proc Roy Soc A 1932; 138:41–8.
- (4) H. P. Grace, Chem. Eng. Commun. 1982; 14, 225
- (5) H. A. Stone, Annu. Rev. Fluid Mech. 1994; 26, 65
- (6) M. Shapira & S. Haber, Int J Multiph Flow 1990;16:305–21
- (7) V. Sibillo, G. Pasquariello, M. Simeone, V. Cristini & S. Guido, Phys Rev Lett 2006; 97, 054502
- (8) A. Vananroye, P. Van Puyvelde & P. Moldenaers, J. Rheol. 2007; 51, 139–153
- (9) P. J. A. Janssen, A. Vananroye, P. Van Puyvelde, P. Moldenaers & P. D. Anderson, J Rheol 2010; 54, 1047–1060
- (10) Y. Renardy, Rheol Acta 2007; 46, 521–529
- (11) Y. Renardy & V. Cristini Phys. Fluids 2007; 13, 7–13
- (12) R. Cardinaels, K. Verhulst & P. Moldenaers, J. Rheol 2009; 53, 1403–1424
- (13) M. Minale, S. Caserta & S. Guido, Langmuir 2010; 26, 126–132
- (14) R. Cardinaels & P. Moldenaers, Microfluid Nanofluid 2011; 10, 1153–1163
- (15) R. B. Bird, R. C Armstrong & O. Hassager, Dynamics of polymeric liquids, J. Wiley & sons (1987).
- (16) J. Onishi, Y. Chen, H. Ohashi, Physica A 2006; 362 84–92
- (17) H. Xi & C. Duncan, Phys. Rev. E. 1999; 59, 3022–3037
- (18) X. Shan & H. Chen, Phys. Rev. E 1993; 47, 1815
- (19) R. Benzi, M. Sbragaglia, S. Succi, M. Bernaschi & S. Chibbaro, J. Chem. Phys. 2009; 131, 104903
- (20) A. E. Komrakovaa, O. Shardta, D. Eskinb & J.J. Derksen , International Journal of Multiphase Flow 2013; 59, 24–43
- (21) R.G.M. Van der Sman & S. van der Graaf, Comput. Phys. Commun. 2008; 178, 492–504
- (22) H. Liu, A. J. Valocchi & Q. Kang, Phys. Rev. E 2012; 85, 046309
- (23) O. Malaspinas, N. Fietier & M. Deville, Jour. Non Newt. Fluid Mech 2010; 165, 1637–1653
- (24) M. Herrchen & H. C. Ottinger, J. Non-Newtonian Fluid Mech. 1997; 68, 17–42
- (25) A. Lindner, J. Vermant, D. Bonn, Physica A 2003; 319, 125–133
- (26) A. Gupta, M. Sbragaglia and A. Scagliarini, arXiv/1406.2686 (2014)
- (27) A. Gupta and M. Sbragaglia, Phys. Rev. E 90, 023305 (2014)