Electrohydrodynamic-induced partial coalescence between a droplet and a liquid–air interface

Abstract

When a droplet coalesces with a flat liquid–air interface, a secondary drop may be left behind resulting in only a partial coalescence rather than complete coalescence. In this paper, we employ an arbitrary Lagrangian–Eulerian method to demonstrate that applying an electric field favours the occurrence of partial coalescence. To understand this phenomenon, we systematically study the effect of an external electric field on the coalescence process between a droplet and a liquid–air interface. In an electric field, the induced electric stresses can overcome the downward flow of the droplet, thus lifting it upwards. As a result, the positive Laplace pressure in the neck region squeezes the droplet towards pinch-off. We observe that both the initial neck expansion and neck shrinkage are suppressed by the electric field. These effects become weaker as the Ohnesorge number $Oh$ increases. Based on the scaling analysis, we report a critical Ohnesorge number $O{h_c} = 14.39{\varGamma ^{3/2}} + 0.029$ to quantify the transition from partial coalescence to complete coalescence in the presence of an electric field, where $\varGamma $ represents the dimensionless electric Bond number. Finally, a relationship between the secondary droplet size and the two key dimensionless numbers of $Oh$ and $\varGamma $ has been developed, which could be useful for producing droplets of desired sizes in microfluidic applications.

1. Introduction

The phenomenon of droplet–interface coalescence has attracted significant attention since it is crucial to a range of industrial technologies, such as liquid–liquid separation (Prabhu 2021), coarsening of emulsions (Sainath & Ghosh 2014), microfluidics (Akartuna et al. 2015) and inkjet printing (Lohse 2022). In these applications, a drop deposited on a liquid pool may result in various outcomes, such as bounce, splashing, full coalescence or partial coalescence (Houssainy, Kabachek & Kavehpour 2020). In partial coalescence, a secondary drop can pinch off from the original drop and deposit on the liquid pool again. This process can repeat several times, which has been termed as a coalescence cascade (Thoroddsen & Takehara 2000). Consequently, the partial coalescence slows down the separation of two phases in an emulsification process, and the generated secondary droplet leads to quality defects in an inkjet printing process. Such partial coalescence can cause unwanted effects. It is therefore necessary to thoroughly understand the mechanism of this phenomenon for a highly efficient and optimal process.

When an injected drop comes in contact with the surface of a miscible liquid, the drop could float on the surface without an instantaneous coalescence owing to the lubrication force generated by a thin layer of surrounding fluid between the drop and surface (Reynolds 1881). The time of this finite period is called the residence time of coalescence (Eow et al. 2001). When the film thickness becomes thin enough, the trapped fluid film drains out and the drop undergoes rapid coalescence with the free surface, where partial coalescence may occur to leave a secondary drop behind.

Over half a century ago, partial coalescence was first observed by Charles & Mason (1960). They found a secondary drop forms from the coalescence of a liquid drop with a liquid–liquid interface, which could be suppressed by adding surfactants or applying an electric field. Using the high-speed video imaging technique, Thoroddsen & Takehara (2000) revealed that partial coalescence occurs several times in a self-similar cascade instead of just once. They also demonstrated that the pinch-off time ${t_{pinch}}$ scales with the drop diameter D as ${t_{pinch}}\sim {D^{2/3}}$. In the previous study, Charles & Mason (1960) suggested that Rayleigh–Plateau instability is the main inducement for partial coalescence, where the original drop becomes a thin column and finally breaks up under the action of surface tension. However, Blanchette & Bigioni (2006) used numerical simulations to confirm that partial coalescence was induced by capillary wave convergence on the top of the droplet rather than the Rayleigh–Plateau instability. Aryafar & Kavehpour (2006) performed a series of experiments to study the effects of various parameters on the coalescence between a drop and a liquid–liquid interface. They noted that the Ohnesorge number is the most important dimensionless parameter in this process. The partial coalescence can occur only when $Oh < 1$ and the size of the secondary droplet mainly depends on this parameter. However, Gilet et al. (2007) identified that the Bond number and the relative density difference between the two fluids are also responsible for partial coalescence. They investigated the effects of capillary wave propagation and damping on partial coalescence and proposed the criterion for partial-to-total coalescence transition. In the further study of Blanchette & Bigioni (2009), they numerically simulated the droplet–interface coalescence with negligible impact velocity. The critical Ohnesorge number for partial coalescence was determined as a function of the Bond number, as well as the density and viscosity differences of the drop and surrounding fluid. However, the occurrence of partial coalescence in the absence of capillary waves was not addressed well by them. In this regard, Ray, Biswas & Sharma (2010) used a coupled level set and volume of fluid method to demonstrate that partial coalescence is affected by the competition between the horizontal and vertical momenta of the drop. A secondary droplet could pinch off once the horizontal momentum exceeds the vertical momentum. Apart from these studies about the effects of physical properties on partial coalescence, researchers have considered some special cases, such as a compound drop on a liquid pool (Deka et al. 2019a), a falling drop on a heated pool (Kirar et al. 2020) or a larger-viscosity pool (Alhareth & Thoroddsen 2020). Additionally, the coalescence processes of two drops or bubbles have been comprehensively studied (Eggers, Lister & Stone 1999; Thoroddsen et al. 2007; Munro et al. 2015; Anthony, Harris & Basaran 2020), where partial coalescence was also observed as reported by Zhang, Li & Thoroddsen (2009) and Deka et al. (2019b).

In the works mentioned above, it is worth noting that no external forces were involved throughout the liquid–air or liquid–liquid system. However, previous experiments showed that the extra Marangoni effect and electric field have significant effects on droplet coalescence. Blanchette, Messio & Bush (2009) studied the coalescence of a drop with a miscible but different fluid reservoir. Under the effect of surface tension gradient, partial coalescence was observed. This phenomenon is primarily dependent on the surface tension ratio of the reservoir to drop. Recently, Constante-Amores et al. (2021) conducted a parametric study of the role of surfactant-induced Marangoni stresses in drop–interface coalescence. Unlike the Marangoni effect caused by the surface tension difference, an external electric field can result in the addition of shear stresses on the two-phase interface and consequently alter the coalescence dynamics. The effect of the electric field on hydrodynamics has been of interest since it was reported by Melcher & Taylor (1969). A curious phenomenon of droplet coalescence under an electric field, namely electrocoalescence, has been an attractive technology for its advantages in promoting or preventing coalescence when it is needed or undesirable. Eow and co-workers reviewed the droplet–droplet and droplet–interface coalescence under an applied electric field (Eow et al. 2001; Eow & Ghadiri 2003a,b). They showed that the residence time before coalescence decreases as the electric field strength increases. The extra shear stress can pull the falling drop upward and promote partial coalescence leaving a secondary drop behind. Aryafar & Kavehpour (2007) found that the introduction of a DC electric field can lead to partial coalescence. However, if the electric field reaches the critical strength for droplet breakup, the droplet will not coalesce but instead form into a fluid jet. Mousavichoubeh and co-workers studied the electrocoalescence of a droplet on an oil–water interface under DC and pulsatile electric fields (Mousavichoubeh, Ghadiri & Shariaty-Niassar 2011a; Mousavichoubeh, Shariaty-Niassar & Ghadiri 2011b; Mousavi, Ghadiri & Buckley 2014). Three patterns were observed: full coalescence, partial coalescence and rebound without coalescence. They found that a pulsatile electric field could promote electrocoalescence and suppress the formation of secondary droplets. Moreover, Pillai et al. (2017) used an electrokinetic model to simulate the charge transfer dynamics during droplet–interface coalescence. They showed that increasing the drop charge can also prevent the secondary droplet. Anand, Juvekar & Thaokar (2020) experimentally confirmed that there exists a critical electrocapillary number determining the transition from full coalescence to partial coalescence. Based on the level-set method, Li et al. (2020) found a different jet-like partial coalescence induced by electric stress. The jet then disintegrated into several progeny droplets, which is undesirable in most cases.

Although many experimental and numerical studies have investigated the effect of an external electric field on droplet–interface coalescence, not much work has been carried out on the coalescence between a droplet and a liquid–air interface, and most work relates to liquid–liquid systems, e.g. water-in-oil systems (Aryafar & Kavehpour 2009; Mousavichoubeh et al. 2011b; Anand et al. 2020). When the surrounding liquid is replaced by air, the density and viscosity differences between the two fluids are not crucial to the occurrence of partial coalescence (Blanchette & Bigioni 2009). The dominant parameters can be summarized as the dimensionless Ohnesorge number (liquid properties), the Bond number (gravity) and, certainly, the electric field. Importantly, we note that previous numerical studies described the two-phase flows by using different modelling techniques, e.g. the level-set method (Li et al. 2020) or the volume-of-fluid method (Deka et al. 2019a), where a so-called distance function (or volume fraction function) is utilized to track the interface. However, for the electrohydrodynamic problem, surface charge conservation requires special attention. For a typical leaky-dielectric fluid (Saville 1997), the density of charge is zero in the bulk but both tangential and normal electric stresses appear at the two-phase interface to balance with the viscous stress. Previous studies usually converted the surface charge to the volume charge as a simplification (Das, Dalal & Tomar 2021; Tian et al. 2022), since it is too difficult to directly solve the jumps of stresses at the interface by the two methods. The electric volumetric forces are then included in the momentum conservation equation as a numerical treatment. Fortunately, these limitations can be addressed by a sharp interface method using an arbitrary Lagrangian–Eulerian (ALE) algorithm (Kjellgren & Hyvärinen 1998), in which all mesh nodes are always moving. The moving interface can be accurately captured by reconstructing the mesh. The leaky-dielectric model and ALE technique have proven to be extremely applicable for simulating electrohydrodynamic problems, such as pinching of charged jets (Nie et al. 2021), tip streaming and emission of charged drops (Collins et al. 2008; Wagoner et al. 2021) and electrospray in cone–jet mode (Herrada et al. 2012; Ponce-Torres et al. 2018).

In this work, we adapt the ALE technique to study the coalescence process of a leaky-dielectric drop on a liquid–air interface under a vertical electric field. The objective is to systematically investigate the influence of the electric field on the characteristic parameters of the droplet–interface coalescence process. The rest of this paper is organized as follows. The problem is formulated in § 2, wherein the physical model, governing equations, initial and boundary conditions and numerical method are presented. Section 3 includes the parameter independence analysis of our model, together with a comparison between our numerical results and experimental observations from the literature. Section 4 discusses the influence of the electric field on the characteristics of droplet–interface coalescence. Finally, concluding remarks are given in § 5.

2. Problem description

2.1. Physical model

The coalescence process of a liquid drop on a liquid–air interface is numerically investigated via axisymmetric ALE simulations using the commercial software COMSOL Multiphysics. The system consists of a spherical liquid drop with radius ${R_0}$ and a large reservoir of identical fluid. The two-phase Navier–Stokes equations are solved in a two-dimensional axisymmetric computational domain $(r,z)$, which is a cylindrical region with a radius of $6{R_0}$ and a height of $15{R_0}$ as shown in figure 1. Here r and z respectively represent the radial and axial coordinates, such that $r = 0$ is the axis of symmetry, and $z = 0$ and $z = 15{R_0}$ are the two parallel positive and negative electrodes, respectively. The width of the domain proves to be large enough to eliminate the influence of boundaries at the time scales in this study (Constante-Amores et al. 2021). The horizontal liquid reservoir has an initially uniform thickness of $3{R_0}$. A vertical electric field is then produced by applying a positive voltage $\phi = {\phi _0}$ at the bottom electrode. The top electrode is grounded and thus the upward electric field is ${E_0} = {\phi _0}/15{R_0}$. The time-dependent computational domains of the two fluid phases are denoted by ${\varOmega _i}(t)$ as depicted in figure 1. Here $i = 1$ represents the drop and the liquid reservoir and $i = 2\; $ represents the surrounding medium, which is air in this study. We assume the two phases are under laminar viscous flow and are an incompressible Newtonian fluid with constant density ${\rho _i}$, dynamic viscosity ${\mu _i}$, electrical permittivity ${\varepsilon _i}$ and conductivity ${\sigma _i}$. The liquid–air interface is denoted by ${S_f}$ with a constant surface tension $\gamma $.

Figure 1. (a) Two-dimensional axisymmetric computational domain used in numerical simulations of the coalescence between a drop and a liquid–air interface under a vertical electric field. Initially, the drop and large fluid reservoir are connected by a neck of radius ${r_{min}} = 0.2{R_0}$ to ensure the continuity of the liquid domain. (b) Three-dimensional schematic diagram of the system. (c) Zoomed-in view of the initial neck region.

We begin our simulations with an initially quiescent flow. Since the ALE technique cannot handle the coalescence of two separate domains, there must exist a neck of finite size between the drop and liquid reservoir for initialization. As shown in figure 1(c), at the initial time $t = 0$, the resting drop is connected to the reservoir by a small neck of radius ${r_{min}} = {r_0} = 0.2{R_0}$. The choice of the initial size of the neck is verified in the next section. From previous studies (Blanchette & Bigioni 2006; Martin & Blanchette 2015; Constante-Amores et al. 2021), it proves convenient to adopt this for a convergent result. Moreover, we note that a small Taylor cone may be formed at the gap between the drop and the interface just before coalescence (Lukyanets & Kavehpour 2008). The locally increased electric field gives rise to the deformation once the gap distance is far smaller than the droplet diameter. Since the time scale is too short to affect the coalescence in this study, we ignore the early deformation and assume that the initial interface is flat. Detailed discussions are shown in the Appendix.

2.2. Governing equations

When the initial drop velocity is not very large, it is reasonable to non-dimensionalize the problem by a characteristic length scale ${l_c} = {R_0}$ and an inertial time scale ${t_c} = \sqrt {{\rho _1}R_0^3/\gamma } $. Then all variables will be made dimensionless (denoted with asterisks) using

(2.1) \begin{equation}\left. {\begin{array}{*{20}{c@{}}} {({r^\ast },{z^\ast }) = \dfrac{{(r,z)}}{{{R_0}}},\quad {t^\ast } = \dfrac{t}{{\sqrt {{\rho_1}R_0^3/\gamma } }},\quad ({u^\ast },{w^\ast }) = \dfrac{{(u,w)}}{{{R_0}/\sqrt {{\rho_1}R_0^3/\gamma } }},\; }\\ {\rho_i^\ast= \dfrac{{{\rho_i}}}{{{\rho_1}}},\quad \mu_i^\ast= \dfrac{{{\mu_i}}}{{{\mu_1}}},\quad \varepsilon_i^\ast= \dfrac{{{\varepsilon_i}}}{{{\varepsilon_0}}},\quad \sigma_i^\ast= \dfrac{{{\sigma_i}}}{{{\sigma_1}}},}\\ {{p^\ast } = \dfrac{p}{{\gamma /{R_0}}},\quad {q^\ast } = \dfrac{q}{{{\varepsilon_0}{E_0}}},\quad {Q^\ast } = \dfrac{Q}{{{\varepsilon_0}{E_0}R_0^2}},\quad {E^\ast } = \dfrac{E}{{{E_0}}},} \end{array}} \right\}\end{equation}

where the flow velocity $\boldsymbol{u} = (u,w)$ represents the radial and vertical components of the velocity vector. The pressure, the surface charge density, the surface charge and the electric strength are denoted by p, q, Q and E, respectively. It is noted that we choose the properties of the drop $(i = 1)$ as characteristic scales, while the dielectric permittivities of both phases are non-dimensionalized by the vacuum permittivity ${\varepsilon _0} = 8.85 \times {10^{ - 12}}\;\textrm{F}\;{\textrm{m}^{ - 1}}$. In what follows, all asterisks have been dropped and the corresponding variables are dimensionless unless stated otherwise.

Since we focus on a liquid–air system, the density and viscosity ratios are fixed at ${\rho _2}/{\rho _1} = 0.001$ and ${\mu _2}/{\mu _1} = 0.0186$ (i.e. air to water) to minimize the influence of the surrounding medium. Then the velocity $\boldsymbol{u}$ and pressure p inside the liquid domain $({\varOmega _1})$ are governed by the Navier–Stokes equations, which in dimensionless form are

(2.2)\begin{gather}\frac{{\partial \boldsymbol{u}}}{{\partial t}} + \boldsymbol{u}\boldsymbol{\cdot }\boldsymbol{\nabla }\boldsymbol{u} =- \boldsymbol{\nabla }p + Bo\,{\boldsymbol{i}_z} + Oh{\nabla ^2}\boldsymbol{u}\;\textrm{in}\;{\varOmega _1},\end{gather}

(2.3)\begin{gather}\boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{u} = 0\;\textrm{in}\;{\varOmega _1},\end{gather}

where ${\boldsymbol{i}_z}$ is the unit vector in the direction of gravity. The two dominant dimensionless numbers in (2.2) are defined as

(2.4a,b)\begin{equation}Bo = \frac{{{\rho _1}gR_0^2}}{\gamma },\quad Oh = \frac{{{\mu _1}}}{{\sqrt {{\rho _1}{R_0}\gamma } }}.\end{equation}

Here, $Bo$ is the bond number that relates gravity ($g$ is the acceleration due to gravity) to the surface tension force and $Oh$ is the Ohnesorge number, which is the ratio of the viscous force to the inertial and surface tension forces. In both domains, the electric potential $\phi $ obeys the axisymmetric form of the electrostatic equation:

(2.5)\begin{equation}{\nabla ^2}\phi = 0\quad \textrm{in}\;{\varOmega _1}\;\textrm{and}\;{\varOmega _2}.\end{equation}

Then the electric field can be computed by the formula $\boldsymbol{E} =- \boldsymbol{\nabla }\phi $. To solve the electrohydrodynamic problem, Laplace's equation (2.5) is coupled with the Navier–Stokes equations by adding the Maxwell stresses into the traction condition (Aris 1990) of tangential and normal stress balance along the interface ${S_f}$:

(2.6)\begin{equation}\boldsymbol{n}\boldsymbol{\cdot }{\bf [}\boldsymbol{T}_i^H + \boldsymbol{T}_i^E{\bf ]}_1^2 = \kappa \boldsymbol{n}\;\textrm{on}\;{S_f},\end{equation}

where $[{x_i}]_1^2 = {x_\textrm{2}} - {x_\textrm{1}}$ represents the jump of a physical quantity x across the two-phase interface and $\boldsymbol{n}$ the outward unit vector of normal direction at the interface pointing to the air phase. On the left-hand side of (2.6), $\boldsymbol{T}_i^H =- {p_i}\boldsymbol{I} + O{h_i}(\boldsymbol{\nabla }\boldsymbol{u} + {(\boldsymbol{\nabla }\boldsymbol{u})^\textrm{T}})$ is the hydrodynamic stress tensor and $\boldsymbol{T}_i^E = 2\varGamma {\varepsilon _i}(\boldsymbol{EE} - {\boldsymbol{E}^2}\boldsymbol{I}/2)$ is the Maxwell stress tensor (Melcher & Taylor 1969). Here, $\boldsymbol{I}$ is the unit tensor and $\varGamma = {\varepsilon _0}E_0^2{R_0}/2\gamma $ denotes the dimensionless electric Bond number, which is the ratio of the electric force to the capillary force. On the right-hand side of (2.6), $\kappa =- \boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{n}$ is the local mean curvature of the interface.

In the leaky-dielectric model (Melcher & Taylor 1969; Saville 1997), the density of charge is zero in the bulk but the charge distribution on ${S_f}$ obeys the surface charge transport equation, which in dimensionless form is

(2.7)\begin{equation}\frac{{\partial q}}{{\partial t}} + {\nabla _s}\boldsymbol{\cdot }(q\boldsymbol{u}) - \frac{1}{{Pe}}\nabla _s^2q = {\varepsilon _1}\alpha \boldsymbol{n}\boldsymbol{\cdot }\left[ {{\boldsymbol{E}_1} - \frac{{{\sigma_2}}}{{{\sigma_1}}}{\boldsymbol{E}_2}} \right]\;\textrm{on}\;{S_f},\end{equation}

where $\; q = \boldsymbol{n}\boldsymbol{\cdot }[{\varepsilon _i}{\boldsymbol{E}_i}]_1^2$ is the surface charge density computed by the Gauss law, ${\nabla _s} = (\boldsymbol{I} - \boldsymbol{nn})\boldsymbol{\cdot }\boldsymbol{\nabla }$ is the surface gradient operator and ${\boldsymbol{E}_1}$ and ${\boldsymbol{E}_2}$ are the electric fields along the interface ${S_f}$ on the liquid and air sides, respectively. The Péclet number, $Pe = R_0^2/{D_c}{t_c}$, is defined as the ratio of the charge diffusion time scale $R_0^2/{D_c}$ to the inertial time scale ${t_c}$, where ${D_c}$ is the surface diffusion coefficient. The dimensionless relaxation parameter $\alpha = {t_c}/{t_e} = \sqrt {{\rho _1}R_0^3\sigma _1^2/\gamma \varepsilon _1^2} $ is introduced for measuring the importance of charge relaxation, where ${t_e} = {\varepsilon _1}/{\sigma _1}$ is the charge relaxation time. It should be noted that here ${\varepsilon _1}$ is the dimensional electrical permittivity of fluid 1, while in (2.7) ${\varepsilon _1}$ is the relative permittivity (asterisk dropped) dimensionalized by ${\varepsilon _0}$. The three terms on the left-hand side of (2.7) represent the charge accumulation, the charge convection and the charge diffusion, respectively. The source term on the right-hand side represents the Ohmic conduction caused by the movement of free electrons. In our simulations, the surrounding air is considered to be free of any electrons. Thus, the conductivity of fluid 2, ${\sigma _2} \equiv 0$, and the second term on the right-hand side of (2.7) has been omitted. When the relaxation parameter $\alpha \to \infty $, the charge is entirely dominated by Ohmic conduction and can accumulate instantaneously to the peaks and troughs of the interface. Then the fluid can be regarded as a perfect conductor. When $\alpha $ is finite (e.g. $\alpha = 5.4$ in the present work), the charge convection and diffusion compete with the Ohmic conduction and thus the charge transport would slow down (Nie et al. 2021).

On the interface, the velocity satisfies the kinematic condition (Scriven 1960):

(2.8)\begin{equation}\boldsymbol{n}\boldsymbol{\cdot }(\boldsymbol{u} - {\boldsymbol{u}_s}) = 0\;\textrm{on}\;{S_f},\end{equation}

where ${\boldsymbol{u}_s}$ is the local velocity on the interface ${S_f}$. For the electric potential, Dirichlet boundary conditions on $\phi $ are imposed on the top $(z = 12)$ and bottom $(z = 0)$ boundaries, and Neumann boundary condition $\partial \phi /\partial r = 0$ is imposed on the right-hand boundary $(r = 6)$. At $r = 0$, symmetry boundary conditions are imposed for both the flow and the electric field. The no-slip and no-penetration conditions are imposed on other boundaries.

From the equations above, the dimensionless groups used to describe the results can be summarized as the Ohnesorge number $Oh$, the electric Bond number $\varGamma $, the Bond number $Bo$ and the Péclet number $Pe$, which represent the roles of fluid viscosity, electric field strength, gravity and charge diffusion, respectively. Moreover, the electrical properties ${\varepsilon _i}$ and $\alpha $ and the density and viscosity ratios ${\rho _2}/{\rho _1}$ and ${\mu _2}/{\mu _1}$ are kept constant in the present work.

2.3. Numerical scheme

The governing equations for mass and momentum conservation (2.2) and (2.3), electrostatic equation (2.5), traction condition (2.6) along with the surface charge transport equation (2.7) are numerically solved by a moving-grid finite-element method (FEM) together with the weak formulation implemented in COMSOL Multiphysics. Following the methodology described by Martínez-Calvo et al. (2020), equations (2.2)–(2.7) are written in weak form by introducing the so-called test functions. First, we derive the weak forms of (2.2) and (2.3) for fluid 1 by defining ${\boldsymbol{v}_{\boldsymbol{u}}}$ and ${\boldsymbol{v}_{\boldsymbol{p}}}$ as the test functions of the velocity and pressure, respectively. Applying the Gaussian theorem and incorporating the boundary conditions, (2.2) and (2.3) become

(2.9) \begin{align} 0 & = \int_{{\varOmega _1}} {{\boldsymbol{v}_{\boldsymbol{u}}}\dfrac{{\partial \boldsymbol{u}}}{{\partial t}}\,\textrm{d}{\varOmega _1}} + \int_{{\varOmega _1}} {{\boldsymbol{v}_{\boldsymbol{u}}}(\boldsymbol{u}\; \boldsymbol{\cdot }\boldsymbol{\nabla }\boldsymbol{u})\,\textrm{d}{\varOmega _1}} + \int_{{\varOmega _1}} {{\boldsymbol{v}_{\boldsymbol{u}}}\boldsymbol{\nabla }p\,\textrm{d}{\varOmega _1}} - Bo\int_{{\varOmega _1}} {{\boldsymbol{v}_{\boldsymbol{u}}}({\boldsymbol{i}_z})\,\textrm{d}{\varOmega _1}} \notag\\ & \quad + Oh\int_{{\varOmega _1}} {(\boldsymbol{\nabla }{\boldsymbol{v}_{\boldsymbol{u}}})(\boldsymbol{\nabla }\boldsymbol{u})\,\textrm{d}{\varOmega _1}} - Oh\int_{{S_1}} {{\boldsymbol{v}_{\boldsymbol{u}}}(\boldsymbol{\nabla }\boldsymbol{u})\boldsymbol{\cdot }\boldsymbol{n}\,\textrm{d}{S_1}} \end{align}

and

(2.10)\begin{equation}0 = \int_{{\varOmega _1}} {{\boldsymbol{v}_{\boldsymbol{p}}}(\boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{u})\,\textrm{d}{\varOmega _1}} .\end{equation}

Here, ${S_1}$ represents the boundary of the domain ${\varOmega _1}$. The solutions and test functions are supposed to belong to the same function space, and we can assume that they satisfy the same boundary conditions, e.g. ${\boldsymbol{v}_{\boldsymbol{u}}}$ satisfies the symmetry boundary condition on the axis of symmetry and the no-slip condition on the other boundaries. Likewise, the weak form of the electrostatic equation (2.5) reads

(2.11)\begin{equation}0 =- \int_{{\Omega _{1,2}}} {\boldsymbol{\nabla }{\boldsymbol{v}_\phi }\boldsymbol{\cdot }\boldsymbol{\nabla }\phi \; \,\textrm{d}{\varOmega _{1,2}}} + \int_S {{\boldsymbol{v}_\phi }(\boldsymbol{\nabla }\phi \boldsymbol{\cdot }\boldsymbol{n})\; \,\textrm{d}S} ,\end{equation}

where ${\boldsymbol{v}_\phi }$ is the test function of the electric potential and S is the boundary of the whole domain (${\varOmega _1}$ and ${\varOmega _2}$). On the interface, the surface charge transport equation (2.7) is given as

(2.12) \begin{align} 0 & = \int_{{S_f}} {{\boldsymbol{v}_q}\dfrac{{\partial q}}{{\partial t}}\,\textrm{d}{S_f}} - \int_{{S_f}} {({\nabla _s}{\boldsymbol{v}_q})\boldsymbol{\cdot }(q\boldsymbol{u})\; \,\textrm{d}{S_f}} + \dfrac{1}{{Pe}}\int_{{S_f}} {({\nabla _s}{\boldsymbol{v}_q})\boldsymbol{\cdot }({\nabla _s}q)\; \,\textrm{d}{S_f}} \notag\\ & \quad - \int_{{S_f}} {{\boldsymbol{v}_q}({\varepsilon _1}\alpha \boldsymbol{n}\boldsymbol{\cdot }{\boldsymbol{E}_1})\; \,\textrm{d}{S_f}} + \int_{{P_f}} {{\boldsymbol{v}_q}(q\boldsymbol{u})\boldsymbol{\cdot }\boldsymbol{n}\,\textrm{d}{P_f}} \; - \dfrac{1}{{Pe}}\int_{{P_f}} {{\boldsymbol{v}_q}({\nabla _s}q)\boldsymbol{\cdot }\boldsymbol{n}\,\textrm{d}{P_f}} , \end{align}

where ${\boldsymbol{v}_q}$ is the test function of the surface charge density and ${P_f}$ is the boundary point of ${S_f}$. Moreover, the traction condition (2.6) along the interface is given as

(2.13)\begin{equation}0 =- \int_{{S_f}} {{\boldsymbol{v}_{\boldsymbol{u}}}\textrm{(}\boldsymbol{n}\boldsymbol{\cdot }[\boldsymbol{T}_i^H + \boldsymbol{T}_i^E]_1^2\textrm{)}\,\textrm{d}{S_f}} + \int_{{S_f}} {{\boldsymbol{v}_{\boldsymbol{u}}}(\kappa \boldsymbol{n})\,\textrm{d}{S_f}} ,\end{equation}

where the effect of the electric field is incorporated into the free-surface flow. The boundary conditions of each equation are imposed to obtain the FEM solutions. Then the computational domain is discretized by free triangular meshes. Importantly, the mesh is heavily concentrated near the neck region while relatively coarser elsewhere. In the ALE technique, the computational mesh inside the domains can change arbitrarily to deform the interface during the coalescence process, which, however, results in an unavoidable large deformation of the mesh. To handle this problem, we adapt the automatic remeshing algorithm to enhance the convergence of the solution. When the mesh quality falls below a specified value, the algorithm will automatically remesh the whole domain and reinitialize the solution. During the calculation, the mesh size is always ensured to be smaller than the neck radius. The time discretization is performed using the built-in backward differentiation formula. All variables are fully coupled by the constant Newton's method (Deuflhard 1974) with a constant damping factor of 0.9. To constrain the mesh displacement, we adapt the Yeoh smoothing method (Yeoh 1993) since it allows the largest displacement before the mesh is inverted. Considering that the ALE technique cannot capture the topology changes when the computational domain is broken into separate parts (Antonopoulou et al. 2020), the calculation will be automatically terminated, and thus only the coalescence dynamics before pinch-off is discussed in this paper.

3. Validation of the numerical method

3.1. Parameter independence analysis

We first validate our numerical simulations by conducting parameter independence analysis. Water and air are chosen as the two working fluids, the density and viscosity ratios of which have been determined as mentioned above. To obtain the electrical properties, we use deionized water (${\varepsilon _1} = 78{\varepsilon _0}$ and ${\sigma _1} = 1 \times {10^{ - 6}}\;\textrm{S}\;{\textrm{m}^{ - 1}}$) as the leaky-dielectric fluid. Thus the charge relaxation time can be calculated as ${t_e} = 0.69\;\textrm{ms}$, which is less than the electrocoalescence characteristic time of ${\sim} 10\;\textrm{ms}$ (Anand et al. 2020). Considering a spherical water drop with ${R_0} = 1\;\textrm{mm}$, ${\rho _1} = 1\; 000\;\textrm{kg}\;{\textrm{m}^{ - 3}}$ and $\gamma = 0.072\;\textrm{N}\;{\textrm{m}^{ - 1}}$, the relaxation parameter $\alpha = 5.4$ is used in the following simulations.

Figure 2(a) shows the grid independence analysis for the variations of minimum neck radii during the coalescence processes. Parameters used are $Oh = 0.02$, $\varGamma = 0.01$ and $Bo = 0$, representing a universal case with the balance of inertial, viscous and electric forces in the absence of gravity (Constante-Amores et al. 2021). Here the selected Péclet number for charge relaxation is set to $Pe = 1\; 000$. Initially the neck radius increases quickly to reach the maximum value and then starts to decrease until pinch-off. We note that the variations of ${r_{min}}$ have no dependence on the grid number N, and thus the choice of $N = 50\;000$ is therefore sufficient. Moreover, comparison with the results ($Oh = 0.02$, $\varGamma = 0$ and $Bo = 0.001$) from Constante-Amores et al. (2021) shows good agreement with our simulations as depicted in figure 2(a).

Figure 2. Temporal evolutions of the minimum neck radius for different (a) grid numbers N, (b) Péclet numbers $Pe$ and (c) initial neck radii ${r_0}$. The open symbols denote the results reported by Constante-Amores et al. (2021) for surfactant-free droplet–interface coalescence. (d) Temporal evolutions of the droplet neck for different initial neck radii ${r_0}$. Here the dimensionless parameters are fixed at $Oh = 0.02$, $\varGamma = 0.01$ and $Bo = 0$.

For the choice of the Péclet number, figure 2(b) shows a universal profile of ${r_{min}}$ versus t for different Péclet numbers. Previous studies have demonstrated that the role of surface charge diffusion is insignificant in charge transport provided that $Pe \gg 1$ (Wagoner et al. 2021). We also test the choice of the initially existing neck radius ${r_0}$. As shown in figure 2(c), decreasing the initial neck radius has a negligible effect on the variations of ${r_{min}}$. Figure 2(d) further confirms that the value of ${r_0}$ is not decisive in the coalescence dynamics. However, a sufficiently small neck brings substantial difficulty in computational efficiency. Therefore, all numerical results are obtained using ${r_0} = 0.2\; {R_0}$. Based on the parameter independence analysis, we use $N = 50\;000$, ${r_0} = 0.2$ and $Pe = 1\; 000$ for all simulations in this study.

3.2. Comparison against experiments

Figure 3 shows snapshots of droplet–interface coalescence under an electric field from the experiments of Aryafar & Kavehpour (2009), together with the present numerical simulations. It should be noted that here the surrounding medium is high-viscosity silicone oil. The corresponding dimensionless numbers are $Oh = 3.39 \times {10^{ - 3}}$, $\varGamma = 0.022$, $Bo = 0.9$, ${\rho _2}/{\rho _1} = 0.971$, ${\mu _2}/{\mu _1} = 9.7 \times {10^3}$, ${\varepsilon _1} = 78$, ${\varepsilon _2} = 2.5$ and $\alpha = 19.6$. A water droplet is initially released and then touches the interface. Under the electric field, the coalescing droplet will form a Taylor-cone-like jet at its apex and then breaks up into a secondary droplet. Our numerical results can well reproduce this phenomenon and show qualitative agreement with the experiments of Aryafar & Kavehpour (2009); however, there are some differences in the shape of the meniscus and sizes of the secondary droplets. We note that this may be attributed to the initial velocity of the deposited droplet in the experiments, which has been neglected in our simulations. Nevertheless, this work confirms that the numerical model can capture the dynamics of droplet–interface coalescence under an electric field.

Figure 3. Comparison of our numerical results (second row) against the experimental results (first row) from Aryafar & Kavehpour (2009) for coalescence of a 4 mm initial diameter drop on a planar interface under an electric field of ${E_0} = 3.27 \times {10^5}\;\textrm{V}\;{\textrm{m}^{ - 1}}$. Here the surrounding medium is silicone oil with viscosity of 10 000 cSt. The dimensionless parameters are $Oh = 3.39 \times {10^{ - 3}}$, $\varGamma = 0.022$, $Bo = 0.9$, ${\rho _2}/{\rho _1} = 0.971$, ${\mu _2}/{\mu _1} = 9.7 \times {10^3}$, ${\varepsilon _1} = 78$, ${\varepsilon _2} = 2.5$ and $\alpha = 19.6$. Reprinted with permission from Aryafar & Kavehpour (2009). Copyright 2009 American Chemical Society.

4. Results and discussion

4.1. Mechanism of partial coalescence under an electric field

We begin our study by presenting the results of droplet–interface coalescence with and without an external electric field. Figure 4(a) shows the time evolutions of droplet shapes during the coalescence processes for $\varGamma = 0$ and $\varGamma = 0.01$. Here the other dimensionless parameters are $Oh = 0.02$ and $Bo = 0$, indicating that gravity is negligible in our study. As predicted by Blanchette & Bigioni (2006), partial coalescence occurs in both cases since the critical Ohnesorge number $O{h_c} = 0.026 \pm 0.003$ has not been reached. Moreover, the time spent from the onset of coalescence to pinch-off, ${\sim} 1.7{t_c}$, is consistent with their experiments. Figure 4(a) shows how partial coalescence occurs: the initial rapid expansion of the neck generates a capillary wave along the whole interface, under which the original droplet starts to oscillate and deform. Once the capillary wave reaches its apex, the droplet will be stretched upward to form a thin column as shown in the sixth profile $(t = 1.25)$. However, the surface tension force will pull the expanded neck inward, and consequently, pinch-off will occur to leave a secondary droplet behind provided that the vertical collapse is not enough to prevent it. The radius of the secondary droplet is nearly half that of the original droplet, ${R_d} = 0.556$ for $\varGamma = 0$ and ${R_d} = 0.585$ for $\varGamma = 0.01$, denoting that $17\,\%$–$20\,\%$ of the volume will be ejected from the original droplet. It is clearly seen that the electric field hardly affects the shapes of the interface during the initial stage of coalescence. However, the electric field gives rise to a slight upward elongation of the droplet as depicted in the last two profiles in figure 4(a). Both the size of the secondary droplet and the pinch-off time increase when applying an external electric field.

Figure 4. (a) Temporal evolutions of the coalescence sequence for $\varGamma = 0$ and $\varGamma = 0.01$. In both cases, $Oh = 0.02$ and $Bo = 0$. Black profiles denote the liquid–air interfaces within $r < 3$ at different moments, where the time interval between the two consecutive profiles is $\Delta t = 0.25$. The times are given in dimensionless units. (b) Temporal evolutions of the three-dimensional interfacial contours of the surface charge density q (left) and the radial interfacial velocity u (right).

Figure 4(b) shows the spatiotemporal contours of the surface charge density and the radial component of the interfacial velocity. It can be seen at the early time $(t \approx 0.5)$ in the velocity field that negative and positive interfacial velocities are created on the interface above and below the neck point $({r_{min}},\;{z_{neck}})$, respectively: $u < 0$ for $z > {z_{neck}}$ and $u > 0$ for $z < {z_{neck}}$. This pulls the neck inward and the higher part of the droplet upward. When pinch-off occurs, the velocity is concentrated near the breaking point and the radial interfacial velocity of the secondary droplet approaches zero. For the case in the presence of an electric field, initially, the surface charge is $q{|_{t = 0}} = 0$ and then is quickly distributed across the interface as seen in figure 4(b). Whether the neck expands out or narrows in, the charge will gradually concentrate at the apex of the secondary droplet, indicating stronger electric stress in this region. Hence the droplet is stretched upward compared with the electric-field-free case.

We now investigate the effect of an electric field on the temporal dynamics of the minimum neck radius ${r_{min}}$ and the maximum height ${h_{max}}$ of the droplet. Figure 5(a) shows the temporal evolutions of neck radius for different electric Bond numbers $\varGamma $. In figure 5(a), the growth of the neck radius ${r_{min}}$ follows a linear relationship with t for a very short time. It can be seen during this expansion stage that these four curves almost overlap, denoting that $\varGamma $ has little or no influence on the initial neck growth. However, at longer times, the rate of neck expansion decreases with $\varGamma $ due to the electric force. The maximum value of the neck radius decreases slightly when the electric field strength becomes stronger. These observations are consistent with the experimental results of Anand et al. (2020). During the subsequent shrinkage stage of the neck, increasing $\varGamma $ leads to a decrease in the neck radius until a transition time. After this time, the decreasing rate of the neck becomes slower as $\varGamma $ increases, and consequently, the bottom pinch-off is delayed due to the presence of an electric field. Unlike neck radius, the maximum height ${h_{max}}$ exhibits a monotonically increasing dependence on $\varGamma $ and the oscillation of the drop apex due to capillary waves is well observed as shown in figure 5(b). Moreover, the evidence for upward elongation of the droplet induced by the electric force can be seen in this figure. For $\varGamma = 0$ and $\varGamma = 0.01$, ${h_{max}}$ starts to decrease when arriving at the final stage of pinch-off. However, for $\varGamma = 0.02$ and $\varGamma = 0.03$, ${h_{max}}$ continuously increases to reach a much higher value than the two former cases. This demonstrates that the electric force has overcome the downward flow, thus lifting the droplet upwards.

Figure 5. Temporal evolutions of (a) the minimum neck radius and (b) the maximum height of the droplet during the coalescence processes for different electric Bond numbers $\varGamma $. (c) Tip streaming phenomenon and the conical tip structure are marked for $\varGamma = 0.03$. Here, $Oh = 0.02$ and $Bo = 0$, which are the same as in figure 4. The surface profile z, the charge density q, the normal electric stress jump $T_n^E$ and the radial interfacial velocity u are shown for (d–g) $t = 1.2$ and (h–k) $t = 1.63$, where the coloured points represent the positions of the neck. The abscissa in (e–g) and (i–k) corresponds to the arclength s, which is zero at the apex of the droplet.

It is noted that for $\varGamma = 0.03$ in figures 5(a) and 5(b), the calculation has been terminated before ${r_{min}}$ decreases to zero since a sharp singularity appears. A conical tip structure is finally formed at the apex of the droplet as shown in figure 5(c), and consequently, a thin fluid jet could be ejected, which is similar to the phenomenon of electrohydrodynamic tip streaming (Collins et al. 2008). Due to the large difference in length scales between the tip streaming $({\sim} \; 0.01)$ and droplet–interface coalescence $({\sim} \; 1)$, we do not simulate the subsequent phenomenon. However, our results show that if the electric field is strong enough, tip streaming will occur at the apex of the droplet before its bottom pinch-off.

Figure 5(d–k) shows the interfacial parameters, including the droplet shape z, the surface charge density q, the normal electric stress jump $T_n^E$ and the radial interfacial velocity u. Here we choose the results of $t = 1.2$ and $t = 1.63$ to study the effect of electric field on the pinch-off dynamics. Droplet shapes for different $\varGamma $ are given in figures 5(d) and 5(h), where the neck positions are all marked. On increasing $\varGamma $, both the droplet height and the neck height are enhanced. At $t = 1.2$, most of the surface charge is accumulated at the droplet apex. When the droplet approaches pinch-off, the surface charge at the droplet apex for $\varGamma = 0.01$ and $\varGamma = 0.02$ increases slowly, whereas this for $\varGamma = 0.03$ rises markedly. Das et al. (2021) demonstrated the importance of charge convection due to fluid flow in the electrohydrodynamic interactions between droplets. However, since ${\varepsilon _1}$ and $\alpha $ are kept constant in our study, we do not discuss the effect of charge convection on the coalescence dynamics. The accumulation of q near the apex modifies the electric field, thereby inducing a local increase of the normal electric stress jump as shown in figure 5(j). Here, $T_n^E$ is defined as the normal component of the electric stress jump across the interface:

(4.1)\begin{equation}T_n^E = \boldsymbol{n}\boldsymbol{\cdot }[{\boldsymbol{T}^E}]_1^2\boldsymbol{\cdot }\boldsymbol{n} = \varGamma [{\varepsilon _2}E_{2,n}^2 - {\varepsilon _1}E_{1,n}^2 - ({\varepsilon _2} - {\varepsilon _1})E_t^2],\end{equation}

where ${E_{i,n}} = \boldsymbol{n}\boldsymbol{\cdot }{\boldsymbol{E}_i}$ and ${E_t} = \boldsymbol{t}\boldsymbol{\cdot }{\boldsymbol{E}_i}$. The normal electric stress is the dominant electric force that acts at the drop interface. It is seen from figures 5(f) and 5(j) that the variation of the normal electric stress is similar to that of the charge density. For $\varGamma = 0.03$, the normal electric stress is much larger at the apex. However, as marked by the coloured points in figure 5(j), the values of the normal electric stress become negative in the neck region for the three electric-field-laden cases. Here $T_n^E \le 0$ indicates that the electric normal stresses are compressive (act inward) along the local interface (Wagoner et al. 2021). Although the inward force favours pinch-off, the coalescence is still delayed under the outward electric normal stress, which is along the rest of the interface and inhibits neck shrinkage. The peak and trough in the curves of the radial interfacial velocity represent the direction of the developing neck as depicted in figures 5(g) and 5(k). It is evident that $u < 0$ near the neck region during the shrinkage stage. The magnitude of u decreases as $\varGamma $ increases and the overall shape in u is similar for all four cases. Further in time, as the curvature of the neck increases, the negative interfacial velocity becomes sharper (large peak or trough) as shown in figure 5(k), denoting the subsequent formation of pinch-off. It is noted in this figure that the curves of u for $\varGamma = 0.02$ and $\varGamma = 0.03$ show little difference in the neck region, while varying greatly in both ends of the interface, especially in the apex region. This reveals that during the later stage of coalescence, the effect of the electric field on neck shrinkage becomes less significant with increasing $\varGamma $, which is consistent with the trend in figure 5(a).

After the above analysis, one remaining problem is how the electric field affects the different stages of the droplet–interface coalescence. From figure 5(a) we conclude that both the initial neck expansion and the eventual pinch-off are suppressed in the presence of an electric field. Previous results show that the outward electric normal stress gives rise to the suppression of pinch-off. However, the initial coalescence is also slowed down owing to the electric force, and therefore the maximum neck radius decreases as $\varGamma $ increases. We plot the distribution of the interface charge density in figure 6. Since the bottom electrode is positive, the surface of the liquid reservoir is positively charged; moreover, the upper and lower parts of the droplet are positively and negatively charged, respectively, as depicted in figure 6(a). The surface charge distribution at different moments is shown in figure 6(b). Here, $Oh = 0.02$ and $\varGamma = 0.03$. At very early time $(t = 0.01)$, positive and negative charges are evenly distributed on both sides of the neck. Over time, the charge starts to accumulate and increase. The peak of q first rises and then falls as shown in figure 6(c). The result is consistent with that of Hamlin, Creasey & Ristenpart (2012). They demonstrated that during the coalescence process, the charge at the droplet bottom is pulled away by convection. However, before coalescence, the top half of the droplet already contains charges of the opposite sign. The remaining charge accumulates on the top and provides the electric driving force for the pinch-off. Additionally, we note that the charge density is always negative at the neck. Then we plot the electric normal stress in figure 6(d), which is found to be compressive (inward) in the neck region and particularly much larger in the early time. For this reason, the electric field has suppressed the initial coalescence and culminated in the reduction of the maximum neck radius. From figures 5 and 6, it is known that the electric normal stress will be compressive in the neck region if there is a singularity, e.g. during the initial coalescence and the eventual pinch-off.

Figure 6. (a) Schematic representation of the charge distribution on the interface during the initial coalescence. (b) Charge density distribution in the neck region at different moments for $Oh = 0.02$ and $\varGamma = 0.03$. (c) The charge density and (d) the normal electric stress jump on the interface with marked neck positions.

Next, we turn our attention towards cases with a larger Ohnesorge number, where the capillary wave is not strong enough to favour the pinch-off due to the highly viscous effect. Blanchette & Bigioni (2009) reported the critical Ohnesorge number for partial coalescence is $O{h_c} \approx 0.026 \pm 0.003$. Therefore, we set the Ohnesorge number to $0.03$ because it allows for the observation of the transition from complete coalescence to partial coalescence. Figure 7(a) demonstrates the time evolutions of droplet shapes during the coalescence processes for $\varGamma = 0$ and $\varGamma = 0.01$. As expected, complete coalescence is observed in the absence of electric field $(\varGamma = 0)$. Upon coalescence, the capillary wave moves towards the apex, leading to the formation of a thin column at $t = 1.25$, which is similar to the results in figure 4. The surface tension pulls the neck inward to form the pinch-off. However, under the highly viscous effect, the droplet starts to drain into the reservoir and the pinch-off fails to occur. The neck will experience a second expansion and shrinkage, resulting in a self-similar behaviour. The inset in figure 7(a) shows that a small bulge forms at the interface centre at $t = 2.5$. The bulge moves downward and finally the droplet completely merges with the liquid reservoir at $t = 4$. The coalescence time is obviously larger than that where pinch-off occurs. Unlike the electric-field-free case, partial coalescence is not prevented by the viscous effect for $\varGamma = 0.01$. The electric force has overcome the interfacial forces, leading to the pinch-off of the droplet at $t = 1.76$. Likewise, the spatiotemporal contours of q and u are given in figure 7(b). The charge is expected to concentrate at the apex of the secondary droplet and the electric force stretches the droplet upward; therefore, the height of the droplet apex for $\varGamma = 0.01$ is higher than that for $\varGamma = 0$. At $t = 2.0$, when the pinch-off does not occur for $\varGamma = 0$, the positive (outward) interfacial velocity is observed on the interface below the neck point, which is responsible for the formation of the small bulge as depicted in figure 7(a).

Figure 7. (a) Temporal evolutions of the coalescence sequence for $\varGamma = 0$ and $\varGamma = 0.01$. In both cases, $Oh = 0.03$ and $Bo = 0$. Black profiles denote the liquid–air interfaces within $r < 3$ at different moments, where the time interval between the two consecutive profiles is $\Delta t = 0.25$. The times are given in dimensionless units. The inset shows the interface shape at $t = 2.5$ for $\varGamma = 0$. (b) Temporal evolutions of the three-dimensional interfacial contours of the surface charge density q (left) and the radial interfacial velocity u (right).

Figure 8(a) shows the temporal evolutions of neck radius for different electric Bond numbers. It can be seen that the initial coalescence is slightly suppressed by increasing $\varGamma $. During the first shrinkage stage, the neck radius decreases to zero, i.e. pinch-off occurs for all electric-field-laden cases. Compared with the curves in figure 5, the effect of the electric field is weakened by increasing the Ohnesorge number. For the electric-field-free case, the neck radius will increase again instead of decreasing to zero, which is consistent with the droplet shapes in figure 7. Note that we define the neck as the location where the radius is minimum; therefore, the neck sometimes disappears during the second shrinkage stage and the final coalescence as shown in figure 7(a). For this reason, the curve of the neck radius will break as marked in figure 8(a). The temporal evolutions of the maximum height are given in figure 8(b). We find a faster vertical collapse velocity for $\varGamma = 0$. Although the dynamics is similar to the case with $Oh = 0.02$, which also exhibits a monotonic dependence on $\varGamma $, there is no emergence of the tip-streaming phenomenon before the final bottom pinch-off for $Oh = 0.03$. This sequence corresponds to the requirement of tip streaming, where the capillary pressure at the droplet apex must be overcome by either electrostatic pressure or dynamic pressure. Figure 8(c) shows a magnified view of the secondary droplet's tail at $t = 1.76$, and also shows a thin thread $(\mathrm{\sim }0.001)$ connecting the about-to-form secondary droplet to the interface. The formation of thin threads was observed for $Oh = 0.03$ in all three electric-field-laden cases, while not for $Oh = 0.02$, denoting that it can only be produced by the electric field when the Ohnesorge number is high enough.

Figure 8. Temporal evolutions of (a) the minimum neck radius and (b) the maximum height of the droplet during the coalescence processes for different electric Bond numbers $\varGamma $. (c) The structure of bottom pinch-off at $t = 1.76$ is marked for $\varGamma = 0.01$. Here, $Oh = 0.03$ and $Bo = 0$, which are the same as in figure 7. The surface profile z, the charge density q, the normal electric stress jump $T_n^E$ and the radial interfacial velocity u are shown for (d–g) $t = 1.3$ and (h–k) $t = 1.75$, where the coloured points represent the positions of the neck. The abscissa in (e–g) and (i–k) corresponds to the arclength s, which is zero at the apex of the droplet.

Likewise, the droplet shape z, the surface charge density q, the normal electric stress jump $T_n^E$ and the radial interfacial velocity u are further provided in figure 8(d–k). Here we choose the results of $t = 1.3$ and $t = 1.75$ to study the effect of the electric field on the coalescence dynamics. At $t = 1.3$, the effect of $\varGamma $ on the instantaneous parameters is similar to that previously observed in figure 5(d–g). The corresponding surface charge density and normal electric stress jump are almost unaltered as compared with the cases of $Oh = 0.02$. Moreover, it is clearly seen that the magnitude of the radial interfacial velocity remains constant as $\varGamma $ increases. This further demonstrates the weakened effect of the electric field on the neck shrinkage process with a higher Ohnesorge number. The droplet shapes at $t = 1.75$ given in figure 8(h) show that the neck radius for $\varGamma = 0$ is much larger than that for $\varGamma > 0$. This arises due to the delayed neck shrinkage velocity caused by the fast vertical collapse. As a result, the negative interfacial velocity is sharpest for $\varGamma = 0.01$ as shown in figure 8(k). Since tip streaming is avoided by the highly viscous effect, the curvature of the apex decreases; therefore, the surface charge density and the normal electric stress jump at the apex for $Oh = 0.03$ decrease slightly as compared with the cases for $Oh = 0.02$.

The results in both figures 5 and 8 comprehensively indicate the effect of the electric field on the partial coalescence, but for the underlying mechanism yet another question has to be answered: Why does applying an electric field favour the partial coalescence between a droplet and a liquid–air interface? It turns out that the electric force overcomes the downward flow of the droplet and stretches the droplet upwards. According to Ray et al. (2010), pinch-off could occur when the horizontal momentum exceeds the vertical momentum. Alhareth & Thoroddsen (2020) suggested that the pinch-off is determined by the interplay between the axial curvature ${\kappa _x} = 1/{r_a}$ and the azimuthal curvature ${\kappa _\theta } = 1/{r_{min}}$. Here ${r_a}$ is the axial curvature radius of the neck point and ${r_{min}}$ is the azimuthal curvature radius (minimum neck radius). If the Laplace pressure $\mathrm{\Delta }p = \gamma ({\kappa _\theta } - {\kappa _x})$ is positive, i.e. ${r_{min}} < {r_c}$, the neck shrinkage will be promoted, finally leading to pinch-off. These studies confirmed that the partial coalescence is related to the vertical and horizontal motion of the neck. To identify this, we first plot the temporal vertical locations of the neck in figure 9 with different $\varGamma $ for $Oh = 0.02$ and $Oh = 0.03$. The figure clearly shows that the initial neck height ${z_{neck}}$ grows in a similar way for all cases. We note that the maximum neck height $\to 4$, and this corresponds to the initial height of the droplet's centre, denoting the upward motion of capillary waves. Then the neck continues to drain into the liquid reservoir, leading to pinch-off as expected for all cases except for $Oh = 0.03$ and $\varGamma = 0$. The electric field slightly increases the neck height and the effect has little dependence on the Ohnesorge number.

Figure 9. Temporal evolutions of the vertical location of the neck for (a) $\; Oh = 0.02$ and (b) $\; Oh = 0.03$.

As shown by figure 9, there is little difference in the vertical motion of the neck in the presence of the electric field. Hence, we turn our attention towards the horizontal motion (expansion and shrinkage) of the neck. The two types of neck curvature radii, ${r_{min}}$ and ${r_a}$, are defined in figure 10(a). Previous studies usually focus on the comparison of the two radii at a certain moment (Zhang et al. 2009; Alhareth & Thoroddsen 2020; Constante-Amores et al. 2021). In contrast to those studies, we plot the temporal evolutions of ${r_{min}}$ and ${r_a}$ during the whole coalescence process. Here the axial curvature radius ${r_a}$ in the neck is calculated by ${r_a} =- 1/(\boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{n})$. For the two cases with $Oh = 0.03$, complete coalescence $(\varGamma = 0)$ and partial coalescence $(\varGamma = 0.01)$ are used for comparison as shown in figures 10(c) and 10(d). During the first expansion stage, ${r_{min}} \gg {r_a}$, and the induced negative Laplace pressure $\mathrm{\Delta }p = \gamma (1/{r_{min}} - 1/{r_a})$ tends to expand the neck as depicted in figure 10(b). Further in time, both ${r_{min}}$ and ${r_a}$ start to decrease during the first shrinkage stage. For $\varGamma = 0$, the two curvature radii are almost equal at $t = 1.85$, but ${r_{min}} > {r_a}$ when $t > 1.5$, and this leads to the second neck expansion as shown in figure 10(b). Finally, at $t = 2.5$, ${r_{min}} \gg {r_a}$ is responsible for the complete coalescence. However, for $\varGamma = 0.01$, we find ${r_{min}} < {r_a}$ during most of the shrinkage stage so that the positive Laplace pressure squeezes the neck towards pinch-off at $t = 1.75$.

Figure 10. (a) Schematic of the neck region during the early coalescence process, with ${r_{min}}$ denoting the radius of curvature in the azimuthal direction and ${r_a}$ denoting the radius of curvature in the axial direction. (b) Schematic of the neck region during the expansion and shrinkage stages for different $\varGamma $. Comparison between the temporal evolutions of ${r_{min}}$ and ${r_a}$ during the whole coalescence processes for $Oh = 0.03$ with (c) $\varGamma = 0$ and (d) $\varGamma = 0.01$.

Figures 10(c) and 10(d) show that the sign of ${r_{min}} - {r_a}$ determines the outcomes of droplet–interface coalescence: ${r_{min}} - {r_a} > 0$ for complete coalescence and ${r_{min}} - {r_a} < 0$ for partial coalescence. However, from the results in figures 5(a) and 8(a) we know that ${r_{min}}$ will be enhanced during the final shrinkage stage as $\varGamma $ increases, and it seems to favour complete coalescence. Actually, the electric field promotes the partial coalescence by increasing the axial curvature radius ${r_a}$, which is attributed to the vertical elongation of the droplet driven by the electric force as depicted in figure 10(b). To verify this, we conduct extra simulations by omitting the electric stresses added on the upper part (for the definition, see figure 10a) of the droplet interface. Figure 11 displays the temporal evolutions of the minimum neck radius and the maximum height for three cases with $Oh = 0.03$. It is clearly seen that when the upper electric stresses are omitted, the electric field does not promote the vertical motion of the droplet, and this gives rise to the final complete coalescence even in the presence of an electric field. Moreover, neck shrinkage is suppressed under the remaining electric force in the neck region. For this reason, we can conclude that the vertical elongation of the droplet is indeed responsible for the occurrence of partial coalescence.

Figure 11. Temporal evolutions of (a) the minimum neck radius and (b) the maximum height of the droplet during the coalescence processes for $Oh = 0.03$. The red dashed line represents the case omitting the upper electric stresses.

4.2. Critical Ohnesorge number for partial coalescence

We now know that the occurrence of partial coalescence is determined by the interplay between the viscous effect and the electric field, i.e. the dimensionless Ohnesorge number $Oh$ and the electric Bond number $\varGamma $. For the electric-field-free case, the critical condition for partial coalescence has been clarified by Blanchette & Bigioni (2009). They determined the critical Ohnesorge number as a linear function of the Bond number if the density and viscosity differences of the two liquids are ignored. In this study, we focus on the liquid–air system and neglect the effect of gravity; thus the critical Ohnesorge number $O{h_c}$ is independent in the absence of an electric field, which is reported to be 0.029 here based on our simulations.

Before determining the critical Ohnesorge number $O{h_c}$ in the presence of an electric field, we first display the final shapes of cases where pinch-off occurs. Figure 12 shows the interface shapes at pinch-off for a wide range of $Oh$ and $\varGamma $. It is clearly seen that the droplet is elongated to deform into a prolate spindle shape with a sharp bottom as the external electric field increases. However, with increasing Ohnesorge number, the secondary droplet tends to be spherical with a thin thread connecting it to the interface. The thin thread is also observed in figure 8 and it may deform into a much smaller satellite droplet after pinch-off. Figure 12 also depicts the drop deformation aspect ratio D as a function of $Oh$ and $\varGamma $. The Taylor deformation parameter is defined as $D = (L - W)/(L + W)$, where L and W are the vertical length and horizontal width of the deformed droplet, respectively. Here $D = 0$ represents a spherical droplet. Figure 12 shows that larger deformation occurs with higher $\varGamma $ and lower $Oh$. Moreover, as $Oh\; $increases to the critical value for partial coalescence, the aspect ratio D decreases to ${\sim} 0.1$ for all electric-field-laden cases.

Figure 12. Interface shapes of the droplet–interface coalescence processes at pinch-off for different $Oh$ and $\varGamma $. The rightmost shape of each row represents the case with critical $Oh$. Inset shows the deformation aspect ratios D as a function of $Oh$ and $\varGamma $.

Based on the analysis above, it is concluded that the electric force promotes the vertical elongation of the droplet, and overcomes the capillary and viscous forces, resulting in the final pinch-off. In the critical state before pinch-off, the characteristic electric stress at the droplet apex can be expressed by

(4.2)\begin{equation}{F_E}\sim {\textstyle{1 \over 2}}{\varepsilon _0}E_0^2.\end{equation}

The results in figure 12 enable us to assume the secondary droplet in the critical state to be a spherical droplet with radius ${R_d}$. Thus the capillary pressure driving the inward flow can be given as

(4.3)\begin{equation}{F_\gamma }\sim \frac{\gamma }{{{R_d}}}\end{equation}

and the resistive viscous stress scales as

(4.4)\begin{equation}{F_v}\sim \frac{{{\mu _1}{U_p}}}{{{R_d}}},\end{equation}

where ${U_p}$ is the stretch velocity required for partial coalescence. We can determine the transition from partial coalescence to complete coalescence by balancing the interfacial stresses:

(4.5)\begin{equation}\frac{{{\mu _1}{U_p}}}{{{R_d}}} - \frac{1}{2}{\varepsilon _0}E_0^2\sim \frac{\gamma }{{{R_d}}}.\end{equation}

Following Hamlin et al. (2012), the apex velocity varies little during the coalescence and is approximately equal to $\sqrt {\gamma /{\rho _1}{R_0}} $, which is the characteristic velocity scale of inertia. Finally, substituting ${U_p} = \sqrt {\gamma /{\rho _1}{R_0}} $ in (4.5) and rearranging, we obtain

(4.6)\begin{equation}Oh - \varGamma \frac{{{R_d}}}{{{R_0}}}\; \sim \; 1.\end{equation}

Equation (4.6) prompts us to find the scale of ${R_d}/{R_0}$ to determine the criterion for partial coalescence. Hamlin et al. (2012) assumed that inertia balances the electric force when the pinch-off occurs: $Q{E_0}\; \sim \; {\rho _1}U_p^2R_d^2$, where Q is the charge remaining on the secondary droplet. Using $Q\; \sim \; {\varepsilon _0}{E_0}R_0^2$, one can obtain the scaling ${R_d}/{R_0}\; \sim \; {\varGamma ^{1/2}}$. Thus, (4.6) becomes

(4.7)\begin{equation}Oh - {\varGamma ^{3/2}}\; \sim \; 1.\end{equation}

Equation (4.7) demonstrates that $Oh$ and ${\varGamma ^{3/2}}$ have opposite effects on the coalescence process. Therefore, it is reasonable to determine the critical Ohnesorge number $O{h_c}$ as a positive correlation of ${\varGamma ^{3/2}}$. We then test this correlation against data from our numerical simulations. In figure 13, $O{h_c}$ is plotted as a function of ${\varGamma ^{3/2}}$ and the straight line $O{h_c} = 14.39{\varGamma ^{3/2}} + 0.029$ is in excellent agreement with the simulations results. It proves that this correlation can accurately capture the transition from partial coalescence to complete coalescence in the presence of an electric field.

Figure 13. Critical Ohnesorge number $O{h_c}$ as a function of the modified electric Bond number with $\varGamma $ ranging from 0 to 0.03. Here the straight line represents the scaling prediction for $O{h_c}$ beyond which partial coalescence does not occur.

4.3. Size of secondary droplet

In this section, we turn our attention to the size of the secondary droplet. We first display the pinch-off times for all cases. Figure 14(a) shows that the pinch-off time increases with increasing $Oh$ and increases sharply when approaching the critical state. Our results in figure 4 demonstrate that the pinch-off time increases as $\varGamma $ increases for $Oh = 0.02$, indicating that the pinch-off has been delayed by the electric field. However, figure 14(a) shows that increasing $\varGamma $ leads to an overall reduction in the pinch-off time when $Oh$ reaches a critical value. We note that the trend becomes inverse once the partial coalescence does not occur without the electric field. This can be attributed to the fact that the effect of the electric field on neck shrinkage becomes much weaker as $Oh$ reaches the critical value; therefore, partial coalescence is only promoted by the upward force caused by the electric field. As a result, the pinch-off time decreases due to the faster upward motion with higher $\varGamma $.

Figure 14. (a) Pinch-off time as a function of the Ohnesorge number for different $\varGamma $. (b) The secondary droplet size ratio ${R_d}/{R_0}$ versus the scaling prediction using the same data as in (a). Here $Oh$ ranges from 0.005 to the critical value; $\varGamma $ ranges from 0 to 0.03. The solid line has a slope equal to one.

Obviously, the ratio of the secondary droplet radius to that of the original droplet depends on both $Oh$ and $\varGamma $. In our simulations, we calculate ${R_d}$ using ${R_d} = {(3{V_d}/4{\rm \pi})^{1/3}}$, where ${V_d} = \int_{{z_{neck}}}^{{h_{max}}\; } {{\rm \pi}{r^2}\,\textrm{d}z} $ represents the volume of the secondary droplet. The studies of Hamlin et al. (2012) indicate that ${R_d}/{R_0}\sim {\varGamma ^{1/2}}$. However, we find this scaling usually holds in cases with high Ohnesorge numbers, where the secondary droplet tends to be spherical. In most cases, the secondary droplet is not spherical except in the critical state as depicted in figure 12.

Following Aryafar & Kavehpour (2006), we consider the volume of the secondary droplet to be the escaped volume of the original droplet during coalescence, resulting in

(4.8)\begin{equation}{\left( {\frac{{{R_d}}}{{{R_0}}}} \right)^3} = \frac{{{V_d}}}{{{V_0}}}\sim \frac{{{V_0} - {Q_v}t}}{{{V_0}}},\end{equation}

where ${Q_v}$ is the average flow rate of the liquid into the reservoir, t is the coalescence time and ${V_0}$ is the volume of the original droplet. We assume the average flow rate is the average cross-sectional area versus the characteristic velocity: ${Q_v} = {\rm \pi}R_0^2\sqrt {\gamma /{\rho _1}{R_0}} $. The coalescence time t can be considered to be the combination of the inertial and viscous time scales, giving

(4.9)\begin{equation}t = a\sqrt {\frac{{{\rho _1}R_0^3}}{\gamma }} + b\frac{{{\mu _1}{R_0}}}{\gamma },\end{equation}

where a and b are scaling constants. Substituting (4.9) into (4.8) produces

(4.10)\begin{equation}\frac{{{R_d}}}{{{R_0}}} = {\left( {1 - \frac{3}{4}a - \frac{3}{4}bOh} \right)^{1/3\; }}.\end{equation}

Equation (4.10) describes the size of the secondary droplet without the electric field. To incorporate the effect of the electric field on the scaling, we can follow the same method as introduced in § 4.2. Replacing the Ohnesorge number in (4.10) with the combination of $Oh$ and ${\varGamma ^{3/2}}$, and using our simulation results to fit the constants, the size ratio finally becomes a function of both the Ohnesorge number and the electric Bond number:

(4.11)\begin{equation}\frac{{{R_d}}}{{{R_0}}} = {[0.22 - 2.5(Oh - 9.53{\varGamma ^{3/2}})]^{1/3}}.\end{equation}

All of the secondary droplet sizes from figure 14(a) are replotted in figure 14(b) versus (4.11). We see that the scaling prediction performs well over a wide range. Examination of figure 14(b) suggests that the size range of the secondary droplet expands with an increase of the electric field. This enables us to produce droplets of desired sizes, which could be useful in microfluidic applications.

5. Conclusion

In the present study, we have numerically investigated the effect of an external electric field on the axisymmetric coalescence between a droplet and a liquid–air interface. A leaky-dielectric model is fully coupled to Navier–Stokes equations in which the Maxwell stresses are directly imposed across the liquid–air interface. The numerical method has been validated against the experimental observation from Aryafar & Kavehpour (2009). We observe that both the initial neck expansion and the eventual neck shrinkage are suppressed by the electric field. The former is due to the contribution of compressive electric stresses in the neck region during the early coalescence, whereas the latter is induced by the outward electric stresses, resulting in the suppression of the final pinch-off. These effects become weaker as the Ohnesorge number $Oh$ increases. For smaller $Oh$, increasing the electric field delays the pinch-off. However, as $Oh$ reaches a critical value above which partial coalescence does not occur, increasing the electric field leads to an overall reduction in the pinch-off time. We have also shown that if the electric field is strong enough, tip streaming will occur at the apex of the droplet before its bottom pinch-off.

Moreover, we find that an increase in electric field strength favours the occurrence of partial coalescence. We demonstrate that this is attributed to the vertical elongation of the droplet caused by the electric force. As a result, the axial curvature radius is increased, and the induced positive Laplace pressure squeezes the neck towards pinch-off. Using the scaling analysis, we report a critical Ohnesorge number $O{h_c} = 14.39{\varGamma ^{3/2}} + 0.029$ to quantify the transition from partial coalescence to complete coalescence in the presence of an electric field, where $\varGamma $ is the dimensionless electric Bond number. Moreover, the ratio of the size of the secondary droplet to that of the original droplet has been determined through our numerical results, which could be useful in microfluidic applications where droplets of desired sizes are needed.