Shear-thinning and shear-thickening emulsions in shear flows

We study the rheology of a two-fluid emulsion in semi-concentrated conditions; the solute is Newtonian while the solvent an inelastic power law fluid. The problem at hand is tackled by means of direct numerical simulations using the volume of fluid method. The analysis is performed for different volume fractions and viscosity ratios under the assumption of negligible inertia and zero buoyancy force. Several carrier fluids are considered encompassing both the shear-thinning and thickening behaviors. We show that the effective viscosity of the system increases for shear-thickening fluids and decreases for the shear-thinning ones for all the viscosity ratios considered. The changes in the emulsion viscosity are mainly due to modifications of the coalescence in the system obtained by changing the carrier fluid property: indeed, local large and low shear rates are found in the regions between two interacting droplets for shear-thickening and thinning fluids, respectively, resulting in increased and reduced local viscosity which ultimately affects the drainage time of the system. This process is independent of the nominal viscosity ratio of the two fluids and we show that it can not be understood by considering only the mean shear rate and viscosity of the two fluids across the domain, but the full spectrum of shear rate must be taken into account.


I. INTRODUCTION
Emulsions are mixtures of two or more liquids that are partially or totally immiscible. They are present in many biological and industrial applications such as waste treatment, oil recovery and pharmaceutical manufacturing and are also relevant in the field of colloidal science where the accuracy and the control of the production process of functional materials rely on the knowledge of the complex microstructure of the suspension 1 . In this work we focus on the rheology of biphasic emulsions by means of direct numerical simulations. The study of twophase flows has recently become of utmost importance related to the ongoing COVID-19 pandemic, caused by airborne transmission of virus-containing saliva droplets transported by human exhalations [2][3][4][5] The study of rheology is motivated by the many fluids in nature and industrial applications which exhibit a non-Newtonian behavior, i.e., a nonlinear relation between the shear stress and the shear rate; the relation between these macroscopic behaviors and the microstructure is often studied assuming suspensions of objects in a Newtonian solvent with dynamic viscosity µ. Einstein 6 was the first to provide a closure for the effective viscosity µ e of a dilute rigid particle suspensions with vanishing inertia, and showed theoretically that µ e is a linear function of the particle volume fraction Φ Non-rigid and deformable objects, such as deformable particles 7 , capsules 8 and droplets 9 behave differently because of their deformation and in the latter case also coalescence and breakup. Taylor 10 extended equation (1) a) Corresponding author: marco.rosti@oist.jp by introducing the viscosity ratio between the two phases λ (defined as the disperse phase over the matrix phase viscosity), thus obtaining which for a unit viscosity ratio reduces to Pal 11,12 derived expressions to evaluate the effective viscosity of infinitely dilute and concentrated emulsions using the effective medium approach. These relations assume limiting cases to model surface tensions effects (either going to zero or infinite) and thus direct numerical simulation is a valuable tool to overcome these limitations. Zhou and Pozrikidis 13 simulated numerically twodimensional emulsions and their work was later extended to three-dimensional flows by Loewenberg and Hinch 14 . Their results revealed the complexity of emulsions rheology, with pronounced shear thinning and large normal stresses associated with an anisotropic microstructure resulting from the alignment of deformed drops in the flow direction. More recently, Srivastava, Malipeddi, and Sarkar 15 also investigated the inertial effects on emulsions and found that, while in the absence of inertia emulsions display positive first and negative second normal stress differences, small amount of inertia alters their signs with the first normal stress difference becoming negative and the second one positive. The sign change is caused by the increasing drop inclination in the presence of inertia, which in turn directly affects the interfacial stresses. These previous studies were focused only on deforming droplets, without taking into account coalescence and breakup, which however plays a key role for the rheological behaviour of emulsions and are the main arXiv:2108.00571v1 [physics.flu-dyn] 2 Aug 2021 elements distinguishing emulsions from particle suspensions.
The effect of coalescence was tackled numerically more recently by Rosti, De Vita, and Brandt 9 and De Vita et al. 16 who showed that the effective viscosity µ e of a droplet suspension can be properly described by extending equation (1) to the second order similarly to what proposed by Batchelor and Green 17 for dilute particle suspensions, i.e., where c 1 and c 2 are fitting parameters depending on the viscosty ratio λ and capillary number Ca. While c 2 is positive in the case of rigid particles, indicating that the effective viscosity always grows with the volume fraction (c 1 is positive consistently with equation (1) and equation (2)), on the other hand c 2 is negative for droplets, and thus µ e is a concave function of Φ, i.e., µ e has a maximum for an intermediate value of Φ and then decreases with the volume fraction. The volume fraction for which the maximum is reached decreases with the viscosity ratio λ and µ e is again a convex function of Φ in the case of droplets when the coalescence of the solute phase is suppressed 16 ; in this case, their behaviour resembles what found for deformable particles 18 , thus confirming that coalescence is a key mechanism that can dramatically change the rheology of emulsions. The role of coalescence in a turbulent channel flow has been recently investigated numerically by Cannon et al. 19 .
In this works we account for the coalescence and study how it affects the rheology of emulsions adding an additional complexity to the system by considering non-Newtonian solvents 20,21 . Few theoretical expressions were proposed in the past to predict the rheology of viscoelstic emulsions [22][23][24] . Palierne 23 derived an expression to predict the shear modulus of emulsions of viscoelastic materials accounting for the mechanical interactions between inclusions. A different expression which however provides similar predictions was later derived by Bousmina 24 who extended Kerner's model 22 for flow of composite elastic media to emulsions of viscoelastic phases with interfacial tension undergoing small deformations. This model is able to predict some general features typical of viscoelastic emulsions, such as shoulders in the storage modulus G at low frequency and a long relaxation time process. However, emulsions with non-Newtonian fluids are characterised by several peculiar behaviours found experimentally which are usually difficult to capture with theoretical models; for example, Mason and Bibette 25 were able to produce monodisperse emulsions using a viscoelastic medium under shear; the authors showed the possibility of controlling the final droplet size by altering the shearing conditions and the emulsion viscoelasticity, suggesting that the monodispersity resulted from the suppression of the capillary instability until the droplet is sufficiently elongated. A similar study was later performed by Zhao and Goveas 26 who investigated the deformation and breakup of a dilute emulsion with a vis-coelastic continuous phase under shear flow. Also these authors found a size selection mechanism of the resulting ruptured drops, not found in Newtonian fluids. To be best of our knowledge, no systematic study of the rheology of concentrated emulsions are available in the literature, with only very few studies focusing on rigid particles suspensions in viscoelastic media [27][28][29][30][31][32][33] , and the present work is aimed to fill in this gap.
The paper is structured as follow: in section II we describe the governing equations and the numerical method used to numerically solve them. After discussing the chosen setup, in section III we present the results of our simulations and discuss the effect of shear-thinning and thickening on the rheology of emulsions. Finally, we summarize the main findings and draw the main conclusions in section IV.

II. METHODOLOGY
We consider a flow with two incompressible viscous fluids separated by an interface in a channel with moving walls, i.e., in a plane Couette geometry. Figure 1(top) shows a view of the geometry considered in the present analysis, where x, y and z (x 1 , x 2 , and x 3 ) denote the streamwise, wall-normal, and spanwise coordinates, while u, v and w (u 1 , u 2 , and u 3 ) denote the corresponding components of the velocity vector field. The lower and upper impermeable moving walls are located at y = −h and y = h, respectively, and move in opposite direction with constant streamwise velocity u = ±U w , providing a nominal shear-rateγ 0 = 2U w /2h.
To identify the two fluid phases, we introduce an indicator (or color) function H. In particular, H = 1 in the regions occupied by the fluid 1 and H = 0 otherwise. In this work, we assume fluid 1 to be the solute (dispersed phase) and fluid 2 the solvent (carrier phase). Since the interface is transported by the fluid velocity, we advect the indicator function as follows where φ is the cell-averaged value of H also called the volume of fluid function. The two fluids are governed by the momentum conservation and the incompressibility constraint and which can be written as one set of equations incorporating the interface conditions, which read Here, ρ is the fluid density assumed to be the same in the two phases, σ ij the Cauchy stress tensors and f i = σκn i δ ≈ σκ∂φ/∂x i is added to account for the interface condition 34 , where σ is the interfacial surface tension coefficient, κ and n i the local curvature and unit normal vector of the interface and δ is the Dirac's delta function at the interface. The Cauchy stress tensor σ ij defined as ij = p ij + 2µ l D ij where p is the presure, D ij = (@u i /@x j + @u j @x i ) /2 the strain rate tensor nd µ l = µ 1 + (1 ) µ 2 the local dynamic viscosity ritten in a mixture form as a simple weighted average ased on the local value of of the dynamic viscosity in he two phases µ 1 and µ 2 . The solute fluid 1 is assumed o be Newtonian, and its dynamic viscosity µ 1 is a contant, while the solvent fluid 2 is a non-Newtonian fluid escribed by a power law model where the local value of 2 is a function of the local shear-rate˙ l .
The equations of motion are solved on an uniform stagered fixed Eulerian grid where the fluid velocity compoents are located on the cell faces and all other variables density, viscosity, pressure, stress and volume of fluid) t the cell centers. The time integration is performed ith a fractional-step method based on the second-order dams-Bashforth scheme. The indicator function H eeded to capture the interface between the two fluids is und following the multidimensional tangent of hyper-

A. Non-Newtonian fluid model
Fluids that exhibit a non-linear behavior between th shear stress and the shear rate are called non-Newtonia and, in particular, fluids whose response does not depen explicitly on time but only on the present shear rate ar denoted generalized Newtonian fluids. When the shea stress increases more than linearly with the shear rat the fluid is called dilatant or shear-thickening, whereas i the case of opposite behavior, i.e., when the shear stres increases less than linearly with the shear rate, the flui is called pseudoplastic or shear-thinning. Several mode have been developed to capture the di↵erent behaviour of various non-Newtonian fluids and in the current stud we focus on the simple inelastic power law models, wher the local viscosity of the fluid is a function of the so local value of shear rate. A relation that can be used t summarize the behaviors previously described for com plex fluids is µ pl = K˙ l n 1 , (7 is defined as σ ij = −pδ ij + 2µ l D ij where p is the pressure, D ij = (∂u i /∂x j + ∂u j ∂x i ) /2 the strain rate tensor and µ l = φµ 1 + (1 − φ) µ 2 the local dynamic viscosity written in a mixture form as a simple weighted average based on the local value of φ of the dynamic viscosity in the two phases µ 1 and µ 2 . The solute fluid 1 is assumed to be Newtonian, and its dynamic viscosity µ 1 is a constant, while the solvent fluid 2 is a non-Newtonian fluid described by a power law model where the local value of µ 2 is a function of the local shear-rateγ l . The equations of motion are solved on an uniform staggered fixed Eulerian grid where the fluid velocity components are located on the cell faces and all other variables (density, viscosity, pressure, stress and volume of fluid) at the cell centers. The time integration is performed with a fractional-step method based on the second-order Adams-Bashforth scheme. The indicator function H needed to capture the interface between the two fluids is found following the multidimensional tangent of hyperbola for interface capturing (MTHINC) method proposed by Ii et al. 35 and transported with the volume of fluid method. More details on the numerical implementation can be found in Rosti, De Vita, and Brandt 9 .

A. Non-Newtonian fluid model
Fluids that exhibit a non-linear behavior between the shear stress and the shear rate are called non-Newtonian and, in particular, fluids whose response does not depend explicitly on time but only on the present shear rate are denoted generalized Newtonian fluids. When the shear stress increases more than linearly with the shear rate, the fluid is called dilatant or shear-thickening, whereas in the case of opposite behavior, i.e., when the shear stress increases less than linearly with the shear rate, the fluid is called pseudoplastic or shear-thinning. Several models have been developed to capture the different behaviours of various non-Newtonian fluids and in the current study, we focus on the simple inelastic power law models, where the local viscosity of the fluid is a function of the sole local value of shear rate. A relation that can be used to summarize the behaviors previously described for complex fluids is where n is the flow index and K the fluid consistency index. A Newtonian behavior is recovered when n = 1 and K = µ, while values of the flow index above and below unity, n > 1 and n < 1, denote shear-thickening and shear-thinning fluids, respectively. The consistency index K measures how strong the fluid responds to the imposed deformation rate but it is not possible to compare different values of K for fluids with different flow indexes n because its dimension is a function of n itself. In general, the local viscosity of the non-Newtonian fluid increases withγ l for shear-thickening fluids, while it reduces for shear-thinning ones, which means that the fluidity of shear-thickening fluids reduces increasing the shear rate, while the opposite is true for shear-thinning fluids. In the previous relation, the local shear rateγ l is the second invariant of the strain-rate tensor D ij and is computed by its dyadic product, i.e.,γ l = 2D ij D ij 36 . The viscosity of a power law shear-thinning fluid becomes infinite for null shear-rate; to overcome this numerical issue, the Carreau fluid model is used instead in which the local viscosity is computed as In this equation, µ ∞ and µ 0 indicate the lower and upper limits of fluid viscosity at infinite and zero shear rates. The flow index n characterizes the behaviour of the fluid: for n < 1 the fluid is shear thinning and the material time constant k represents the degree of shear-thinning. More details on the Carreau and power-low models can be found in 37 .

B. Setup
We consider a two-phase system in a Couette flow consisting of a suspension of droplets, see figure 1. The droplets are initially spherical with radius R and are randomly distributed in the computational domain, which is a rectangular box of size 16R × 10R × 16R in the stream, wall-normal and spanwise directions. The computational domain is discretized with a uniform Cartesian grid with 16 grid points per initial droplet radius R (i.e., 160 grid points per channel height 2h). Note that, the actual resolution of the droplets dynamic improves as droplets coalescence and grow in size. The suspended fluid is Newtonian with uniform constant viscosity µ 1 , while three different kinds of carrier fluids are studied: Newtonian fluids with viscosity µ 2 = µ, shear-thickening fluids with viscosity µ 2 = µ pl defined by equation (7) and shear-thinning fluids with viscosity µ 2 = µ ca defined by equation (8). The flow index n is varied in the range 0.25 ÷ 1.75 and we choose to define the three fluids such that the viscosity at the applied reference shear rateγ 0 is the same in all of them, i.e., µ. Thus, the parameter K in the power law fluid is determined by imposing µ = µ pl (γ 0 ) (yielding K = µγ 1−n 0 ), and similarly the Carreau fluid model parameter k is found by imposing µ = µ ca (γ 0 ). Note that, in our simulations with the Carreau model we have fixed the ratio µ 0 /µ ∞ to 10 4 . Figure 2 shows the rheology of all the carrier fluids considered in the present study. We vary the volume fractions of the solute phase Φ (defined as φ xyz with the brackets indicating the average operator) in the range 0 ÷ 0.3 and we consider two different viscosity ratios λ = µ 1 /µ equal to 1 and 0.01. The surface tension coefficient σ is fixed, such that the droplets have an initial capillary numbers Ca = µγ 0 R/σ equal to 0.2. Finally, the Reynolds number Re = ργ 0 R 2 /µ is fixed to 0.1, so that we can consider the inertial effects negligible.
Topological changes of droplets such as coalescence and breakup are treated here naturally without any additional model, since with the volume of fluid method these happen naturally when the distance between two interface falls within a grid cell. However, this means that the coalescence and breakup processes can be strongly influenced by the grid size resolution; a detailed discussion of this matter is reported in the appendix. Note also that, the parameters are chosen similar to that of previous works that can be found in the literature to ease comparisons 7,9,16,18 . Finally, the numerical code used in the present work has been tested and validated in the past in several works for laminar and turbulent flows of single and multiphase systems 18,33, .

III. RESULTS
We study the rheology of a droplet suspension immersed in power law fluids and compare the results with those obtained in a Newtonian fluid. To do so, we start our rheological analysis by considering the effective viscosity µ e which is the viscosity of a Newtonian fluid that gives the same shear stress at the nominal shear rate, and is thus defined as where σ w 12 is the shear component of the stress tensor evaluated at the walls. In general, we expect the rheology of a two fluid system to be a function of the Reynolds number Re, the capillary number Ca, the viscosity ratio λ, the solute volume fraction Φ, the confinement ratio 2R/2h and the non-Newtonian property of the carrier fluid here summarised by the power law index n. In the present work, we limit our analysis to inertialess flows, i.e., Re ≪ 1, the capillary number Ca is not varied and fixed to a value such that a single droplet of size R subject to the applied shear rateγ 0 does not breakup, and the domain size is chosen such that confinement effects are negligible. Thus, we can simplify the analysis to µ e ≈ F (Φ, λ, n). Figure 3 shows the effective viscosity µ e as a function of the droplet volume fraction Φ for two different viscosity ratio λ and for three different power law index n. In the Newtonian carrier fluid (black color), the behaviour is the same observed by De Vita et al. 16 : when the viscosity of the two fluids is the same, i.e., λ = 1, the effective viscosity µ e first grows with the volume fraction similarly to a rigid particle suspension, then reaches a maximum value for an intermediate Φ and then starts decreasing again. When the viscosity of the dispersed phase is reduced, i.e., λ decreases, the volume fraction for which the maximum effective viscosity is reached reduces, and for a sufficiently low λ the maximum is reached at Φ ≈ 0 and the effective viscosity curve decreases with Φ. Note that, the effective viscosity µ e can be even smaller than the carrier fluid one µ in the case of λ < 1. The non-monotonic behaviour of µ e with Φ is a direct consequence of the limiting behavior for Φ → 0 and Φ → 1 where µ e is equal to µ 2 and µ 1 by definition 9 . When the carrier fluid is non-Newtonian the situation is further modified. In general, for all the power law index n and viscosity ratio λ considered in the present study, we observe that the non-monotonic behaviour is preserved and that the effective viscosity µ e is larger for shear-thickening fluids with n > 1 and smaller for shear-thinning fluids with n < 1 than their Newtonian counterparts with n = 1. The difference in µ e between the Newtonian and non-Newtonian fluids grows with the volume fraction Φ and larger difference are evident between the Newtonian and shear-thinning fluids than what observed with the shear-thickening fluid for both the viscosity ratio λ considered. These changes of µ e with n are qualitatively similiar to what observed for a rigid particle suspension immersed in a power law carrier fluid 32 , except for the non-monotonicity of the viscosity curve which is due to the coalescence mechanism in the emulsions and absent in particle suspensions.
To understand the change in the effective viscosity curve with n, we compute the mean shear rateγ m = φγ l xyz,t / φ xyz,t and viscosity µ m = φµ l xyz,t / φ xyz,t in the carrier phase which are shown in figure 4. In the figure we observe that the mean shear rate exhibits a behaviour similar to the effective viscosity µ e :γ m grows with the volume fraction Φ, then reaches a maximum value and finally decreases for large Φ for λ = 1, while for λ < 1 the maximum is reached at lower volume fraction and for λ = 0.01 the mean shear rate deceases for all Φ. Also the effect of the flow index is analogous to what previously observed for µ e , and indeedγ m is larger for n > 1 and smaller for n < 1 than the Newtonian case with n = 1, with larger differences evident for the shearthinning fluids. Less trivial is the behaviour of the mean viscosity µ m ; for λ ≥ 1, being the shear rate larger than the reference oneγ 0 , the mean viscosity in the shearthickening fluid is larger than the Newtonian one while it is smaller in the shear-thinning fluid, as expected from the rheology of the fluid previously reported in figure 2. On the other hand, when λ < 1, the mean shear rate is actually smaller thanγ 0 and thus the behaviour of the non-Newtonian fluids is opposite, with the mean viscosity of the shear-thinning fluid being larger than the Newtonian one and smaller for the shear-thickening fluid. Thus, while for λ ≥ 1 the mean viscosity and shear rate are both contributing to increase or decrease the effective viscosity of the suspension, when λ < 1 because of this peculiar behaviour, the mean viscosity is actually counteracting and reducing the effect of the change in mean shear rate. This reverse trend obviously never happens in a particle suspension, where the shear rate is always larger than the single-phase value in the presence of suspended particles, and is thus a prerogative of emulsions.
Notwithstanding the counter-effect of the mean vis- cosity for λ < 1, the effective viscosity µ e in the shearthinning and thickening fluids is always smaller and larger than in the Newtonian fluid, respectively. In order to understand what is really affecting the emulsion viscosity, we need to study the value of the local shear rate not in terms of its mean value but by considering its full range of values assumed. Figure 5(top) reports the probability density function pdf of the shear rateγ l in the carrier phase for all the fluids studied with different n and λ and a selected volume fraction Φ equal to 10%. The graphs show that a large spectrum of shear rate exists for all the fluids, with the pdf characterised by a strong peak at its mean value (see figure 4(top)) around the nominal shear rateγ 0 and with a positive skewness, i.e., with more high shear rates than low ones. The latter result indicates an higher probability of finding low local viscosity in a shear-thinning fluid and of high viscosity in a shear-thickening one, independently of the viscosity ratio considered. The origin of these high shear rate events is the droplet interaction. To prove this, we show in figure 6 some instantaneous visualisation of the local viscosity in the solvent and solute phases. As we can clearly see from the figure, when two droplets interact in a shear-thickening fluid (panel 1 and 3) the local viscosity in the interstitial region between the two is large independently of the value of viscosity ratio λ considered, thus indicating a local high shear rate, while on the contrary when the two interacting droplets are suspended in a shear-thinning fluid (panel 2 and 4) the local viscosity in between is low, thus indicating a local low shear rate. Since the occurrence of the high shear rate events is linked to the droplet interactions, it is a function of the volume fraction Φ and grows with it, as shown in figure 5(bottom). The different viscosity found for shear-thickening and work to study the interaction of droplets. He proposed that the dynamics of two colliding droplets is an interplay of an external flow -the driving flow -responsible of the collision, and an internal flow -the drainage film between the two droplets -responsible of the interface deformation and rupture. The former is described in terms of a collision duration ⌧ c , roughly proportional to the inverse of the shear rate, i.e., ⌧ c˙ ⇡ 1 43 , the latter in terms of a drainage time ⌧ d which by scaling analysis can be shown to be proportional to the capillary number Ca of the droplets, i.e., ⌧ d˙ ⇡ Ca m where m = 4/3 if the drainage film is assumed to be flat 42 or m = 1 for a dimpled-film shape 44 . When the collision duration is larger than the drainage time ⌧ c > ⌧ d , droplets coalesce whereas in the opposite case ⌧ c < ⌧ d they repel. In the former condition the emulsion is often called attractive, while in the latter it is defined repulsive 45 . These ideas are based on several simplifying hypothesis, with the main ones being almost spherical droplets and head on collision; notwithstanding this, the theory proved able to estimate the general coalescence behaviour of emulsions and we will use it here as well. For all the viscosity ratio considered, the viscosity in the drainage film between two droplets is lower in the shear-thinning fluid and larger in the shear-thickening one than their Newtonian counterpart; because of this change, the capillary number Ca relevant for the collision is smaller/larger in the shear-thinning/thickening fluid than in the Newtonian case with a consequent decrease/increase of the drainage time ⌧ d . Thus, Chester's theory predicts an increase of the coalescence in shearthinning fluids and a decrease in the shear-thickening ones. Note however that, if one considers only the mean value µ m reported in figure 4, than this results would be opposite for the cases with < 1. To disentangle this opposite result, we quantify the droplets coalescence in our simulations by measuring the mean total surface area S in all the considered fluids, reported in figure 7(top). Note that, a reduction of surface area indicates droplet merging and thus a reduction of the number of droplets. As shown in the figure, the surface area decreases with the droplet volume fraction , thus indicating the merging of small droplets into large ones, and as expected is larger in the cases with unit viscosity ratio = 1 than in those with = 0.01 16 . When considering di↵erent flu- thinning fluids in the interstitial region between two interacting droplets affects the coalescence mechanism. Chesters 42 was the first to develop a theoretical framework to study the interaction of droplets. He proposed that the dynamics of two colliding droplets is an interplay of an external flow -the driving flow -responsible of the collision, and an internal flow -the drainage film between the two droplets -responsible of the interface deformation and rupture. The former is described in terms of a collision duration τ c , roughly proportional to the inverse of the shear rate, i.e., τ cγ ≈ 1 43 , the latter in terms of a drainage time τ d which by scaling analysis can be shown to be proportional to the capillary number Ca of the droplets, i.e., τ dγ ≈ Ca m where m = 4/3 if the drainage film is assumed to be flat 42 or m = 1 for a dimpled-film shape 44 . When the collision duration is larger than the drainage time τ c > τ d , droplets coalesce whereas in the opposite case τ c < τ d they repel. In the former condition the emulsion is often called attractive, while in the latter it is defined repulsive 45 . These ideas are based on several simplifying hypothesis, with the main ones being almost spherical droplets and head on collision; notwithstanding this, the theory proved able to estimate the general coalescence behaviour of emulsions and we will use it here as well. For all the viscosity ratio considered, the viscosity in the drainage film between two droplets is lower in the shear-thinning fluid and larger in the shear-thickening one than their Newtonian counterpart; because of this change, the capillary number Ca relevant for the collision is smaller/larger in the shear-thinning/thickening fluid than in the Newtonian case with a consequent decrease/increase of the drainage time τ d . Thus, Chester's theory predicts an increase of the coalescence in shearthinning fluids and a decrease in the shear-thickening ones. Note however that, if one considers only the mean value µ m reported in figure 4, than this results would be opposite for the cases with λ < 1. To disentangle this opposite result, we quantify the droplets coalescence in our simulations by measuring the mean total surface area S in all the considered fluids, reported in figure 7(top). Note that, a reduction of surface area indicates droplet merging and thus a reduction of the number of droplets. As shown in the figure, the surface area decreases with the droplet volume fraction Φ, thus indicating the merg-ing of small droplets into large ones, and as expected is larger in the cases with unit viscosity ratio λ = 1 than in those with λ = 0.01 16 . When considering different fluids, i.e., different power law index n, we clearly observe that when the fluid is shear-thickening the surface area S is always larger than in the Newtonian case, while on the contrary when the fluid is shear-thinning S is smaller than in the Newtonian case. This result is consistent with the prediction based on Chester's theory when using the local viscosity in the region in between droplets, i.e., the right tale of the pdf of the shear rate in figure 5(top), while it does not agree when considering only the mean shear (viscosity) value in figure 4. It is thus important to consider the full-spectrum of available shear rates in the domain when studying power-law fluids and theory and models uniquely based on mean values would result in wrong predictions. Figure 7(bottom) shows how the droplets redistribute across the domain for the different fluids analysed. Two representative volume fractions are considered, a dilute case with Φ = 0.1 and the most concentrated case studied with Φ = 0.3. In general, we observe that droplets tend to concentrate in the middle of the domain away from the walls. The migration of a droplet away from the wall in a shear flow was predicted theoretically by Magnaudet, Takagi, and Legendre 46 who provided an analytical expression for the the migration velocity; in our case, coalescence further enhance the concentration at the channel centerline in agreement with previous numerical studies in the literature 9,16,47,48 . Indeed, the concentration in the center of the domain is larger when the viscosity ratio λ is small, due to the tendency of droplets to merge and form large conglomerates, while this is reduced for large values of λ where the coalescence process is slowed down. In the latter case, we also find a reminiscence of the formation of wall-parallel layers as typically found in rigid particle suspensions which tend to form particle layers. When changing the carrier fluid property, we note that for shear-thickening fluids the concentration is always less peaked than in the shear-thinning case, thus further confirming the reduced level of coalescence and a more rigid-like behaviour of the droplets. The presence of the suspended droplets and their distribution in the channel affect the shear-rate, which differently from a single phase flow where it is constant, shows a strong inhomogeneity across the wall-normal direction, as reported in the inset of figure 7(bottom).
Finally, we analyse in more details how the results changes with the power index n. The top left panel of figure 8 reports the effective viscosity µ e as a function of n for the two viscosity ratio λ investigated and for a fixed volume fraction Φ = 0.2. The figure shows that µ e grows monotonically with n for both the viscosity ratio tested. Furthermore, the figure also suggests that the rate of change of the effective viscosity is approximately independent of λ. This is proved in the top right panel of the figure, where we report the effective viscosity µ e divided by the corresponding effective viscosity of the Newtonian case µ N ewt e (i.e., n = 1). The two curves obtained for different viscosity ratio λ well collapse into a single curve, whose expression can be obtained with a simple numerical fit as µ e = µ N ewt e 0.801 + 0.256n − 0.058n 2 .
The fit is also reported in the figure with a black line showing a good level of agreement with the results of our simulations. This result provides the important conclusion that the effective viscosity of an emulsion in a power law fluid can be well predicted with equation (10) knowing only the effective viscosity for the Newtonian configuration. When we study the mean values of the shear rateγ m and of the viscosity µ m the result is quite different. In particular, while the mean shear rateγ m still grows monotonically with n, its rate of change is enhanced by lower values of the viscosity ratio λ. This is a consequence of the enhanced coalesce of the droplets arising from the combined effect of the low viscosity ratio and of the shear-thinning property of the fluid. The mean viscosity µ m increases with n for the unit viscosity ratio case λ = 1 while it decreases for the case with λ = 0.01. This results further shows that the largest variations of mean viscosity are found in the case with low viscosity ratio λ while much smaller variations are evident as λ grows. We conclude our investigation by showing in figure 9 the first and second normal stress differences of the system, defined as From the figure we observe that the magnitude of both the normal stress differences grows with the flow index n. The first normal stress difference N 1 is greater then zero, i.e., σ 11 > σ 22 , and thus droplets are elongated in the direction of the flow and compressed in the wallnormal direction. On the other hand, the second normal difference N 2 is always negative. While N 1 substantially reduces for λ = 0.01, N 2 is mostly unaffected by the viscosity ratio.

IV. CONCLUSIONS
We study the rheology of a two-fluid emulsion in semiconcentrated conditions with the solute being a Newtonian fluid and the solvent an inelastic power law fluid. The problem is studied numerically by the volume of fluid method. Different volume fractions and viscosity ratios are considered together with several carrier fluids exhibiting both shear-thinning and thickening behaviours. We do not vary the capillary number which is fixed to a value such that a single droplet subject to the applied shear rate does not breakup, in order to focus our attention on the droplet coalescence mechanism.
First, we study how the effective viscosity of the system changes with the different carrier fluid properties. We show that the effective viscosity grows with the volume fraction of the dilute phase, reaches a maximum and then decreases with it beyond a certain concentration; this critical concentration corresponding to the maximum effective viscosity reduces when decreasing the viscosity ratio. The general behaviour remains mostly unaltered when considering power law fluids, except for an increase in effective viscosity for shear-thickening fluids and a decrease for the shear-thinning ones.
We show that the mean shear rate exhibit a behaviour similar to the effective viscosity, while the mean viscosity does not. Indeed, the mean viscosity is less than the nominal one for shear-thinning fluids and larger than the nominal one for shear-thickening fluids with unit viscosity ratios and the opposite when the viscosity ratio is reduced. This is consistent with the values of local shear rates measured in the domain but in contrast with the value of the effective viscosity. We explain this contradiction by studying the coalescence efficiency of the system which is modified by the nature of the carrier fluid. To support our claim, we show that the local shear rate assumes a very wide range of values; although its probability density function is peaked at the mean value, the distribution is strongly skewed exhibiting strong tails especially for large shear rates. This local high shear rate corresponds to the regions in between two interacting droplets, thus increases with the solute volume fraction. As a consequence, local high and low viscosity arises in these regions for shear-thickening and thinning fluids, respectively, which ultimately affect the coalescence of the droplets. Indeed, we relate the changes in the emulsion viscosity mainly to modifications of the coalescence in the system obtained by changing the carrier fluid property: coalescence is enhanced for shear-thinning fluids and reduced for shear-thickening ones due to modifications of the drainage time of the system, caused by modifications of the viscosity in the system. This process is mainly dominated by the power law fluid index n, with shearthinning and thickening fluids exhibiting two different behaviours which are (qualitatively) independent of the nominal viscosity ratio of the two fluids; also, this can not be understood by considering only the mean shear rate and viscosity of the two fluids across the domain.
Finally, we provide an expression able to successfully estimate the effective viscosity of the emulsion, as a function of power law index n and of the effective viscosity of the emulsion in a Newtonian fluid at the same volume fraction and viscosity ratio. This is a consequence of the fact that µ e grows with n with a rate of change which is approximately independent of λ.
To conclude, our results show that the coalescence efficiency which strongly affects the system rheology can be controlled by properly choosing the non-Newtonian property of the carrier fluids. Fully neglecting droplet merging can lead to erroneous or incomplete predictions in several flow conditions. In this work we have investigated a limiting cases where the nominal coalescence efficiency is close to unity; in a real scenario, the coalescence efficiency is likely to have intermediate values and future works may therefore be devoted to handling both coalescence and collisions in order to simulate a more wide range of emulsions. Anyway, the present choice helped to highlight the importance of coalescence.

Appendix: Verification of the numerical method and of the results
When two droplets approach each other, they may coalesce or not, with the coalescence happening without the addition of any non-hydrodynamic force due to the finite (although small) inertia. Topological changes of droplets such as coalescence and breakup are treated in the present methodology automatically without any additional model, since with the volume of fluid method these happen automatically when the distance between two interface falls within a grid cell. However, this does not mean that the coalescence process is correctly captured. On the other hand, fully resolving the coalesce is numerically challenging and not feasible, due to the very small scales involved, which are order of magnitude smaller than the characteristic size of the droplets. Notwithstanding the fact that coalescing events are here numerically driven, due to the finite size of our grid which will cause droplets to merge, we can verify that the present simulations still provides meaningful macroscopic results. To do that, we have proceeded as follows: first, we use the same setup employed in a previous work where the rheological data of a Newtonian droplet suspension have been compared with available experimental measurements 16 ; next, we have verified that our grid resolutions is fine enough to provide results which are weakly dependent on the grid size resolution, by evaluating their changes when doubling the grid size, and found variations less than 5%. Note that, while this ensure that the bulk results are accurate, it does not mean that we are fully capturing the very small dynamics corresponding to the actual drainage process. To prove this, figure 10 reports the time evolution of the total surface area for a case with Φ = 10% and viscosity ratio λ = 0.01. This test was chosen because the surface area is the quantity that is mostly affected by coalescence, being its changes directly related to it, and is thus a good indicator to show if coalescence is properly resolved; also the case chosen has a low viscosity ratio which is the one where coalescence is promoted, thus being the most demanding in terms of resolution. We can observe that after an initial transient phase where the surface area decreases, i.e., droplets are coalescing, the system reaches a statistically steady state where the surface area slightly oscillates around a mean value, in-dicated by the dashed line and extracted by time averaging over the second half of the time signal. The figure report the same quantity obtained using three different grid sizes, 8, 16 and 32 grid points per initial droplet radius R (corresponding to 80, 160 and 320 grid points per channel height 2h). We can observe that, although the time-histories are different, with the coarse grid leading to a faster coalescence than the finer ones, the steady state value obtained with the two finer grids are comparable, while a slightly lower value is obtained for the coarsest grid. The figure suggests that the grid used in the present simulations (16 grid points per initial droplet radius R) is appropriate to properly describe the macroscopic phenomenon under investigation.
We have also similarly tested that the results are independent of the size of the homogeneous directions x and z and that the wall-to-wall distance guarantees a low level of confinement similarly to what found in previous works 49,50 .