2.1. Governing equations

Several applications discussed in § 1, such as cross-flow microemulsification (Charcosset et al. Reference Charcosset, Limayem and Fessi2004) and oil recovery (Eckmann et al. Reference Eckmann, Cavanagh and Branger2001), involve long and narrow channels, which corresponds to $\epsilon = W/L \ll 1$ for our model. Consistent with lubrication theory, the maximum droplet height $h_{max}$ is assumed to be much smaller than its diameter $D$, with $h_{max} \sim W$ and $D \sim L$ at most. In principle, lubrication theory requires the topography slope to be small, but in practice lubrication theory can yield qualitatively accurate predictions for two-phase problems even outside of this regime, as demonstrated by comparison with numerical simulations of the full governing equations (Zhou & Kumar Reference Zhou and Kumar2012).

The vertical and horizontal distances are non-dimensionalized with $W$ and $L$, respectively. All stresses are scaled with a characteristic capillary pressure $W\sigma /L^2$, where $\sigma$ is the interfacial tension (Schwartz & Eley Reference Schwartz and Eley1998). The horizontal velocity is scaled with a capillary spreading speed $u^* = \epsilon ^3\sigma /3\mu ^d$, where $\mu ^d$ is the droplet viscosity, the vertical velocity is scaled with $\epsilon u^*$ and time is scaled with $L/u^*$.

At leading order, the mass and momentum conservation equations in each phase are

(2.1)$$\begin{gather} u^i_{x} + v^i_{y} = 0, \end{gather}$$

(2.2)$$\begin{gather}(\mu^i/\mu^d) u^i_{yy} = 3P^i_{x}, \end{gather}$$

(2.3)$$\begin{gather}P^i_y + (\rho^i/\rho^d) Bo = 0, \end{gather}$$

where the superscript $i$ denotes the phase ($i = d$ for droplet, and $i = s$ for surrounding fluid), $u$ is the horizontal velocity, $v$ is the vertical velocity, $\mu$ is the viscosity, $P$ is the pressure and $\rho$ is the density. Three dimensionless parameters arise here: (i) the Bond number $Bo = \rho ^d g L^2 / \sigma$, where $g$ is the magnitude of the gravitational acceleration, which represents a ratio of gravitational forces to surface-tension forces, (ii) the viscosity ratio $\mu _r = \mu ^s / \mu ^d$ and (iii) the density ratio $\rho _r = \rho ^s / \rho ^d$. In this work, we assume $Bo \ll 1$ to isolate the influence of the imposed pressure gradient on the droplet dynamics.

Motivated by prior work on droplet motion on a solid substrate, we use a precursor-film/disjoining-pressure approach to model contact-line motion (Schwartz Reference Schwartz1998; Schwartz & Eley Reference Schwartz and Eley1998; Espín & Kumar Reference Espín and Kumar2015, Reference Espín and Kumar2017; Charitatos, Pham & Kumar Reference Charitatos, Pham and Kumar2021; Mhatre & Kumar Reference Mhatre and Kumar2024). As opposed to approaches that impose a slip law on the substrate, here, the contact-line position is extracted from the droplet height profile (see § 2.2) and is not an extra variable, which makes the resulting equations less cumbersome to solve (Savva & Kalliadasis Reference Savva and Kalliadasis2011). This approach assumes the presence of a precursor film of thickness $b$ along the entire substrate (figure 1). The total pressure within the droplet is $P = P^{'} + \varPi$, where $P^{'}$ is the hydrodynamic pressure and $\varPi$ is the disjoining pressure. A two-term disjoining pressure is used

(2.4)\begin{equation} \varPi = A \left[\left(\frac{b}{h}\right)^n - \left(\frac{b}{h}\right)^m\right], \end{equation}

where $A$ is the dimensionless Hamaker constant and $h$ is the droplet thickness. The term with exponent $n$ represents contributions from repulsive intermolecular forces, and the term with exponent $m$ represents contributions from attractive intermolecular forces. We choose $n=3$ and $m=2$, as prior studies show that these values provide qualitatively accurate predictions of contact-line motion at a reasonable computational cost (Espín & Kumar Reference Espín and Kumar2015, Reference Espín and Kumar2017; Charitatos et al. Reference Charitatos, Pham and Kumar2021). It should be noted that our model can be extended to incorporate chemical heterogeneity on the substrate by spatially varying $A$ (Schwartz Reference Schwartz1998; Schwartz & Eley Reference Schwartz and Eley1998; Kalpathy, Francis & Kumar Reference Kalpathy, Francis and Kumar2012).

The two-term disjoining pressure determines an equilibrium contact angle of magnitude $\theta _{eq}$, which the droplet attains in the absence of an imposed pressure gradient (Schwartz & Eley Reference Schwartz and Eley1998)

(2.5)\begin{equation} \theta_{eq} = \sqrt{\frac{2(n-m)A b}{(n-1)(m-1)}}. \end{equation}

Here, $\theta _{eq}$ is the scaled equilibrium contact angle and is related to the actual laboratory-frame equilibrium contact angle by $\theta _{eq, lab} = \epsilon \theta _{eq}$. All the angles discussed in the remainder of this paper are scaled angles.

The interface between the droplet and the surrounding fluid is described as $H(x,t) = h(x,t) + \eta (x)$, and we impose the following conditions:

(2.6)$$\begin{gather} u^d = u^s, \end{gather}$$

(2.7)$$\begin{gather}v^d = v^s, \end{gather}$$

(2.8)$$\begin{gather}P^d - P^s ={-}H_{xx}, \end{gather}$$

(2.9)$$\begin{gather}u^d_y = \mu_r u^s_y, \end{gather}$$

(2.10)$$\begin{gather}h_t + \left( \int_{\eta}^{H} u^d \,{{\rm d} y} \right)_x = 0, \end{gather}$$

where (2.6) and (2.7) are continuity of velocity, (2.8) is the normal stress balance, (2.9) is the tangential stress balance and (2.10) is the kinematic condition. No-slip and no-penetration conditions are imposed at the substrate ($y = \eta (x)$) and the top boundary ($y = 1$)

(2.11a,b)$$\begin{gather} u^d (y = \eta, t) = 0, \quad u^s (y = 1, t) = 0, \end{gather}$$

(2.12a,b)$$\begin{gather}v^d (y = \eta, t) = 0, \quad v^s (y = 1, t) = 0. \end{gather}$$

An imposed pressure gradient in the channel generates a mean flow rate in the $x$-direction, and we represent its dimensionless form by $Q_p$. Due to mass conservation, $Q_p$ does not vary along the $x$-direction, and we calculate it by considering a portion of the channel that only contains the precursor film, as it is a steady two-layer flow with a flat interface of height $b$

(2.13)\begin{equation} Q_p = \frac{1}{4} \Delta P \left({-}b^3 + \frac{3b(b-1)}{1+(\mu_r-1)b} + \frac{(b-1)^3}{\mu_r} \right). \end{equation}

Here, $\Delta P$ denotes the magnitude of the dimensionless imposed pressure gradient. Since the precursor-film thickness is much smaller than the channel width ($b \ll 1$), the flow rate simplifies to $Q_p = \Delta P / 4\mu _r$, which is used in the flow-rate condition

(2.14)\begin{equation} \int_{\eta}^{H} u^d \,{{\rm d} y} + \int_{H}^{1} u^s \,{{\rm d} y} = Q_p. \end{equation}

We find that the qualitative nature of our results does not change if we use (2.13) in the flow-rate condition, and the quantitative difference is approximately $5$ %–$10\,\%$. We use the simplified flow rate here to eliminate the dependence on $b$.

Integrating (2.2) and (2.3) with respect to $y$ and using conditions (2.6)–(2.14) gives the droplet height evolution equation

(2.15) \begin{align} \frac{\partial h}{\partial t} &= \frac{\partial}{\partial x} \left[ \frac{(H-1)^3h^3 f(H,\mu_r, \eta)}{g(H,\mu_r, \eta)} \left( \frac{\partial H^3}{\partial x^3} + \frac{\partial \varPi}{\partial x}\right) \right.\nonumber\\ & \quad + \left. \Delta P \frac{h^3 (3+H^2(\mu_r-1))-2H(1+(\mu_r-2)\eta + \eta(\mu_r \eta - 4))}{4 f(H,\mu_r, \eta)} \right], \end{align}

where

(2.16)\begin{gather} f(H,\mu_r,\eta) = 1+H(\mu_r-1)-\mu_r\eta, \end{gather}

(2.17)$$\begin{gather} g(H,\mu_r,\eta) = (H-\eta)^4 \mu_r^2 - 2\mu_r(H-1)(H^3-H^2(1+2\eta) \nonumber\\ +\eta({-}2+3\eta-2\eta^2)+H(2-2\eta+3\eta^2)) +(H-1)^4.\end{gather}$$

We set the length of the computational domain to $L$, and require the interface between the precursor film and the surrounding fluid to be flat at the left and right ends

(2.18a,b)$$\begin{gather} h(x = 0, t) = b, \quad \hfill h_x (x = 0, t) = 0, \end{gather}$$

(2.19a,b)$$\begin{gather}h(x=L, t) = b, \quad h_x (x = L, t) = 0. \end{gather}$$

We use a fourth-order polynomial which satisfies (2.18a,b) and (2.19a,b), and specify a dimensionless droplet volume $v_0$, to define the initial droplet shape. Evolution equation (2.15) is discretized using a fully implicit centred fourth-order finite-difference scheme with 6000–8000 spatial nodes, and the MATLAB built-in solver ode15s is used for time integration. All results were checked to ensure that they are converged with respect to spatial and temporal discretization.

In the simulations presented here, we set $L =6$ or $L = 9$ to obtain results that are independent of the length of the computational domain. A value of $b=0.001$ is used as it allows us to recover Tanner's spreading law for a two-dimensional droplet (§ 3) while allowing for reasonable computation times. For all calculations with rough surfaces, the precursor-film thickness is much smaller than the defect amplitude ($b/h_d \sim$ $10^{-2}$–$10^{-1}$), indicating that the bump accounts for surface roughness in a physically realistic manner. Smaller values of $b$ lead to slower droplet spreading and a more computationally stiff problem, but do not change the qualitative behaviour of the results presented here. For many of our results, we use representative values of $v_0=0.2$ and $\theta _{eq} = 10^{\circ }$ (which corresponds to $A=10^5$ using (2.5)), although we also study the effect of changing these parameters. In this paper, the largest value of $\theta _{eq}$ is $16^\circ$, which corresponds to $\tan {\theta _{eq}} \sim 10^{-1}$ and small interfacial slopes, consistent with the lubrication approximation.

2.2. Contact angles and contact lines

For a smooth substrate, we define apparent advancing ($\theta _{acl}$) and receding ($\theta _{rcl}$) contact angles (figure 1a) as the largest angles between the substrate and the tangents to the interface between the droplet and the surrounding fluid, on the advancing and receding sides of the droplet. Advancing and receding contact-line locations, $x_{acl}$ and $x_{rcl}$, are defined as the points where the tangents corresponding to the apparent contact angles and the substrate intersect. In the absence of an external pressure gradient $\theta _{acl} = \theta _{rcl}$. If $\theta _{acl} > \theta _{eq}$, the droplet spreads until $\theta _{acl} = \theta _{eq}$, and if $\theta _{acl} < \theta _{eq}$ the droplet retracts until $\theta _{acl} = \theta _{eq}$.

For substrates with topographical defects, a mesoscopic contact angle $\theta _m$ is defined (Espín & Kumar Reference Espín and Kumar2015; Park & Kumar Reference Park and Kumar2017)

(2.20)\begin{equation} \tan(\theta_m) = \frac{h_x}{1+(h_x + \eta_x)h_x}. \end{equation}

The apparent contact angles $\theta _{acl}$ and $\theta _{rcl}$ are obtained by finding the points on the interface between the droplet and the surrounding fluid where $\theta _m$ is the largest, on the advancing ($\theta _{ma}$ in figure 2a) and receding ($\theta _{mr}$ in figure 2b) sides of the droplet, respectively, and then extrapolating the tangents at these points to the substrate. The contact-line positions $x_{acl}$ and $x_{rcl}$ are defined as the points where the tangents corresponding to the apparent contact angles and the substrate intersect.

Figure 2. (a) Contact-line region when the advancing contact line moves over a defect, where $\eta$ represents the substrate shape, $H$ represents the interface shape, $\theta _{acl}$ is the apparent advancing contact angle and $\theta _{ma}$ is the largest mesoscopic angle on the advancing half of the droplet interface. (b) Contact-line region when the receding contact line moves over a defect, where $\theta _{rcl}$ is the apparent receding contact angle and $\theta _{mr}$ is the largest mesoscopic angle on the receding half of the droplet interface.

4.1. Pinning–depinning transition

Figures 6(a) and 6(b) show the variation of $\theta _{acl}$ and $x_{acl}$ with $t$ for $\Delta P = 0.05$ and $\Delta P = 0.07$. For $\Delta P = 0.05$, $x_{acl}$ does not vary with $t$, and $\theta _{acl}$ increases with $t$ to eventually reach a steady value. This behaviour is indicative of droplet pinning. For $\Delta P = 0.07$, $\theta _{acl}$ increases with $t$ and $x_{acl}$ remains constant until $t = 266$. Then, for $t > 266$, $\theta _{acl}$ decreases with $t$ and eventually reaches a steady value, and $x_{acl}$ increases with $t$. This behaviour is indicative of droplet depinning. Overall, these results suggest that there is a critical pressure gradient $\Delta P_{crit}$ above which droplet depinning occurs.

Figure 6. Droplet dynamics on a rough substrate. (a) Plot of $\theta _{acl}$ vs $t$. (b) Plot of $x_{acl}$ vs $t$. The kinks in $x_{acl}$ at $t = 266$ indicate that the contact line rapidly depins from the defect at that time and shifts to the right. The smaller kink to the left of the larger kink arises while numerically resolving $x_{acl}$ from droplet profiles. The other parameters are $L=9$, $\mu _r = 0.01$, $A = 10^5$ ($\theta _{eq} = 10^{\circ }$), $b=0.001$, $v_0 = 0.2$, $h_d = 0.02h_{max}$ and $w_d = 2h_d$.

Figure 7(a) shows the variation of the steady $\theta _{acl}$ value with $\Delta P$, and it is seen that $\theta _{acl}$ increases with $\Delta P$ for $\Delta P < 0.063$. For $\Delta P < 0.063$, it is found that $x_{acl}$ remains equal to its initial value, similar to what is seen for $\Delta P = 0.05$ (blue line in figure 6b). As $\Delta P$ increases from zero, the droplet remains pinned but becomes increasingly deformed, as shown in figure 7(b).

Figure 7. (a) Steady value of $\theta _{acl}$ vs $\Delta P$. (b) Steady droplet profiles for different $\Delta P$. (c) Droplet profiles at different $t$ for $\Delta P = 0.07$. The solid black lines show substrate topography and the other lines show droplet profiles. All other parameters are the same as in figure 6.

At $\Delta P = 0.063$, a significant decrease in $\theta _{acl}$ is observed, indicating droplet depinning, and for $\Delta P > 0.063$, there is a slight increase in $\theta _{acl}$ with $\Delta P$. Furthermore, for $\Delta P > 0.063$, the behaviour of $x_{acl}$ as a function of $t$ is qualitatively similar to that for $\Delta P = 0.07$ (red line in figure 6b), where $x_{acl}$ remains equal to the initial value for a period of time and then gradually increases with $t$. Thus, after the droplet depins it slides on the substrate with a steady shape, as shown in figure 7(c) (see supplementary material for videos of droplet pinning and depinning).

4.2. Droplet pinning mechanism

As discussed in § 1, the pinning–depinning transition has been rationalized in prior work using simple force-balance models (Fan et al. Reference Fan, Wilson and Kapur2011; Lu et al. Reference Lu, Xu, Xie, Wang, Sun and Yang2019). We apply those ideas here to see how well they describe our results. In those models, the drag force acting on the droplet due to the surrounding fluid flow drives depinning, and the surface-tension force acting along the droplet contact line resists depinning. The total drag force acting on the droplet is predominantly composed of skin drag and the pressure drag is negligible. Also, the total drag force is nearly equal to the shear force due to the droplet being thin (figure 2 of supplementary material). We note that pressure drag is expected to become significant for droplets not well described by the lubrication approximation. The balance between the shear force and the surface-tension force can be represented in a dimensionless form as $(D^2/h_{max})(\cos {\theta _{rcl}} - \cos {\theta _{acl}}) \sim \int _{s} \boldsymbol {n}\boldsymbol {\cdot }\boldsymbol {T}\boldsymbol {\cdot }\boldsymbol {t} \,{\rm d}s$ (see Appendix A). Here, $D = x_{acl} - x_{rcl}$ is the droplet width, $h_{max}$ is the maximum droplet height, $\boldsymbol {n}$ is the unit normal vector at the interface between the droplet and the surrounding fluid that points into the surrounding fluid, $\boldsymbol {T}$ is the droplet stress tensor, $\boldsymbol {t}$ is the unit tangent vector at the interface and $s$ is the interface arclength coordinate such that $s = 0$ at the receding contact line and $s = 1$ at the advancing contact line.

We extract the values of $D$, $h_{max}$, $\theta _{acl}$, $\theta _{rcl}$, $\boldsymbol {n}$, $\boldsymbol {T}$ and $\boldsymbol {t}$ from the steady droplet shapes in numerical simulations to calculate the terms in the force balance. Figure 8(a) shows these terms for different $\Delta P < \Delta P_{crit}$, with the open red circles presenting the simulation results. As can be seen, there is indeed a linear relationship, consistent with the force-balance model. It should be noted that inertial forces may become important for larger droplets or larger $\Delta P$ values (Mortazavi & Jung Reference Mortazavi and Jung2023), but they are neglected in this work.

Figure 8. (a) Terms in the force-balance model for pinned droplets at different $\Delta P$ values. The open red circles are results from numerical simulations and the straight blue line has a slope of unity. (b) Enlarged view of the pinned advancing contact line for $\Delta P = 0.03$, where the blue line shows the droplet and the red line shows the defect. (c) Enlarged view of the pinned receding contact line for $\Delta P = 0.03$, where the blue line shows the droplet and the red line shows the defect. The other parameters are $L=9$, $\mu _r = 0.01$, $A = 10^5$ ($\theta _{eq} = 10^{\circ }$), $b=0.001$, $v_0 = 0.2$, $h_d = 0.02h_{max}$ and $w_d = 2h_d$.

Next, we discuss the pinning locations of the two contact lines. It can be seen from figure 8(b,c) that the following geometric relations hold at the pinned advancing and receding contact lines:

(4.1)$$\begin{gather} \theta_{acl} = \theta_{ma} + \gamma_{a}, \end{gather}$$

(4.2)$$\begin{gather}\theta_{rcl} = \theta_{mr} - \gamma_{r}, \end{gather}$$

where $\gamma _a$ is the magnitude of the slope at the point on the defect on the right that coincides with $x_{acl}$, and $\gamma _r$ is the magnitude of the slope at the point on the defect on the left that coincides with $x_{rcl}$. Our calculations show that, for all the cases of a pinned droplet, $\theta _{ma} \approx \theta _{mr} \approx \theta _{eq}$, so the values of $\theta _{rcl}$ and $\theta _{acl}$ depend primarily on $\gamma _r$ and $\gamma _a$, respectively.

For all the cases where the droplet is pinned, it is found that $x_{rcl}$ and $x_{acl}$ coincide with the points on the defects where there is a maximum negative slope, whose magnitude we denote as $\gamma _{max}$. Using (4.1) and (4.2), it can be seen that, as a consequence of these pinning locations, $\theta _{rcl} \approx \theta _{eq} - \gamma _{max}$ is minimized, and $\theta _{acl} \approx \theta _{eq} + \gamma _{max}$ is maximized, which maximizes the surface-tension term in the force-balance model $(D^2/h_{max})(\cos {\theta _{rcl}} - \cos {\theta _{acl}})$. These findings are consistent with predictions from a lubrication-theory-based model for a droplet sliding down an inclined substrate with a single Gaussian-shaped bump, where the advancing contact line pins at the point on the bump that has the maximum negative slope, as it maximizes the surface-tension force (Park & Kumar Reference Park and Kumar2017). The present work shows that similar behaviour occurs even in more complicated situations, where the geometry and driving force for flow are different, and the dynamics of both phases are accounted for.

The pinning–depinning transition can also be understood through an analysis of the pressure gradients in the droplet (Espín & Kumar Reference Espín and Kumar2015; Charitatos et al. Reference Charitatos, Pham and Kumar2021). Figure 9 shows the variation of the pressure in a pinned droplet with $x$ for two $\Delta P$ values. The total pressure includes contributions from the capillary pressure and the disjoining pressure, and figure 9 shows the former quantity as well. For each case, the capillary pressure has a negative gradient and drives flow from the droplet interior toward the advancing contact line (located at $x = 5.56$), which promotes depinning. However, the disjoining pressure increases closer to the contact line, and causes the total pressure gradient within the droplet to become positive, which drives flow from the advancing contact line toward the droplet interior. Near the receding contact line, the total pressure gradient within the droplet is negative (figure 3 of supplementary material), which drives flow from the receding contact line toward the droplet interior. The flows near the two contact lines oppose each other and result in a pinned droplet.

Figure 9. Droplet pressure vs $x$ for pinned droplets. The dashed blue line shows the contribution from the capillary pressure and the solid red line shows the total pressure (capillary and disjoining). All the other parameters are the same as in figure 8.

As $\Delta P$ increases, the pinned droplet becomes increasingly deformed (see figure 7b), which increases the droplet thickness near the advancing contact line and decreases the disjoining pressure there. As a result, the total pressure gradient near the advancing contact line decreases, as shown in figure 9. At $\Delta P_{crit}$, the total pressure gradient near the advancing contact line becomes negative, and the total pressure gradient near the receding contact line remains negative. This results in a net negative pressure gradient within the droplet, which causes the two contact lines to depin simultaneously, and the droplet slides to the right.

4.3. Droplet depinning time

The time it takes the droplet to depin when $\Delta P > \Delta P_{crit}$ is of both fundamental and practical interest. The droplet depinning time, $t_{depinning}$, is calculated as the time at which the receding contact-line position, $x_{rcl}$, crosses the defect on the right. A scaling relation can be obtained for $t_{depinning}$ by dividing the droplet width, $D$, with the droplet velocity (3.1) obtained from the flat-interface approximation (§ 3). For $\mu _r \ll 1$, this yields $t_{depinning} \sim D / {\Delta P h_{max}}$.

Figure 10 shows $\log {t_{depinning}}$ vs $\log {D / ({\Delta P h_{max}})}$, where the open blue circles are results from numerical simulations and the dashed black line has a slope of $1$. The intercept for the dashed black line is chosen such that it overlaps with the simulation result at $\log {D / (\Delta P h_{max})} = 4.1814$. Fitting a straight line to the numerical calculations yields a slope of $1.423$. Thus, the scaling law derived above predicts a faster depinning time (${\sim D / (\Delta P h_{max})}$) than what is observed from numerical simulations (${\sim (D/(\Delta P h_{max}))^{1.423}}$). This is expected based on the results in § 3, as the flat-interface approximation does not account for viscous dissipation and surface-tension forces near the contact lines, which slow down droplet motion.

Figure 10. Plot of $\log {t_{depinning}}$ vs $\log {D/(\Delta P h_{max})}$, where the open blue circles show results from numerical simulations, and the dashed black line shows a slope of $1$. The other parameters are $L=9$, $\mu _r = 0.01$, $A = 10^5$ ($\theta _{eq} = 10^{\circ }$), $b=0.001$, $v_0 = 0.2$, $h_d = 0.02h_{max}$ and $w_d = 2h_d$.

5.1. Influence of viscosity ratio

Figure 11(a) shows $\Delta P_{crit}$ vs $\mu _r$. For $\mu _r < 1$, $\Delta P_{crit}$ decreases slightly with $\mu _r$, and this result is qualitatively consistent with the experimental observations of Fan et al. (Reference Fan, Wilson and Kapur2011) where, as the viscosity of the glycerol–water droplet increases (by increasing the glycerol concentration), a larger critical flow rate of the surrounding air is required to depin it on a treated glass surface. For $\mu _r > 1$, $\Delta P_{crit}$ decreases with $\mu _r$ to reach a minimum value at $\mu _r \approx 10$, and then increases with $\mu _r$. These results can be rationalized by examining the shear force acting on the droplet. Our calculations show that the surface-tension force does not vary significantly with $\mu _r$ ($\approx 6\,\%$), so we do not consider it here.

Figure 11. (a) Plot of $\Delta P_{crit}$ vs $\mu _r$. (b) Shear force vs $\mu _r$ for $\Delta P = 0.03$. (c) Plot of $\partial u^s / \partial y$ at the maximum height of pinned droplet vs $\mu _r$ for $\Delta P = 0.03$. The other parameters are $L=9$, $\mu _r = 0.01$, $A = 10^5$ ($\theta _{eq} = 10^{\circ }$), $b=0.001$, $v_0 = 0.2$, $h_d = 0.02h_{max}$ and $w_d = 2h_d$.

Figure 11(b) shows shear force vs $\mu _r$ for a fixed $\Delta P < \Delta P_{crit}$. Using lubrication theory, the shear force $\int _{s} \boldsymbol {n}\boldsymbol {\cdot }\boldsymbol {T}\boldsymbol {\cdot }\boldsymbol {t} \,{\rm d}s$ simplifies to $\int _{s} \mu _r \partial u^s / \partial y \,{\rm d}s$, where $\partial u^s /\partial y$ is the horizontal velocity gradient with respect to $y$ in the surrounding fluid at the interface. For a constant $\Delta P$, as $\mu _r$ increases, the surrounding fluid becomes more viscous and exerts a larger shear force on the droplet. However, for $\mu _r \gg 1$, the surrounding fluid becomes so viscous that the net flow rate along the $x$-direction decreases. Consequently, $\partial u^s / \partial y$ decreases with $\mu _r$, as exemplified in figure 11(c), where $\partial u^s / \partial y$ at the maximum height of the pinned droplet is shown. As a result, the shear force decreases with $\mu _r$ for $\mu _r \gg 1$, and the droplet experiences a maximum shear force, and consequently a minimum $\Delta P_{crit}$ for $\mu _r \approx 10$. This observation is potentially important for applications such as cross-flow microemulsification and oil recovery, where the viscosity of the surrounding fluid could be tuned to minimize $\Delta P_{crit}$.

5.2. Influence of droplet volume

Figure 12(a) shows $\Delta P_{crit}$ vs $v_0$ for two different values of $\mu _r < 1$ (droplet more viscous than the surrounding fluid). It is seen that $\Delta P_{crit}$ decreases with $v_0$ for each case, and a larger $\Delta P_{crit}$ is required for a smaller $\mu _r$. This result is qualitatively consistent with the experimental findings of Fan et al., where a lower critical flow rate of surrounding air is required to depin a larger glycerol–water droplet on a treated glass surface, and a larger flow rate is required to depin a more viscous droplet of the same volume (Fan et al. Reference Fan, Wilson and Kapur2011). For $\mu _r > 1$, $\Delta P_{crit}$ decreases with $v_0$ as well (figure 4 of supplementary material).

Figure 12. (a) Plot of $\Delta P_{crit}$ vs $v_0$ for $\mu _r = 0.001$ and $\mu _r = 0.9$. (b) Plot of $\log {\Delta P_{crit}}$ vs $\log {v_0}$ for $\mu _r = 0.001$, where the open blue circles show results from droplet simulations, and the dashed black line shows a linear fit with a slope of $-0.5$. The other parameters are $L=9$, $A = 10^5$ ($\theta _{eq} = 10^{\circ }$), $b=0.001$, $h_d = 0.02h_{max}$ and $w_d = 2h_d$.

These results can be understood in terms of the force-balance model discussed in § 4.2. An analytical expression is obtained for the shear force acting on the droplet using the flat-interface approximation (see § 3), where the shear stress acting on a flat interface of height $H$ is multiplied by the length of the horizontal domain $L$

(5.1)\begin{equation} F_s^{flat} = 3\Delta P L \left( H - \frac{1+H^2 (\mu_r - 1)}{2(1+H(\mu_r - 1))} \right). \end{equation}

Using (5.1) to calculate the shear force by setting $H = h_{max}$ obtained from the steady droplet profiles, and assuming $L \sim D$, the dimensionless force balance at the point of depinning can be expressed as $3\Delta P_{crit} (h_{max}-1)/2 \sim {(}D/h_{max}(\cos {\theta _{rcl}} - \cos {\theta _{acl}})$ for $\mu _r \ll 1$, and $3\Delta P_{crit} h_{max}/2 \sim (D/h_{max})(\cos {\theta _{rcl}} - \cos {\theta _{acl}})$ for $\mu _r \gg 1$. As discussed in § 3, both $D$ and $h_{max}$ $\sim v_0^{0.5}$. Thus, $\Delta P_{crit} \sim v_0^{-0.5}$ for both $\mu _r \ll 1$ and $\mu _r \gg 1$. Figure 12(b) shows $\log {\Delta P_{crit}}$ vs $\log {v_0}$ for $\mu _r = 0.001$, where the open blue circles are results from numerical simulations and the dashed black line shows a slope of $-0.5$. The close agreement indicates that the flat-interface approximation accurately captures the influence of $v_0$ on $\Delta P_{crit}$. Results from numerical simulations for $\mu _r = 100$ follow the scaling relation as well (figure 4 of supplementary material).

6.1. Influence of substrate wettability

Figure 13(a) shows $\Delta P_{crit}$ vs $\theta _{eq}$ for $\mu _r = 0.01$, and it is seen that $\Delta P_{crit}$ increases with $\theta _{eq}$. This result is qualitatively consistent with the experimental findings of Madani & Amirfazli (Reference Madani and Amirfazli2014), where a lower surrounding water flow rate is required to depin an oil droplet ($\mu _r < 1$) on a more wettable substrate.

Figure 13. (a) Plot of $\Delta P_{crit}$ vs $\theta _{eq}$. (b) Plot of $(D^2/h_{max}) (\cos {\theta _{rcl}} - \cos {\theta _{acl}})$ vs $\theta _{eq}$ for $\Delta P = 0.03$. The other parameters are $L=9$, $\mu _r = 0.01$, $v_0 = 0.2$, $b=0.001$, $h_d = 0.02h_{max}$ and $w_d = 2h_d$.

These results can be understood by examining the surface-tension term in the force-balance model, $(D^2/h_{max}) (\cos {\theta _{rcl}} - \cos {\theta _{acl}})$ (see § 4.2), which increases with $\theta _{eq}$ as shown in figure 13(b). The shear force acting on the droplet is not considered as it does not vary significantly with $\theta _{eq}$ ($\approx 7\,\%$). This can be understood from the flat-interface approximation. According to (5.1), the shear force depends on $h_{max}$, which does not vary significantly with $\theta _{eq}$ ($\approx 8\,\%$).

As $\theta _{eq}$ increases, $\theta _{rcl} \approx \theta _{eq} - \gamma _{max}$ increases, which tends to decrease the surface-tension force. Additionally, $\theta _{acl} \approx \theta _{eq} + \gamma _{max}$ increases, which tends to increase the surface-tension force. Overall, a net increase is observed in the surface-tension force (figure 13b), and a larger $\Delta P_{crit}$ is required for droplet depinning on a less wettable substrate.

6.2. Influence of defect geometry

The defect shape is governed by its maximum height $h_d$ and maximum width $w_d$ (see § 2). It can be seen from figure 1(b) that the maximum slope along the defect $\gamma _{max}$ increases with $h_d$ (constant $w_d$), and decreases with $w_d$ (constant $h_d$). As will be discussed below, $\Delta P_{crit}$ depends only on $\gamma _{max}$, and increasing $h_d$ has the same qualitative effect on the droplet dynamics as decreasing $w_d$. Here, we vary $h_d$ while keeping $w_d$ constant.

Figure 14(a) shows $\Delta P_{crit}$ vs the ratio of $h_d$ to the initial maximum droplet height, and it is seen that $\Delta P_{crit}$ increases with $h_d$. This result can be understood by examining the surface-tension force acting along the contact line. As discussed in § 6.1, the shear force acting on the droplet for a fixed $\Delta P$ depends on $h_{max}$, which does not vary significantly with $h_d$ (${\approx }4\,\%$), and $\mu _r$, which does not vary with $h_d$. As a consequence, the shear force does not vary significantly with $h_d$ (${\approx }6\,\%$).

Figure 14. (a) Plot of $\Delta P_{crit}$ vs $h_d$. (b) Plot of $(D^2/h_{max})(\cos {\theta _{rcl}} - \cos {\theta _{acl}})$ vs $h_d$ for $\Delta P = 0.015$. The other parameters are $L=9$, $\mu _r = 0.01$, $b=0.001$, $A = 10^5$ ($\theta _{eq} = 10^{\circ }$), $v_0 = 0.2$ and $w_d = 0.008$.

As discussed in § 4.2, $\theta _{acl} \approx \theta _{eq} + \gamma _{max}$ and $\theta _{rcl} \approx \theta _{eq} - \gamma _{max}$. It can be seen that for a fixed $\theta _{eq}$, as $h_d$, and consequently $\gamma _{max}$, increases, $\theta _{acl}$ increases and $\theta _{rcl}$ decreases. As a result, the surface-tension term in the force-balance model $(D^2/h_{max})(\cos {\theta _{rcl}} - \cos {\theta _{acl}})$ increases with $h_d$, as shown in figure 14(b). Thus, a larger $\Delta P_{crit}$ is required for droplet depinning if taller defects are present on the substrate. Similarly, a larger $\Delta P_{crit}$ is required for narrower defects (figure 5 of supplementary material).

We also find that when the ratio of the defect height to the initial maximum droplet height is greater than $0.025$ (see dashed line in figure 14a), smaller residual droplets are formed behind and in front of the defect on the right, after the original droplet depins. (The value of 0.025 is specific to the parameters used in the calculations.) We discuss this phenomenon in the next section.

7.1. Geometric mechanism

Figure 15 shows droplet profiles as time progresses after the droplet depins, and it is seen that a residual droplet is formed behind the defect on the right as the droplet slides over it. This happens when $\theta _{rcl} < \gamma _{max}$ as the receding contact line approaches the defect. The disjoining pressure above the defect increases due to droplet thinning, which drives flow away from this region and a smaller residual droplet is pinched off behind the defect. However, if $\theta _{rcl} > \gamma _{max}$, the droplet slides over the defect without leaving behind a residual droplet. It should be noted that $\theta _{rcl}$ is obtained as an output from our calculations (see § 2.2) and generally depends on the viscosity ratio, droplet volume, substrate wettability and defect geometry. These findings are consistent with a prior lubrication-theory-based model for a droplet sliding down an inclined substrate with a single Gaussian-shaped bump, where a residual droplet is formed behind the bump, if $\theta _{rcl} < \gamma _{max}$ as the receding contact line approaches the bump (Park & Kumar Reference Park and Kumar2017).

Figure 15. Droplet profiles at (a) $t = 500$, (b) $t = 520$ and (c) $t = 540$. Enlarged view of the contact line at (d) $t = 500,$ (e) $t = 520$ and (f) $t = 540$. The solid black lines show substrate topography and the blue lines show droplet profiles. The parameters are $L=9$, $A = 10^5$ ($\theta _{eq} = 10^{\circ }$), $v_0 = 0.2$, $b=0.001$, $\Delta P = 0.1$, $h_d = 0.02h_{max}$ and $w_d = 2h_d$.

It is also found that an additional residual droplet forms in front of the defect at later times if $\theta _{rcl} < \gamma _{max}$. Figure 16 shows droplet profiles as time progresses after the formation of the residual droplet behind the defect. As seen in figure 15, the receding contact line remains in the vicinity of the defect during the pinch-off of the residual droplet behind the defect, while the advancing contact line slides to the right. As a result, the droplet elongates and attains a slender concave upward shape on its receding side to the right of the defect as shown in figure 16(a,e). Consequently, the disjoining pressure in this region increases, which drives flow away from the thinnest point of the interface, and a smaller residual droplet is pinched off in front of the defect (figures 16b–d and 16f–h). These findings are also consistent with the prior work discussed above on droplets sliding down an inclined substrate (Park & Kumar Reference Park and Kumar2017), indicating that the condition for residual droplet formation ($\theta _{rcl} < \gamma _{max}$) holds more generally. Indeed, one would expect this criterion to hold for a broad range of geometries and flows as long as the topography amplitude is much smaller than the length scales characterizing the overall flow. We note that increasing the droplet volume increases the value of $\theta _{rcl}$ when the receding contact line encounters the defect, eventually leading to $\theta _{rcl} > \gamma _{max}$ and suppression of residual droplet formation.

Figure 16. Droplet profiles at (a) $t = 560$, (b) $t = 580$, (c) $t = 600$ and (d) $t = 620$. Enlarged view of the contact line at (e) $t=560$, (f) $t=580$, (g) $t=600$ and (h) $t=620$. The solid black lines show substrate topography and the blue lines show droplet profiles. The parameters are the same as in figure 15.

The formation of residual droplets has important implications for practical applications. Our results indicate that the presence of topographical defects on the substrate may make it difficult to recover the entire volume of the droplet in applications such as surface cleaning and oil recovery. However, for applications such as microemulsification, our results suggest that the size of the dispersed-phase droplets could be controlled by designing defects to break-up larger droplets into smaller ones.

7.2. Influence of pressure gradient on residual droplet volume

We now consider the influence of $\Delta P$ on the volume of the residual droplets formed at the defect. The volume of a residual droplet is calculated as $\int _{x_l}^{x_r} h \, \mathrm {d}\kern0.06em x$, where $x_l$ and $x_r$ are the left and right end points of the residual droplet. These values are calculated such that $h>b$ $\forall$ $x_l < x < x_r$ in the region where the residual droplet is present. Figure 17(a) shows residual droplet volume vs $\Delta P$, where the open blue circles are for droplets formed behind the defect, and the open red triangles are for droplets formed in front of the defect. It can be seen that the volume of the droplet formed behind the defect decreases with $\Delta P$. This result can be rationalized by examining the droplet profiles when the receding contact line approaches the defect. Figure 17(b) shows these droplet profiles for different $\Delta P$ values. For each case, the receding contact line slows down as it approaches the defect. As $\Delta P$ increases, the depinned advancing contact line slides faster, leading to a more elongated droplet shape (figure 17b). This results in a larger portion of the droplet being trapped behind the defect just before the residual droplet is pinched off. Thus, a larger residual droplet is formed behind the defect for a larger $\Delta P$ value.

Figure 17. (a) Residual droplet volume vs $\Delta P$. (b) Droplet profiles at $t = 300$ for different $\Delta P$ values. The solid black line shows substrate topography and the other lines show droplet profiles. The parameters are $L=9$, $A = 10^5$ ($\theta _{eq} = 10^{\circ }$), $v_0 = 0.2$, $b=0.001$, $\Delta P = 0.1$, $h_d = 0.02h_{max}$ and $w_d = 2h_d$.

Figure 17(a) also shows that the volume of the residual droplet in front of the defect decreases with $\Delta P$. This can be rationalized by noting that because the receding contact line slides faster after depinning from the defect for larger $\Delta P$ values, there is a shorter time for droplet pinch-off and a smaller residual droplet is formed in front of the defect.