1. Introduction
Coalescence of drops is one of the most common capillarity-driven phenomena which can be observed in multiphase fluid dynamics. The early-time dynamics of coalescence is dependent on both the viscosity of the drops and their geometry. Different power laws for the growth of the connecting structure (referred to as neck or bridge) have been found for viscous and inviscid freely suspended drops (Eggers, Lister & Stone Reference Eggers, Lister and Stone1999; Duchemin, Eggers & Josserand Reference Duchemin, Eggers and Josserand2003; Aarts et al. Reference Aarts, Lekkerkerker, Guo, Wegdam and Bonn2005; Thoroddsen et al. Reference Thoroddsen, Qian, Etoh and Takehara2007; Paulsen, Burton & Nagel Reference Paulsen, Burton and Nagel2011), as well as for sessile drops in the viscous and inviscid limits (Ristenpart et al. Reference Ristenpart, McCalla, Roy and Stone2006; Narge, Beysens & Pomeau Reference Narge, Beysens and Pomeau2008; Hernández-Sánchez et al. Reference Hernández-Sánchez, Lubbers, Eddi and Snoeijer2012; Lee et al. Reference Lee, Kang, Yoon and Yarin2012; Eddi, Winkels & Snoeijer Reference Eddi, Winkels and Snoeijer2013). The study of coalescence phenomena is also relevant for many applications where the underlying substrate of the coalescing drops is a liquid. Some examples are wet-on-wet printing (Hack et al. Reference Hack, Costalonga, Segers, Karpitschka, Wijshoff and Snoeijer2018), emulsions (Shaw Reference Shaw2003; Kamp, Villwock & Kraume Reference Kamp, Villwock and Kraume2016) and lubricant impregnated substrates (Anand et al. Reference Anand, Paxson, Dhiman, Smith and Varanasi2012; Smith et al. Reference Smith, Dhiman, Anand, Reza-Garduno, Cohen, McKinley and Varanasi2013).
Here, we focus on liquid lenses (de Gennes, Brochard-Wyart & Quéré Reference de Gennes, Brochard-Wyart and Quéré2004), consisting of liquid drops floating on a quiescent pool of another liquid. This case was studied for Newtonian drops (Burton & Taborek Reference Burton and Taborek2007) and liquid crystals (Delabre & Cazabat Reference Delabre and Cazabat2010), where the authors analysed the growth of the bridge in top-view experiments. Recent work considered the coalescence of lenses using side-view experiments (Hack et al. Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020). This perspective is sketched in figure 1, providing a quasi-two-dimensional view of the problem. The experiments revealed a self-similar dynamics of the bridge profiles, with scaling laws for the bridge height 
$h_0$ with time 
$t$ that depend on the viscosity of the lenses (Hack et al. Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020). At low viscosity, the dominant balance during coalescence is between surface tension and inertia, and it was found that 
$h_0\sim t^{2/3}$. At high viscosity, the dominant balance between surface tension and viscosity leads to 
$h_0 \sim t$. These scaling laws are the same as those described in the merging of liquid wedges (Billingham & King Reference Billingham and King2005). Owing to the slender geometry of the drops – typically the contact angle 
$\theta$ in figure 1 is small – the coalescence of liquid lenses can be analysed using the thin-sheet equations (Ting & Keller Reference Ting and Keller1990; Erneux & Davis Reference Erneux and Davis1993). Using a similarity analysis, the experimentally observed inertial and viscous scaling laws are recovered (Hack et al. Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020). For example, the viscous coalescence speed was found to be 
$\mathrm {d}\bar {h}_0/\mathrm {d}\bar {t} \approx \bar {V}_0 = 0.5525 \gamma \theta ^2/\eta$, where 
$\gamma$ and 
$\eta$ respectively are the drop surface tension and viscosity (for consistency we here use overbars to indicate dimensional variables). This prediction was found to be in very good agreement with experiments.
Side-view sketch of two coalescing lenses initial radius 
$R$ and equilibrium contact angle 
$\theta$. The drop profile is described by 
$\bar h(\bar x,\bar t)$, while the minimal height in the neck region 
$\bar h_0(t)=\bar h(0,\bar t)$. The coalescence velocity 
$\bar V$ is defined as 
$\bar {V}=\mathrm {d}\bar {h}_0/\mathrm {d}\bar {t}$.
The viscous similarity analysis, however, contains a salient feature that remains to be explained: the obtained self-similar velocity profile does not decay at large distance from the thin bridge region, but reaches a finite value. This is rather unusual for problems involving coalescence (or drop breakup, cf. Eggers & Fontelos Reference Eggers and Fontelos2015). Namely, the similarity analysis is typically based on the assumption that the flow remains confined to the scale of the bridge – at large scale, i.e. the scale of the drop, the flow is usually assumed to vanish. Such is the case for the coalescence of sessile drops as illustrated in figure 2(a). It shows an experimental top-view sequence of merging drops that are in contact with a solid substrate. During the initial growth of the bridge the global features of the drops appear nearly stationary – away from the bridge region one observes only a minute spreading of the drops. Figure 2(b) shows the equivalent top-view sequence for liquid lenses, for which the situation is manifestly different. Clearly, these floating drops do not remain stationary, but their centres of mass exhibit an inward motion as soon as the drops establish contact. This inward motion is not a small effect. Figure 2(c) compares the centre of mass velocity 
$\bar {U}$ (measured in top view) with the bridge coalescence velocity 
$\bar {V} = \mathrm {d}\bar {h}_0/\mathrm {d}\bar {t}$ (measured in side view), showing that the two velocities are proportional. These velocities were obtained using the experimental method described in Hack et al. (Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020). We remark that this inward motion does not at all arise for liquid lenses of very low viscosity – this is in line with the inviscid similarity solutions, whose velocity rapidly decays away from the bridge (Hack et al. Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020).
Top-view coalescence of drops in the viscous regime. (a) Sessile drops in contact with a substrate, (b) liquid lenses floating on a water bath. During lens coalescence we see a clear inward motion of the drops that does not occur for sessile drops. (c) Comparison of the centre of mass velocity 
$\bar {U}$ and the ‘bridge’ coalescence velocity 
$\bar {V}=\mathrm {d}\bar {h}_0/\mathrm {d}\bar {t}$ for liquid lenses, taken for three oils of different viscosity and nearly identical contact angles (blue: 
$\eta =115\ \textrm {Pa}\,\textrm {s}$ and 
$\theta =27^\circ$, green: 
$\eta =33\ \textrm {Pa}\,\textrm {s}$ and 
$\theta =32^\circ$, red: 
$\eta =9\ \textrm {Pa}\,\textrm {s}$ and 
$\theta =31^\circ$, see Hack et al. (Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020) for experimental details). The centre of mass velocity 
$\bar {U}$ is not a small effect, as it is comparable in magnitude to the coalescence velocity 
$\bar {V}$.
In this paper, we provide a detailed analysis of the coalescence of highly viscous lenses and elucidate the coupling between the inner ‘bridge’ solution and the global dynamics of the drops. We treat the problem using the two-dimensional thin-sheet equations, reflecting an analysis along the cross-section shown in figure 1. Such a two-dimensional approximation turned out successful for the geometry of spherical caps (Ristenpart et al. Reference Ristenpart, McCalla, Roy and Stone2006; Hernández-Sánchez et al. Reference Hernández-Sánchez, Lubbers, Eddi and Snoeijer2012; Hack et al. Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020), since close to the bridge the length scale in the third dimension, 
$O((R \bar {h}_0/\theta )^{1/2})$, is much larger than the horizontal and vertical length scales of the bridge, 
$O(\bar {h}_0/\theta )$ and 
$O(\bar {h}_0)$.
We demonstrate that, in general, the coalescence velocity exhibits significant logarithmic corrections, 
$O(1/ \ln {t})$, where 
$t\ll 1$ is the (dimensionless) time after coalescence. On the other hand, the thin-sheet equations admit an outer solution where the drop's centre of mass can migrate freely, closely resembling the motion observed in figure 2(b). In this latter case, the corrections to the leading-order result are much smaller. The analysis is confirmed in detail by comparison with time-dependent numerical simulations of the thin-sheet equations.
The article is organised as follows: the governing equations for coalescing lenses are presented in § 2. In § 3, we study the inner scale dynamics, near the point of coalescence, followed in § 4 by an analysis of the outer region, which is where the two different boundary conditions manifest themselves. We end with our conclusion and outlook in § 5.
2. The viscous thin-sheet equations
2.1. Formulation
Following the approach of Hack et al. (Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020), the process of coalescence is modelled by the two-dimensional viscous thin-sheet equations (Ting & Keller Reference Ting and Keller1990; Erneux & Davis Reference Erneux and Davis1993). The underlying approximations are the following: (i) similar to sessile drops (Ristenpart et al. Reference Ristenpart, McCalla, Roy and Stone2006; Hernández-Sánchez et al. Reference Hernández-Sánchez, Lubbers, Eddi and Snoeijer2012) the flow in the bridge region is quasi-two-dimensional in the early stage of coalescence; (ii) the equilibrium contact angle 
$\theta$ is small, such that a slender body approximation can be employed; (iii) the influence of the bath on the dynamics is negligible, i.e. free slip boundary conditions can be employed at both interfaces of the two-dimensional lenses; (iv) due to negligible differences in surface tension between the bath and liquid lens and the liquid lens and air, the liquid lenses are assumed to be symmetrical with respect to the bath–air interface – although asymmetric surface tensions can actually be mapped to an ‘effective’ symmetric surface tension (cf. Hack et al. Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020, Supplementary Material). We note that these assumptions are in accordance to the previous experiments of Hack et al. (Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020) where the ratio of viscosities of the lenses and the bath is more than a thousand and the corresponding surface tension asymmetry is only approximately ten per cent. The former justifies (iii) and the latter (iv).
The resulting model equations for negligible inertia of the flow, in dimensionless variables, take the form
\begin{equation} \dot h + (hu)' = 0, \quad h h''' + 4 (h u')' = 0. \end{equation}
Here, 
$h(x,t)$ and 
$u(x,t)$, respectively, are the dimensionless interface height and horizontal velocity, and dots and primes indicate derivatives with respect to 
$t$ and 
$x$, respectively. Dimensional variables 
$(\bar {x},\bar {h},\bar {u},\bar {t})$ are scaled as
\begin{equation} \bar{x}=R x, \quad \bar{h}=\theta R h, \quad \bar{u}=\frac{\gamma}{\eta}\theta u, \quad \bar{t}= \frac{R \eta}{\theta\gamma}t, \end{equation}
where 
$R$ is the initial lens radius and 
$\theta$ is the contact angle (cf. figure 1). The surface tension is denoted as 
$\gamma$ and 
$\eta$ is the viscosity of the liquid inside the lenses. The thin-sheet equations (2.1) correspond to mass and (horizontal) momentum conservation. The latter gives the balance between capillary forces (first term) and viscous forces (second term), while inertia has here been neglected.
The momentum equation (2.1b) can be readily integrated to give the horizontal force balance
\begin{equation} h h'' - \tfrac{1}{2} h^{\prime 2} + 4 u'h ={-}\tfrac{1}{2} + F( t). \end{equation}
The terms on the left-hand side represent the horizontal force transmitted through the thin drop by pressure, surface tension and viscous stresses, respectively. We have introduced the constant 
$-1/2$ on the right-hand side corresponding to the total force in the static solution, see (2.6) below. Therefore 
$F(t)$ measures any additional horizontal force that arises during the coalescence dynamics.
2.2. Two-dimensional numerical simulations
In order to illustrate the interplay between the dynamics of the small bridge region and the large bulk of the drops, we simulate the coalescence of two-dimensional lenses numerically. Due to symmetry about the point of coalescence 
$x=0$, we can impose the boundary conditions 
$h' = u = 0$ at 
$x = 0$ and focus on the domain 
$x \geq 0$. We consider two different sets of outer boundary conditions.
The first case corresponds to the experimentally realised setup of two freely floating lenses. In this case, the length 
$L(t)$ of each lens decreases with time (as can be seen in figure 2) from its initial value 
$L(t=0) = 2$, and at the edge of the lens we impose the thickness 
$h(L)$, the contact angle and conservation of mass
\begin{equation} h(L) = 0,\quad h'(L) ={-}1, \quad u(L) = \dot{L}. \end{equation}
Substitution into the momentum equation (2.3) then yields 
$F(t) \equiv 0$, i.e. the horizontal force transmitted through the lens is equal to its initial value at all times. With 
$F=0$ in (2.3), the condition 
$h'(L) = -1$ is redundant as it follows from 
$h(L) = 0$, so we are left with three boundary conditions 
$h'(0) = u(0) = h(L) = 0$ for the third-order governing equations and a fourth condition 
$\dot L = u(L)$ that determines the evolution of 
$L(t)$.
The second case we consider is one where the lenses do not move, due to a symmetry condition being imposed about the centres of the lenses, such as in a periodic array of simultaneously coalescing lenses,
\begin{equation} h'(1) = u(1) = 0. \end{equation}
This yields three boundary conditions 
$h'(0) = u(0) = h'(1) = 0$ on the governing equations and a fourth condition 
$u(1) = 0$ that determines the unknown 
$F(t)$.
The key finding is that the different outer boundary conditions will lead to different spatial and temporal dynamics also in the inner region, at the scale of the bridge. As we shall see, the leading-order solution in the two cases remains the same – however, the coalescence velocity can exhibit significant corrections at the next order.
The theoretical initial condition is taken to be the static solution of (2.1) with non-dimensional radius 
$1$ and contact angle 
$1$ (corresponding to the dimensional values 
$R$ and 
$\theta$, respectively),
\begin{equation} h(x,t=0) = h_s(x) = \tfrac{1}{2} x (2 - x), \end{equation}
but in the numerical simulations a small-scale perturbation 
$(h_{0i}-x/2)\exp (-x/2h_{0i})$ is added to initialise the coalescence. The perturbation profile is chosen to satisfy the symmetry boundary condition 
$h'(0) = 0$ and have an initial bridge height 
$h_{0i} = h_0(t=0)$, which is taken to be 
$10^{-10}$ unless stated otherwise.
We solve the governing equations (2.1a) and (2.3) numerically in Matlab using a finite difference method with a Crank–Nicolson time-stepping scheme. Both the spatial and temporal grids are non-uniform with a relative resolution of 
$1\,\%$, and the discretisation is second-order accurate in both space and time. For the periodic lenses (2.5), the computational domain is 
$0 \leq x \leq 1$. For the free lenses (2.4), where the domain 
$0 \leq x \leq L(t)$ changes with time, we use a rescaled position variable 
$\hat {x} = 2x/L(t)$ and solve the rescaled equations on the fixed domain 
$0 \leq \hat {x} \leq 2$ instead. As explained above, four boundary conditions are imposed on the third-order system of equations, in order to determine the additional unknown 
$\dot L(t)$ or 
$F(t)$.
Some snapshots of height and velocity profiles from the simulations are shown in figure 3. The free-floating lenses immediately begin to move towards each other at a constant velocity, and eventually merge into one larger lens with double the volume (and hence 
$\sqrt {2}$ times the radius and height). The periodic lenses instead flatten out towards a uniform height of 
$1/3$, determined by the initial volume in the lens. Importantly, and this is one of the central points of the paper, the velocity profiles are significantly different between the two cases, even at very early times, due to the different outer boundary conditions. We will show how this affects the coalescence velocity, 
$\dot {h}_0$.
3. The inner region
In order to study the dynamics of the region close to the point of coalescence, we rescale the variables as
\begin{equation} \xi = \dfrac{x}{h_0 (t)},\quad h(x, t) = h_0 ( t )\ \mathcal{H}( \xi , \tau
), \quad u(x, t) = \dot{h}_0 ( t )\ \mathcal{U}(\xi,
\tau ),
\end{equation}
with 
$\tau = t$. In choosing the scaling used above, we are motivated by the fact that in the inner region the length scale of importance is the bridge height 
$h_0$. The resulting governing equations become
$$\begin{gather} \frac{h_0}{V}\dot{\mathcal{H}} + \mathcal{H} - \xi \mathcal{H}' + (\mathcal{H}\mathcal{U})' = 0, \end{gather}$$
$$\begin{gather}\mathcal{H}\mathcal{H}'' - \tfrac{1}{2} \mathcal{H}'^2 + 4 V \mathcal{H} \mathcal{U}' + \tfrac{1}{2} - F(\tau) = 0, \end{gather}$$
where 
$V = \dot h_0$ is the unknown coalescence velocity that we wish to determine. At 
$\xi = 0$, the definition 
$h(0,t) = h_0(t)$ and the symmetry about 
$\xi = 0$ yield
\begin{equation} \mathcal{H} (0) = 1, \quad \mathcal{H}'(0) = \mathcal{U}(0) = 0. \end{equation} The system of equations (3.2) has three spatial derivatives and accordingly three boundary conditions (3.3). However, the equations contain two unknowns 
$V$ and 
$F$ that are determined on matching to the outer solution. The matching necessitates that we evaluate 
$\mathcal {H}$ and 
$\mathcal {U}$ as 
$\xi \rightarrow \infty$, and from (3.2) we obtain, on neglecting the time-derivative term,
$$\begin{gather} \mathcal{U} = \frac{F}{4\alpha V}\ln \xi - \beta + O(\ln\xi/\xi), \quad \mathcal{H} = \alpha \xi + O(\ln \xi), \quad \text{as } \xi\to\infty. \end{gather}$$
This asymptotics contain two degrees of freedom 
$\alpha$ and 
$\beta$ that are determined during the solution process. (Given that the system is third order in space, there should be a further degree of freedom as 
$\xi \to \infty$. However, it can be shown by perturbation analysis of (3.2) (cf. Eggers & Fontelos Reference Eggers and Fontelos2015), that this degree of freedom appears as a prefactor of a term decaying through an exponential in 
$\xi$ as long as 
$\alpha >0$, which is always the case in this coalescence problem. The degree of freedom does thus not appear in our leading-order asymptotics.) This matching will be performed explicitly in § 4, but below we already anticipate some of the results necessary to evaluate the inner solutions.
3.1. The leading-order similarity solution
The leading-order solution to the foregoing equations is a steady self-similar solution, which was previously obtained by Hack et al. (Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020). The leading-order quantities 
$\mathcal {H}_0(\xi )$, 
$\mathcal {U}_0(\xi )$ and 
$V_0$ are independent of 
$\tau$, so that (3.2) reduces to
$$\begin{gather} \mathcal{H}_0 - \xi \mathcal{H}_0' + (\mathcal{H}_0\mathcal{U}_0)' = 0, \end{gather}$$
$$\begin{gather}\mathcal{H}_0\mathcal{H}_0'' - \tfrac12 \mathcal{H}_0^{'2} + 4 V_0 \mathcal{H}_0 \mathcal{U}_0^{'} + \tfrac{1}{2} = 0, \end{gather}$$
where we have anticipated that 
$F(\tau ) \ll 1$ – this indeed is true as we see later. The boundary conditions required to evaluate these quantities are given by (3.3), complemented by 
$\mathcal {H}'_0 = 1$ as 
$\xi \rightarrow \infty$ (
$\alpha = 1$) set by the leading-order outer solution given in (2.6). These equations (and also the higher-order equations) are solved numerically in Mathematica using a shooting method on the domain 
$0 \leq \xi \leq 10^5$. We find the leading-order coalescence and far-field translational velocities to be (Hack et al. Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020)
\begin{equation} V_0 = 0.5525, \quad \mathcal{U}_0 \to -\beta_0 ={-}0.5734 \text{ as }\ \xi\to\infty. \end{equation}3.2. Next-order corrections
The solutions 
$\mathcal {H}_0$ and 
$\mathcal {U}_0$ are plotted in figure 4(
$a$,
$b$), where they are compared with the numerical solutions. The two columns correspond to the two types of boundary conditions, free-floating lenses and periodic lenses. For both cases an excellent agreement is found at early times after coalescence. However, the velocity profiles in figure 3(a,b) revealed large differences between these two situations. While these differences do not turn up in the leading-order inner solutions 
$\mathcal {H}_0$ and 
$\mathcal {U}_0$, they impact the next-order corrections. As we will see, this also has a strong effect on the evolution of the coalescence velocity 
$V$. We therefore now address the inner solution beyond the leading order.
Similarity solution profiles – leading order (3.6) and first correction ((3.9) and (3.12)) – for (
$a,c$) free-floating lenses and (
$b{,}d$) periodic lenses. Rescaled numerical results are shown for comparison, evaluated at (
$a{,}c$) 
$t = 10^{-3},\ 10^{-2},\ 10^{-1}$ and (
$b{,}d$) 
$t = 10^{-4},\ 10^{-3},\ 10^{-2}$. The numerical height profile is transformed as 
$\mathcal {H}_{0,{num}} = h/h_0$ and 
$\mathcal {H}_{1, num} = (h/h_0 - \mathcal {H}_{0, asy})/\epsilon$ where 
$\epsilon = h_0$ for free lenses and 
$\epsilon = F$ for periodic lenses, and the velocity is transformed similarly. (In order to resolve the 
$O(10^{-3})$ corrections in (
$c$), the numerical simulation was performed with a resolution of 
$0.1\,\%$, which restricted 
$h_{0i}$ to the larger value 
$10^{-6}$.) The arrows indicate increasing time.
3.2.1. The case 
$F(\tau ) \equiv 0$
In the case of free-floating lenses, for which 
$F(\tau )\equiv 0$, it turns out that we can set up a consistent expansion based on the expansion parameter 
$h_0(\tau )$. On substituting the expansions
$$\begin{gather} \{\mathcal{H}(\tau),\mathcal{U}(\tau),V(\tau), \beta(\tau)\} = \{\mathcal{H}_0,\mathcal{U}_0,V_0, \beta_0\} + h_0(\tau)\,\{\mathcal{H}_1,\mathcal{U}_1,V_1, \beta_1\} + O(h_0^2), \end{gather}$$
(where the dependence of 
$\mathcal {H}$ and 
$\mathcal {U}$ on 
$\xi$ is understood) into (3.2)–(3.3) we obtain at 
$O(h_0)$,
$$\begin{gather} 2\mathcal{H}_1 - \xi\mathcal{H}_1' + (\mathcal{H}_0\mathcal{U}_1 + \mathcal{H}_1\mathcal{U}_0)' = 0, \end{gather}$$
$$\begin{gather}\mathcal{H}_0\mathcal{H}_1'' + \mathcal{H}_1 \mathcal{H}_0'' - \mathcal{H}_0'\mathcal{H}_1' + 4(V_0\mathcal{H}_0 \mathcal{U}_1' + V_0\mathcal{H}_1\mathcal{U}_0' + V_1\mathcal{H}_0\mathcal{U}_0') = 0, \end{gather}$$
$$\begin{gather}\mathcal{H}_1(0) = \mathcal{H}_1'(0) = \mathcal{U}_1(0) = 0. \end{gather}$$ The differential equations for 
$\mathcal {H}_1$ and 
$\mathcal {U}_1$ are of third order with three boundary conditions, so one further condition is required in order to determine the unknown 
$V_1$. This is obtained from a matching with the outer region. As 
$\xi \to \infty$, the generic solutions to (3.8) are quadratic, 
$\mathcal {H}_1 \propto \xi ^2$, and so will need to match the second derivative 
$h_s''(0) = -1$ (an 
$O (h_0)$ quantity when expressed in inner variables) of the outer solution. We thus impose 
$\mathcal {H}_1'' \to -1$ as 
$\xi \to \infty$, and obtain the numerical solutions plotted in figure 4(
$c$), with the coefficients
\begin{equation} V_1 ={-}0.7625, \quad \mathcal{U}_1 \to -\beta_1 ={-}0.3492, \quad \text{as } \xi \to\infty. \end{equation}3.2.2. The case of non-zero $F(\tau )$
When the horizontal force 
$F(\tau )$ is non-zero, the next-order corrections are at 
$O(F)$. As we shall see in the next section, these corrections dominate the 
$O( h_0)$ corrections that arise due to time evolution of the outer solution, since 
$F$, though unknown yet, is evaluated to be 
$O(1/\ln (h_0))$. For this reason, we now expand the solution in terms of 
$F$, i.e.
\begin{align} \left\{ \mathcal{H}(\tau), \mathcal{U}(\tau), V(\tau), \beta(\tau) \right\} = \left\{ \mathcal{H}_0, \mathcal{U}_0, V_0, \beta_0 \right\} + F(\tau) \left\{ \mathcal{H}_1, \mathcal{U}_1, V_1, \beta_1 \right\} + O(F^2). \end{align}
After substituting into the governing equations, we obtain, at 
$O(F)$,
$$\begin{gather} \mathcal{H}_1 - \xi \mathcal{H}_1' + (\mathcal{H}_0 \mathcal{U}_1 + \mathcal{H}_1 \mathcal{U}_0)' = 0, \end{gather}$$
$$\begin{gather}\mathcal{H}_0 \mathcal{H}_1'' + \mathcal{H}_1 \mathcal{H}_0'' - \mathcal{H}_0' \mathcal{H}_1' + 4 (V_0 \mathcal{H}_0 \mathcal{U}_1' + V_0 \mathcal{H}_1 \mathcal{U}_0' + V_1 \mathcal{H}_0 \mathcal{U}_0') = 1, \end{gather}$$
$$\begin{gather}\mathcal{H}_1'(0) = \mathcal{U}_1(0) = \mathcal{H}_1(0) = 0. \end{gather}$$
Since the corrections to the outer solution at 
$O(h_0)$ are to be neglected when matching the inner and outer solutions at 
$O(F)$, we now impose the condition 
$\mathcal {H}_1' \to 0$ as 
$\xi \to \infty$. On solving these equations numerically, we obtain the solution profiles plotted in figure 4(
$d$) and the coefficients
\begin{equation} V_1 ={-}0.6227, \quad \beta_1 = 0.7079. \end{equation} As 
$F = O(1/\ln h_0)$ will typically not be very small, it is useful to calculate the second-order 
$O(F^2)$ correction, which includes contributions proportional to 
$F^2$ as well as contributions proportional to 
$(h_0/\dot h_0)\dot F$. We perform this calculation in Appendix A.
4. The outer region and matching
We now turn to the outer solution, and calculate the corrections to the initial static profile 
$h_s$ (2.6) that are generated by the coalescence in the bridge region. Time integration of the evolution equation (2.1a) from this initial condition yields
\begin{equation} h(x,t) = h_s(x) - \int_0^t \left(u(x,\hat t) h(x,\hat t) \right)'\,\mathrm{d}\hat t. \end{equation}
The outer solution, where 
$h,u,x = O(1)$, is thus given by the initial profile 
$h_s(x)$ with an 
$O(t) = O(h_0)$ correction.
Integration of the momentum equation (2.3) using 
$h = h_s + O(h_0)$ then yields the outer-region velocity profile
\begin{equation} u_o = \dfrac{F(t)}{4} \ln{\dfrac{x}{2-x}} + B(t) + O(h_0). \end{equation}
The coefficients 
$B(t)$ and 
$F(t)$ are determined by boundary conditions and by matching to the inner solution.
4.1. The case of freely floating lenses
For coalescing free lenses, the boundary conditions (2.4) result in no additional horizontal force, 
$F(t) \equiv 0$, and hence the leading-order velocity profile (4.2) is a uniform translation. The translation velocity is given by the far-field behaviour (3.6) of the inner solution, so at leading order 
$u_o \approx -V_0\beta _0 = -0.3168$, which is independent of time.
The first corrections in the outer solution then come in at 
$O(t) = O(h_0)$. Time integration of the mass conservation equation (2.1a) yields the result
\begin{equation} h(x,t) = h_s(x) - t u_o h_s'(x) + O(h_0^2) = h_s(x - u_o t) + O(h_0^2), \end{equation}
which reveals that the first-order correction simply represents a translation of the initial profile 
$h_s(x)$ by the steady leading-order velocity 
$u_o < 0$. (In fact, it can be shown that to all orders in 
$h_0$, the velocity profile is spatially uniform and hence the drop is undergoing pure translation, with deformation only occurring in the bridge region 
$x = O(h_0)$.)
In order to match with the bridge region, we substitute 
$x = h_0 \xi$ into (4.3) and expand in powers of 
$h_0$, making use of 
$t u_o \approx -V_0 \beta _0 t \approx -\beta _0 h_0$. This yields
\begin{equation} h = h_0 \left[ \left(\xi - h_0\frac{\xi^2}{2}\right) + \beta_0 \left(1 - h_0 \xi\right) + \cdots \right]. \end{equation}Comparing this with the definition of the inner variables (3.1), we now identify the appropriate far-field behaviour of the inner solution to be
\begin{equation} \mathcal{H}_1(\xi) \sim{-}\frac{\xi^2}{2} - \beta_0 \xi + \cdots \quad \text{as} \quad \xi \to \infty. \end{equation}
This yields the condition 
$\mathcal {H}_1'' \to -1$ as 
$\xi \to \infty$ that we anticipated in § 3.2.1, which was imposed to obtain the solution (3.9).
The resulting expression for the coalescence velocity is, from (3.6a) and (3.9a),
\begin{equation} \dot h_0( t) = 0.5525 - 0.7625 h_0 + O(h_0^2), \end{equation}
predicting a linear correction to the coalescence velocity. This prediction is tested quantitatively in comparison with numerical simulations in figure 5. The result (4.6) is shown as a dashed line labelled ‘free’, while the solid lines are numerics that were initialised at three different initial bridge heights 
$h_{0i}$. After a short time, the numerical curves for free lenses rapidly converge to the predicted asymptotics (4.6). Note that in typical experimental conditions where 
$h_0 \sim 10^{-3}\cdots 10^{-2}$, the coalescence velocity is very close to the leading-order value 
$V_0$.
The dependence of the coalescence velocity 
$V=\dot h_0$ on the bridge height 
$h_0$, for both free-floating (2.4) and periodic (2.5) lenses. Numerical results for three different values of the initial bridge height 
$h_{0i}$ are shown (solid lines). The asymptotic results for the free-floating and periodic lenses are given (dashed lines), respectively, by (4.6) and (4.10).
4.2. The case of periodic lenses
The periodic lenses come with the symmetry boundary condition (2.5) on the velocity profile (4.2). This yields 
$B(t)\equiv 0$ and hence
\begin{equation} u_o = \dfrac{F(t)}{4} \ln{\dfrac{x}{2-x}}. \end{equation}
This logarithmic velocity at small 
$x$ can indeed be matched to the velocity at large 
$\xi$ from the inner region calculated in § 3.2.2. Specifically, we equate (4.7) (taking the limit 
$x \ll 1$) with the far-field inner velocity profile (3.4a) (using 
$u = V\mathcal {U}$ and 
$\xi = x/h_0$),
\begin{equation} \dfrac{F(t)}{4} \ln{\dfrac{x}{2}} = \dfrac{F(t)}{4} \ln{\dfrac{x}{h_0}} - V \beta, \end{equation}
where it is noted that 
$V$ and 
$\beta$ are expanded in 
$F$ themselves.
Equation (4.8) finally enables us to express 
$F$ in terms of 
$h_0$, which to leading order gives the result
\begin{equation} F = \frac{4V\beta}{\ln(2/h_0)} = \frac{4V_0\beta_0}{\ln(2/h_0)} + O(\ln(2/h_0)^{{-}2}). \end{equation}This confirms the slow, logarithmic decay of the force term, suggesting strong corrections with respect to the leading-order similarity solution. Most importantly, this leads us to the sought-after coalescence velocity for periodic drops
\begin{equation} \dot h_0 = V_0 + V_1 F + O(F^2) = V_0 + \frac{\hat{V}_1}{\ln(2/h_0)} + O(\ln(2/h_0)^{{-}2}), \end{equation}
where 
$V_0 = 0.5525$, 
$V_1 = -0.6227$ and 
$\hat{V}_1 = 4V_0\beta _0V_1 = 0.7892$. The 
$O(\ln (2/h_0)^{-2})$ correction to the coalescence velocity is calculated in Appendix A.
Once again, the asymptotic prediction is tested quantitatively in comparison with numerical simulations in figure 5. The result (4.10) is shown as a dashed line labelled ‘periodic’, while the solid lines are numerics that were initialised at three different initial heights 
$h_{0i}$; again the numerics are in excellent agreement with the asymptotics. We notice a dramatic difference between the coalescence velocities for periodic drops as compared with the free drops, owing to the logarithmic decay of 
$F$. At values where 
$h_0 \sim 10^{-3}\cdots 10^{-2}$, the coalescence velocity for periodic drops is significantly below the leading-order value 
$V_0$.
5. Conclusion
In the present work, we have analysed the coalescence dynamics of viscous liquid lenses using both matched asymptotics and numerical simulations. We restricted our attention to the case where the flow in the bath is negligible, which is a consistent approximation for lenses of high viscosity, and where the contact angles are small to enable a slender thin-sheet description. In addition, following a previously used assumption in coalescence, we treated the problem as quasi-two-dimensional. The common scenario in coalescence and pinch-off is that the flow remains localised into the narrow bridge region, while the far field remains stationary. Here, we have found that this is not the case for viscous lens coalescence, and demonstrated that the bridge region affects the global dynamics, using both matched asymptotics and numerical simulations.
For freely floating two-dimensional viscous lenses, as soon as the lenses start coalescing there is a motion of the outer contact line towards the point of coalescence. In dimensional form, the growth rate of the bridge height 
$\bar {h}_0(\bar t) = \bar {h}(\bar x = 0, \bar {t})$ (4.6), and the horizontal translation velocity of the drop's centre of mass 
$\bar {u}_o$ ((3.6) and (3.9)) can be written as
$$\begin{gather} \bar V = \dot{\bar{h}}_0 = \left(0.5525 + 0.7625 \frac{\bar{h}_s''(0)\bar{h}_0}{\theta^2} + \cdots\right) \frac{\gamma\theta^2}{\eta}, \end{gather}$$
$$\begin{gather}\bar U = \bar{u}_o ={-}\left(0.3168 - 0.2443 \frac{\bar{h}_s''(0)\bar{h}_0}{\theta^2} + \cdots\right) \frac{\gamma\theta}{\eta}. \end{gather}$$
Here, we have eliminated the lens radius 
$R$ in favour of the second derivative 
$\bar {h}_s''(0) = -\theta /R$ of the initial condition (2.6), in order to highlight that the coalescence process depends only on the local shape of the lens near the bridge.
For a periodic array of two-dimensional viscous lenses, the centres of mass of the drops do not move due to symmetry. Instead, a horizontal force is generated that resists the inward motion towards the bridge, which results in a logarithmic velocity profile in the lens. In dimensional form, the growth rate of the bridge (4.10) and the additional force generated (4.9) are then given by
\begin{align} \bar V = \dot{\bar{h}}_0 = \left(0.5525 + \frac{0.7892}{\ln(2R\theta/\bar{h}_0)} + \cdots\right) \frac{\gamma\theta^2}{\eta}, \quad \bar{F} = \left(\frac{1.267}{\ln(2R\theta/\bar{h}_0)} + \cdots\right)\gamma \theta^2. \end{align}
A further, second-order, correction is calculated in Appendix A. It is important to note that the corrections are logarithmic in time, so these can be significant even at the very early stages of coalescence. In fact, logarithmic corrections also arise for viscous coalescence of freely suspended drops (Hopper Reference Hopper1990; Eggers et al. Reference Eggers, Lister and Stone1999). The structure of the problem is, however, different. In the freely suspended case the dynamics is 
$h\sim t \ln t$, which in contrast to (5.2a) does not exhibit a finite velocity at early times.
It is of interest to compare these two-dimensional, slender predictions with the experiments shown in figure 2(c). In Hack et al. (Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020) it was already shown that the experimental vertical coalescence velocity 
$\bar V$ was in very good quantitative agreement with the leading-order prediction of (5.1). Indeed, for freely floating drops, the higher-order terms are expected to be negligible in the experimental range. The current theory predicts a ratio of horizontal to vertical velocity 
$\bar U/\bar V \approx 0.57/\theta$ for the freely moving drops, which for the experimental contact angles amounts to 
$\bar U/\bar V \approx 1.1$. The experimentally measured ratio in figure 2(c) was found 
$\bar U/\bar V \approx 1.7$, which implies an even stronger centre of mass motion than predicted. This quantitative disparity could be due to three-dimensional effects, or due to the fact that the contact angles are not very small. Still, our theory offers an explanation for the appearance of a centre of mass motion for viscous lenses, and provides the relevant scaling laws. This centre of mass motion is not observed for inertial drops, which is in accordance with the decaying velocity in the inertial similarity solution (Hack et al. Reference Hack, Tewes, Xie, Datt, Harth, Harting and Snoeijer2020). For future experiments, it might be of interest to study coalescence while the centres of the drop are prevented from translating towards each other (e.g. by attaching the drops to fixed capillaries). In that case, we expect the appearance of significant (logarithmic) corrections to the coalescence velocity.
Besides three-dimensional effects, it is also of interest to discuss the influence of gravity. Gravity becomes important if the (dimensional) radius 
$R$ of the lenses is no longer small compared with the capillary length 
$\ell _c = \sqrt {\gamma /\varDelta \rho \,g}$, where 
$\varDelta \rho$ is the difference in density between the lenses and either of the external fluids. In this case, we expect gravity to flatten the static lens profile 
$h_s(x)$ (2.6), but near the rim of the lens, when 
$x \ll \ell _c$, gravity becomes negligible and the capillary thin-sheet equations (2.1) and contact angle 
$\theta$ are recovered, and hence the leading-order inner solution (3.6) will still hold. For (two-dimensional) free lenses, which can undergo uniform translation even in the presence of gravity, the 
$O(h_0)$ correction (5.1) also holds, but with a modified value of 
$h_s''(0)$ that simply follows from the static drop. For periodic lenses, the correction becomes more involved as the outer velocity profile (4.2) is affected by the change in 
$h_s$.
More generally, the dynamical structure of the problem bears a strong similarity with drop spreading on a rigid substrate, where the motion of the contact line also induces a weak flow on the scale of the drop (Bonn et al. Reference Bonn, Eggers, Indekeu, Meunier and Rolley2009). The spreading velocity exhibits a logarithmic dependence on the scale separation between drop size and the characteristic scale of the contact line in that case too. An important difference, however, is that for drop spreading the universal leading-order similarity solution for the inner problem captures the phenomenon, as it has only algebraically small corrections. Here we found that for viscous lens coalescence, the corrections are themselves only logarithmically small, and therefore can be significant.
Acknowledgements
We are grateful to J. Eggers for discussions. We would also like to thank the referees for all their very valuable comments and suggestions, and in particular for highlighting the different large-scale boundary conditions.
Funding
W.T. acknowledges support by the European Union's Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreement No 722497 – LubISS. M.A.H. acknowledges support from an Industrial Partnership Programme of the Netherlands Organisation for Scientific Research (NWO), cofinanced by Canon Production Printing Netherlands B.V., University of Twente, and Eindhoven University of Technology. C.D. and J.H.S. acknowledge support from NWO through VICI Grant No. 680-47-632.
Declaration of interests
The authors report no conflict of interest.
Appendix A. Details of the inner solution for non-zero $F$
We calculate both the 
$O(F)$ and 
$O(F^2)$ corrections to the leading-order result (3.6) by expanding
\begin{equation} \mathcal{H} = \mathcal{H}_0 + F \mathcal{H}_1 + F^2 \mathcal{H}_2 + \frac{h_0}{\dot h_0}\dot F \mathcal{H}_T + O(F^3), \end{equation}
and similarly for 
$\mathcal {U}$, 
$V$ and 
$\beta$. Note that, since we anticipate that 
$F = O(1/\ln t)$, the terms with subscript 
$2$ and 
$T$ are both 
$O(F^2)$.
We substitute the expansion into the (3.2)–(3.3), together with the matching condition 
$\mathcal {H}'(\infty ) = 1$ which is accurate to all orders in 
$F$, and identify coefficients. The resulting governing equations for the corrections 
$\mathcal {H}_\varDelta$, 
$\mathcal {U}_\varDelta$ and 
$V_\varDelta$, where 
$\varDelta = 1, 2, T$, are given by
$$\begin{gather} \begin{pmatrix} \mathcal{H}_\varDelta - \xi \mathcal{H}_\varDelta' + (\mathcal{H}_0 \mathcal{U}_\varDelta + \mathcal{H}_\varDelta \mathcal{U}_0)' \\ \mathcal{H}_0 \mathcal{H}_\varDelta'' + \mathcal{H}_\varDelta \mathcal{H}_0'' - \mathcal{H}_0' \mathcal{H}_\varDelta' + 4 (V_0 \mathcal{H}_0 \mathcal{U}_\varDelta' + V_0 \mathcal{H}_\varDelta \mathcal{U}_0' + V_\varDelta \mathcal{H}_0 \mathcal{U}_0') \end{pmatrix} {} \end{gather}$$
$$\begin{gather}{} = \underbrace{\begin{pmatrix} 0 \\ 1 \end{pmatrix}}_{\varDelta = 1}, \underbrace{\begin{pmatrix} -(\mathcal{H}_1 \mathcal{U}_1)' \\ -\mathcal{H}_1 \mathcal{H}_1'' + \tfrac12 \mathcal{H}_1'^2 - 4(V_0 \mathcal{H}_1 \mathcal{U}_1' + V_1 \mathcal{H}_0 \mathcal{U}_1' + V_1 \mathcal{H}_1 \mathcal{U}_0') \end{pmatrix}}_{\varDelta = 2}, \underbrace{\begin{pmatrix} - \mathcal{H}_1 \\ 0 \end{pmatrix}}_{\varDelta = T}, \end{gather}$$
and the boundary conditions are 
$\mathcal {H}_\varDelta '(0) = \mathcal {U}_\varDelta (0) = \mathcal {H}_\varDelta (0) = \mathcal {H}_\varDelta '(\infty ) = 0$.
Solving these equations numerically using Mathematica yields the results
\begin{gather} V_0 = 0.5525, \quad V_1 ={-}0.6227, \quad V_2 ={-}0.1267, \quad V_T = 0.3028, \end{gather}
\begin{gather} \beta_0 = 0.5734, \quad \beta_1 = 0.7079,\quad \beta_2 = 0.8728, \quad \beta_T ={-}0.5277, \end{gather}for the coefficients in the far-field behaviour (3.4).
We can then use the matching (4.9) to obtain the results
\begin{align} F &= \frac{4(V_0 + V_1 F + O(F^2))(\beta_0 + \beta_1 F + O(F^2))}{\ln(2/h_0)} \end{align}
\begin{align} &= \frac{4V_0\beta_0}{\ln(2/h_0)} + \frac{16 V_0\beta_0(V_0 \beta_1 + V_1 \beta_0)}{\ln(2/h_0)^2} + O(\ln(2/h_0)^{{-}3}), \end{align}and
\begin{align} \dot h_0 &= V_0 + V_1 F + V_2 F^2 + V_T (h_0/\dot h_0) \dot F + O(F^3) \end{align}
\begin{align} &= V_0 + \frac{\hat{V}_1}{\ln(2/h_0)} + \frac{\hat{V}_2}{\ln(2/h_0)^2} + O(\ln(2/h_0)^{{-}3}) \end{align}
\begin{align} &= V_0 + \frac{\hat{V}_1}{\ln(c/h_0)} + O(\ln(2/h_0)^{{-}3}), \end{align}
where 
$V_0 = 0.5525$, 
$\hat{V}_1 = 4V_0\beta _0 V_1 = 0.7892$ and
\begin{align} \hat{V}_2 &= 16V_0^2\beta_0^2 V_2 + 16V_0\beta_0(V_0\beta_1 + V_1\beta_0)V_1 + 4V_0\beta_0 V_T = 0.07272, \end{align}
\begin{align} c &= 2\exp(-\hat{V}_2/\hat{V}_1) = 2.19. \end{align}
The second-order coefficient 
$\hat{V} _2$ is coincidentally quite small when 
$\dot h_0$ is expanded in terms of 
$\ln (2/h_0)$, so including the second-order correction in figure 5 would not have a noticeable effect. However, when 
$\dot h_0$ is expanded in terms of e.g. 
$\ln (1/h_0)$ or 
$\ln (1/t)$, the second-order correction yields a significant improvement.
Open access
