1. Introduction
Understanding the deformation and breakup of a droplet in a viscous flow plays an essential role in designing stable emulsions and foams for diverse industrial applications (Langevin Reference Langevin2000; Dan et al. Reference Dan, Gochev, Krägel, Aksenenko, Fainerman and Miller2013; van Kempen et al. Reference van Kempen, Schols, van der Linden and Sagis2013; Lam, Velikov & Velev Reference Lam, Velikov and Velev2014; Hamed, Schenck & Fiegel Reference Hamed, Schenck and Fiegel2020). Droplet breakup analysis can guide the generation of dispersed-phase droplets with a size distribution aptly tuned for chemical reactors, droplet microencapsulation, inkjet printing and spraying systems (Frankel & Acrivos Reference Frankel and Acrivos1970; Kennedy, Pozrikidis & Skalak Reference Kennedy, Pozrikidis and Skalak1994; Belyaeva et al. Reference Belyaeva, Della Valle, Neufeld and Poncelet2004; Lemenand et al. Reference Lemenand, Dupont, Della Valle and Peerhossaini2013; Daly et al. Reference Daly, Harrington, Martin and Hutchings2015; Narsimhan, Wang & Xiang Reference Narsimhan, Wang and Xiang2019; He et al. Reference He, Yang, Qin, Wen and Zhang2017). There is extensive research on the deformation and breakup of a clean droplet under different domains of external flow types. Barthes-Biesel & Acrivos (Reference Barthes-Biesel and Acrivos1973) first developed perturbation theories for examining the critical conditions for the breakup of a droplet under a general linear flow field in the limit of small droplet deformations. Some theoretical studies have examined the deformation and breakup of long slender droplets in shear and extensional flows where the internal fluid has a very low viscosity compared with the suspending fluid (Acrivos & Lo Reference Acrivos and Lo1978; Hinch & Acrivos Reference Hinch and Acrivos1979, Reference Hinch and Acrivos1980). Some studies have experimentally investigated the transient effects in the deformation and burst of viscous droplets (Bentley & Leal Reference Bentley and Leal1986; Stone, Bentley & Leal Reference Stone, Bentley and Leal1986). Rallison & Acrivos (Reference Rallison and Acrivos1978) and Rallison (Reference Rallison1981) have numerically explored the deformation and breakup of a viscous droplet in extensional and shear flows.
The role of simple surfactants in the deformation and breakup of droplets has also been widely investigated. Here, the role of surfactants is to change the surface tension of the interface as well as introduce Marangoni flows (Stone & Leal Reference Stone and Leal1990; Milliken, Stone & Leal Reference Milliken, Stone and Leal1993; Pawar & Stebe Reference Pawar and Stebe1996; Li & Pozrikidis Reference Li and Pozrikidis1997; Vlahovska, Bławzdziewicz & Loewenberg Reference Vlahovska, Bławzdziewicz and Loewenberg2009; Kamat et al. Reference Kamat, Wagoner, Thete and Basaran2018). Some studies have also focused on the deformation and burst of elastic capsules in extensional and simple shear flows using theories and experiments (Li, Barthes-Biesel & Helmy Reference Li, Barthes-Biesel and Helmy1988; Chang & Olbricht Reference Chang and Olbricht1993). The mechanical behaviour of artificial capsules and vesicles and their different modelling strategies are reviewed in Barthès-Biesel (Reference Barthès-Biesel2009).
Suspensions of droplets and bubbles laden with lipid bilayers, polymers, proteins and solid particles exhibit complex interfacial rheology. The interplay among surface viscoelastic effects, Marangoni stresses, surface pressure changes, surfactant solubility and diffusivity governs the mechanics of these complex interfaces (Jaensson, Anderson & Vermant Reference Jaensson, Anderson and Vermant2021; Pourali et al. Reference Pourali, Kröger, Vermant, Anderson and Jaensson2021). The role of interfacial rheology in droplet mechanics remains underexplored. Droplet interfaces embedded with certain low-molecular-weight surfactants (eicosanol, hexadecanol) and proteins (
$\beta$-casein) form highly viscous membranes (Brooks et al. Reference Brooks, Fuller, Frank and Robertson1999; Fuller & Vermant Reference Fuller and Vermant2012; Gunning et al. Reference Gunning, Kirby, Wilde, Penfold, Woodward and Morris2013; Verwijlen, Moldenaers & Vermant Reference Verwijlen, Moldenaers and Vermant2013). In such systems, the in-plane friction that arises when the molecules/particles within the interface slide past each other results in viscous dissipation, i.e. interfacial viscosities. The surface rheology of viscous interfaces can be modelled using the Boussinesq–Scriven law that provides an expression for the interfacial traction that arises from the surface shear and dilational viscosities (Boussinesq Reference Boussinesq1913; Scriven Reference Scriven1960). Erni (Reference Erni2011) reviews the influence of interfacial viscosity in the shear and dilational modes of deformation in systems with complex interfaces. Several experimental studies have also explored the different techniques to measure the interfacial shear and dilational rheological properties of complex fluid interfaces (Miller et al. Reference Miller, Ferri, Javadi, Krägel, Mucic and Wüstneck2010; Choi et al. Reference Choi, Steltenkamp, Zasadzinski and Squires2011; Fuller & Vermant Reference Fuller and Vermant2012; Balemans, Hulsen & Anderson Reference Balemans, Hulsen and Anderson2016; Jaensson et al. Reference Jaensson, Anderson and Vermant2021).
In recent years there have been advances in understanding the effect of surface rheology on droplet deformation and breakup using theories and numerical simulations. Flumerfelt (Reference Flumerfelt1980) first examined the impact of surface viscosities on droplet deformation in shear and extensional flows by extending the first-order perturbation analysis developed in Cox (Reference Cox1969) for a clean droplet. Narsimhan (Reference Narsimhan2019) developed second-order perturbation theories to describe the shape and rheology of a droplet with interfacial shear and dilational viscosity. The condition for the breakup of a droplet with surface viscosity in a linear flow field is presented in Singh & Narsimhan (Reference Singh and Narsimhan2020) in the limit of small droplet deformation. The deformation and breakup of a droplet in shear flow using boundary-element simulations have been explored in Gounley et al. (Reference Gounley, Boedec, Jaeger and Leonetti2016) and Luo, Shang & Bai (Reference Luo, Shang and Bai2019). Herrada et al. (Reference Herrada, Ponce-Torres, Rubio, Eggers and Montanero2022) also examined the steady-state solutions and stability of a gaseous droplet with interfacial viscosity using linear stability analysis and numerical computations in extensional flow. These previous studies have shown that the surface shear viscosity reduces the droplet deformation and increases the critical capillary number for droplet breakup compared with a clean droplet at the same viscosity ratio. In contrast, surface dilational viscosity enhances the droplet deformation and reduces the critical capillary number for droplet breakup compared with a clean droplet at the same viscosity ratio. In this work, we perform boundary-element simulations to investigate the role of surface viscosity in droplet deformation and breakup under axisymmetric extensional flow for a wide range of viscosity ratio values, interfacial viscosity and surfactant parameters. We present the mechanism behind the effect of surface shear and dilational viscosity on droplet deformation using traction arguments from interfacial viscosities. We also explore the combined impact of surface viscosity and surfactant transport on droplet relaxation in a quiescent fluid for different viscosity contrasts between the inner and outer fluid and surfactant transport parameters. We compare the results from our numerical simulations with those of second-order analytical theories in the limit of small capillary number (Narsimhan Reference Narsimhan2019; Singh & Narsimhan Reference Singh and Narsimhan2020). Numerical investigation of this problem allows us to examine stronger deformations, the effects of Marangoni stresses and the effects of pressure-thickening/thinning surface viscosities, the latter of which have not been addressed in previous studies.
2. Problem set-up and methodology
2.1. Problem overview
We numerically investigate how interfacial rheology alters the dynamics of a droplet suspended in an unbounded immiscible fluid and subject to an external flow. Figure 1 shows a schematic of the problem. We consider an initially spherical droplet of radius 
$R^{\prime }$ placed in a uniaxial extensional flow symmetric about the 
$z^{\prime }$ axis. We model the problem in polar-cylindrical coordinates with origin at the centre of the droplet. The internal fluid viscosity is 
$\lambda \eta ^{\prime }$, and the external fluid viscosity is 
$\eta ^{\prime }$. The interface of the droplet has an insoluble monolayer of a surface-active agent. We assume that the surfactant monolayer is primarily viscous. In this work, we neglect surface elasticity effects and employ the Boussinesq–Scriven constitutive relationship to describe the rheology of the viscous Newtonian interface (Boussinesq Reference Boussinesq1913; Scriven Reference Scriven1960). The surfactant interface has a surface shear viscosity 
$\eta _\mu ^{\prime }$, a surface dilational viscosity 
$\eta _\kappa ^{\prime }$ and a surface tension 
$\sigma ^{\prime }$. As the droplet deforms under external flow, the flow around the droplet can affect the local surface tension 
$\sigma ^{\prime }$, surfactant concentration 
$\varGamma ^{\prime }$ and surface viscosity (
$\eta _\kappa ^{\prime }$ and 
$\eta _\mu ^{\prime }$) on the interface. We discuss the formulations of variable surface tension and surface viscosity in the next subsection.
Figure 1. Problem overview.
The steady-state deformation of the droplet under flow is characterized by the Taylor deformation parameter 
$D_{Taylor}=(L^{\prime }-B^{\prime })/(L^{\prime }+B^{\prime })$, where 
$L^{\prime }$ and 
$B^{\prime }$ are the lengths of the major and minor axes, respectively, as shown in figure 1. The undisturbed far-field velocity vector 
$\boldsymbol {u}_{\boldsymbol {\infty }}^{\boldsymbol {\prime }}$ in dimensional form is expressed as
\begin{equation}
\boldsymbol{u}_{\boldsymbol{\infty}}^{\boldsymbol{\prime}}
= G^{\prime} \boldsymbol{\cdot} \begin{bmatrix} -1 & 0 &
0\\0 & -1 & 0\\0 & 0 & 2\end{bmatrix} \boldsymbol{\cdot}
\boldsymbol{x}^{\boldsymbol{\prime}}.
\end{equation}
In the above equation, 
$G^{\prime }$ represents the external strain rate and 
$\boldsymbol {x}^{\boldsymbol {\prime }}$ represents the position vector. Throughout the paper, primed variables are in dimensional form while unprimed variables are in dimensionless form.
2.2. Governing equations
We non-dimensionalize all lengths by radius 
$R^{\prime }$ of the initial spherical droplet, viscosities by the outer fluid viscosity 
$\eta ^{\prime }$, times by the inverse strain rate 
${G^{\prime }}^{-1}$, velocities by 
$G^{\prime } R^{\prime }$, bulk stresses by 
$\eta ^{\prime } G^{\prime }$ and surface stresses by 
$R^{\prime } \eta ^{\prime } G^{\prime }$. The interfacial surfactant concentration and surface tension are normalized by their initial equilibrium values 
$\varGamma _{eq}^{\prime }$ and 
$\sigma _{eq}^{\prime }$ at time 
$t=0$.
The fluid flow inside and outside of the droplet is described by the Stokes and continuity equations. The dimensionless forms of these equations are
$$\begin{gather} \lambda\nabla^{2}\boldsymbol{u}^{in}= \boldsymbol{\nabla} p^{in}, \quad \boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{u}^{in}=0, \end{gather}$$
$$\begin{gather}\nabla^{2}\boldsymbol{u}^{out}= \boldsymbol{\nabla} p^{out}, \quad \boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{u}^{out}=0. \end{gather}$$The boundary conditions at the interface are as follows.
(i) Continuity of velocity:
(2.3)where
\begin{equation} \boldsymbol{u}^{\boldsymbol{in}} = \boldsymbol{u}^{{\boldsymbol{out}}} = \boldsymbol{u}_{{\boldsymbol{s}}}, \end{equation}$\boldsymbol {u}_{\boldsymbol {s}}$ represents the velocity of the interface.(ii) Force balance:
(2.4)\begin{equation} (\boldsymbol{\tau}^{in} - \boldsymbol{\tau}^{out}) \boldsymbol{\cdot} \boldsymbol{n} = \boldsymbol{f}_{\boldsymbol{\mu}} + \boldsymbol{f}_{\boldsymbol{\kappa}} + \frac{1}{Ca} \left( \boldsymbol{\nabla}_{\boldsymbol{s}}\sigma - \sigma\boldsymbol{n}\boldsymbol{\nabla}\boldsymbol{\cdot}\boldsymbol{n} \right). \end{equation}
In (2.4), 
$\boldsymbol {\tau }^{in} - \boldsymbol {\tau }^{out}$ is the traction jump across the interface, 
$Ca=G^{\prime }\eta ^{\prime } R^{\prime }/\sigma _{eq}^{\prime }$ represents the capillary number, 
$\sigma$ is the dimensionless surface tension and 
$\boldsymbol {n}$ is the outward-pointing normal. The expressions for the interfacial shear traction 
$\boldsymbol {f}_{\boldsymbol {\mu }}$ and the interfacial dilational traction 
$\boldsymbol {f}_{\boldsymbol {\kappa }}$ in (2.4) are given by the Boussinesq–Scriven constitutive relationship (Boussinesq Reference Boussinesq1913; Scriven Reference Scriven1960):
$$\begin{gather} \boldsymbol{f}_{\boldsymbol{\kappa}} = \boldsymbol{\nabla}_{\boldsymbol{s}} \boldsymbol{\cdot}[Bq_\kappa\boldsymbol{\mathsf{P}}(\boldsymbol{\nabla}_{\boldsymbol{s}}\boldsymbol{\cdot} \boldsymbol{u}_{\boldsymbol{s}})], \end{gather}$$
$$\begin{gather}\boldsymbol{f}_{\boldsymbol{\mu}} ={-} \boldsymbol{\nabla}_{\boldsymbol{s}}\boldsymbol{\cdot}[Bq_\mu \boldsymbol{\mathsf{P}}(\boldsymbol{\nabla}_{\boldsymbol{s}}\boldsymbol{\cdot}\boldsymbol{u}_{\boldsymbol{s}}) - \boldsymbol{\mathsf{P}}\boldsymbol{\cdot}(\boldsymbol{\nabla}_{\boldsymbol{s}}\boldsymbol{u}_{\boldsymbol{s}} + {\boldsymbol{\nabla}_{\boldsymbol{s}}\boldsymbol{u}_{\boldsymbol{s}}}^{T})\boldsymbol{\cdot} \boldsymbol{\mathsf{P}}]. \end{gather}$$ In the above expressions, 
$\boldsymbol{\mathsf{P}} = \boldsymbol{\mathsf{I}}-\boldsymbol {nn}$ is the projection operator on the surface, 
$\boldsymbol {\nabla }_{\boldsymbol {s}} = \boldsymbol{\mathsf{P}} \boldsymbol {\cdot} \boldsymbol {\nabla }$ is the surface gradient, 
$Bq_{\mu }={\eta _\mu ^{\prime }}/{R^{\prime }\eta ^{\prime }}$ is the Boussinesq number for the surface shear viscosity and 
$Bq_{\kappa }={\eta _\kappa ^{\prime }}/{R^{\prime }\eta ^{\prime }}$ is the Boussinesq number for the surface dilational viscosity.
As an alternative way of quantifying the surface viscous effects, we introduce additional dimensionless parameters 
$Bq=Bq_\mu +Bq_\kappa$ and interfacial viscosity ratio 
$\lambda _{ds}={\eta _\kappa ^{\prime }}/{\eta _\mu ^{\prime }}$. The former represents the ratio of combined interfacial forces to the bulk viscous forces, while the latter represents the relative ratio of surface dilational to shear viscosity.
The changes in surfactant concentration 
$\varGamma$ along the droplet's surface are taken into account using the time-dependent convection–diffusion equation (Stone Reference Stone1990; Wong, Rumschitzki & Maldarelli Reference Wong, Rumschitzki and Maldarelli1996):
\begin{equation} \frac{\partial\varGamma}{\partial t} + \boldsymbol{\nabla}_{\boldsymbol{s}} \boldsymbol{\cdot}(\varGamma \boldsymbol{u}_{\boldsymbol{t}}) + \varGamma (\boldsymbol{\nabla}_s \boldsymbol{\cdot}\boldsymbol{n})(\boldsymbol{u}_{\boldsymbol{s}}\boldsymbol{\cdot} \boldsymbol{n}) = \frac{1}{Pe_s}\boldsymbol{\nabla}_s^{2}\varGamma, \end{equation}
where 
$Pe_s = {Ca R^{\prime }\sigma _{eq}^{\prime } }/{ D_s^{\prime }\eta ^{\prime }}$ is the surface Péclet number, 
$\boldsymbol {u}_{\boldsymbol {t}}$ is the tangential component of the interfacial velocity 
$\boldsymbol {u}_{\boldsymbol {s}}$ and 
$D_s^{\prime }$ is the surface diffusivity of the surfactant molecules. The surface Péclet number is the ratio of surface convection rate to the surface diffusion rate of the surfactant. In (2.7), the interpretation of the time derivative is such that the surface coordinates are fixed. In our numerical code, we displace the grid points in the direction normal to the interface. Therefore, the term 
$-\dot {\boldsymbol {x}} \boldsymbol {\cdot} \boldsymbol {\nabla }_s \varGamma$ (that arises because of the displacement of the grid points along the surface as the droplet deforms) is zero in the evolution equation given by Wong et al. (Reference Wong, Rumschitzki and Maldarelli1996).
We consider the Langmuir equation of state to relate the surface tension to the surfactant concentration at the interface:
\begin{equation} \sigma = \frac{\sigma_c^{\prime}}{\sigma_{eq}^{\prime}} + E \ln{\left(1-\frac{\varGamma}{\varGamma_\infty}\right)}. \end{equation} In the above equation, 
$\sigma _c^{\prime }$ is the surface tension of a clean interface without surfactant. The dimensionless parameter 
$E={R_G^{\prime } T^{\prime } \varGamma _\infty ^{\prime }}/{\sigma _{eq}^{\prime }}$ is the surface elasticity number, where 
$R_G^{\prime }$ is the ideal gas constant, 
$T^{\prime }$ is the absolute temperature and 
$\varGamma _\infty ^{\prime }$ is the maximum packing density for a given surfactant system.
At initial time 
$t=0$, we assume a homogeneous distribution of surfactant over the interface of the droplet (
$\varGamma =1$) and (2.8) can be rewritten as
\begin{equation} \sigma = 1 + E\ln{\left( \frac{\varGamma_\infty-\varGamma}{\varGamma_\infty-1}\right)}. \end{equation} In the case of non-homogeneous distribution of surfactant over the droplet's surface, the surface viscosities can vary strongly with surface pressure 
$\varPi ^{\prime } = \sigma _c^{\prime } - \sigma ^{\prime }$, as has been explored in several experimental studies (Kurtz, Lange & Fuller Reference Kurtz, Lange and Fuller2006; Kim et al. Reference Kim, Choi, Zasadzinski and Squires2011; Fuller & Vermant Reference Fuller and Vermant2012; Kim et al. Reference Kim, Choi, Zell, Squires and Zasadzinski2013; Samaniuk & Vermant Reference Samaniuk and Vermant2014). The dependence of surface viscosity on surface tension is taken into account using the following expressions (Manikantan & Squires Reference Manikantan and Squires2017):
$$\begin{gather} {\eta_{\kappa}^{\prime}}^{{\pm}} = {\eta_{\kappa, eq}^{\prime}} \exp\left({\pm\frac{\varPi^{\prime}-\varPi_{eq}^{\prime}}{\varPi_c^{\prime}}}\right), \end{gather}$$
$$\begin{gather}{\eta_{\mu}^{\prime}}^{{\pm}} = {\eta_{\mu, eq}^{\prime}}\exp\left({\pm \frac{\varPi^{\prime}-\varPi_{eq}^{\prime}}{\varPi_c^{\prime}}}\right), \end{gather}$$
where 
$+$ and 
$-$ signs represent a pressure-thickening and pressure-thinning surfactant, respectively; 
${\eta _{\mu, eq}^{\prime }}$ and 
${\eta _{\kappa, eq}^{\prime }}$ denote the equilibrium shear and dilational surface viscosity at the initial value of surface pressure 
$\varPi _{eq}^{\prime }$; and 
$\varPi _c^{\prime }$ is the pressure scale over which significant surface viscosity changes occur (Manikantan & Squires Reference Manikantan and Squires2017).
The dimensionless forms of (2.10) and (2.11) are
$$\begin{gather} Bq_{\kappa}^{{\pm}} = Bq_{\kappa, eq}\exp\left({\pm \frac{1-\sigma}{\varPi_c}}\right), \end{gather}$$
$$\begin{gather}Bq_{\mu}^{{\pm}} = Bq_{\mu, eq}\exp\left({\pm \frac{1-\sigma}{\varPi_c}}\right). \end{gather}$$
In the above equations, 
$\varPi _c=\varPi _c^{\prime }/\sigma _{eq}^{\prime }$ is the dimensionless surface pressure scale.
The dimensionless parameters are summarized in table 1. The ranges of surface parameters (
$Bq_\mu, Bq_\kappa, \varPi _c, E$) explored in this study are based on their experimental values reported in previous literature (Brooks et al. Reference Brooks, Fuller, Frank and Robertson1999; Kurtz et al. Reference Kurtz, Lange and Fuller2006; Georgieva et al. Reference Georgieva, Schmitt, Leal-Calderon and Langevin2009; Kim et al. Reference Kim, Choi, Zasadzinski and Squires2011; Fuller & Vermant Reference Fuller and Vermant2012; Kim et al. Reference Kim, Choi, Zell, Squires and Zasadzinski2013; Verwijlen et al. Reference Verwijlen, Moldenaers and Vermant2013; Hermans & Vermant Reference Hermans and Vermant2014; Samaniuk & Vermant Reference Samaniuk and Vermant2014; Zell et al. Reference Zell, Nowbahar, Mansard, Leal, Deshmukh, Mecca, Tucker and Squires2014; Manikantan & Squires Reference Manikantan and Squires2017; Singh & Narsimhan Reference Singh and Narsimhan2021). The surface Péclet number is based on surface diffusivity measurements listed in Shmyrov & Mizev (Reference Shmyrov and Mizev2019) that are collected from various sources. The surface diffusivity can vary widely depending on the phase behaviour of the adsorbed surfactant. A rule of thumb is that many macroscopic droplets (millimetre or above) will have 
$Pe_s \gg 1$ while those that are small (e.g. nanoemulsions or microemulsions) will have 
$Pe_s \sim O(1)$ or smaller.
Table 1. Dimensionless parameters.
2.3. Numerical implementation and validation
We implement axisymmetric boundary-element simulations to numerically compute the droplet's interfacial velocity. The velocity at the interfacial location 
$\boldsymbol {x}_{\boldsymbol {0}}$ can be expressed as
\begin{align} \boldsymbol{u}_{\boldsymbol{s}}(\boldsymbol{x}_{\boldsymbol{0}}) ={-}\frac{1}{4{\rm \pi}(1+\lambda)} \int_C \boldsymbol{\mathsf{M}}(\boldsymbol{x}, \boldsymbol{x}_{\boldsymbol{0}}) \boldsymbol{\cdot} [[\boldsymbol{\tau} \boldsymbol{\cdot}\boldsymbol{n}]] \,{\rm d}l(\boldsymbol{x}) + \frac{1}{4{\rm \pi}}\frac{1-\lambda}{1+\lambda}\int_C^{PV }\boldsymbol{\mathsf{q}}(\boldsymbol{x}, \boldsymbol{x}_{\boldsymbol{0}}) \boldsymbol{\cdot} \boldsymbol{u}_{\boldsymbol{s}}(\boldsymbol{x}) \,{\rm d}l(\boldsymbol{x}). \end{align}
In the above equation, 
$\boldsymbol{\mathsf{M}}$ and 
$\boldsymbol{\mathsf{q}}$ represent the single-layer and the double-layer potentials (Pozrikidis Reference Pozrikidis1990, Reference Pozrikidis1992), 
$[[\boldsymbol {\tau } \boldsymbol {\cdot}\boldsymbol {n}]]$ is the traction jump across the interface and 
${\rm d}l$ is the differential arclength along the droplet contour.
The droplet surface is discretized into 
$N$ elements connecting 
$N+1$ nodes. The detailed numerical implementation to solve the boundary-integral equation (2.14) and the convection–diffusion equation (2.7) can be found in the axisymmetric code developed in Singh & Narsimhan (Reference Singh and Narsimhan2021). A typical simulation starts with 
$N=100$ nodes on a spherical droplet interface. The time step 
$\Delta t$ is set as 
$(2\unicode{x2013}6)\times 10^{-3} Ca$ to update the droplet shape using the numerical scheme and the re-meshing procedure described in Singh & Narsimhan (Reference Singh and Narsimhan2021). To ensure that the total surfactant is conserved on the droplet's surface, we rescale the surfactant concentration at each time step. Therefore, the surfactant concentration integrated over the interface is constant over time. In the absence of such rescaling, we find that the error in total surfactant concentration integrated over the interface is less than 0.6 % for the simulations shown in the paper. The droplet is considered to have attained a steady-state deformation when the interfacial normal velocity 
$|\boldsymbol {u}_{\boldsymbol {n}}|<0.01$ at all collocation points. The critical capillary number 
$Ca_C$ denotes the largest value of 
$Ca$ below which the droplet attains a steady-state deformation under flow.
To validate our code, we compare our results with boundary-element simulations from Pawar & Stebe (Reference Pawar and Stebe1996). Figure 2 shows the variation of Taylor deformation parameter 
$D_{Taylor}$ with capillary number 
$Ca$ for a clean droplet with viscosity ratio 
$\lambda =1$. In the same plot, we also compare our numerical results with those of second-order perturbation theories developed in previous studies (Barthes-Biesel & Acrivos Reference Barthes-Biesel and Acrivos1973; Vlahovska et al. Reference Vlahovska, Bławzdziewicz and Loewenberg2009; Narsimhan Reference Narsimhan2019; Singh & Narsimhan Reference Singh and Narsimhan2020) and the lowest-order perturbation analysis from Flumerfelt (Reference Flumerfelt1980). The lowest-order perturbation theories solve the droplet shape up to 
${O}{(Ca)}$ for 
$Ca \ll 1$, while the second-order theories solve the shape up to 
${O} {(Ca^{2})}$ (Barthes-Biesel & Acrivos Reference Barthes-Biesel and Acrivos1973; Vlahovska et al. Reference Vlahovska, Bławzdziewicz and Loewenberg2009; Narsimhan Reference Narsimhan2019). At 
${O} {(Ca)}$, the analytical expression of the Taylor deformation parameter for a droplet with surface shear and dilational viscosity under uniaxial extensional flow is given as (Taylor Reference Taylor1934; Flumerfelt Reference Flumerfelt1980)
\begin{equation} D_{Taylor} = \frac{15}{32}\, \frac{19\lambda +16 + 24 Bq_\kappa + 8 Bq_\mu}{5\lambda + 6 Bq_\kappa + 4 Bq_\mu + 5} \, Ca . \end{equation}
Our results compare well with the numerical simulation from Pawar & Stebe (Reference Pawar and Stebe1996) and the second-order perturbation theories at all 
$Ca$ values. At 
$Ca>0.01$, we see that the lowest-order perturbation analysis shows deviations from the data predicted from simulations and underpredicts the 
$D_{Taylor}$ values.
Figure 2. Plot of 
$D_{Taylor}$ versus 
$Ca$ for a clean droplet without surfactant and viscosity ratio 
$\lambda =1.0$. The blue dots represent our numerical results and the red triangles represent boundary-element simulations from Pawar & Stebe (Reference Pawar and Stebe1996). The dashed blue curve is from second-order perturbation analysis (Barthes-Biesel & Acrivos Reference Barthes-Biesel and Acrivos1973; Vlahovska et al. Reference Vlahovska, Bławzdziewicz and Loewenberg2009; Singh & Narsimhan Reference Singh and Narsimhan2020), and the black solid curve is the lowest-order perturbation analysis from Flumerfelt (Reference Flumerfelt1980).
3. Results: droplet deformation and breakup
In §§ 3.1–3.3, we explore the effect of surface viscosity and surfactant transport on the deformation and breakup of a droplet placed in an extensional flow.
3.1. Comparison with small deformation perturbation theories for constant surface viscosity
In this section, we discuss the effect of constant surface viscosity on the steady-state deformation and the critical capillary number 
$Ca_C$ for droplet breakup. This analysis assumes that surface tension is constant and that surfactant is homogeneously distributed on the interface. Surfactant transport effects are neglected to analyse the isolated impact of interfacial viscosity on droplet dynamics. Sections 3.2 and 3.3 examine the coupling between interfacial viscosity, Marangoni stresses and surface dilution.
Figure 3 shows the variation of 
$D_{Taylor}$ with 
$Bq=Bq_\kappa +Bq_\mu$ for a droplet with interfacial viscosity ratio 
$\lambda _{ds}= Bq_{\kappa }/Bq_{\mu }=1$. We observe that at a given capillary number, the steady-state deformation 
$D_{Taylor}$ of the droplet decreases upon increasing 
$Bq$. Our numerical results compare well with those of the small deformation perturbation theories at 
${O} {(Ca^{2})}$ (Narsimhan Reference Narsimhan2019; Singh & Narsimhan Reference Singh and Narsimhan2020). Again, we see that the lowest-order perturbation analysis using (2.15) underpredicts the 
$D_{Taylor}$ values, and this effect is more significant at 
$Ca=0.05$ than at 
$Ca=0.01$.
Figure 3. Plots of 
$D_{Taylor}$ versus 
$Bq$ for a droplet with interfacial viscosity ratio 
$\lambda _{ds}=1$ at capillary number values 0.01 and 0.05. The dots represent our numerical results, dashed curves are from second-order perturbation theories (Narsimhan Reference Narsimhan2019; Singh & Narsimhan Reference Singh and Narsimhan2020) and bold curves are from using first-order perturbation theory (2.15).
Figure 4 shows the variation of 
$D_{Taylor}$ with interfacial viscosity ratio 
$\lambda _{ds}$ for a droplet with viscosity ratio 
$\lambda =1$, 
$Ca=0.05$ at 
$Bq$ values of 2, 5 and 10. Ratio 
$\lambda _{ds}=0$ represents a droplet with pure surface shear viscosity, and 
$\lambda _{ds}=\infty$ represents a droplet with pure surface dilational viscosity. The horizontal dashed line shows 
$D_{Taylor}$ for a clean droplet without surfactant at 
$Ca=0.05$ and 
$\lambda =1$. We observe that at 
$\lambda _{ds} \approx 2.0$, 
$D_{Taylor}$ of a droplet with surface viscosity is the same as that of a clean droplet for different 
$Bq$ values. Droplets with 
$\lambda _{ds}$ below this critical value exhibit smaller deformations than a clean droplet, and the deformation decreases as the interfacial effects (
$Bq$) increase. Droplets with 
$\lambda _{ds}$ above the critical value demonstrate the opposite effect – i.e. they exhibit larger deformations than clean droplets, and the deformation increases as the interfacial effects (
$Bq$) increase. We observe qualitatively similar behaviour if we examine droplets at different viscosity ratios 
$\lambda$. The critical 
$\lambda _{ds}$ where the droplet deformation is the same as that for a clean droplet shifts to 
$\lambda _{ds}=0.9$ for a droplet with viscosity ratio 
$\lambda =0.1$ (see figure 5), while the critical value shifts to 
$\lambda _{ds} = 40$ for a droplet with viscosity ratio 
$\lambda = 10$ (see figure 6).
Figure 4. Plots of 
$D_{Taylor}$ versus 
$\lambda _{ds}$ for a droplet with viscosity ratio 
$\lambda =1$, 
$Ca=0.05$ and 
$Bq$ values of 2, 5 and 10. Here 
$D_{Taylor}$ for a clean droplet without surfactant at 
$Ca=0.05$ and 
$\lambda =1$ is shown using a horizontal black dashed line.
Figure 5. Plots of 
$D_{Taylor}$ versus 
$\lambda _{ds}$ for a droplet with viscosity ratio 
$\lambda =0.1$, 
$Ca=0.05$ and at 
$Bq$ values of 2, 5 and 10. Here 
$D_{Taylor}$ for a clean droplet without surfactant at 
$Ca=0.05$ and 
$\lambda =0.1$ is shown using a horizontal black dashed line.
Figure 6. Plots of 
$D_{Taylor}$ versus 
$\lambda _{ds}$ for a droplet with viscosity ratio 
$\lambda =10$, 
$Ca=0.05$ and at 
$Bq$ values of 2, 5 and 10. Here 
$D_{Taylor}$ for a clean droplet without surfactant at 
$Ca=0.05$ and 
$\lambda =10$ is shown using a horizontal black dashed line.
The experimentally reported values of surface dilational viscosity are found to be much larger than surface shear viscosity as 
$\lambda _{ds}>{O} {(1)}$ for most interfacial monolayers. For instance, the interfacial viscosity ratio estimates of lipid monolayer dipalmitoylphosphatidylcholine and dimyristoylphosphatidylethanolamine are reported to be 
$\lambda _{ds} \sim {O} (10^{5})$ (Krägel et al. Reference Krägel, Kretzschmar, Li, Loglio, Miller and Möhwald1996; Kim et al. Reference Kim, Choi, Zasadzinski and Squires2011). Surface rheological studies of various interfacial monolayers have shown 
$\lambda _{ds} \sim {O} (10^{3})$ for poly(
$t$-butyl methacrylate) system (Krägel et al. Reference Krägel, Kretzschmar, Li, Loglio, Miller and Möhwald1996; Samaniuk & Vermant Reference Samaniuk and Vermant2014), 
$\lambda _{ds} \sim {O} (10^{4})$ for hexadecanol interface (Verwijlen et al. Reference Verwijlen, Moldenaers and Vermant2013; Samaniuk & Vermant Reference Samaniuk and Vermant2014), 
$\lambda _{ds} \sim {O} (1)$ for protein 
$\beta$-casein monolayer (Freer et al. Reference Freer, Yim, Fuller and Radke2004; Erni, Windhab & Fischer Reference Erni, Windhab and Fischer2011) and 
$\lambda _{ds} \sim {O} (10^{6})$ for viscoelastic globular protein ovalbumin interface (Erni et al. Reference Erni, Fischer, Windhab, Kusnezov, Stettin and Läuger2003; Xiong et al. Reference Xiong, Ren, Tian, Yang, Li and Li2018). The results shown in figures 4, 5 and 6 suggest that the combined influence of surface shear and dilational viscosity will increase the droplet deformation compared with a clean droplet at the same capillary number for most surfactant systems.
The next part of this section examines how interfacial viscosity alters the critical capillary number 
$Ca_C$ for drop breakup. The critical capillary number 
$Ca_C$ is defined as the largest value of 
$Ca$ below which the droplet will attain a steady-state deformation under flow (
$|\boldsymbol {u_n}|<0.01$ at all collocation points). Above this capillary number, the droplet will keep deforming and become unstable. Figure 7 compares the variation of 
$Ca_C$ with viscosity ratio 
$\lambda$ from the small deformation theories (represented by dashed lines) and numerical simulations (represented by bold lines). From figure 7(a), we can see that the surface dilational viscosity reduces 
$Ca_C$ of the droplet compared with a clean droplet at the same viscosity ratio 
$\lambda$. When 
$\lambda >4$, we observe that 
$Ca_C$ of a droplet with dilational viscosity is nearly the same as 
$Ca_C$ of a clean droplet from both theory and simulations. In contrast, we observe that the surface shear viscosity increases 
$Ca_C$ compared with a clean droplet at the same 
$\lambda$ (figure 7b). From figure 7(c), we can see the effect of equal surface shear and dilational viscosities on 
$Ca_C$. In this case, the 
$Ca_C$ values are lower than that of the clean droplet for viscosity ratio 
$\lambda <0.4$ and are higher than that of the clean droplet for viscosity ratio 
$\lambda >0.4$. The numerical results from boundary-element simulations are in good qualitative agreement with the small deformation theories for the three cases shown in figure 7. Previous studies examining clean droplets under extensional flow have also reported the 
$Ca_C$ values from analytical theories to be within 
$15\,\%$ agreement with the numerical results (Rallison & Acrivos Reference Rallison and Acrivos1978), which is what we observe here.
Figure 7. Critical capillary number 
$Ca_C$ versus viscosity ratio 
$\lambda$ for a clean droplet and a droplet with surface viscosity. The solid lines correspond to boundary-element simulations, while the dashed lines correspond to second-order perturbation theories (Narsimhan Reference Narsimhan2019; Singh & Narsimhan Reference Singh and Narsimhan2020). (a) Effect of surface dilational viscosity. (b) Effect of surface shear viscosity. (c) Effect of equal surface shear and dilational viscosity.
Figure 8 shows steady-state droplet shapes at 
$Ca=0.04$ for a droplet with pure surface shear viscosity 
$Bq_\mu =2, Bq_\kappa =0$ (figure 8a,c) and a droplet with pure surface dilational viscosity (figure 8b,d) at viscosity ratio values 
$\lambda = 0.1$ and 
$\lambda = 10$. The vectors shown in the droplet's upper half represent the computed traction from interfacial viscosity (2.5) and (2.6). As can be seen from figure 8, the droplet shapes from boundary-element simulations (drawn using dashed lines) match well with the droplet shapes from second-order perturbation theories (drawn using red dots). For a clean droplet, 
$D_{Taylor}=0.137$ at 
$\lambda =0.1$ and 
$D_{Taylor}=0.187$ at 
$\lambda =10$. In the case of a droplet with pure surface shear viscosity (figure 8a,c), the shear traction opposes the surface velocity gradients and reduces the droplet deformation compared with the clean droplet. For a droplet with pure dilational viscosity, the expression of dilational traction (2.5) mimics the formulation of variable surface tension. The influence of surface dilational viscosity can be seen as an ‘effective’ surface tension (Singh & Narsimhan Reference Singh and Narsimhan2021):
\begin{equation} \sigma^{effective} = \sigma + Bq_\kappa(\boldsymbol{\nabla}_{\boldsymbol{s}}\boldsymbol{\cdot} \boldsymbol{u}_{\boldsymbol{s}}). \end{equation}
The negative divergence of velocity in the droplet's pole reduces the effective surface tension in the pole region, and we observe that the droplet with surface dilational viscosity is more deformed than a clean droplet at the same 
$Ca$ (figure 8b,d). The effect of surface dilational viscosity appears similar to the influence of surfactant convection. When a droplet is placed under extensional flow, the surfactant gets swept towards the droplet's pole. The higher surfactant concentration in the pole region reduces the effective surface tension and increases the droplet deformation (Pawar & Stebe Reference Pawar and Stebe1996; Feigl et al. Reference Feigl, Megias-Alguacil, Fischer and Windhab2007).
Figure 8. Steady-state droplet shapes at 
$Ca=0.04$ for four different cases: (a) 
$Bq_\mu =2, Bq_\kappa =0$ and 
$\lambda =0.1$; (b) 
$Bq_\mu =0, Bq_\kappa =2$ and 
$\lambda =0.1$; (c) 
$Bq_\mu =2, Bq_\kappa =0$ and 
$\lambda =10$; (d) 
$Bq_\mu =0, Bq_\kappa =2$ and 
$\lambda =10$. Shapes drawn using red dots are from second-order perturbation theories (Narsimhan Reference Narsimhan2019) and shapes drawn with dashed lines are from boundary-element simulations. The vectors shown in the droplet's upper half represent the computed traction from interfacial viscosity (2.5) and (2.6). For a clean droplet, 
$D_{Taylor}=0.137$ at 
$\lambda =0.1$ and 
$D_{Taylor}=0.187$ at 
$\lambda =10$.
From figure 8, we can see that the interfacial traction acting on the droplet's surface is more significant at 
$\lambda =0.1$ than at 
$\lambda =10$ for both surface shear and dilational viscosity (see scale for heat map). As a result, the surface viscosity influences the droplet deformation and 
$Ca_C$ more strongly at lower values of viscosity ratio than higher values, as can be seen from figure 7.
3.2. Combined influence of constant surface viscosity and surfactant transport
As the droplet deforms under extensional flow, the surfactant can get swept from the droplet's equator towards the droplet's poles leading to non-homogeneous surfactant distribution. In this section, we include the effect of Marangoni stresses on droplet deformation and breakup, assuming the surface viscosity remains constant on the interface.
Figure 9 shows the variation of 
$Ca_C$ with 
$\lambda _{ds}$ for a droplet with 
$Bq=2$ at surface Péclet number values 
$0.01$, 
$1$, 
$10$ and 
$1000$. The values of dimensionless parameters are: 
$\lambda =1$, 
$E=0.2$ and 
$\varGamma _\infty ^{-1}=0.5$. The dashed horizontal lines in the figure denote 
$Ca_C$ of a droplet without surface viscosity, i.e. 
$Bq=0$. We observe that the critical capillary number 
$Ca_C$ decreases upon increasing 
$\lambda _{ds}$ and approaches 
$Ca_C$ of a droplet with pure surface dilational viscosity at 
$\lambda _{ds}=\infty$. The surface viscosity does not significantly impact droplet breakup at higher values of surface Péclet number (
$Pe_s=10$ and 
$Pe_s=1000$). We find the 
$Ca_C$ values for a droplet with interfacial viscosity to be the same as that without interfacial viscosity at 
$Pe_s=1000$.
Figure 9. Variation of 
$Ca_C$ with 
$\lambda _{ds}$ for a droplet with 
$Bq=2$ at 
$Pe_s$ values 0.01, 1, 10 and 1000. The dashed horizontal lines denote 
$Ca_C$ of the droplet in the absence of surface viscosity, i.e. 
$Bq=0$. The dimensionless parameters are: 
$\lambda =1$, 
$E=0.2$ and 
$\varGamma _\infty ^{-1}=0.5$.
Figure 10(a) shows the droplet shapes at 
$Pe_s=0.01$ and 
$Pe_s=1000$ for a droplet with interfacial viscosity 
$Bq=2$ and 
$\lambda _{ds}=1$. The values of dimensionless parameters are: 
$Ca=0.05$, 
$\lambda =1$, 
$E=0.4$ and 
$\varGamma _\infty ^{-1}=0.5$. The droplet at surface Péclet number 
$Pe=0.01$ appears less deformed than the droplet at 
$Pe=1000$. Figure 10(b) shows the variation of tangential interfacial traction and Marangoni forces acting along the droplet's surface with normalized arclength 
$s$ for the shapes shown in figure 10(a). At 
$Pe=0.01$, we observe that the Marangoni effects are negligible compared with the surface viscosity effects. In contrast, at 
$Pe=1000$, the traction contribution from interfacial viscosity is much smaller than the Marangoni forces. As a result, the interfacial viscosity has a negligible impact on droplet deformation and breakup at large 
$Pe_s$, and the 
$Ca_C$ values approach the results for a droplet with only surfactant transport (as was observed in figure 9).
Figure 10. (a) Steady-state droplet shapes at 
$Pe_s=0.01$ and 
$Pe_s=1000$ for a droplet with 
$Bq=\!2$ and 
$\lambda _{ds}=\!1$. (b) Variation of tangential component of interfacial traction (shown by bold curves) and Marangoni traction (shown by dashed curves) with normalized arclength 
$s$ for the droplet shapes shown in (a). The values of dimensionless parameters are: 
$Ca=0.050$, 
$\lambda =1$, 
$E=0.4$ and 
$\varGamma _\infty ^{-1}=0.5$.
Figure 11 shows the variation of 
$Ca_C$ with 
$Bq$ for a droplet with equal surface shear and dilational viscosity 
$\lambda _{ds}=1$. The values of dimensionless parameters are: 
$\lambda =1$, 
$E=0.2$ and 
$\varGamma _\infty ^{-1}=0.5$. Similar to the case of constant surface viscosity (figure 7c), for a droplet with surfactant transport, we observe that the droplet with equal interfacial viscosities has a higher 
$Ca_C$ than that of a droplet without interfacial viscosity. We observe that increasing the interfacial viscosity has a higher impact on 
$Ca_C$ at 
$Pe_s=0.01$ than at 
$Pe_s=10$. The interfacial viscosity has a negligible impact on 
$Ca_C$ at 
$Pe_s=1000$.
Figure 11. Variation of 
$Ca_C$ with 
$Bq=Bq_\mu +Bq_\kappa$ for a droplet with interfacial viscosity ratio 
$\lambda _{ds}=1$ at surface Péclet number values 
$Pe_s=0.01$, 
$1$, 
$10$ and 
$1000$. The values of dimensionless parameters are: 
$\lambda =1$, 
$E=0.2$ and 
$\varGamma _\infty =2$.
3.3. Influence of pressure-dependent surface viscosity
Previous literature has shown that the surface viscosity can change sharply with surface pressure (Kurtz et al. Reference Kurtz, Lange and Fuller2006; Kim et al. Reference Kim, Choi, Zasadzinski and Squires2011; Fuller & Vermant Reference Fuller and Vermant2012; Kim et al. Reference Kim, Choi, Zell, Squires and Zasadzinski2013; Samaniuk & Vermant Reference Samaniuk and Vermant2014). Here we discuss how the pressure-dependent surface shear and dilational viscosity enhance/reduce the droplet deformation and critical capillary number compared with the case of 
$\varPi$-independent surface viscosity. We include the effect of surfactant transport as in § 3.2, but now also let the surface viscosity vary with surfactant concentration.
We first discuss the impact of pressure-dependent surface shear viscosity on droplet deformation and breakup. Figure 12 shows the evolution of droplet shapes with time for a pressure-thinning surfactant with 
$\varPi _c=0.25$ (shown in figure 12a) and a pressure-thickening surfactant with 
$\varPi _c=0.25$ (shown in figure 12b) at capillary number 
$Ca=0.08$. The values of dimensionless parameters are: 
$\lambda =1$, 
$Pe_s=0.01$, 
$E=0.4$, 
$\varGamma _\infty =2$, 
$Bq_{\mu, {eq}}=2$ and 
$Bq_{\kappa,{eq}}=0$. The droplet with 
$\varPi$-thinning surfactant eventually attains a stable shape under flow at time 
$t=5.4$. In contrast, at the same capillary number, the droplet with 
$\varPi$-thickening surfactant deforms continuously with time and becomes unstable. In figure 13, we plot the variation of 
$Bq_\mu$ and surfactant concentration 
$\varGamma$ with normalized arclength 
$s$ (
$s=0$ and 
$s=1$ at the droplet's two poles) for the droplet shapes shown in figures 12(a) and 12(b) at 
$t=5.4$. At 
$Pe=0.01$, the diffusion of the surfactant along the droplet's interface dominates over the surfactant convection. As a result, we observe a homogeneous distribution of surfactant on the droplet's surface and 
$\varGamma <1$ over the entire droplet interface for 
$\varPi$-thickening and 
$\varPi$-thinning surfactant. The lowered surfactant concentration over the droplet's surface also lowers the droplet's 
$Bq_\mu$ for a 
$\varPi$-thickening surfactant and increases the droplet's 
$Bq_\mu$ for a 
$\varPi$-thinning surfactant. As a result, the droplet with 
$\varPi$-thinning surfactant breaks at a higher 
$Ca$ than the droplet with 
$\varPi$-thickening surfactant. Tables 2 and 3 compare 
$Ca_C$ of a droplet with a 
$\varPi$-thickening, 
$\varPi$-thinning and 
$\varPi$-independent surfactant at surface elasticity numbers 
$E=0.4$ and 
$E=2$, respectively.
Figure 12. Time evolution of a droplet with pressure-dependent interfacial shear viscosity 
$Bq_{\mu,{eq}}=2.0$: (a) 
$\varPi _c=0.25$, pressure-thinning; (b) 
$\varPi _c=0.25$, pressure-thickening. The values of dimensionless parameters are: 
$\lambda =1$, 
$E=0.4$, 
$\varGamma _\infty =2$, 
$Ca = 0.08$ and 
$Pe_s=0.01$.
Figure 13. Variation of 
$Bq_\mu$ and surface concentration 
$\varGamma$ with the normalized arclength 
$s$ corresponding to the droplet shapes shown in figure 12 at time 
$t = 5.4$.
Table 2. Values of 
$Ca_C$ of a droplet with 
$Bq_{eq}=2$, 
$\lambda =1.0$, 
$E=0.4$, 
$\varGamma _\infty =2$ and 
$Pe_s=0.01$.
Table 3. Values of 
$Ca_C$ of a droplet with 
$Bq_{eq}=2$, 
$\lambda =1.0$, 
$E=2$, 
$\varGamma _\infty =2$ and 
$Pe_s=0.01$.
Figure 14 shows the evolution of droplet shapes with time for a pressure-thinning surfactant with 
$\varPi _c=0.25$ (shown in figure 14a) and a pressure-thickening surfactant with 
$\varPi _c=0.25$ (shown in figure 14b) at 
$Ca=0.056$, but now the Péclet number is large (
$Pe_s = 10$). Here also we observe that the droplet with 
$\varPi$-thinning surfactant attains a stable shape under flow while the droplet with 
$\varPi$-thickening surfactant becomes unstable. Figure 15 shows the variation of 
$Bq_\mu$ and surfactant concentration 
$\varGamma$ with normalized arclength 
$s$ for the droplet shapes shown in figures 14(a) and 14(b) at 
$t=8.4$. At 
$Pe_s=10$, we observe that the surfactant gets convected towards the droplet's pole leading to 
$\varGamma >1$ in the pole region and 
$\varGamma <1$ in the equator region. The lowered surfactant concentration in the equator region for a droplet with 
$\varPi$-thinning surfactant increases the droplet's 
$Bq_\mu$ and surface shear traction acting on the droplet's equator. For a droplet with only surface shear viscosity, what matters is traction near the droplet's equator and not traction near the droplet's poles, as was discussed for figure 8. As a result, the droplet with pressure-thinning surfactant breaks at a larger 
$Ca$ than the droplet with pressure-thickening surfactant (
$Ca_C$ of a droplet with pressure-thinning surfactant is 0.056, while 
$Ca_C$ of a droplet with pressure-thickening surfactant is 0.054). We observe similar behaviour of pressure-dependent surface viscosity at a higher value of surface elasticity number. At 
$E=2$, 
$Ca_C$ of a droplet with pressure-thinning surfactant is 0.067, while 
$Ca_C$ of a droplet with pressure-thickening surfactant is 0.064.
Figure 14. Time evolution of a droplet with pressure-dependent interfacial shear viscosity 
$Bq_{\mu,{eq}}=2.0$: (a) 
$\varPi _c=0.25$, pressure-thinning; (b) 
$\varPi _c=0.25$, pressure-thickening. The values of dimensionless parameters are: 
$\lambda =1$, 
$E=0.4$, 
$\varGamma _\infty =2$, 
$Ca = 0.056$ and 
$Pe_s=10.0$.
Figure 15. Variation of 
$Bq_\mu$ and surface concentration 
$\varGamma$ with the normalized arclength 
$s$ corresponding to the droplet shapes shown in figure 14 at time 
$t = 8.4$.
Next we discuss the impact of a 
$\varPi$-thickening and 
$\varPi$-thinning surfactant on a droplet with only surface dilational viscosity. Figure 16 shows the droplet shape evolution for a pressure-thinning and pressure-thickening surfactant (figure 16a,b). Unlike surface shear viscosity, here we see that the droplet with 
$\varPi$-thickening surfactant attains a stable shape and the droplet with 
$\varPi$-thinning surfactant becomes unstable. Here, for a droplet with 
$\varPi$-thickening surfactant, the lowered surfactant concentration due to strong diffusive effects lowers the droplet's 
$Bq_\kappa$ and increases 
$Ca_C$ compared with a 
$\varPi$-thinning surfactant. Table 2 compares 
$Ca_C$ of a droplet with a 
$\varPi$-thickening, 
$\varPi$-thinning and 
$\varPi$-independent surfactant at 
$Pe_s=0.01$. From table 2 we can see that the effects of 
$\varPi$-thickening and 
$\varPi$-thinning surfactant are not very pronounced for a droplet with only dilational viscosity compared to a droplet with only surface shear viscosity.
Figure 16. Time evolution of a droplet with pressure-dependent interfacial dilational viscosity 
$Bq_{\kappa,{eq}}=2.0$: (a) 
$\varPi _c=0.25$, pressure-thinning; (b) 
$\varPi _c=0.25$, pressure-thickening. The values of dimensionless parameters are: 
$\lambda =1$, 
$E=0.4$, 
$\varGamma _\infty =2$, 
$Ca = 0.059$ and 
$Pe_s=0.01$.
Figure 17. Variation of 
$Bq_\kappa$ and surface concentration 
$\varGamma$ with the normalized arclength 
$s$ corresponding to the droplet shapes shown in figure 16 at time 
$t = 6.7$.
In the case of a droplet with only surface dilational viscosity at 
$Pe_s=10$, we do not observe a difference in 
$Ca_C$ for 
$\varPi$-thickening and 
$\varPi$-thinning surfactant. The convection of surfactant towards the droplet's pole at 
$Pe_s=10$ increases/decreases the droplet's 
$Bq_\kappa$ in the pole region for 
$\varPi$-thickening/thinning surfactant, and we observe the droplet with 
$\varPi$-thickening surfactant is slightly more deformed compared with the droplet with 
$\varPi$-thinning surfactant. Because we do not observe a significant difference in deformation at large 
$Pe_s$ for a droplet with only surface dilational viscosity, the results for the 
$\varPi$-thickening and 
$\varPi$-thinning surfactant for this case are omitted in the paper.
4. Results: droplet relaxation in quiescent fluid
In this section, we discuss the combined influence of surface viscosity and surfactant transport on the relaxation of an initially extended droplet in a quiescent flow. For the results presented in this section, a spherical droplet of radius 
$R^{\prime }$ at a given set of non-dimensional parameters is first stretched under extensional flow to a specified initial non-dimensional length 
$L=L^{\prime }/R^{\prime }$, where 
$L^{\prime }$ is the length of the droplet's major axis. The flow is then stopped to examine the relaxation of the extended droplet in the now quiescent external fluid.
In the absence of external flow, we non-dimensionalize all lengths by 
$R^{\prime }$, viscosities by the outer fluid viscosity 
$\eta ^{\prime }$, velocities by 
$\sigma _{eq}^{\prime }/\eta ^{\prime }$, time by 
$R^{\prime } \eta ^{\prime }/ \sigma _{eq}^{\prime }$, bulk stresses by 
$\sigma _{eq}^{\prime }/R^{\prime }$ and surface stresses by 
$\sigma _{eq}^{\prime }$. The interfacial surfactant concentration and surface tension are normalized by their equilibrium values 
$\varGamma _{eq}^{\prime }$ and 
$\sigma _{eq}^{\prime }$ on a spherical droplet of radius 
$R^{\prime }$ placed under no external flow.
The non-dimensionalized time-dependent convection–diffusion equation in the absence of flow can be written as
\begin{equation} \frac{\partial\varGamma}{\partial t} + \boldsymbol{\nabla}_{\boldsymbol{s}} \boldsymbol{\cdot}(\varGamma \boldsymbol{u}_{\boldsymbol{t}}) + \varGamma (\boldsymbol{\nabla}_s \boldsymbol{\cdot}\boldsymbol{n})(\boldsymbol{u}_{\boldsymbol{s}}\boldsymbol{\cdot} \boldsymbol{n}) = \frac{1}{\gamma}\boldsymbol{\nabla}_s^{2}\varGamma. \end{equation}
The non-dimensional parameter 
$\gamma = {R^{\prime }\sigma _{eq}^{\prime } }/{ D_s^{\prime }\eta ^{\prime }}$ in the above equation denotes the ratio of diffusive time scale 
${R^{\prime }}^{2}/D_s^{\prime }$ to capillary time scale 
$R^{\prime } \eta ^{\prime }/ \sigma _{eq}^{\prime }$.
Figure 18 shows the droplet relaxation of an initially extended droplet with initial length 
$L=16$ for three cases: droplet with only surface shear viscosity (
$Bq_\mu =1$, 
$Bq_\kappa =0$), clean droplet (
$Bq_\mu =0$, 
$Bq_\kappa =0$) and droplet with only surface dilational viscosity (
$Bq_\mu =0$, 
$Bq_\kappa =1$). The values of dimensionless parameters are: 
$\lambda =1$, 
$E=0.2$, 
$\varGamma _\infty =2$ and 
$\gamma =1.0$. As we can see, both clean droplet and droplet with surface dilational viscosity eventually break up by an end-pinching mechanism. In the case of the droplet with surface shear viscosity, the droplet becomes stable and eventually returns to a spherical shape.
Figure 18. Relaxation of an initially extended droplet in a stagnant fluid for three cases having the same initial shape (
$L=16.0$) over time: (a) droplet with surface shear viscosity (
$Bq_\mu =1$, 
$Bq_\kappa =0$), (b) clean droplet (
$Bq_\mu =0$, 
$Bq_\kappa =0$) and (c) droplet with surface dilational viscosity (
$Bq_\mu =0$, 
$Bq_\kappa =1$). The values of dimensionless parameters are: 
$\lambda =1$, 
$E=0.2$, 
$\varGamma _\infty =2$ and 
$\gamma =1.0$.
In the next example shown in figure 19, we keep all the parameters the same as in figure 18 but increase the 
$E$ from 
$0.2$ to 
$2$. At 
$E=2$, the higher sensitivity of surface tension to surfactant concentration enhances the Marangoni stresses and reduces the effective surface tension. As a result of these two effects, the initial shapes at 
$E=2$ have less pronounced bulging ends than the initial droplet shapes shown at 
$E=0.2$ in figure 18. Here also, we observe that the droplet with surface shear viscosity eventually returns to a spherical shape. In the case of clean droplet and droplet with only surface dilational viscosity, the droplet breaks into two smaller droplets.
Figure 19. Relaxation of an initially extended droplet in a stagnant fluid for three cases having the same initial shape (
$L=16.0$) over time: (a) droplet with surface shear viscosity (
$Bq_\mu =1$, 
$Bq_\kappa =0$), (b) clean droplet (
$Bq_\mu =0$, 
$Bq_\kappa =0$) and (c) droplet with surface dilational viscosity (
$Bq_\mu =0$, 
$Bq_\kappa =1$). The values of dimensionless parameters are: 
$\lambda =1$, 
$E=2$, 
$\varGamma _\infty =2$ and 
$\gamma =1$.
In figure 20, we examine the relaxation of a droplet with very low viscosity ratio 
$\lambda =0.01$ and initial extension 
$L=10$ for three cases: a clean droplet (
$Bq_\mu =0$, 
$Bq_\kappa =0$), a droplet with equal surface viscosity (
$Bq=0.4$, 
$\lambda _{ds}=1$) and a droplet with only surface shear viscosity (
$Bq=0.4$, 
$\lambda _{ds}=0$). The values of dimensionless parameters are: 
$E=2$, 
$\varGamma _\infty =2$ and 
$\gamma =1.0$. We observe that the clean droplet becomes unstable with the eventual breaking of the thin filaments from the droplet's poles. In the case of a droplet with surface viscosity (figure 20b,c), the middle portion of the droplet gradually reduces to an almost spherical shape. The droplet's poles form very small droplets that appear to pinch off from the parent droplet at times 
$t=3.9$ and 
$t=6.2$, respectively. We observe that the droplet with surface viscosity has a higher breakup time than the clean droplet. Figure 21 shows droplet shapes at different values of 
$\lambda _{ds}$ having the same axial extension 
$L=10$. The values of dimensionless parameters are: 
$\lambda =0.01$, 
$Bq=0.4$, 
$E=2$, 
$\varGamma _\infty =2$ and 
$\gamma =1$. We see that droplets with a higher value of interfacial viscosity ratio 
$\lambda _{ds}$ have more pronounced pointed ends than droplets with a lower value of 
$\lambda _{ds}$. At higher values of interfacial viscosity ratio 
$\lambda _{ds}>7$, the droplet breakup under relaxation will be similar to that of a clean droplet.
Figure 20. Relaxation of an initially extended droplet in a stagnant fluid for three cases having the same initial shape (
$L=10.0$) over time: (a) clean droplet (
$Bq_\mu =0$, 
$Bq_\kappa =0$), (b) droplet with equal surface viscosity (
$Bq_\mu =0.2$, 
$Bq_\kappa =0.2$) and (c) droplet with only surface shear viscosity (
$Bq_\mu =0.4$, 
$Bq_\kappa =0.0$). The values of dimensionless parameters are: 
$\lambda =0.01$, 
$E=2$, 
$\varGamma _\infty =2$ and 
$\gamma =1$.
Figure 21. Droplet shapes at different values of 
$\lambda _{ds}$ having the same axial extension 
$L=10$. The values of dimensionless parameters are: 
$\lambda =0.01$, 
$Bq=0.4$, 
$E=2$, 
$\varGamma _\infty =2$ and 
$\gamma =1$.
In figure 22, we examine the relaxation of a droplet with viscosity ratio 
$\lambda =0.05$ and initial extension 
$L=10$ for four cases: droplet with surface shear viscosity (
$\lambda _{ds}=0$) and 
$E=0.4$, droplet with equal surface viscosity (
$\lambda _{ds}=1$) and 
$E=0.4$, droplet with surface shear viscosity (
$\lambda _{ds}=0$) and 
$E=2.0$ and droplet with equal surface viscosity (
$\lambda _{ds}=1$) and 
$E=2.0$. The values of dimensionless parameters are: 
$\lambda =0.05$, 
$Bq=0.4$, 
$\varGamma _\infty =2$ and 
$\gamma =100$. At 
$E=0.4$, the droplet with only surface shear viscosity (
$\lambda _{ds}=0$) and the droplet with equal viscosity (
$\lambda _{ds}=1$) form very small droplets at the poles that appear to pinch off from the parent droplet at times 
$t=14$ and 
$t=11.2$, respectively. At 
$E=2$, the droplet with only surface shear viscosity gradually relaxes into a spherical shape while the droplet with equal surface viscosity breaks and forms daughter droplets at the poles.
Figure 22. Relaxation of an initially extended droplet in a stagnant fluid for four cases having the same initial shape (
$L=10.0$) over time: 
$(a)$ droplet with surface shear viscosity (
$\lambda _{ds}=0$) and 
$E=0.4$, 
$(b)$ droplet with equal surface viscosity (
$\lambda _{ds}=1$) and 
$E=0.4$, 
$(c)$ droplet with surface shear viscosity (
$\lambda _{ds}=0$) and 
$E=2.0$ and 
$(d)$ droplet with equal surface viscosity (
$\lambda _{ds}=1$) and 
$E=2.0$. The values of dimensionless parameters are: 
$\lambda =0.05$, 
$Bq=0.4$, 
$\varGamma _\infty =2$ and 
$\gamma =100$.
Figure 23 shows droplet shapes at different combinations of 
$\lambda _{ds}$ and 
$\gamma$ with the same axial extension 
$L=10$. The values of dimensionless parameters are: 
$\lambda =10$, 
$Bq=0.4$, 
$E=0.4$ and 
$\varGamma _\infty =2$. At a higher value of viscosity ratio 
$\lambda =10$, the droplet shapes at different values of 
$\lambda _{ds}$ are similar for the same initial axial extension and 
$\gamma$ value. Droplet shapes at 
$\gamma =100$ form pointed ends due to stronger convection of surfactant towards the droplet's ends. We see that irrespective of the interfacial viscosity ratio 
$\lambda _{ds}$ and 
$\gamma$ values, the droplet eventually returns to a spherical shape. The relaxation of a droplet with 
$\lambda _{ds}=\infty$ and 
$\gamma =100$ is shown in figure 24.
Figure 23. Droplet shapes at different combinations of 
$\lambda _{ds}$ and 
$\gamma$ with the same axial extension 
$L=10$. The values of dimensionless parameters are: 
$\lambda =10$, 
$Bq=0.4$, 
$E=0.4$ and 
$\varGamma _\infty =2$.
Figure 24. Relaxation of an initially extended droplet with surface dilational viscosity (
$\lambda _{ds}=\infty, Bq=0.4$) and viscosity ratio 
$\lambda =10$ in a stagnant fluid with initial shape 
$L=10.0$ over time. The values of dimensionless parameters are: 
$E=0.4$, 
$\varGamma _\infty =2$ and 
$\gamma =100$.
5. Conclusion
In this paper, we numerically investigated the effect of surface shear and dilational viscosity on single-droplet dynamics in an axisymmetric domain under the Stokes flow regime. We modelled the droplet rheology using the Boussinesq–Scriven constitutive relationship. We examined how the interfacial viscosity alters the steady-state deformation and the critical capillary number for droplet breakup under uniaxial extensional flow at different values of the droplet's viscosity ratio and surfactant transport parameters. We also discussed the impact of surface viscosity on the relaxation dynamics of an initially extended droplet in a stagnant external fluid.
In the case of a droplet with constant surface viscosity, we found that our numerical results for the Taylor deformation parameter are in good quantitative agreement with second-order perturbation theories up to 
${O} {(Ca^{2})}$. The critical capillary numbers from our boundary-element analysis are found to be within 
$15\,\%$ of the theoretical results. Similar correlations have been reported in previous literature examining clean droplet deformation and breakup under extensional flow (Rallison & Acrivos Reference Rallison and Acrivos1978). We observed that the surface shear/dilational viscosity increases/decreases the critical capillary number for droplet breakup by reducing/increasing the droplet deformation at a given capillary number and viscosity ratio compared with a clean droplet. We found that the isolated effect of surface viscosity increases the steady-state droplet deformation compared with a clean droplet for droplets with 
$\lambda _{ds}>{O} {(1)}$. This will be the case for the majority of surfactant systems as the experimentally reported values of surface dilational viscosity are reported to be several orders of magnitude higher than the surface shear viscosity (Krägel et al. Reference Krägel, Kretzschmar, Li, Loglio, Miller and Möhwald1996; Erni et al. Reference Erni, Fischer, Windhab, Kusnezov, Stettin and Läuger2003; Freer et al. Reference Freer, Yim, Fuller and Radke2004; Erni et al. Reference Erni, Windhab and Fischer2011; Kim et al. Reference Kim, Choi, Zasadzinski and Squires2011; Verwijlen et al. Reference Verwijlen, Moldenaers and Vermant2013; Samaniuk & Vermant Reference Samaniuk and Vermant2014; Xiong et al. Reference Xiong, Ren, Tian, Yang, Li and Li2018).
We also discussed the combined influence of surface viscosity and surfactant transport on droplet deformation and breakup. We incorporated the effect of surfactant transport by assuming the Langmuir equation of state to correlate interfacial tension with surfactant concentration changes. We observed that the critical capillary number decreases upon increasing the interfacial viscosity ratio 
$\lambda _{ds}$ at different values of surface Péclet number. At high values of surface Péclet number, the Marangoni effects dominate over the effect of interfacial viscosities and 
$Ca_C$ of a droplet with surface viscosity is found to be the same as that of a clean droplet. We investigated the effect of 
$\varPi$-thickening and 
$\varPi$-thinning surfactants on droplet deformation and breakup by assuming an exponential dependence of surface viscosity on surface pressure. For a droplet with surface shear viscosity, at both low and high values of surface Péclet number, we observed that droplets with pressure-thinning/thickening surfactant have higher/lower 
$Ca_C$ than droplets with pressure-independent surfactant. However, the effect of surface viscosity on 
$Ca_C$ is not very pronounced at 
$Pe_s=10$. In contrast, for a droplet with surface dilational viscosity, at a low value of 
$Pe_s$, a droplet with pressure-thinning/thickening surfactant has lower/higher 
$Ca_C$ compared with a droplet with pressure-independent surfactant. For a droplet with surface dilational viscosity, at a high value of 
$Pe_s$, a droplet with 
$\varPi$-thickening surfactant is observed to be only slightly more deformed than a droplet with 
$\varPi$-thinning surfactant, and we do not observe a significant difference in 
$Ca_C$. Tables 4 and 5 summarize how the Taylor deformation parameter results for a droplet with pressure-dependent surface shear and dilational viscosity differ from those for the pressure-independent case, respectively.
Table 4. Changes in Taylor deformation parameter 
$D_{Taylor}$ of a droplet with pressure-dependent surface shear viscosity relative to the pressure-independent case.
Table 5. Changes in Taylor deformation parameter 
$D_{Taylor}$ of a droplet with pressure-dependent surface dilational viscosity relative to the pressure-independent case.
We also discussed the influence of interfacial viscosity on the relaxation of an initially extended droplet in a quiescent fluid at different values of droplet viscosity ratio 
$\lambda$ (0.01, 0.05, 1 and 10), elasticity number 
$E$ (0.2, 0.4 and 2) and surfactant parameter 
$\gamma$ (1 and 100) that denotes the ratio of diffusive time scale to capillary time scale. At a low value of 
$\lambda =0.01$, a droplet with surface viscosity eventually breaks and forms daughter droplets irrespective of the interfacial viscosity ratio. In contrast, at a high value of 
$\lambda$, a droplet with surface viscosity relaxes into a stable spherical shape irrespective of the interfacial viscosity ratio. At a moderate value of 
$\lambda =1$, we observe that surface dilational viscosity aids in end pinching while surface shear viscosity suppresses the effect, and the droplet relaxes into a stable spherical shape. In cases where a droplet with surface viscosity breaks under relaxation, we observe that the droplet takes a longer time to pinch off compared with a clean droplet. We observe that a high value of elasticity number 
$E=2$ suppresses/delays the pinch-off compared to low values of elasticity number (
$E=0.2$ and 
$E=0.4$).
Acknowledgements
The authors would like to thank the Michael and Carolyn Ott Endowment of Purdue University.
Declaration of interests
The authors report no conflict of interest.
Open access
