6.1. Viscous flow

From Appendices B and C, viscoplastic flows with $R\lesssim \kappa$ are viscous-dominated in both inner and outer regions. With $\ell =0.2$ (upper bound of the thick-channel limit), this corresponds to approximately $R\lesssim 50$ . In this range, the creeping-flow semianalytical model should provide accurate predictions. Since $d\lt 1$ throughout, the outer flow sets the upper Reynolds number limit. Accordingly, results in this section are restricted to $R\leqslant 50$ , i.e. the viscous-dominated regime.

The variations in slip velocity with slip area fraction, shown in figure 2(a), highlight the role of groove width in modulating the slippery behaviour. Overall, the larger slip area fractions lead to greater slip velocities. The semianalytical model also accurately reproduces the axial velocity profiles (figure 2 b) and predicts the position of the yield surfaces. Meanwhile, the asymmetry index $I_s$ (figure 2 c) provides additional insight into the deviation of the slip distribution from symmetric behaviour induced by the interface deflection. Figure 2(c) further shows that, as the slip area fraction increases, the asymmetry of the slip velocity profile becomes more pronounced, i.e. implying an enhanced interface deflection. As shown in figure 2, the semianalytical results show an excellent agreement with the CFD predictions, while the approximate solution indicate slightly deviated results.

In figure 2, a small Bingham number ( $B=0.1$ ) is used because the semianalytical solution applies only prior to SH plug formation. For the parameters of figure 2, an SH plug is expected at $B=1$ ; thus, the admissible Bingham-number range depends on the remaining flow parameters. There, $d=0.05$ , but for smaller $d$ the model remains valid at larger $B$ , and subsequent figures include higher- $B$ cases computed semianalytically. For very small groove aspect ratios ( $d\lt 0.05$ ), CFD becomes computationally demanding due to the fine mesh required to resolve the air flow in the groove, whereas the semianalytical model remains efficient. A combined semianalytical and CFD approach is therefore necessary to capture the flow over a broad parameter range and across regimes.

Figure 2.

The CFD and semianalytical results: (a) slip velocity, (b) axial velocity profiles, and (c) local asymmetry index ( $I_s$ ) from the semianalytical solution (4.32), the approximate solution ((4.33)–(4.35)), and numerical simulations (symbols: squares, triangles, circles for $\varphi =0.3$ , $0.5$ , $0.9$ ). Parameters: $d=0.05$ , $B=0.1$ , $\textit{Ca}=0.03$ , $R=0.01$ , $\varTheta =95^\circ$ , $\ell =0.2$ . In (b), grey shading denotes the numerically obtained centre plug; dashed black lines indicate the lower ( $h$ ) and upper ( $h^u$ ) yield surfaces from the semianalytical solution (4.32).

As discussed, for very small groove aspect ratios, CFD becomes impractical due to mesh-resolution constraints and high computational cost; we therefore use the semianalytical solution in this regime. Figure 3 shows that both meniscus shape and slip velocity depend on groove geometry ( $d$ ), yield stress ( $B$ ) and surface tension (via $\textit{Ca}$ ). Here, $\eta =\hat \eta / \hat d$ denotes interface deflection relative to groove depth. A flatter meniscus occurs for larger $d$ (figure 3 a) and smaller $B$ and $\textit{Ca}$ (figure 3 b,c). At fixed $\textit{Ca}$ , decreasing $d$ increases interfacial shear stress and thus deflection. At larger $B$ , maintaining a fixed flow rate requires a higher pressure gradient, increasing interfacial shear and deflection. Increasing $\textit{Ca}$ (weaker surface tension) likewise enhances deflection as surface tension cannot counteract shear forces.

Figure 3.

Semianalytical results. (a–c) Liquid–air interface profiles for varying (a) groove aspect ratio ( $B=0.1$ , $\textit{Ca}=0.03$ ), (b) Bingham number ( $d=10^{-2}$ , $\textit{Ca}=0.03$ ) and (c) capillary number ( $B=0.1$ , $\varTheta =95^\circ$ , $d=10^{-2}$ ). All cases: $R=0.01$ , $\ell =0.2$ , $\varphi =0.5$ . (d–f) Corresponding slip velocity profiles for varying (d) groove aspect ratio, (e) Bingham number and (f) capillary number, under the same conditions. In (a–c), $\eta =-1$ denotes the groove bottom wall. Insets in (d–f): asymmetry index of the slip velocity profile ( $I_s$ ).

Figure 3. Long description

The image contains six line graphs arranged in two rows and three columns. Each graph shows the relationship between a variable on the x-axis and another variable on the y-axis. Panel A: The line graph shows the liquid-air interface profiles for varying groove aspect ratios. The x-axis is labeled ‘x’ with units of ℓ/4, and the y-axis is labeled ‘η’. Three lines represent different groove aspect ratios: d = 1 x 10^-1, 5 x 10^-2, and 1 x 10^-2. Panel B: The line graph shows the liquid-air interface profiles for varying Bingham numbers. The x-axis is labeled ‘x’ with units of ℓ/4, and the y-axis is labeled ‘η’. Three lines represent different Bingham numbers: B = 1 x 10^-2, 1 x 10^-1, and 1. Panel C: The line graph shows the liquid-air interface profiles for varying capillary numbers. The x-axis is labeled ‘x’ with units of ℓ/4, and the y-axis is labeled ‘η’. Three lines represent different capillary numbers: Ca = 1 x 10^-3, 1 x 10^-2, and 3 x 10^-2. Panel D: The line graph shows the corresponding slip velocity profiles for varying groove aspect ratios. The x-axis is labeled ‘x’ with units of ℓ/4, and the y-axis is labeled ‘u_s’. Three lines represent different groove aspect ratios: d = 1 x 10^-1, 5 x 10^-2, and 1 x 10^-2. An inset shows the asymmetry index of the slip velocity profile. Panel E: The line graph shows the corresponding slip velocity profiles for varying Bingham numbers. The x-axis is labeled ‘x’ with units of ℓ/4, and the y-axis is labeled ‘u_s’. Three lines represent different Bingham numbers: B = 1 x 10^-2, 1 x 10^-1, and 1. An inset shows the asymmetry index of the slip velocity profile. Panel F: The line graph shows the corresponding slip velocity profiles for varying capillary numbers. The x-axis is labeled ‘x’ with units of ℓ/4, and the y-axis is labeled ‘u_s’. Three lines represent different capillary numbers: Ca = 1 x 10^-3, 1 x 10^-2, and 3 x 10^-2. An inset shows the asymmetry index of the slip velocity profile.

The corresponding slip velocity profiles in figure 3(d–f) illustrate the effects of interface deflection. Increased meniscus deformation, occurring for smaller $d$ or larger $B$ and $\textit{Ca}$ , leads to a more asymmetric slip velocity distribution along the liquid–air interface, as also quantified by the local asymmetry index ( $I_s$ ) shown in the insets. As $d$ decreases, the slip number is reduced and the average slip correspondingly decreases (figure 3 d), whereas increasing $B$ results in a larger average slip (figure 3 e). An increase in $\textit{Ca}$ , associated with reduced surface tension, primarily amplifies the asymmetry of the slip profile.

In figure 4, the average slip velocity, $\langle u_s\rangle$ , is calculated using the semianalytical and CFD models over a wide range of groove aspect ratios, Bingham numbers and capillary numbers. In figure 4(a), the coloured curve is obtained from the semianalytical approach and indicates the onset of SH plug formation along the liquid–air interface. This line represents the limit up to which the semianalytical solution is valid. Above this line (i.e. the region highlighted in grey), the CFD model is the only modelling option, and its results are shown by markers. Figure 4(a) shows that the average slip increases with $d$ and $B$ . Also, there is good consistency between the semianalytical and CFD predictions. In figure 4(b), similar results are presented but for varying $\textit{Ca}$ and at $B=0.1$ . Again, there is a limit for the semianalytical solution, shown by the vertical line on the right-hand side. However, this figure also reveals the limitation of the CFD simulation: for very small groove aspect ratios, numerical results could not be obtained due to mesh-resolution constraints and the prohibitive computational cost of resolving the highly confined geometry. The variation of $\langle u_s\rangle$ with the capillary number shown in figure 4(b) further indicates that, although increasing $\textit{Ca}$ promotes deformation of the liquid–air interface, the resulting change in the average slip velocity remains modest. Furthermore, figure 4 implies that with increasing groove depth, the average slip eventually converges to its saturated value.

Figure 4.

Semianalytical and CFD results: average slip velocity along the liquid–air interface, $\langle u_s \rangle$ , versus groove aspect ratio $d$ : (a) $B=0.1$ , $1$ , $10$ at $\textit{Ca}=1.37\times 10^{-4}$ ; (b) $\textit{Ca}=1.37\times 10^{-4}$ , $0.01$ , $0.03$ at $B=0.1$ . In both panels $\ell =0.2$ , $\varphi =0.5$ , $\varTheta =95^\circ$ , $R=0.01$ . Coloured dashed lines: semianalytical solution. Symbols: CFD. Gradient-coloured line (see colour bar): $\langle u_s \rangle$ at the critical groove aspect ratio $d_{\mathit{cr}}$ for SH plug onset (semianalytical). In (a), grey shading is SH plug formation domain; in (b), blue and green shading is regions where CFD and semianalytical solutions are not feasible, respectively.

Figure 5, which illustrates the CFD results, elucidates the influence of groove depth and the fluid’s yield stress on the airflow dynamics and interface deflection. For both Newtonian (figure 5 a–d) and Bingham (figure 5 e–h) outer-flow scenarios, the airflow inside the groove sustains a primary recirculation vortex whose morphology is governed by the aspect ratio, $d$ : for the shallower groove, the vortex becomes more elongated. In addition, the liquid–air interface is more deflected for the shallower groove. Indeed, for a shallower groove, the airflow shear gradient becomes more significant, which can markedly affect the interface dynamics (as suggested by the second term on the left-hand side of (4.15)).

Comparing figure 5(a–d) and figure 5 (e–h), we observe that the yield stress affects the flow dynamics in two ways: (i) the interface is more deflected due to enhanced shear in the outer flow and at the interface (as suggested by the first term on the left-hand side of (4.15)), which is associated with the increased pressure gradient required to sustain the flow rate; (ii) an SH plug forms on the interface, which exhibits an asymmetric shape (with respect to the $y$ -axis). Although the results at $R=50$ are still viscous-dominated, inertia may still affect the morphology of the SH plug. Moreover, the deflected interface should also contribute to this asymmetric SH plug profile. Note that, inside the SH plug, the norm of the strain-rate is very small and falls below a critical value calculated at the yield surface location for a regularized Bingham model.

Figure 5.

The CFD results: streamwise velocity contours within and near the groove for $d=0.05$ (a,e), $0.1$ (b, f), $0.2$ (c,g), $0.5$ (d,h), at $B=0$ (a–d) and $B=10$ (e–h); $R=50$ ( $\textit{Ca}\approx 0.69$ ), $\ell =0.2$ , $\varTheta =95^\circ$ , $\varphi =0.5$ . Solid black lines are liquid–air interface. For $B=10$ (e–h), dashed white lines are boundary of the SH plug region. Solid and dashed contours denote positive and negative streamwise velocity, respectively.

Figure 5. Long description

Panel A: A heat map showing streamwise velocity contours within and near the groove. The x-axis ranges from -ℓ/4 to ℓ/4 and the y-axis ranges from -5 to 5 times 10^-3. Solid black lines denote the liquid-air interface. Panel B: A heat map showing streamwise velocity contours within and near the groove. The x-axis ranges from -ℓ/4 to ℓ/4 and the y-axis ranges from -10 to 0. Solid black lines denote the liquid-air interface. Panel C: A heat map showing streamwise velocity contours within and near the groove. The x-axis ranges from -ℓ/4 to ℓ/4 and the y-axis ranges from -20 to 0. Solid black lines denote the liquid-air interface. Panel D: A heat map showing streamwise velocity contours within and near the groove. The x-axis ranges from -ℓ/4 to ℓ/4 and the y-axis ranges from -50 to 0. Solid black lines denote the liquid-air interface, and dashed white lines denote the boundary of the SH plug region. Panel E: A heat map showing streamwise velocity contours within and near the groove. The x-axis ranges from -ℓ/4 to ℓ/4 and the y-axis ranges from -5 to 5 times 10^-3. Solid black lines denote the liquid-air interface. Panel F: A heat map showing streamwise velocity contours within and near the groove. The x-axis ranges from -ℓ/4 to ℓ/4 and the y-axis ranges from -10 to 0. Solid black lines denote the liquid-air interface. Panel G: A heat map showing streamwise velocity contours within and near the groove. The x-axis ranges from -ℓ/4 to ℓ/4 and the y-axis ranges from -20 to 0. Solid black lines denote the liquid-air interface. Panel H: A heat map showing streamwise velocity contours within and near the groove. The x-axis ranges from -ℓ/4 to ℓ/4 and the y-axis ranges from -50 to 0. Solid black lines denote the liquid-air interface, and dashed white lines denote the boundary of the SH plug region.

6.2. Inertial flow

As discussed in the previous sections (see also Appendices B and C), for $R \lesssim 50$ the semianalytical solution is reliable and captures the creeping-flow dynamics. This conclusion follows from the negligible advection terms and the relatively weak deformation of the liquid–air interface. However, for $R\gt 50$ , inertial effects become important, for a typical thick channel with $\ell =0.2$ . In addition, as the Reynolds number increases, the interface deformation also intensifies. Collectively, these effects push the system beyond the validity range of the semianalytical model. Therefore, to resolve the flow field and characterize the dynamics at higher Reynolds numbers, we employ full CFD simulations, which can accommodate both the nonlinear inertial contributions and the enhanced interface deformation.

Figure 6 illustrates how rheology, inertia and groove aspect ratio influence the local slip-velocity distribution along the liquid–air interface. Figure 6(a,b) present the variations of $u_s(x)$ for a range of Bingham numbers, groove aspect ratios, and Reynolds numbers, while figure 6(c) shows the corresponding distribution of the asymmetry index along the interface for a broad range of Reynolds numbers. Figure 6(a) indicates that the slip-velocity profile $u_s(x)$ remains nearly symmetric with respect to $y$ axis across all conditions considered. For Newtonian fluids ( $B=0$ ), the slip velocity reaches its maximum at the groove centre and decreases monotonically towards the groove edges, forming a characteristic dome-shaped profile. Introducing yield stress ( $B\gt 0$ ) fundamentally alters this behaviour by generating a central plateau where $u_s$ remains constant. This plateau corresponds to the unyielded plug (i.e. the SH plug) observed in figure 5. Both the plateau width and the maximum slip velocity increase with the Bingham number. Figure 6(b) examines the influence of Reynolds number at $B=10$ and $d=0.5$ . Increasing the Reynolds number leads to a reduction in the dimensionless slip velocity, while simultaneously enlarging the central plateau associated with the SH plug region. However, as $R$ increases, the profiles become increasingly asymmetric. This behaviour is quantified in figure 6(c), which shows that the asymmetry grows with Reynolds number, particularly near the groove edges. In contrast, the central region where the unyielded SH plug appears indicates symmetry.

Figure 6.

The CFD results. (a) Slip velocity profiles at the liquid–air interface for varying groove aspect ratio $d$ and $B=0$ (solid), $1$ (dashed), $10$ (dash–dotted), at fixed $R=50$ . (b) Slip velocity profiles for varying $R$ at $d=0.5$ , $B=10$ . (c) Slip velocity asymmetry index ( $I_s$ ) versus $R$ for $B=10$ , $d=0.5$ . All panels: $\ell =0.2$ , $\varphi =0.5$ , $\varTheta =95^\circ$ .

To elucidate the influence of inertia on the flow field within and adjacent to the groove, figure 7 shows contours of the streamwise velocity over a range of Reynolds ( $R$ ) and Bingham numbers ( $B$ ). The results reveal a complex evolution of the recirculation vortex with increasing inertia. In the creeping flow regime ( $R=0.01$ ), the vortex structure is symmetric. As the Reynolds number increases into the intermediate range ( $50 \leqslant R \leqslant 100$ ), the vortex core shifts towards the downstream edge of the groove. Interestingly, at higher Reynolds number ( $R=500$ ), this trend reverses and the vortex migrates back towards the upstream edge. This vortex behaviour is connected to the deflection mode of the interface. Under creeping flow condition, $R=0.01$ , the interface remains nearly flat and pinned to both groove edges (figure 7 a,e), leading to a symmetric vortex structure within the groove. For $50 \leqslant R \leqslant 100$ , where the shear increases but the viscous effects are still significant, the interface shape changes from concave (penetrating into the groove) to convex (penetrating into the viscoplastic flow) from the upstream to the downstream edge. Consistent with the discussion in § 4.1.1, this is the typical interface shape captured by the semianalytical model (see figure 3). The convex interface near the downstream edge explains the displacement of the vortex core towards that edge. In contrast, at $R=500$ , inertial effects become significant and control the interface dynamics, such that the interface becomes convex near the upstream edge. This causes the vortex to shift back towards the upstream edge. Generally, we may argue that the concave–convex mode of the interface shape is a signature of Stokes flow, where viscous effects are dominant. However, when inertia becomes significant, the interface deformation becomes more complex and typically shows a convex shape near the upstream edge. For example, in figure 7(d,h), where $R=500$ , we observe a convex–concave–convex mode for the interface shape. The SH plug also shows a complex evolution. It is evident that the SH plug expands with increasing $B$ , which is associated with a larger yield stress. On the other hand, with increasing $R$ , the SH plug becomes asymmetric, which is attributed to the interface shape and the asymmetric flow field.

Figure 7.

The CFD results: streamwise velocity contours within and near the groove for $R=0.01$ (a,e), $50$ (b, f), $100$ (c,g), $500$ (d,h), at $B=1$ (a–d) and $B=10$ (e–h); $d=0.05$ , $\varTheta =95^\circ$ , $\varphi =0.5$ , $\ell =0.2$ . Solid black lines are liquid–air interface; dashed white lines are SH plug boundary. Solid and dashed contours denote positive and negative streamwise velocity, respectively.

The decrease of the dimensionless slip velocity, $u_s$ , with increasing Reynolds number should be interpreted as a feature of the present transverse-groove Bingham configuration, rather than as a general inertial effect in SH flows. In fact, in this work, we use the average velocity as the velocity scale, so that the characteristic stress becomes a viscous stress defined based on the average velocity. This choice generally leads to a decrease in the dimensionless slip velocity with increasing Reynolds number. Apart from the choice of scaling, another aspect is the possible effect of the Reynolds number on the interface dynamics and hence on the slip velocity, which may differ between transverse and longitudinal flow configurations. For example, for Newtonian flows over longitudinal ridges, inertial effects may either increase or decrease the flow-rate gain and the slip velocity depending on the meniscus geometry and channel configuration, as reported in Game et al. (Reference Game, Hodes and Papageorgiou2019). For our transverse flow configuration, as also reported in Rahmani & Taghavi (Reference Rahmani and Taghavi2022) and Rahmani et al. (Reference Rahmani, Kumar, Greener and Taghavi2023), the dimensionless slip velocities generally decrease with increasing Reynolds number, and the slip profile becomes more asymmetric. Similar effects were also reported for Newtonian transverse flow, as reported in Teo & Khoo (Reference Teo and Khoo2014).

Figure 8.

The CFD results showing the liquid–air interface morphology for different Reynolds and Bingham numbers: $B=0$ in (a,d), $1$ in (b,e), $10$ in (c, f). Panels (a–c) correspond to $\varTheta =95^\circ$ , while (d–f) correspond to $\varTheta =120^\circ$ . In all cases, $d=0.05$ , $\varphi =0.5$ , $\ell =0.2$ . The line $\eta =-1$ denotes the groove bottom wall. Arrows indicate the locations where the liquid–air interface depins from the groove edges.

Figure 8. Long description

Panel A: A line graph shows the liquid-air interface morphology for Reynolds number (R) values of 0.01, 35, 50, 500, and 1000. The x-axis ranges from -ℓ/4 to ℓ/4, and the y-axis (η) ranges from -1.0 to 1.0. Panel B: Another line graph displays the interface morphology for the same R values. Panel C: A third line graph illustrates the interface morphology for the same R values. Panel D: A line graph shows the liquid-air interface morphology for Bingham numbers (Bn) of 0.01, 50, 100, and 500. The x-axis ranges from -ℓ/4 to ℓ/4, and the y-axis (η) ranges from -1.0 to 1.0. Panel E: Another line graph displays the interface morphology for the same Bn values. Panel F: A third line graph illustrates the interface morphology for the same Bn values. The line denotes the groove bottom wall, and arrows indicate where the liquid-air interface depins from the groove edges.

We take a closer look at the interface dynamics in figure 8. In this figure, the interface position is plotted for several Reynolds numbers from $R=0.01$ to $R=1000$ , for Newtonian flow (figures 8 a and 8 d) and for viscoplastic flows with $B=1$ (figures 8 b and 8 e) and $B=10$ (figures 8 c and 8 f). The contact angle is fixed at $\varTheta =95^\circ$ and $120^\circ$ for the figure 8(a–c) and figure 8(d–f), respectively. To quantify the onset of depinning in our simulations, we define a depinning criterion based on the displacement of the liquid–air interface from the groove edge: if the interface, located at $\alpha _a=0.5$ , remains within the first computational cell adjacent to the edge, it is considered pinned otherwise it is depinned. In typical VOF solvers, such as the OpenFOAM solver used in this work, at each time step the solver modifies the interface curvature at the wall using the imposed contact angle condition, which eventually leads to a change in the surface tension source term in the momentum equation. This change in the source term modifies the velocity and pressure fields near the wall. Then, the updated velocity field advects the interface (the $\alpha$ phase), and eventually the contact line moves because of the interface advection. Based on the above description, numerically, for a given mesh, depinning occurs when the contact line displaces exactly one cell from the groove edge along the groove sidewall. Of course, the absolute value of the displacement and the onset of depinning may depend on the cell size adjacent to the groove edge; however, we believe this is a decent way to characterize the depinning onset, because such a criterion is developed based on the numerical approach that the VOF solver uses to advect the interface. However, a mesh study is necessary to ensure that the depinning onset, characterized by critical dimensionless numbers (e.g. Reynolds number), is independent of the cell size. A mesh-sensitivity study confirmed that further refinement does not significantly change the detected depinning threshold; only minor changes in the detailed interface morphology were observed, which became negligible for the finer meshes (see Appendix G). Let us first discuss the results for $\varTheta =95^\circ$ . Considering the Newtonian flow, the interface is again nearly flat and pinned at the edges for $R=0.01$ . At $R=35$ and $50$ , where the flow is still viscous-dominated, we observe the concave–convex mode of the interface. In addition, for both Reynolds numbers $R = 35$ and $50$ , the interface first depins from the upstream edge. Interface state at the downstream edge, however, differs between the two cases: at $R = 35$ the interface remains pinned, whereas at $R = 50$ it depins, although the resulting displacement into the groove is less pronounced than at the upstream edge. This behaviour that the depinning first occurs at the upstream edge is expected for viscous-dominated flow, as discussed in § 4.1.1. When the Reynolds number increases to $R=500$ and $R=1000$ , the interface dynamics change completely: the interface shape mode becomes convex–concave–convex, and depinning occurs at the downstream edge. For viscoplastic flows, as the Bingham number increases, the interface deflection is amplified while the interface shape modes are preserved. In fact, increasing $B$ magnifies the interface deflection and causes greater penetration of the depinned interface into the groove from the downstream edge. This is due to the increased shear at the interface as $B$ increases.

For $\varTheta = 120^\circ$ (figure 8 d–f), at $R=50$ , depinning occurs only at the upstream edge for all values of Bingham number, in contrast to the $\varTheta = 95^\circ$ case, where depinning takes place at both groove edges. Another notable difference concerns the extent of interface penetration into the groove. Under identical conditions at $B=1$ and $R=50$ , the interface depins from the upstream edge and then penetrates into the groove up to an $\eta \approx -0.05$ for $\varTheta = 120^\circ$ , whereas such a penetration extends to $\eta \approx -0.17$ for $\varTheta = 95^\circ$ . These observations highlight the effective role of a larger contact angle in delaying depinning, suppressing interface penetration and thus sustaining the Cassie state.

The CFD regime boundaries in figure 8 should be interpreted as the detected onset of depinning from the initially pinned state during the imposed parameter sweep. They do not resolve possible hysteresis loops or metastable intermediate wetting states after depinning, which would require a dedicated history-dependent or continuation analysis (Varchanis et al. Reference Varchanis, Dimakopoulos and Tsamopoulos2017; Pettas et al. Reference Pettas, Karapetsas, Dimakopoulos and Tsamopoulos2022).

The results in figure 8 were obtained by prescribing a single equilibrium contact angle. This treatment should be interpreted as representing the limiting case of weak contact-angle hysteresis. Such an approximation is relevant for low-pinning hydrophobic or SH surfaces, for which the advancing and receding contact angles may differ only slightly. For example, liquid-like hydrophobic coatings can exhibit contact-angle hysteresis as small as $1^\circ$ – $3^\circ$ (Barrio-Zhang et al. Reference Barrio-Zhang, Ruiz-Gutiérrez, Armstrong, McHale, Wells and Ledesma-Aguilar2020). In addition, ultrahydrophobic surfaces with advancing/receding water contact angles of $177^\circ /177^\circ$ and $172^\circ /169^\circ$ have also been reported (Chen et al. Reference Chen, Fadeev, Hsieh, Öner, Youngblood and McCarthy1999; Öner & McCarthy Reference Öner and McCarthy2000). Thus, the calculations in figure 8, based on the equilibrium contact angle, provide a useful reference for the response of the liquid–air interface when contact-angle hysteresis is weak, particularly while the contact line remains pinned.

Figure 9.

The CFD results: liquid–air interface morphology for varying Reynolds and Bingham numbers using the dynamic contact-angle condition. Panels correspond to $B=0$ (a), $B=1$ (b) and $B=10$ (c), while the Reynolds number is varied as indicated in the legend. In all cases, $d=0.05$ , $\varphi =0.5$ , $\ell =0.2$ , $\varTheta =95^\circ$ , $\theta _{{adv}}=120^\circ$ and $\theta _{\textit{rec}}=70^\circ$ . The line $\eta =-1$ denotes the groove bottom wall. Arrows indicate the detected liquid–air interface depinning locations.

To assess the influence of larger contact-angle hysteresis, we also performed an additional set of CFD simulations using a dynamic contact-angle condition. In these simulations, the equilibrium contact angle is set to $\varTheta =95^\circ$ , while the advancing and receding limits are prescribed as $\theta _{{adv}}=120^\circ$ and $\theta _{\textit{rec}}=70^\circ$ , respectively. This wider advancing–receding interval reflects the fact that hydrophobic surfaces can exhibit substantial hysteresis depending on surface chemistry, roughness and contact-line pinning. Previous dynamic-wetting simulations and measurements on hydrophobic substrates support advancing angles in the range $\theta _{{adv}}\simeq 105^\circ$ – $120^\circ$ and receding angles in the approximate range $\theta _{\textit{ rec}}\simeq 70^\circ$ – $100^\circ$ . For example, Yokoi et al. (Reference Yokoi, Vadillo, Hinch and Hutchings2009) used measured static advancing and receding angles of $107^\circ$ and $77^\circ$ , respectively, while Kim et al. (Reference Kim, Kavehpour and Rothstein2015) reported static advancing and receding angles of $98^\circ$ and $81^\circ$ on smooth polytetrafluoroethylene (PTFE). In addition, Malgarinos et al. (Reference Malgarinos, Nikolopoulos, Marengo, Antonini and Gavaises2014) considered hydrophobic VOF cases with $\theta _{{adv}}\simeq 105^\circ$ – $120^\circ$ , and polymeric hydrophobic surfaces such as PTFE may exhibit even larger hysteresis, with receding angles near $100^\circ$ despite much larger advancing angles (Pericet-Camara et al. Reference Pericet-Camara, Moraila-Martinez, Rodriguez-Valverde and Cabrerizo-Vilchez2009). The selected interval, $\theta _{\textit{ rec}}=70^\circ$ to $\theta _{{adv}}=120^\circ$ , is therefore used as a representative hysteretic hydrophobic wetting condition. The corresponding interface morphologies obtained with this dynamic wetting model are shown in figure 9(a–c) for Newtonian and viscoplastic fluid with $B=1$ , $10$ , respectively.

As shown in figure 9, the dynamic contact-angle condition preserves the main inertia-driven evolution of the interface morphology observed in figure 8. At low Reynolds number, $R=0.01$ , the interface remains nearly flat and pinned at both groove edges for all values of $B$ , indicating that the dynamic contact-angle condition has little influence while the contact line is not moving. As $R$ increases, the interface progressively develops an asymmetric shape and the local contact-angle condition is first violated at the upstream groove edges. Notably, the effect of the Bingham number remains qualitatively similar to that observed in the constant- $\varTheta$ simulations (cases modelled with equilibrium contact angle; see figure 8). Increasing $B$ amplifies the interface deformation and enhances the penetration of the liquid–air interface into the groove, particularly at larger Reynolds numbers.

A comparison between figure 9 and the constant- $\varTheta =95^\circ$ cases in figure 8(a–c) shows that the dynamic contact-angle model mainly affects the depinning behaviour near the groove edges, while the overall inertia-induced change in interface shape is preserved. In particular, allowing the apparent contact angle to vary between $\theta _{\textit{rec}}=70^\circ$ and $\theta _{{adv}}=120^\circ$ modifies the extent of contact-line displacement after depinning. For several cases, the dynamic contact-angle condition reduces excessive penetration compared with the constant- $\varTheta =95^\circ$ case, because the interface can adjust within the imposed advancing–receding interval before a stronger displacement into the groove develops. Thus, the dynamic wetting model does not change the main physical mechanism of interface evolution, namely inertia-driven asymmetric interface shapes, but it provides a more accurate prediction of the interface depinning onset and postdepinning contact line motion in the groove.

The SH plug morphology is examined in greater detail in figure 10, which presents the plug geometry across a wide range of groove aspect ratios and Reynolds numbers. The results clearly show that the plug remains symmetric in the creeping-flow regime but becomes increasingly asymmetric as $R$ increases. At $R=50$ and $R=100$ , the SH plug is oriented towards the downstream edge; however, at $R=500$ , it is oriented towards the upstream edge. This behaviour is related to the interface deflection, as discussed previously. Furthermore, increasing the groove aspect ratio enlarges the SH plug area for all Reynolds numbers considered. This is because, as the groove aspect ratio increases, the slip number increases, which promotes SH plug growth.

Figure 10.

The CFD results: SH plug morphology versus groove aspect ratio $d$ for (a) $R=0.01$ , (b) $R=50$ , (c) $R=100$ , (d) $R=500$ , at $B=10$ ; $\varphi =0.5$ , $\varTheta =95^\circ$ , $\ell =0.2$ .

6.3. Effective slip length

The overall performance of SH surfaces can be characterized by their effective slip length, which represents the surface’s capacity to reduce drag and enhance flow. Previous analyses of the flow fields have revealed the dynamics of slip velocity, interface deflection and SH plug formation. We now evaluate how the geometrical and flow parameters influence the effective slip length in figure 11. This analysis examines the dependence of the effective slip length, $\chi _e$ , on the groove aspect ratio $d$ , Bingham number $B$ , Reynolds number $R$ , slip-area fraction $\varphi$ and groove periodicity length $\ell$ .

In figures 11(a) and 11(b), the effective slip length is presented over a wide range of groove aspect ratios for the creeping and inertial flow regimes, respectively. In figure 11(a), for cases without an SH plug, the effective slip length is obtained from the semianalytical, approximate (5.3), and explicit (5.4) solutions, whereas for cases with an SH plug, the CFD results are used. The regions with and without an SH plug are separated by the colour-gradient line, which represents the onset of SH plug formation for different Bingham numbers, as indicated by the colourmap. The figure shows that the effective slip length increases with the groove aspect ratio $d$ and the Bingham number $B$ . Furthermore, comparison of the semianalytical, approximate and explicit solutions in figure 11(a) shows that, when the SH plug is absent (i.e. at small values of $d$ ), the approximate and explicit approaches agree well with the semianalytical solution, particularly the semianalytical and explicit solutions, whose results nearly match. We may also use the explicit and approximate solutions for the effective slip length as useful correlations when an SH plug is present. Extending these solutions to the SH-plug regime, we observe a reasonable agreement between the CFD and explicit results, whereas a more considerable discrepancy exists between the CFD results and the approximate solution. This analysis reveals that the explicit solution for the effective slip length serves as an accurate prediction before SH plug formation and as a reasonably accurate correlation in the presence of an SH plug.

Figure 11(b) further demonstrates that, under inertial flow conditions, $\chi _e$ increases monotonically with $d$ across all examined Bingham numbers. However, for a fixed groove aspect ratio and Bingham number, increasing the Reynolds number leads to a reduction in the effective slip length. This can be mainly attributed to the decrease in the dimensionless slip velocity, as shown in figure 6. Treating the approximate and explicit solutions as useful correlations for the regime with an SH plug (typically at large $B$ ), figures 11(c) and 11(d) present the scaled effective slip length ( $\chi _e/\ell$ ) versus the Bingham number for different slip-area fractions and groove periodicity lengths, respectively. These comparisons reveal that $\chi _e$ increases with increasing $B$ , and this enhancement becomes more pronounced for larger values of $\varphi$ and $\ell$ . This implies that, as the grooves become wider and the channel becomes thinner, the effective slip length increases. Figures 11(c) and 11(d) also indicate that the effective slip length ( $\chi _e$ ) is typically one to two orders of magnitude smaller than the groove periodicity length $\ell$ .

Figure 11.

Semianalytical and CFD results: effective slip length variation with geometric and flow parameters at $\ell =0.2$ , $\varphi =0.5$ . (a,b) Effective slip length versus groove aspect ratio $d$ for $B=0.1$ , $1$ , $10$ , at (a) $R=0.01$ and (b) $R=50$ , $500$ . (c,d) Normalized effective slip length versus $B$ at $d=0.05$ : (c) $\varphi =0.3$ , $0.5$ , $0.7$ with $\ell =0.2$ ; (d) $\ell =0.1$ , $0.2$ with $\varphi =0.5$ . In (a), (c) and (d), solid lines denote the explicit solution (5.4) in the regime where the shallow-groove semianalytical model provides an exact prediction; however, dotted lines denote the use of the same expression only as a correlation-type estimate. Colour-gradient line in (a) (see colourbar): effective slip length at $d_{\textit{cr}}$ from semianalytical solution. In (a), (c) and (d), dashed line, approximate solution (5.3). In (a) black and coloured symbols denote semianalytical and CFD results, respectively; in (b)–(d) symbols denote CFD. In (a), grey shading marks SH plug formation regime. All panels $\varTheta =95^\circ$ .

It would be insightful to compare the pressure drop provided by the SH surface for our SH channel with the pressure drop of an equivalent no-slip channel with increased height due to the addition of the groove depth, i.e. $\hat H+\hat d$ . This could demonstrate how the SH surface strategy could be effective in reducing the pressure drop compared with simply increasing the channel height for a no-slip channel. For $B=10$ , $R=0.01$ , $\hat L=0.2 \times 10^{-4}$ and $d=1$ , the pressure drop over one groove period is $0.275\,\mathrm{Pa}$ , for the SH surface. For a no-slip reference case in which the channel half-height is increased to account for the groove depth, i.e. $\hat H+\hat d$ , the pressure drop over the same streamwise period is found to be $0.302\,\mathrm{Pa}$ , which is still larger than the SH value of $0.275\,\mathrm{Pa}$ . Thus, even after accounting for the increased channel height, the SH surface provides an additional pressure-drop reduction of approximately 9.8 %.

Figure 12.

Semianalytical, explicit and CFD results: comparison of effective slip length from the explicit solution, $\chi _{e,2}$ (5.4), with (a) semianalytical and (b) CFD predictions. Symbol face colour $B$ – log scale in (a), linear in (b); edge colour: slip area fraction $\varphi$ ; size increases with groove aspect ratio $d$ . Triangles, $\ell =0.2$ ; circles, $\ell =0.1$ . In (a) $\textit{Ca}=1.37\times 10^{-4}$ ; in (b) $R=50$ . Both panels $\varTheta =95^\circ$ .

Figure 12. Long description

Panel A: A scatter plot compares semianalytical predictions of effective slip length. The x-axis is labeled ‘χe (Semianalytical)’ and the y-axis is labeled ‘χe,2’. The color of the symbols represents the log scale of B, ranging from -3 to 1. The size of the symbols increases with the groove aspect ratio, and the edge color indicates the slip area fraction φ. Triangles represent ℓ = 0.2, and circles represent ℓ = 0.1. The plot includes a dashed line indicating a reference or expected trend. Panel B: Another scatter plot compares CFD predictions of effective slip length. The x-axis is labeled ‘χe (CFD)’ and the y-axis is labeled ‘χe,2’. The color of the symbols represents the linear scale of B, ranging from 0.1 to 50. The size and edge color of the symbols follow the same conventions as in Panel A. The plot also includes a dashed line indicating a reference or expected trend.

To better compare the effective slip length obtained from the explicit solution ( $\chi _{e,2}$ ) with the semianalytical and CFD results, these values are plotted against one another in figure 12. In figure 12(a), where the semianalytical results are presented, the SH plug has not formed; however, in figure 12(b), the SH plug is present for most of the data, and the comparison is therefore made with the CFD results. The strong agreement observed in both comparisons, indicated by data points lying close to the line of perfect agreement, confirms the high predictive accuracy of the explicit model. This consistency is maintained across the entire range of Bingham numbers, groove aspect ratios and slip-area fractions examined. In addition, the close correspondence between the explicit solution and the CFD results show that the explicit formulation remains reliable even when an SH plug forms along the liquid–air interface. A small deviation appears only at large Bingham numbers ( $B \geqslant 10$ ) and at a groove periodicity length of $\ell = 0.2$ , where the SH plug becomes more pronounced. However, this difference remains minor, and the explicit model can still be treated as an accurate correlation for predicting the effective slip length in the presence of SH plug.

Figure 13.

The CFD results showing the regime classification of the liquid–air interface morphology in the groove-aspect-ratio–Reynolds number ( $d$ – $R$ ) plane. (a,b), (c,d) and (e, f) correspond to $B=0$ , $B=1$ and $B=50$ , respectively; (a,c,e) correspond to $\varTheta =95^\circ$ , while (b,d, f) correspond to $\varTheta =120^\circ$ . Green, purple and red symbols denote Regimes I, II-a and II-b, respectively; Regime I corresponds to an interface pinned at both groove edges; Regime II-a corresponds to depinning from the upstream edge; Regime II-b corresponds to depinning from the downstream edge. Black cross symbols ( $\times$ ) denote CFD cases for which depinning was observed, but the simulations did not allow us to clearly determine whether depinning first occurred at the upstream or downstream groove edge.

Figure 13. Long description

Panel A: A scatter plot shows the regime classification for a specific groove aspect ratio. The x-axis represents the Reynolds number (R) on a logarithmic scale from 10^1 to 10^3, and the y-axis represents the dimensionless parameter d ranging from 0.05 to 1.00. Green circles denote Regime I, purple stars denote Regime II-a, and red squares denote Regime II-b. Black crosses indicate cases where depinning was observed but the exact edge could not be determined. Panel B: Another scatter plot for a different groove aspect ratio with the same axes and symbols. Panel C: A scatter plot for another groove aspect ratio with the same axes and symbols. Panel D: A scatter plot for yet another groove aspect ratio with the same axes and symbols. Panel E: A scatter plot for a different groove aspect ratio with the same axes and symbols. Panel F: A scatter plot for another groove aspect ratio with the same axes and symbols.

6.4. Regime mapping

Based on the detailed investigation of the liquid–air interface morphology presented in figures 5, 7 and 8, two distinct interfacial regimes are identified: a pinned interface (Regime I) and a depinned interface (Regime II), where the latter may involve depinning from the upstream edge (Regime II-a) or from the downstream edge (Regime II-b). Figure 13(a–f) classify these regimes in the $d$ – $R$ parameter space for Newtonian ( $B=0$ ) and viscoplastic fluids with $B=1$ and $B=50$ , at two contact angles $\varTheta =95^\circ$ and $\varTheta =120^\circ$ . For the Newtonian case (figures 13 a and 13 b), at low Reynolds numbers, the interface remains pinned (Regime I, green symbols) across all groove aspect ratios. As inertia increases ( $R \gtrsim 70$ for $\varTheta =95^\circ$ and $R \gtrsim 200$ for $\varTheta =120^\circ$ ), the interface transitions to Regime II-b (red symbols), in which detachment occurs only at the downstream edge, predominantly for deeper grooves ( $d \gtrsim 0.2$ ). At moderate Reynolds numbers ( $50 \lesssim R \lesssim 70$ for $\varTheta =95^\circ$ and $80 \lesssim R \lesssim 200$ for $\varTheta =120^\circ$ ), this transition follows a sequential pattern: the interface first depins from the upstream edge (Regime II-a, purple stars) and eventually detaches only at the downstream edge (Regime II-b, red squares). Depinning from the upstream edge also occur in shallow grooves ( $d \lesssim 0.2$ ), but at wider ranges of Reynolds number ( $30 \lesssim R \lesssim 60$ for $\varTheta =95^\circ$ and $40 \lesssim R \lesssim 80$ for $\varTheta =120^\circ$ ). The introduction of yield stress substantially suppresses Regime II-a. For $B=1$ (figures 13 c and 13 d), the region corresponding to depinning from upstream edge (Regime II-a, purple symbols) shrinks considerably. As the yield stress becomes large ( $B=50$ , figures 13 e and 13 f), the influence is even more pronounced.

A direct comparison between the regime classification in figure 13 for $\varTheta = 95^\circ$ (figure 13 a,c,e) and $\varTheta = 120^\circ$ (figure 13 b,d, f) further highlights the stabilizing role of a larger contact angle. Increasing $\varTheta$ not only delays the onset of depinning from the upstream edge but also shifts the entire transition sequence (Regimes II-a $\rightarrow$ II-b) towards higher Reynolds numbers. Consequently, the depinned regimes require stronger inertial forcing to emerge. This trend is consistent across all Bingham numbers and groove aspect ratios, confirming that a larger contact angle systematically enhances the interface stability and promotes a sustained Cassie state. To further assess the influence of contact-angle hysteresis, figure 14 shows the regime classification obtained with the dynamic contact-angle condition, $\varTheta =95^\circ$ , $\theta _{{adv}}=120^\circ$ and $\theta _{\textit{rec}}=70^\circ$ . This case provides a direct comparison with the two constant-contact-angle limits in figure 13: its equilibrium value is identical to the lower case, $\varTheta =95^\circ$ , whereas its advancing limit coincides with the larger case, $\varTheta =120^\circ$ . Compared with the constant- $\varTheta =95^\circ$ results, the onset of upstream edge depinning is shifted to larger $R$ , because the contact line can sustain apparent angles up to $\theta _{{adv}}=120^\circ$ before losing pinning. However, the dynamic case does not reduce simply to the constant- $\varTheta =120^\circ$ limit, since the downstream edge is governed by the lower receding limit, $\theta _{\textit{rec}}=70^\circ$ , which controls the transition to downstream/receding-edge depinning at larger $R$ . Thus, the dynamic contact-angle condition preserves the main inertia-driven transition mechanism from pinned to depinned states, while providing a modified transition onset and regimen edge-specific classification: Regime II-a corresponds to depinning from the upstream/advancing edge (controlled by the advancing contact angle), and Regime II-b to depinning from the downstream/receding edge (controlled by the receding contact angle).

Figure 14.

The CFD regime classification of the liquid–air interface morphology in the groove aspect ratio–Reynolds number $(d-R)$ plane for $B=0$ (a), $B=1$ (b) and $B=50$ (c). In all cases, $\varphi =0.5$ , $\ell =0.2$ and the dynamic contact-angle condition is prescribed with $\varTheta =95^\circ$ , $\theta _{{adv}}=120^\circ$ and $\theta _{\textit{rec}}=70^\circ$ . Symbols denote the identified interface state: green diamonds, Regime I, pinned at both groove edges; magenta stars; Regime II-a, depinning from the upstream/advancing edge at $x=-\varphi \ell /2$ ; red squares; Regime II-b, depinning from the downstream/receding edge at $x=\varphi \ell /2$ . The regime classification highlights the edge-specific nature of depinning under transverse flow, where the first loss of pinning occurs when either the advancing or receding contact-angle limit is reached.

A closer examination of figure 13 reveals that for grooves with $d \gtrsim 0.2$ , the onset of depinning exhibits a threshold Reynolds number. In this range, there exists a critical Reynolds number beyond which the interface transitions from a pinned state (Regime I) to partial depinning at the upstream edge (Regime II-a). To evaluate this transition, we introduce the effective Reynolds number $R^* = R/(1+B)$ , which is defined based on the effective viscosity of viscoplastic flow and thus incorporates the increase in effective viscosity due to yield stress. Using this definition, the $R^*$ – $B$ plane is constructed to classify Regime I and Regime II for all cases with $d \gt 0.2$ .

As shown in figures 15(a) and 15(b) for $\varTheta =95^\circ$ and $\varTheta =120^\circ$ , respectively, the pinned regime (Regime I) occupies the lower portion of the $R^*$ – $B$ plane, while the depinned regimes (Regime II) emerges only above a well-defined threshold. The black dashed line represents a correlation extracted from the CFD data, given by $R^* = 30 B^{-1}$ and $R^* = 50 B^{-1}$ for $\varTheta =95^\circ$ and $\varTheta =120^\circ$ , respectively, which provides a reasonably accurate demarcation of the transition boundary for $B\gt 0$ and $d \gt 0.2$ . For the dynamic contact-angle case, shown in figure 15(c), the equilibrium contact angle is kept at $\varTheta =95^\circ$ , while the advancing and receding limits are prescribed as $\theta _{{adv}}=120^\circ$ and $\theta _{\textit{rec}}=70^\circ$ , respectively. The corresponding CFD data are described by the correlation $R^*=32B^{-1}$ . This value is close to the constant- $\varTheta =95^\circ$ correlation, but slightly shifted because the advancing limit allows the upstream contact line to sustain larger apparent angles before depinning, while the lower receding limit controls the downstream-edge depinning at higher Reynolds numbers. Comparing the correlation $R^* = 30 B^{-1}$ and $R^* = 50 B^{-1}$ , together with the dynamic-contact-angle correlation $R^*=32B^{-1}$ , with the effective Reynolds number relation $R^* = R/(1+B)$ , at sufficiently large $B$ (i.e. $B\gt 10$ ), we can also argue that the transition Reynolds number $R$ is approximately fixed around $R\approx 40$ and $R\approx 70$ for $\varTheta =95^\circ$ and $\varTheta =120^\circ$ , respectively, and around $R\approx 32$ for the dynamic contact-angle case, which can also be observed in figures 13 and 14. This is an interesting observation and may be justified as follows. Depinning is mainly sensitive to the flow dynamics near the groove edges, where the shear rate ( $\dot \gamma$ ) is considerably large. This means that near the edges, the effective viscosity $\mu _e = 1 + B/\dot \gamma$ approaches the plastic viscosity $\mu _{\!p} = 1$ (viscosities are in dimensionless form). This also implies that near the edges, local definition of the Reynolds numbers based on local viscosities may converge (i.e. $R^* \to R$ locally near the edges). Therefore, we may argue that the yield stress contribution to the regime transition for deep grooves is less significant, and a fixed critical $R \approx 40$ and $R\approx 70$ for $\varTheta =95^\circ$ and $\varTheta =120^\circ$ , respectively, with $R\approx 32$ for the dynamic contact-angle condition, can effectively distinguish between the pinned and depinned regimes. We emphasize that the transition Reynolds numbers reported here are specific to $\ell =0.2$ and $\varphi =0.5$ , as well as to the prescribed contact-angle model, and that for other geometric parameters or wetting conditions the transition value of Reynolds numbers may differ.

Figure 15.

The CFD regime classification of liquid–air interface morphology for $d\gt 0.2$ in the $R^{*}$ – $B$ plane, where $R^{*}=R/(1+B)$ . In all cases, $\varphi =0.5$ and $\ell =0.2$ . Panels (a) and (b) correspond to the equilibrium contact-angle cases $\varTheta =95^\circ$ and $\varTheta =120^\circ$ , respectively. Panel (c) corresponds to the dynamic contact-angle condition with $\varTheta =95^\circ$ , $\theta _{{adv}}=120^\circ$ , and $\theta _{\textit{rec}}=70^\circ$ . Green symbols denote Regime I, where the interface remains pinned at both groove edges. Magenta stars denote Regime II-a, corresponding to depinning from the upstream edge at $x=-\varphi \ell /2$ . Red squares denote Regime II-b, corresponding to depinning from the downstream edge at $x=\varphi \ell /2$ . The black dashed lines in (a) and (b) show the CFD-based correlations $R^{*}=30B^{-1}$ and $R^{*}=50B^{-1}$ , respectively, while the solid black curve in (c) shows the updated CFD-based correlation $R^{*}=32B^{-1}$ , valid for $B\gt 0$ and $d\gt 0.2$ .

We now examine the influence of capillarity on the interfacial morphology using the semianalytical model. Figure 16 presents the regime classification in the $d$ – $\textit{Ca}$ plane for two representative Bingham numbers, $B=0.1$ and $B=1$ . These results correspond to the low-Reynolds-number limit ( $R=0.01$ ), where inertial contributions are negligible and the interface dynamics is governed solely by the balance between capillary forces, imposed pressure gradients and the viscoplastic response. For both Bingham numbers, the regime map focuses on Regime I, in which the liquid–air interface remains pinned at the groove edges. Within this regime, the interface may be either flat or weakly deflected depending on the local interfacial slope: it is classified as flat when $\max (|\eta |) \leqslant 0.01$ , and as deflected otherwise. Figure 16 also delineates the thresholds for the onset of depinning from the upstream edge and for interface collapse, the latter corresponding to contact of the interface with the groove bottom wall. Both transition boundaries are predicted by the semianalytical model. Accordingly, at small capillary numbers ( $\textit{Ca} \lesssim 10^{-3}$ ) and for a broad range of groove aspect ratios, the interface remains flat and pinned at the groove edges. However, with increase in $\textit{Ca}$ , the interface, while still pinned, undergoes deflection for both values of $B$ considered (the onset of such deflection is shown by the red dashed line). Further increase in $\textit{Ca}$ lead to the onset of depinning, predicted by the semianalytical criterion defined in (4.17), and shown by the solid black line ( $d \gtrsim 10^{-2}$ ). Beyond this threshold, the interface enters Regime II-a. For very shallow grooves ( $d \lesssim 10^{-2}$ ), the semianalytical solution is also able to predict the onset of interface collapse, as discussed in § 4.1.1. This transition, marked by the dashed black line. As expected, the collapsed interface regime becomes more prominent for small values of $d$ , where even moderate deflection results in contact with the groove bottom wall. Furthermore, the onset of interface depinning occurs at a higher $\textit{Ca}$ for a larger contact angle $\varTheta$ .

Figure 16.

Semianalytical results: regime classification of liquid–air interface morphology in the groove aspect ratio–capillary number ( $d$ – $\textit{Ca}$ ) plane for (a) $B=0.1$ and (b) $B=1$ ; $\varphi =0.5$ , $\ell =0.2$ , $R=0.01$ . Green shading is pinned morphology (Regime I). In regime I, circles indicate flat pinned interfaces ( $max(|\eta |) \leqslant 0.01$ ), while diamonds indicate pinned but deflected interfaces ( $max(|\eta |) \gt 0.01$ ). The dashed red line separates the two subregimes. Solid black line is upstream depinning threshold from (4.17) with static contact angle $\varTheta =110^\circ$ and $120^\circ$ ; dashed black line is interface collapse criterion.

Figure 16. Long description

Two scatter plots labeled (a) and (b) depict the relationship between capillary number (Ca) on the horizontal axis and dimensionless distance (d) on the vertical axis. Both plots feature a logarithmic scale for both axes. Panel (a) shows data for a contact angle of 110 degrees and 120 degrees, with green shading indicating pinned morphology (Regime I). Within Regime I, circles represent flat pinned interfaces, while diamonds represent pinned but deflected interfaces. A dashed red line separates these two subregimes. Solid and dashed black lines indicate upstream depinning thresholds and interface collapse criteria, respectively. Panel (b) similarly shows data for a contact angle of 120 degrees, with the same color and symbol scheme as Panel (a).

Based on the observations in figures 5, 7 and 10, the viscoplastic flow exhibits an unyielded region, referred to as the SH plug, located above the liquid–air interface. To systematically determine the conditions under which this SH plug forms, figure 17 is constructed using the relevant geometrical parameters, namely the groove aspect ratio, slip area fraction and groove periodicity length. This figure is inspired by figure 19 of Rahmani & Taghavi (Reference Rahmani and Taghavi2024a ), where the parameters shown on the $x$ - and $y$ -axes were obtained from analytical relations for the critical slip number. Here, we instead substitute the slip number using its definition, i.e. $\bar b \approx \varphi \ell d/4\mu _a$ . The figure identifies the onset of SH plug formation using the semianalytical, approximate and explicit solutions for two representative Bingham numbers. The comparison between the semianalytical prediction (solid line) and the CFD results for $B = 0.1$ demonstrates that the semianalytical formulation accurately captures the onset of SH plug formation across the entire range of geometrical parameters. Further examination of the semianalytical, explicit and approximate criteria shows that the explicit relation provides an excellent estimate of the critical groove aspect ratio $d_{cr}$ , with only a minor deviation from the semianalytical boundary. In contrast, the approximate solution exhibits a somewhat larger offset relative to the semianalytical curve; however, it still reproduces the correct qualitative dependence of $\textrm{d}_{cr}$ on key geometrical parameters such as the groove aspect ratio, slip area fraction and groove periodicity length.

Figure 17.

Semianalytical and CFD results: regime map of SH wall plug presence at the liquid–air interface; light blue, no SH plug; red, SH plug. Blue and red markers: CFD cases for $B=0.1$ without and with SH plug, respectively. Transition boundaries from semianalytical solution (solid), explicit expression ((5.8), dash–dotted) and approximate relation ((5.7), dashed). Boundaries for $B=1$ also shown (thicker lines). Here $\textit{Ca}=1.37\times 10^{-4}$ , $\varTheta =95^\circ$ .