3.1. Numerical methods

The computations presented in this study were performed using our recently developed numerical method, which was initially described and validated in our previous work on the flow past bubbles (Fang et al. Reference Fang, Zhang, Xiong, Xu and Ni2021). Previously, this numerical tool was limited to handling oblate inviscid bubbles ( $\mu ^* = 0$ ) and rigid spheroids ( $\mu ^* = \infty$ ). However, recent advancements in the numerical method, as thoroughly described by Wei et al. (Reference Wei, Meng, Ni and Zhang2023), have extended its applicability to droplets with irregular shapes and arbitrary viscosities. In this section, together with Appendix A, we will briefly introduce the key features of the method, while more comprehensive details can be found in a separate publication (Wei et al. Reference Wei, Meng, Ni and Zhang2023).

As described in § 2, the main challenge in solving (2.1) is to accurately impose the interfacial boundary conditions (2.2) at the droplet interface. The boundary conditions (2.2) can be reformulated as follows:

(3.1) \begin{eqnarray} \left . \begin{array}{ll} u_{e,n} = u_{i,n} = 0 \\[0.1cm] \mu _e \left (\frac {\partial {u}_{e,\tau _1}}{\partial n} - \kappa _1{u}_{e,\tau _1} \right ) = \mu _i \left (\frac {\partial {u}_{i,\tau _1}}{\partial n} - \kappa _1{u}_{i,\tau _1} \right ) \\[10pt] \mu _e \left (\frac {\partial {u}_{e,\tau _2}}{\partial n} - \kappa _2{u}_{e,\tau _2} \right ) = \mu _i \left (\frac {\partial {u}_{i,\tau _2}}{\partial n} - \kappa _2{u}_{i,\tau _2} \right ) \end{array} \qquad \right \} \quad \mathrm {at} \quad \Gamma , \end{eqnarray}

with $n$ , $\tau _1$ and $\tau _2$ denoting the normal component and the two principal tangential components along $\boldsymbol {n}$ , $\boldsymbol {t}_1$ and $\boldsymbol {t}_2$ , respectively, while $\kappa _{1}$ and $\kappa _{2}$ are the corresponding interface curvatures in these two directions. This set of equations represents the kinematic jump conditions across the droplet interface. Specifically, when $\mu _i = 0$ , these conditions simplify to the free-slip boundary condition and when $\mu _i = \infty$ , they simplify to the no-slip boundary condition.

Moreover, such kinematic jumps can be easily imposed on body-fitted grids, where the grid boundaries perfectly coincide with the droplet interface. However, while body-fitted grids can accurately impose these jump conditions on a spherical droplet, they are less effective for oblate droplets (Legendre et al. Reference Legendre, Rachih, Souilliez, Charton and Climent2019; Rachih et al. Reference Rachih, Legendre, Climent and Charton2020). To address this challenge, we have developed a novel numerical methodology within the framework of Cartesian grids (Wei et al. Reference Wei, Meng, Ni and Zhang2023). This approach is non-trivial because the grid boundaries do not align with the droplet interface, leading to the partitioning of interfacial grids into interior and exterior regions of the droplet. Accurately imposing jump conditions on such meshes is particularly challenging when the droplet has an irregular shape. This task heavily relies on the precise estimation of local principal curvatures within the interfacial cells, denoted as $\kappa _1$ and $\kappa _2$ . The numerical method proposed in our recent work incorporates two key techniques to address this issue. First, a height function method is employed to determine the principal curvatures and their directions in the interfacial cells. Second, the internal and external flows (2.1) are solved simultaneously using an embedded boundary method, which ensures the correct implementation of the jump conditions described in (3.1). Detailed descriptions of these techniques can be found from Wei et al. (Reference Wei, Meng, Ni and Zhang2023) and additional information about the numerical algorithms is provided in Appendix A. For simplicity, we will not elaborate further here.

We implemented this numerical method using the open-source library Basilisk (Popinet Reference Popinet2003, Reference Popinet2015), which employs the finite volume method to discretise the equations, thereby preserving their conservative properties. In general, Basilisk uses an approximate projection method to decouple the velocity and pressure equations (2.1). The well-known Bell–Colella–Glaz scheme (Bell et al. Reference Bell, Colella and Glaz1989) is employed for discretising the convection term, while a fully implicit scheme is used for the viscous term. Additionally, Basilisk incorporates an adaptive mesh refinement technique to enhance computational efficiency by refining the mesh in regions close to the interface and those with concentrated vorticity. Detailed descriptions of these algorithms can be found in the references by Popinet (Reference Popinet2003, Reference Popinet2015) and we do not reiterate them here. It is worth noting that the original Basilisk library was designed to simulate flow past rigid particles with no-slip boundary conditions. Extending these algorithms to handle droplets with arbitrary viscosity constitutes our novel contribution in developing the numerical method.

Appendix A includes validation tests, specifically focusing on benchmark problems involving shear flow past a frozen body. These tests have been conducted to verify the accuracy and validity of the aforementioned numerical methods for simulating three-dimensional flows past a viscous droplet.

3.2. Numerical settings

Figure 2.

Grid distribution within the computational domain for the parameters $(Re_e, Re^*, \mu ^*, \chi ) = (200, 1.0, 100, 2.0)$ . (a) Global view of the computational grid. (b) Enlarged view near the droplet surface. The grid resolution reaches a minimum size of $\varDelta = R/80$ in the vicinity of the droplet, with 10 layers extending both inside and outside of the interface.

The simulations were conducted in a three-dimensional cubic domain with dimensions $L \times L \times L$ . The droplet is positioned at the origin $(x, y, z) = (0, 0, 0)$ . The inlet (left boundary) is located $L/3$ to the left of the droplet, while the outlet (right boundary) is situated $2L/3$ to the right. The droplet is centred along the other four walls of the domain. Figure 2 illustrates the non-uniform mesh employed in the simulations for the case with $Re_e = Re_i = 200$ , $\mu ^* = 100$ , $\chi = 2.0$ and $Sr = 0.2$ . Figure 2(a) provides a global view of the mesh system, while figure 2(b) shows a close-up view near the droplet. The mesh features a smallest grid size of $\varDelta = R/80$ near the droplet interface, with 10 layers of grids both interior and exterior to the interface. This configuration allows for accurate resolution of the flow boundary layer, which has a thickness $\delta \propto 2R/\sqrt {Re} \simeq R/7$ . Consequently, there are at least 11 grid points within the boundary layer. Further, grids with a slightly larger size of $\varDelta = R/20$ extend to the upstream and downstream regions of the droplet, covering distances of $2R$ and $5R$ , respectively. Grids with a size of $\varDelta = R/10$ then expand to a broader range of $[-5R, 15R]$ . This strategy of spatial resolution ensures precise capture of the shear flow around the droplet and accurate depiction of the wake region. Additional numerical tests for grid independence, detailed in Appendix A, confirm that this spatial resolution is sufficient for converging results on the forces experienced by the droplet.

Another important consideration is determining the optimal size for the computational domain. This is crucial when investigating flows around a blunt body, particularly at low Reynolds numbers, where confinement effects from a smaller computational domain can significantly influence the results (Legendre & Magnaudet Reference Legendre and Magnaudet1998). To address this, a size-independent study was conducted for two different Reynolds numbers on a bubble: $Re_e = 0.1$ and $Re_e = 200$ , representing low- and moderately high- $Re$ regimes, respectively. The domain size was varied from $L = 50R$ to $L = 800R$ and results from different domain sizes are detailed in Appendix A. It is evident that a domain size of $L = 400R$ provides convergence for the low-Reynolds-number case, while $L = 100R$ already achieves convergence for the high-Reynolds-number case. Consequently, $L = 400R$ and $200R$ were selected for the simulations in these two $Re$ regimes reported in the manuscript. Note that simulations were not performed at Reynolds numbers lower than $0.1$ due to the excessive domain size requirements. Furthermore, the case with $Re_e = 0.1$ sufficiently captures the asymptotic behaviour for the low-Reynolds-number limit.

In addition to the interfacial boundary conditions at the droplet, appropriate boundary conditions must be specified for the computational domain. A straightforward approach would be to set an inflow velocity $\boldsymbol {U} = (U_\infty + \alpha y)\boldsymbol {e}_x$ at the left boundary ( $x = -L/3$ ) and a parabolic outflow condition at the right boundary ( $x = 2L/3$ ). However, these settings can lead to numerical instability for this shear flow problem. Specifically, the unperturbed velocity $U_x = \boldsymbol {U} \cdot \boldsymbol {e}_x$ becomes negative in the region where $y \lt -U_\infty /\alpha$ , resulting in a reversal of the inlet and outlet boundaries in this negative $U_x$ region. Consequently, the unique inflow (outflow) condition at the left (right) boundary becomes non-physical, potentially leading to numerical disturbances that compromise the accuracy of the simulations. To address this issue, Legendre & Magnaudet (Reference Legendre and Magnaudet1998) proposed a hybrid boundary condition at the downstream boundary. This method involves imposing the standard outflow condition for the region within the disk $(y^2 + z^2)^{1/2} \lt U_\infty /\alpha$ , while applying the inflow condition $(U_\infty + \alpha y)$ for the region outside this disk $(y^2 + z^2)^{1/2} \gt U_\infty /\alpha$ . In our study, we tested several boundary condition configurations to minimise numerical contamination in the flow field. Ultimately, we found that applying periodic boundary conditions at both the upstream and downstream boundaries yielded the most reliable results. Validation tests presented in Appendix A also support the effectiveness of this approach. Additionally, no-slip conditions, $\boldsymbol {U} = (U_\infty + \alpha y)\boldsymbol {e}_x$ , are applied to the four lateral walls (top, bottom, front and back) of the domain.

As previously noted, the problem under investigation is characterised by five dimensionless parameters $(Re_e, Sr, Re^*, \mu ^*, \chi )$ . However, our focus is not on creating a comprehensive phase map of flow features based on these parameters. Instead, we concentrate primarily on the phenomenon of lift reversal as the droplet transitions from an inviscid oblate bubble to a rigid spheroid. Therefore, in §§ 4 and 5, we focus on the shear rate of $Sr = 0.2$ , while the effects of varying the shear rate are discussed in § 6. Throughout this study, we investigate only one Reynolds number ratio, $Re^* = 1$ . Although exploring additional values of $Re^*$ could be insightful, given that some studies (Edelmann et al. Reference Edelmann, Clercq, Patrick and Noll2017; Rachih Reference Rachih2019; Shi et al. Reference Shi, Climent and Legendre2024) have identified secondary flow instabilities at higher $Re^*$ values, a detailed examination of $Re^*$ is beyond the scope of the current work. A brief discussion of such secondary instabilities is provided in Appendix E , while a more comprehensive discussion of this scenario can be found in the very recent work by Shi et al. (Reference Shi, Climent and Legendre2024), almost in synchronisation with the present study. For simplicity, we will use $Re$ throughout the manuscript, as $Re_e$ and $Re_i$ are equivalent in this study. In summary, our primary focus is on the lift response of the droplet by varying $Re$ , $\mu ^*$ and $\chi$ , within the ranges $Re \in [0.1, 300]$ , $\mu ^* \in [0.01, 100]$ and $\chi \in [1.0, 2.5]$ . Particularly, cases with $\mu ^* = 0.01$ and $\mu ^* = 100$ are used to model an inviscid bubble and a rigid spheroid, respectively, in the numerical simulations.

4.1. Spherical droplets

Figure 3.

Lift coefficients $C_L$ experienced by spherical droplets in the low-but-finite $Re$ regime, with variations in the Reynolds number and viscosity ratio. (a) $C_L$ as a function of $\mu ^*$ at different $Re$ . (b) Comparison between numerical results and analytical correlations for spherical droplets at $Re = 0.1$ ; the solid line represents (4.1), (4.2a ) from Mei (Reference Mei1992), and the dashed line is (4.1), (4.2b ) from Legendre & Magnaudet (Reference Legendre and Magnaudet1998), while for both lines, we use $C (\mathcal {R}) \propto \mathcal {R} ^2$ in (4.1) predicted by Legendre & Magnaudet (Reference Legendre and Magnaudet1997). The grey triangles and squares in panel (b) indicate the pressure and viscous components of the lift coefficients, respectively.

To begin with, we focus on spherical droplets. Figure 3(a) presents the lift forces experienced by these droplets, which have varying Reynolds numbers and viscosity ratios but a constant shear rate of $Sr = 0.2$ . The results show that the lift coefficient $C_L (Re, \mu ^*, \chi = 1.0)$ increases with $\mu ^*$ but decreases with $Re$ . To understand this trend, it is important to review the state-of-the-art of this problem in the low-but-finite Reynolds number regime. It was Saffman (Reference Saffman1965) who first used the asymptotic method to investigate the lift force acting on a particle in the low-Reynolds-number regime with a very high shear rate, under the assumption that the flow in the far field is equivalent to that produced by a point force (whose magnitude is the Stokes drag on the particle) located at the particle’s position. This point force, combined with advective effects induced by the ambient shear, results in a small correction to the far-field flow, which takes the form of a uniform flow at leading order. This uniform flow is not aligned with the upstream flow, due to the asymmetry introduced by the shear. The force resulting from this small correction is a small Stokes drag, whose direction is collinear to the flow correction, i.e. it is not collinear to the upstream flow. The force component perpendicular to the upstream flow is the Saffman lift force, also known as the Saffman mechanism. In addition, the origin of the primary Stokes drag is the azimuthal vorticity at the particle surface: in the Stokes limit, the force and the vorticity share the same pre-factor. Hence, increasing the magnitude of the vorticity, by changing the particle shape or the boundary condition, results in an increase of the drag and vice versa. Advection in the wake is key to obtaining a non-zero lift force, as reversibility in the Stokes limit prevents the existence of a lift component for a sphere (Bretherton Reference Bretherton1962). For the wake to produce a lift force, some asymmetry is required; otherwise, the wake only provides the Oseen drag correction. Clearly, shear is responsible for this asymmetry and, thus, the Saffman mechanism can also be understood as stemming from the asymmetric convection of azimuthal vorticities produced at the particle’s interface in the presence of shear.

Several works have been devoted to extending Saffman’s asymptotic scheme to obtain the lift force on a spherical object moving in shear flow in the low-Reynolds-number regime. These efforts, such as those by McLaughlin (Reference McLaughlin1991) and Legendre & Magnaudet (Reference Legendre and Magnaudet1997, Reference Legendre and Magnaudet1998), refine the asymptotic approach to account for more complex flow conditions. In summary, the analytical solution for the lift coefficient, to the first order, has the following formulation:

(4.1) \begin{eqnarray} C_L (Re \ll 1) &=& \frac {C (\mathcal {R})}{\pi ^2}\varepsilon J_L(\varepsilon ) - C_L^{\prime(2)}, \end{eqnarray}

where $\varepsilon = (\mid Sr \mid / Re)^{1/2}$ represents the ratio of the Oseen to the Saffman length scales. Here, $C_L^{\prime(2)}$ is the second-order inertial contribution to lift, originally derived by Saffman (Reference Saffman1965) as $11Sr/8$ for rigid sphere, with more recent computations by Candelier et al. (Reference Candelier, Mehaddi, Mehlig and Magnaudet2023) yielding smaller corrections. The constant $C$ , which will be elaborated in more detail later, depends on the dimensionless parameter $\mathcal {R} = (3\mu ^* + 2) / (\mu ^* + 1)$ that is the strength of the Stokeslet in the unbounded solution. The term $J_L(\varepsilon )$ denotes a three-dimensional integral that is independent of the droplet’s interfacial condition. In the limit as $\varepsilon \rightarrow \infty$ , analytical solutions for $J_L(\varepsilon )$ were obtained by Saffman (Reference Saffman1965), with $J_L(\infty ) = 2.254$ . However, the exact manner in which $J_L(\varepsilon )$ approaches $J_L(\infty )$ as $\varepsilon$ increases remains analytically unresolved. To solve this, semi-empirical approximations for $J_L(\varepsilon )$ at limited value of $\varepsilon$ have been proposed in some studies, given by

(4.2) \begin{eqnarray} J_L (\varepsilon ) &=& J_L(\infty )(1+0.2\varepsilon ^{-2})^{-3/2} \\ \nonumber \mathrm {or} \qquad J_L (\varepsilon ) &=& 0.6765\left [1+\mathrm {tanh}[2.5(\mathrm {lg}\varepsilon + 0.191)]\right ]\left [ 0.667 + \mathrm {tan}[6(\varepsilon -0.32)] \right ], \end{eqnarray}

where the first correlation was derived through numerical data fitting by Legendre & Magnaudet (Reference Legendre and Magnaudet1998) and the second correlation was proposed by Mei (Reference Mei1992) based on earlier numerical results provided by McLaughlin (Reference McLaughlin1991). Therefore, both of these approximations for $J_L(\varepsilon )$ are known to be highly dependent on their respective numerical results. Additionally, in the low-but-finite $Re$ regime, it is important to note that the numerical solutions can be significantly affected by the confinement of the computational domain, as reviewed by Shi & Rzehak (Reference Shi and Rzehak2019).

The value of $C (\mathcal {R})$ warrants further discussion. Saffman (Reference Saffman1965) derived an analytical value of $C (\mu ^* \rightarrow \infty ) = 18$ for a rigid sphere. In contrast, Mei & Klausner (Reference Mei and Klausner1994) reported that $C (\mu ^* \rightarrow 0) = C (\mu ^* \rightarrow \infty ) \cdot \mathcal {R}(\mu ^* \rightarrow 0) / \mathcal {R}(\mu ^* \rightarrow \infty ) = 12$ for an inviscid bubble. Subsequently, Legendre & Magnaudet (Reference Legendre and Magnaudet1997) challenged this result and provided a more rigorous analysis of $C$ for arbitrary values of $\mu ^*$ (or $\mathcal {R}$ ). They argued that $C (\mathcal {R})$ is proportional to $C_{D0}^2$ or $\mathcal {R}^2$ , where $C_{D0}$ is the drag coefficient experienced by the droplet in homogeneous flow at the same Reynolds number. Accordingly, for an inviscid bubble, this implies $C (\mu ^* \rightarrow 0) = C (\mu ^* \rightarrow \infty ) \cdot \mathcal {R}^2(\mu ^* \rightarrow 0) / \mathcal {R}^2(\mu ^* \rightarrow \infty ) = 8$ . However, the correlation proposed by Legendre & Magnaudet (Reference Legendre and Magnaudet1997) for this low-but-finite Reynolds number regime has not yet been validated for droplets with arbitrary $\mu ^*$ or $\mathcal {R}$ .

Our initial approach involves comparing the present numerical solutions with the aforementioned formulae for varying viscosity ratios $\mu ^*$ , including (4.1), (4.2), and the proportionality $C(\mathcal {R}) \propto \mathcal {R}^2$ . We keep the Reynolds number fixed at $Re = 0.1$ . Figure 3(b) presents the lift coefficients as a function of the viscosity ratio, while $J_L (\varepsilon )$ is estimated using either (4.2a ) (Legendre & Magnaudet Reference Legendre and Magnaudet1998) or (4.2b ) (Mei Reference Mei1992). The results show excellent agreement with the theoretical predictions, although they are slightly lower than those predicted by (4.2a ), with a discrepancy of less than $5\,\%$ . This deviation is considered acceptable, as both of the correlations (4.2a ) and (4.2b ) were derived from numerical data specific to the referenced studies. Similar findings are reported in the recent work by Shi et al. (Reference Shi, Climent and Legendre2024) (figure 2 therein), where the ratio of the inner to outer Reynolds number was 5, in contrast to 1 in the present work. Together, these results align with the low-Re analysis conducted by Magnaudet et al. (Reference Magnaudet, Takagi and Legendre2003), which showed that the hydrodynamic force in the terminal state does not depend on the density ratio (and hence, the $Re^*$ ). Furthermore, figure 3(b) illustrates both the pressure and viscous components of the lift coefficients, demonstrating that the viscous component is more dominant than the pressure component in this low-but-finite Reynolds number regime.

Figure 4.

Comparison between numerical results and normalised correlations for $C_L$ on spherical droplets at low-but-finite $Re$ . The numerical data collapse onto a single curve expressed by (4.3) after normalisation.

Moreover, we try to normalise the lift coefficients with the formulation of $C_{L}^* = [C_{L}(Re,\mu ^*) - C_{L}(Re, \mu ^*\rightarrow 0)]/[C_{L}(Re, \mu ^*\rightarrow \infty ) - C_{L}(Re, \mu ^*\rightarrow 0)]$ . If we assume that the correlation $C_L (Re \rightarrow 0, \mu ^*) \propto C_{D0}^2 (Re \rightarrow 0, \mu ^*)$ (Legendre & Magnaudet Reference Legendre and Magnaudet1997) holds for spherical droplets in the low-but-finite $Re$ regime, and given that the Oseen formulation can be used to estimate the drag coefficient as $C_{D0}(Re \rightarrow 0, \mu ^*) \simeq 8\mathcal {R}/Re$ , we derive the following normalised lift in this regime:

(4.3) \begin{eqnarray} C_{L}^* = \frac {C_L (Re,\mu ^*) - C_{L}(Re,\mu ^*\rightarrow 0)}{C_{L}(Re,\mu ^*\rightarrow \infty ) - C_{L}(Re,\mu ^*\rightarrow 0)} &=& \frac {C_{D0}^2(Re,\mu ^*) - C_{D0}^2(Re, \mu ^*\rightarrow 0)}{C_{D0}^2(Re,\mu ^*\rightarrow \infty ) - C_{D0}^2(Re, \mu ^*\rightarrow 0)} \nonumber \\ &\approx & \frac {\left (\frac {8}{Re}\right )^2\left [\mathcal {R}(\mu ^*)^2 - 2^2\right ]}{\left (\frac {8}{Re}\right )^2\left [(3^2 - 2^2) \right ]} \nonumber \\ &\approx & \frac {(4+5\mu ^*)\mu ^*}{5(1 + \mu ^*)^2}. \end{eqnarray}

This normalised lift coefficient, $C_{L}^*$ , is compared with our numerical results, as shown in figure 4 for $0.1 \leqslant Re \leqslant 1$ . Excellent agreement is observed between the theoretical correlation and the numerical data. Thus, this normalised correlation (4.3) effectively predicts the lift force experienced by a spherical droplet in the low-but-finite $Re$ regime. In summary, our numerical results validate the conclusion by Legendre & Magnaudet (Reference Legendre and Magnaudet1997) that, for a spherical droplet in this regime, the lift force is proportional to the square of the drag force acting on the droplet.

Figure 5.

Streamlines internal and external of the droplets at $Re = 0.1$ , with the viscosity ratio varying in the range of $0.01 \leqslant \mu ^* \leqslant 100$ . As $\mu ^*$ increases, a stronger recirculating zone develops on the side with higher velocity past the droplet, thereby generating a more pronounced Saffman mechanism.

The next step is to explore why the lift force increases with larger $\mu ^*$ in the low-but-finite $Re$ regime. It is known that the Saffman mechanism is dominant in this regime (Legendre & Magnaudet Reference Legendre and Magnaudet1997, Reference Legendre and Magnaudet1998). This mechanism is associated with the asymmetric advection of azimuthal vorticities generated at the droplet in the presence of shear in the $XOY$ plane. Consequently, it is pertinent to investigate the asymmetric distribution of $\omega _z$ (i.e. the projection of $\omega _\phi$ onto the $XOY$ plane) around the spherical droplet as $\mu ^*$ varies. Figure 5 displays the internal and external streamlines on the $XOY$ plane for different values of $\mu ^*$ , while maintaining a Reynolds number of $Re = 0.1$ . It is evident that a larger $\mu ^*$ results in a more pronounced recirculating zone on the side with higher velocity past the droplet. This flow structure is consistent with those reported by Ervin & Tryggvason (Reference Ervin and Tryggvason1997) and Hidman et al. (Reference Hidman, Ström, Sasic and Sardina2022), who investigated the rise motion of bubbles in low- $Re$ shear flows. Under these conditions, the negative $\omega _z$ observed in the figure is expected to generate a positive lift force $\textbf {F}_{L}\cdot \boldsymbol {e}_y$ , given that $\textbf {F}_L \propto \boldsymbol {U} \times \boldsymbol {\omega }$ . Additionally, a stronger recirculation for larger $\mu ^*$ leads to a greater lift force. Figure 6 illustrates the distribution of $\omega _z$ on the top ( $\phi = 0$ ) and bottom ( $\phi = \pi$ ) surfaces of the droplet, showing the dependence of $\omega _{\phi = 0}$ and $\omega _{\phi = \pi }$ on $\theta$ . The figure reveals that a more viscous droplet exhibits a larger $\Delta \omega _\phi = \omega _{\phi = 0} - \omega _{\phi = \pi }$ , confirming that the asymmetric advection of $\omega _\phi$ is responsible for the lift force experienced by the droplet in this low-but-finite $Re$ regime.

Figure 6.

Distribution of $\omega _z$ at the (a) top and (b) bottom surfaces of the droplet. A more viscous droplet exhibits a larger $\Delta \omega _\phi = \omega _\phi (\phi = 0) - \omega _\phi (\phi = \pi )$ at the equator, thereby enhancing the Saffman mechanism.

4.2. Oblate droplets

We now shift our focus to oblate droplets in the low-but-finite $Re$ regime ( $0.1 \leqslant Re \leqslant 1$ ), with aspect ratios varying in the range of $1 \leqslant \chi \leqslant 2.5$ , which encompasses moderate to high shape deformations. As with spherical droplets discussed in § 4.1, our primary concerns are to understand how and why the lift force evolves with the changing parameters. It is well known that deformable drops and bubbles moving in shear flows can experience lateral migration, even in the zero-Reynolds-number limit, as demonstrated by Leal (Reference Leal1980). Even in the low-Reynolds-number regime, the problem becomes nonlinear due to the fact that the matching of velocities and stresses across the droplet surface must be satisfied on a deformed interface, whose shape is a part of the solution to the problem. Notable contributions to this phenomenon include the works by Chan & Leal (Reference Chan and Leal1979) and Magnaudet et al. (Reference Magnaudet, Takagi and Legendre2003), who analytically solved for the deformation rate of the droplet and the corresponding deformation-induced lift using the reciprocal theorem. Their studies focused on deformation caused by the presence of a wall and shear flow. They found that droplets migrated in a wall-bounded shear flow, with the direction of deformation-induced migration changing with the viscosity ratio, and the migration was directed towards the wall when the viscosity ratio was $\mu ^*{-} O (1)$ . In contrast, for an almost inviscid bubble rising in a quiescent liquid, the deformation-induced migration was directed away from the wall, in agreement with experimental results by Takemura et al. (Reference Takemura, Takagi, Magnaudet and Matsumoto2002). Sugiyama & Takemura (Reference Sugiyama and Takemura2010) also analysed the deformation-induced lift when the bubble was extremely close to the wall. However, the deformation-induced lift force in these studies was derived for a linear shear flow near a wall, where the deformation is somewhat asymmetric compared with the case of an unbounded flow. This differs from the present study, where we assume the deformable droplet maintains an axi-symmetrical shape. The asymmetric deformation of bubbles has been shown to produce a negative lift force in unbounded flow, referred to as the $A\hbox{-}$ mechanism (Ervin & Tryggvason Reference Ervin and Tryggvason1997; Hidman et al. Reference Hidman, Ström, Sasic and Sardina2022), which acts in the opposite direction to the Saffman mechanism. This effect, however, is beyond the scope of the present study.

Figure 7.

Interpolation scheme for predicting $C_L$ on oblate droplets in the low-but-finite $Re$ regime for (a) $Re = 0.1$ , (b) $Re = 0.5$ and (c) $Re = 1$ . Points represent numerical results, while the solid lines correspond to the correlation proposed as (4.4).

Figure 7 illustrates the dependence of $C_L (Re, \mu ^*, \chi )$ on $\mu ^*$ at $Re = 0.1$ , $0.5$ and $1$ , with aspect ratios varying as $\chi = 1.0$ , $\chi = 1.5$ , $2.0$ and $2.5$ in each panel. The first two pictures demonstrate that larger aspect ratios correspond to higher lift coefficients, regardless of the viscosity ratio under consideration. This effect arises because the oblateness of the droplet influences the curvature of the interface and the vorticities generated at the interface are positively correlated with the curvature of the equator. For example, $\omega _{\phi , {max}} \sim \chi ^{8/3}$ for an inviscid oblate bubble (Magnaudet & Mougin Reference Magnaudet and Mougin2007). Consequently, the Saffman mechanism is more pronounced for more oblate droplets in this $Re$ regime. In contrast, figure 7(c) shows that at $Re = 1$ , the earlier trend – that larger aspect ratios result in higher lift coefficients – largely persists. However, an anomalous reversal is observed for $\mu ^* \geqslant 10$ at $\chi = 1.5$ , where the lift force is smaller than that for $\chi = 1.0$ . This is not unexpected, as at $(Re, \mu ^*) = (5, 100)$ , the lift force was already found to decrease monotonically with $\chi$ (see figure 12 c). This phenomenon can be attributed to the dominance of $S\hbox{-}$ mechanism , which generates negative lift on highly viscous droplets. Consequently, the $S\hbox{-}$ mechanism may already be present at $Re = 1$ , partially offsetting the Saffman mechanism. To substantiate the observation of an intensified Saffman mechanism for oblate droplets, a detailed examination of the vorticity distribution around the droplet and its evolution in the flow field is essential. Figure 8 presents the distribution of $\Delta \omega _\phi = (\omega _{\phi = 0} - \omega _{\phi = \pi })$ along the droplet interface for different aspect ratios, with $Re$ fixed at $0.1$ and viscosity ratio maintained at $\mu ^* = 100$ . The results confirm that a larger $\chi$ enhances the asymmetry of the vorticity at the top ( $\phi = 0$ ) and bottom ( $\phi = \pi$ ) of the droplet, thereby intensifying the Saffman mechanism. Additionally, we aim to develop a semi-empirical correlation to predict $C_L (Re, \mu ^*, \chi )$ in this low-but-finite Reynolds number regime for oblate droplets. Although no theoretical work is available for inviscid oblate bubbles or rigid spheroids in this context, we employ an interpolation technique inspired by the correlation established for spherical droplets in (4.3). The following formula is proposed based on this approach:

(4.4) \begin{equation} C_{L} (Re, \mu ^*, \chi ) = \frac {a\mu ^* + b^2}{(\mu ^* + b)^2}{C_{L} (Re, \mu ^*\rightarrow 0, \chi )} + \frac {{(\mu ^*)}^2 + a\mu ^*}{(\mu ^* + b)^2}{C_{L} (Re, \mu ^*\rightarrow \infty , \chi )}, \end{equation}

where $a$ and $b$ are coefficients that can be determined empirically. After fitting the numerical results, we find that the coefficients are $a = - 0.1954{\chi ^2} + 0.7833\chi + 0.2730$ and $b = - 0.1403{\chi ^2} + 0.5355\chi + 0.5483$ . The validity of this correlation is confirmed by the results presented in figure 7. In this figure, the solid lines corresponding to (4.4) align very well with the numerical data, demonstrating that the proposed correlation effectively captures the variation of lift force with both the viscosity ratio and aspect ratio for spherical to oblate droplets in the low-but-finite Reynolds number regime.

Figure 8.

Same plot as in figure 6, with the Reynolds number fixed at $Re = 0.1$ and the viscosity ratio fixed at $\mu ^* = 100$ , while the aspect ratio of the droplet is varied.

Additionally, it is important to note that in real-world scenarios, a deformable droplet elongates in the extensional direction of the background strain rate and rotates in the rotational direction of the background vorticity. In the low-Reynolds-number regime, the level of asymmetric deformation scales with the capillary number $Ca$ ( $Ca = \mu U_\infty / \sigma$ , with $\sigma$ the surface tension coefficient). This deformation-induced lift force scales as $Ca \mu R U_\infty$ (Magnaudet et al. Reference Magnaudet, Takagi and Legendre2003). In contrast, the inertial lift force in the presence of shear scales as $(Re \, Sr)^{1/2} \mu R U_\infty$ implied by (4.1). For the deformation-induced lift to be negligible compared with the inertial lift, the condition $Ca \ll (Re \, Sr)^{1/2}$ must hold. Since $Sr = 0.2$ is fixed for the majority of the present work and $Re$ values down to 0.1 are considered, the above constraint implies that $Ca \ll 0.1$ is required for the inertial lift to dominate in real-world scenario.

5.1. Spherical droplet

Figure 9.

Lift coefficients experienced by spherical droplets in the moderate- to high- $Re$ regime, with variations in Reynolds number and viscosity ratio: (a) $C_L$ versus $\mu ^*$ for different $Re$ ; (b) $C_L$ versus $\mu ^*$ with specified $Re = 100$ . In panel (b), the dash-dotted line represents the value predicted by correlation (5.1) (Legendre & Magnaudet Reference Legendre and Magnaudet1998), the dotted line indicates $C_L = 0$ , and the grey triangles and squares denote the pressure and viscous components of the lift coefficients.

Maintaining the droplet in a spherical shape ( $\chi = 1.0$ ), figure 9(a) illustrates the variation of lift coefficient $C_L$ with viscosity ratio $\mu ^*$ across different Reynolds numbers $Re$ . For highly viscous droplets, the lift coefficient may exhibit unsteady behaviour due to vortex shedding, therefore, the time-averaged $C_L$ is presented in the figure for those unsteady scenarios. It is observed that $C_L (Re, \mu ^*)$ decreases with increasing $\mu ^*$ and can even turn negative for high-Reynolds-number droplets ( $Re \geqslant 100$ ) at viscosity ratio exceeding a certain critical value. The same trend has been reported by Shi et al. (Reference Shi, Climent and Legendre2024) (figures 4 and 5 therein), where the critical viscosity ratio depends only weakly on the Reynolds number ratio. This behaviour suggests a fundamental shift in the lift generation mechanism between the moderate- to high- and low-but-finite Reynolds number regimes. Focusing on spherical droplets at $Re = 100$ , figure 9(b) shows the dependence of lift coefficients on $\mu ^*$ . A sign reversal in the lift coefficient is noted within the range $10 \lt \mu ^* \lt 50$ . The figure also displays the pressure and viscous components of the lift force. It is evident that the viscous component is negative across all viscosity ratios, while the positive pressure component decreases with $\mu ^*$ until it becomes smaller than the magnitude of the negative viscous component. This trend indicates that viscous effects dominate and lead to a reversal of the net lift force for viscous droplets in the moderate- to high-Reynolds-number regime. This negative contribution of the viscous component to lift has been reported previously by Legendre & Magnaudet (Reference Legendre and Magnaudet1998) and Kurose et al. (Reference Kurose, Misumi and Komori2001) for inviscid bubbles, as well as by Sugioka & Komori (Reference Sugioka and Komori2007) and Bagchi & Balachandar (Reference Bagchi and Balachandar2002) for water/rigid spheres. Additionally, Legendre & Magnaudet (Reference Legendre and Magnaudet1998) proposed an empirical correlation for predicting the lift coefficient of an inviscid spherical bubble in this moderate- to high-Reynolds-number regime, given by

(5.1) \begin{eqnarray} C_L (Re, \mu ^*\rightarrow 0, \chi = 1.0) = \frac {2Sr}{3}\frac {1 + 16 Re^{-1}}{1 + 29 Re^{-1}}. \end{eqnarray}

According to this correlation, the predicted lift coefficient for $Re = 100$ and $\mu ^* \rightarrow 0$ is approximately $0.12$ , as indicated by the dash-dotted line in figure 9(b). This value is in excellent agreement with our numerical results at $\mu ^* = 0.01$ . Conversely, we have not included reference data for the case of a rigid sphere ( $\mu ^* = 100$ ) due to the variability in available direct numerical simulation results for this Reynolds number (Shi & Rzehak Reference Shi and Rzehak2019). The negative lift coefficient observed for $Re = 100$ and $\mu ^* \rightarrow \infty$ is notably sensitive to the shear rate and the specific numerical settings employed in different simulations. The reported values in different numerical simulations and the corresponding numerical settings are summarised in Appendix B, where the information on the parameter choices, grid resolutions and boundary conditions is provided.

To investigate the mechanism underlying the lift reversal observed in figure 9 and the dependency of $C_L$ on $\mu ^*$ , we first consider the lift generation mechanisms in inviscid shear flow. In such a flow, the lift force arises from the advection and stretching of the basic vorticity, represented by $\omega _z = \alpha \boldsymbol {e}_z$ in the unperturbed shear flow. The asymmetric stretching of $\omega _z$ , characterised by $\omega _z \partial U_x / \partial z$ due to the presence of a spherical droplet, results in a pair of counter-rotating streamwise vortices in the wake (Lighthill Reference Lighthill1956; Auton Reference Auton1987; Auton et al. Reference Auton, Hunt and Prud’Homme1988). These vortices induce a positive lift force directed towards the side with the higher velocity, a mechanism known as the $L$ -mechanism (Adoua et al. Reference Adoua, Legendre and Magnaudet2009). This mechanism has been confirmed by Legendre & Magnaudet (Reference Legendre and Magnaudet1998) and Adoua et al. (Reference Adoua, Legendre and Magnaudet2009) for inviscid bubbles in the moderate- to high-Reynolds-number regime.

Additionally, it is important to consider another source of counter-rotating streamwise vortices, which arises from the asymmetric advection of the azimuthal vorticity $\omega _\phi$ generated at the body surface. In the presence of shear, the stretching and tilting of $\omega _\phi$ , represented by $\omega _\phi \partial U_x / \partial \phi$ , lead to the formation of counter-rotating streamwise vortices that produce a negative lift force directed towards the side with the lower velocity. This mechanism, referred to as the $S$ -mechanism (Adoua et al. Reference Adoua, Legendre and Magnaudet2009), becomes dominant for highly deformed bubbles where $\omega _{\phi , max}$ increases proportionally to $\chi ^{8/3}$ (Magnaudet & Mougin Reference Magnaudet and Mougin2007), as discussed by Adoua et al. (Reference Adoua, Legendre and Magnaudet2009) and Hidman et al. (Reference Hidman, Ström, Sasic and Sardina2022). Furthermore, increasing the viscosity of the droplet, represented by a higher $\mu ^*$ , may lead to the formation of attached eddies in the droplet’s wake due to flow separation. Rachih et al. (Reference Rachih, Legendre, Climent and Charton2020) demonstrate that such attached eddies form at Reynolds numbers of approximately $25$ for homogeneous flow past a droplet with $\mu ^* \gt 10$ . In the presence of shear, the inhomogeneous flow disrupts the axisymmetry of these attached eddies, causing them to tilt and generate counter-rotating streamwise vortices through the stretching and tilting effects in the vorticity equations. These streamwise vortices contribute to a negative lift force on the droplet and represent a specific case of the $S$ -mechanism. Thus, in the following study, we will also attribute such tilting eddy-induced lift to the $S$ -mechanism.

Consequently, the fundamental distinction between the $L$ -mechanism and the $S$ -mechanism lies in their origins of streamwise vorticity. In the $L$ -mechanism, the streamwise vorticity is generated from $\omega _z$ present in the upstream shear flow, while in the $S$ -mechanism, it originates from $\omega _\phi$ produced at the droplet’s interface. Referring to figure 9(b), it is evident that an increase in the viscosity ratio $\mu ^*$ leads to a significant rise in $\omega _\phi$ at the interface, thereby making a transition of dominance from the $L$ - to $S$ -mechanism. This transition results in a decrease in lift, potentially reversing it from positive to negative values. Additionally, in the scenario dominated by the tilting of the azimuthal vorticity, the lift force is influenced not only by the counter-rotating streamwise vortices (associated with the previously mentioned $S\hbox{-}$ mechanism) but also by the asymmetric advection of azimuthal vorticities. This latter effect is referred to as the ‘extended Saffman mechanism’ in the moderate- to high- $Re$ regime because it shares some similarities with the Saffman mechanism dominating in the low- $Re$ regime. However, it is important to emphasise that the Saffman and extended Saffman mechanisms are fundamentally different from a mathematical perspective. as we point out in § 4.1 that the Saffman mechanism is rooted in an asymptotic approach, which no longer applies in the moderate- to high-Reynolds-number regime. In this regime, the radius of the Stokes region becomes small, making it impossible to distinguish between inner and outer regions. Since these two mechanisms, arising from azimuthal vorticity, produce opposing effects on the lift – positive lift from the extended Saffman mechanism and negative lift from the $S$ -mechanism – they tend to counter balance each other, particularly for highly viscous droplets.

Figure 10.

Flow patterns past spherical droplets of $Re = 100$ with varying viscosity ratios: (a) $\mu ^* = 0.01$ ; (b) $\mu ^* = 1$ ; (c) $\mu ^* = 10$ ; (d) $\mu ^* = 100$ . For each viscosity ratio, the left panel displays the iso-surfaces of streamwise vorticity $\omega _x = \pm 0.3$ in the wake of the droplet, while the right panel illustrates the three-dimensional streamlines, with colour indicating $\omega _x$ values ranging from $-0.3$ (blue) to $0.3$ (red). In the right panels, attached eddies appear and are tilted at $\mu ^* \geqslant 10$ .

To investigate the transition from the $L$ -mechanism to the $S$ -mechanism as $\mu ^*$ increases in the moderate- to high-Reynolds-number regime, we examine the streamwise vortices $\omega _x = \pm 0.3$ in the wake of droplets, with the Reynolds number maintained at $Re = 100$ . The left panels of figure 10 illustrate these vortices. For weakly to moderately viscous droplets at $\mu ^* = 0.01$ and $1$ , where the $L$ -mechanism is predominant, the produced $\omega _x$ exhibits positive values at $z \lt 0$ and negative values at $z \gt 0$ . Consequently, the fluid situated between these vortices is entrained inwards ( $y \lt 0$ ), generating a positive counteracting lift on the droplet, consistent with the observations in figure 9(b). In contrast, for moderately to highly viscous droplets $\mu ^* = 10$ and $100$ , the signs of $\omega _x$ reverse compared with the cases with lower viscosity ratios, indicating a shift to the $S$ -mechanism. The right panels of figure 10 show that at these higher viscosity ratios, attached eddies form in the wake of the droplet and become tilted due to the ambient shear. This tilting generates streamwise vortices consistent with the $S$ -mechanism. Unlike in homogeneous flow, where the direction of the lift force can be influenced by random disturbances (Johnson & Patel Reference Johnson and Patel1999; Thompson et al. Reference Thompson, Leweke and Provansal2001), the shear in the present case acts as a controlled disturbance that consistently tilts the eddies in the $y$ -direction, thus producing a lift in the same direction. Additionally, Rachih et al. (Reference Rachih, Legendre, Climent and Charton2020) demonstrated that at $Re = 100$ , the formation of attached eddies behind the droplet occurs only when the viscosity ratio exceeds $\mu ^* \gt 5$ . This finding is in agreement with our observations presented in figure 10.

It is noteworthy that although the counter-rotating streamwise vortices at $\mu ^* = 10$ and $100$ share the same sign, the corresponding lift forces on the droplets exhibit opposite signs – specifically, a positive $C_L$ for $\mu ^* = 10$ and a negative $C_L$ for $\mu ^* = 100$ (see figure 9 b). This observation suggests that the reversal of lift force lags behind the reversal of streamwise vortices, indicating that the lift in this transitional regime is not solely determined by the sign of the streamwise vorticities. This scenario contrasts with the well-understood behaviour in homogeneous flow past blunt bodies, where the streamwise vortices in the wake predominantly dictate the direction of lift, whether for an oblate bubble (Magnaudet & Mougin Reference Magnaudet and Mougin2007; Zhang & Ni Reference Zhang and Ni2017) or a rigid sphere (Johnson & Patel Reference Johnson and Patel1999; Thompson et al. Reference Thompson, Leweke and Provansal2001). Hence, an additional mechanism must be at play to account for the positive lift observed at $\mu ^* = 10$ . This positive lift could counterbalance the negative lift induced by the $S$ -mechanism, resulting in a net positive lift. As conjectured earlier, this ‘positive’ lift may be attributed to the asymmetric advection of azimuthal vorticities, which is termed the extended Saffman mechanism and becomes significant at higher viscosity ratios. We will explore this hypothesis further in § 6, employing a force decomposition method grounded in the theory of vorticity moments to quantitatively analyse this effect.

Figure 11.

Variation of the lift coefficient $C_L$ with Reynolds number ranging from low to moderately high values for spherical droplets. The solid line represents the correlation $C_L (Re, Sr) = ([C_L^{\mathrm {low}\,Re} (Re, Sr)]^2 + [C_L^{\mathrm {high}\,Re} (Re, Sr)]^2)^{1/2}$ proposed by Legendre & Magnaudet (Reference Legendre and Magnaudet1998) for spherical bubbles. The dashed line denotes the correlation $C_L (Re, Sr) = 2.518 J_L (\varepsilon ) (\alpha /Re)^{1/2}$ suggested by McLaughlin (Reference McLaughlin1991) for rigid spheres at low-but-finite $Re$ , with $J_L (\varepsilon )$ being estimated by (4.2b ).

Next, we examine the influence of Reynolds number on the lift coefficient across the entire low- to moderately high- $Re$ regime. Figure 11 illustrates the dependency of $C_L$ on $Re$ . For clarity, the graph is presented on a logarithmic scale, with negative values of $C_L$ not shown but indicated by an inset arrow. The figure includes the correlation proposed by Legendre & Magnaudet (Reference Legendre and Magnaudet1998) (black solid line), given by $C_L (Re, Sr) = ([C_L^{\mathrm {low}\,Re} (Re, Sr)]^2 + [C_L^{\mathrm {high}\,Re} (Re, Sr)]^2)^{1/2}$ , where $C_L^{\mathrm {low}\,Re} (Re, Sr)$ is predicted using (4.1) and (4.2a ), and $C_L^{\mathrm {high}\,Re} (Re, Sr)$ is computed from (5.1). This correlation effectively predicts $C_L$ for inviscid bubbles and shows excellent agreement with the present results across the full parameter space, from the Saffman mechanism to the regime dominated by the $L$ -mechanism. Specifically, $C_L$ initially decreases with $Re$ until a critical value of $Re \approx 5$ , due to the weakening of the Saffman mechanism, and then increases as the $L$ -mechanism becomes more influential with increasing inertial effects. This trend is consistent with the numerical findings of Legendre & Magnaudet (Reference Legendre and Magnaudet1998) and Kurose et al. (Reference Kurose, Misumi and Komori2001), and is also applicable to weakly to moderately viscous droplets with $\mu ^* \leqslant 10$ . Additionally, the correlation proposed by McLaughlin (Reference McLaughlin1991) (red dashed line), given by $C_L (Re, Sr) = 2.518 J_L (\varepsilon ) (\alpha /Re)^{1/2}$ where $J_L (\varepsilon )$ is estimated by (4.2b ), is used to model the lift coefficient for a rigid sphere. This model aligns with our numerical data up to $Re \approx 5$ . For moderately to highly viscous droplets ( $\mu ^* \geqslant 10$ ), we observe a monotonic decline in $C_L$ with increasing $Re$ , attributed to the dominance of the $S$ -mechanism which reduces and can even reverse the lift. However, Rachih et al. (Reference Rachih, Legendre, Climent and Charton2020) indicates that attached eddies form in the wake of the droplet at $(Re, \mu ^*) \approx (50, 50)$ . Beyond this point, the shear tends to tilt the eddies, potentially reversing the streamwise vorticities and making negative lift more probable. Despite this, positive lift is still observed from our data displayed in figure 11, again suggesting the presence of an additional mechanism contributing extra positive lift.

Figure 12.

Evolution of lift coefficients with aspect ratio in the moderate- to high-Reynolds-number regime: (a) inviscid bubble with $\mu ^* = 0.01$ ; (b) droplet with $\mu ^* = 1$ ; (c) rigid spheroid with $\mu ^* = 100$ .

5.2. Oblate droplet

We now examine the behaviour of oblate droplets in the moderate- to high-Reynolds-number regime. Figure 12 presents the lift coefficient $C_L$ for oblate droplets with aspect ratios $1.0 \leqslant \chi \leqslant 2.5$ , with viscosity ratios $\mu ^* = 0.01$ , $1$ and $100$ depicted in separate panels. For cases where the lift coefficient is unsteady, such as for oblate bubble at $(Re, \mu ^*, \chi ) = (300, 0.01, 2.5)$ and for rigid sphere at $(Re, \mu ^*, \chi ) = (300, 100, 1.0)$ , we still use time-averaged values of $C_L$ . For inviscid bubbles ( $\mu ^* = 0.01$ ), figure 12(a) shows that at $Re = 5$ and $10$ , $C_L$ decreases monotonically with increasing $\chi$ . This trend is attributed to the fact that the $L$ -mechanism does not dominate at this moderate Reynolds number and the increase in oblateness strengthens the $S$ -mechanism, thereby reducing the positive value of the lift. As the Reynolds number increases to $Re \geqslant 50$ , $C_L$ initially increases with $\chi$ until it reaches a Reynolds-number-dependent threshold, $\chi _{cr}$ . For $Re = 50$ and $100$ , $\chi _{cr}$ is approximately $1.5$ , while for $Re = 300$ , $\chi _{cr}$ is approximately $1.8$ . Beyond this critical aspect ratio, $C_L$ begins to decrease with increasing $\chi$ . To understand this non-monotonic dependence of $C_L$ on $\chi$ , we must consider the interplay between the $S$ -mechanism and the $L$ -mechanism for oblate bubbles. According to Naciri (Reference Naciri1992), for inviscid shear flow around a body, such as a potential flow past an oblate bubble, $C_L (Re \rightarrow \infty , \mu ^* \rightarrow 0, \chi )$ closely mirrors the added mass coefficient and thus increases with aspect ratio. This theoretical prediction was numerically validated by Adoua et al. (Reference Adoua, Legendre and Magnaudet2009) at very high Reynolds numbers ( $Re \approx 4\times 10^{3}$ ) and is also supported by our numerical data presented in Appendix A. However, at moderate to high Reynolds numbers ( $Re \leqslant 15\times 10^{3}$ ), Adoua et al. (Reference Adoua, Legendre and Magnaudet2009) also observed that $C_L (Re, \mu ^* \rightarrow 0, \chi )$ does not always increase with $\chi$ , but can decline after reaching a critical threshold $\chi _{cr}$ , consistent with our numerical results. This behaviour arises because the enhanced $S$ -mechanism grows faster than the $L$ -mechanism when $\chi \gt \chi _{cr}$ . Figure 12(b) presents results for droplets with $\mu ^* = 1$ , where the variation in $C_L$ with $\chi$ follows a similar trend to that observed for $\mu ^* = 0.01$ . This suggests that the $L$ - and $S$ -mechanisms exhibit comparable dominating regimes for droplets with $\mu ^* = 0.01$ and $\mu ^* = 1$ . For droplets with higher viscosity ratios, corresponding to rigid spheroids ( $\mu ^* = 100$ ), figure 12(c) shows that $C_L$ consistently decreases with $\chi$ and may even become negative. This outcome indicates that the $S$ -mechanism predominates for such highly viscous droplets, leading to a reduction in $C_L$ as the droplet deformation increases.

We now explore the dependency of $C_L$ on $\mu ^*$ for droplets with different oblateness. Figure 13 presents the numerical results for different Reynolds numbers: $Re = 10$ , $100$ and $300$ . In figure 13(a), at $Re = 10$ , it is evident that the lift coefficient $C_L$ decreases both with $\chi$ and with $\mu ^*$ . This behaviour is attributed to the enhanced $S$ -mechanism associated with increasing $\mu ^*$ , as previously discussed. At higher Reynolds numbers, specifically $Re = 100$ and $300$ , figures 13(b) and 13(c) show that $C_L$ consistently declines with $\mu ^*$ for a given aspect ratio, similar to observations for spherical droplets. However, the relationship between $C_L$ and aspect ratio $\chi$ varies significantly depending on the viscosity ratio $\mu ^*$ . For weakly to moderately viscous droplets ( $0.01 \leqslant \mu ^* \leqslant 1$ ), the variation of $C_L$ with $\chi$ is non-monotonic. In contrast, for moderately to highly viscous droplets ( $1 \leqslant \mu ^* \leqslant 100$ ), $C_L$ decreases monotonically with $\chi$ . The differing impacts of $\chi$ in these two viscosity regimes can be explained by the roles of the $L$ - and $S$ -mechanisms. While the $S$ -mechanism not only strengthens with increasing $\mu ^*$ but also intensifies with $\chi$ , the $L$ -mechanism grows with $\chi$ up to a certain threshold for weakly to moderately viscous droplets, but always decreases with $\chi$ for moderately to highly viscous droplets. Consequently, for droplets with $\mu ^* \leqslant 1$ , the effect of $\chi$ on $C_L$ is complex due to the concurrent influences of both $L$ - and $S$ -mechanisms. Conversely, for droplets with $\mu ^* \gt 1$ , where the $S$ -mechanism predominates, $C_L$ is expected to have a purely negative relationship with $\chi$ .

Figure 13.

Evolution of lift coefficients with viscosity ratio and aspect ratio in the moderate- to high-Reynolds-number regime: (a) $Re = 10$ ; (b) $Re = 100$ ; (c) $Re = 300$ .

Figure 14.

Counter-rotating streamwise vortices behind the droplet at different aspect ratios, with iso-surfaces corresponding to $\omega _x = \pm 0.3$ : (a) inviscid oblate bubbles with $(Re, \mu ^*) = (300, 0.01)$ ; (b) rigid spheroids with $(Re, \mu ^*) = (50, 100)$ .

As previously discussed in § 5.1, the $L$ - and $S$ -mechanisms generate opposing streamwise vortices, resulting in reversed lift forces. To further elucidate these vortex structures behind oblate droplets in different regimes, we examine the iso-surfaces of $\omega _x = \pm 0.3$ in the wake of oblate droplets. Figure 14(a) illustrates these iso-surfaces for bubbles with $(Re, \mu ^*) = (300, 0.01)$ across varying aspect ratios, $1 \leqslant \chi \leqslant 2.5$ . For weakly to moderately deformed bubbles ( $1 \leqslant \chi \leqslant 1.5$ ), which fall under the regime dominated by the $L$ -mechanism, a pair of counter-rotating streamwise vortices produces a positive lift. However, as the aspect ratio increases to $\chi = 2.0$ , a second pair of vortices with opposite signs begins to emerge on the inner side of the wake. This additional pair of vortices is associated with the $S$ -mechanism, specifically the tilting of attached eddies that develop beyond a critical aspect ratio of $\chi \approx 1.75$ (Blanco & Magnaudet Reference Blanco and Magnaudet1995; Magnaudet & Mougin Reference Magnaudet and Mougin2007). As $\chi$ increases to $2.5$ , the $S$ -mechanism intensifies, leading to negative lift that partially offsets the positive lift generated by the $L$ -mechanism. This explains the observed decrease in $C_L$ for $\chi \gt 1.8$ , as shown in figure 12(a). For rigid spheroids at $(Re, \mu ^*) = (50, 100)$ , with aspect ratios varying in the range of $1 \leqslant \chi \leqslant 2.5$ , figure 14(b) depicts the iso-surfaces of the streamwise vortices. In this case, $\omega _x$ does not change sign with increasing $\chi$ , rather, it is consistently produced by the tilting of the attached eddies. Notably, figure 12(c) reveals that negative $C_L$ is observed only for $\chi \geqslant 1.2$ , while $C_L$ remains positive for $\chi = 1.0$ . This discrepancy, again, indicates that the reversal of the lift force lags behind the reversal of the streamwise vortices. Thus, there must be an additional mechanism contributing to the positive lift that counteracts part of the negative lift induced by the $S$ -mechanism. We will address this further in § 6, where we demonstrate that the asymmetric distribution of azimuthal vorticities is responsible for this additional positive lift.

6.1. Lift dependency on the counter-rotating streamwise vortices

From the wake structures presented in figures 10 and 14, we observe that in the moderate- to high-Reynolds-number regime, the lift force on a droplet in shear flow is closely related to the signs of the counter-rotating streamwise vortices. These vortices can arise from either the $L\hbox{-}$ or $S\hbox{-}$ mechanism. However, the reversal of lift, induced by increasing $\mu ^*$ or $\chi$ , occurs with a delay relative to the sign reversal of these vortices. This discrepancy is attributed to an undetermined mechanism, which contributes a positive lift that partially offsets the negative lift generated by the $S\hbox{-}$ mechanism.

To confirm this, we conduct a quantitative analysis based on the theory of vorticity dynamics (Saffman Reference Saffman1995; Wu et al. Reference Wu, Ma and Zhou2007b ). This framework indicates that the total force acting on a body moving through a viscous flow arises from body-induced vorticity moments. This theoretical approach has been successfully applied to study solid–fluid interactions (Wu et al. Reference Wu, Pan and Lu2005, Reference Wu, Lu and Zhuang2007a ; Wang et al. Reference Wang, He and Liu2019; Tong et al. Reference Tong, Yang and Wang2021). More recently, Hidman et al. (Reference Hidman, Ström, Sasic and Sardina2022) employed this vorticity-based approach to evaluate the lift force on a freely rising bubble in a shear flow, finding strong agreement with numerical results from a force-balanced model based on bubble motion dynamics. In this section, we apply a rigorous force decomposition using vorticity moments to elucidate the precise origins of the lift experienced by the droplet.

Figure 15.

Decomposition of vorticity forces using (6.1), with droplets having fixed parameters of $(Re, \chi ) = (100, 1.0)$ . (a) Evolution of vorticity forces with respect to the viscosity ratio $\mu ^*$ . Here, $C_{L, z}^\omega$ is less significant for weakly to moderately viscous droplets, but becomes more important for moderately to highly viscous droplets. (b) Integration domain for computing $F_{L, \omega _z}$ using (6.1 a), and the iso-contours represent $-0.3 \leqslant \omega _z \leqslant 0.3$ , with $r_e = 4R$ indicating the outer radius of $\Omega _e$ used in the integration. (c) Integration domain for computing $F_{L, \omega _x}$ using (6.1 b), and the iso-surfaces correspond to $\omega _x = \pm 0.3$ , with $l_A = 14R$ being the distance from $A_{wake}$ to the droplet centre and $r_A = 4.4R$ representing the radius for area integration.

The lift forces $F_L$ parallel to $\boldsymbol {e}_y$ can be decomposed into two components: one arising from the $z$ -component azimuthal vorticity $\omega _z$ and the other from the $x$ -component streamwise vorticity $\omega _x$ . The detailed derivations of these vorticity forces, based on the theory of vorticity moments, are presented in Appendix C. Note that these derivations are based on the Cartesian coordinate system rather than the cylindrical system, as both the flow shear and the induced lift force are aligned with the $y$ -direction, making the use of the Cartesian coordinate system more convenient for these investigations. As depicted in figure 15(b), the approximate expressions for the lift forces due to $\omega _z$ and $\omega _x$ can be formulated as

(6.1) \begin{eqnarray} F_{L, z}^\omega &=& -\frac {1}{2} \left ( \rho _l U_x \int _{\Omega _i + \Omega _e} \omega _z {\rm d} \Omega + \int _{S_e} \omega _z x \boldsymbol {u} \cdot \boldsymbol {n} {\rm d}S \right ),\\ \nonumber F_{L, x}^\omega &=& \frac {1}{2} \rho _l U_x \int _{A_{wake}} \omega _x z {\rm d} A. \end{eqnarray}

Here, $\Omega _i$ represents the volume of the internal fluid occupied by the droplet, $\Omega _e$ is the volume of the external fluid surrounding the droplet and $S_e$ denotes the outer boundary of $\Omega _e$ . The term $\boldsymbol {u} \cdot \boldsymbol {n}$ refers to the normal velocity on $S_e$ . Additionally, $U_x$ is the imposed velocity in the $x$ -direction, and $A_{wake}$ is the cross-sectional area downstream of the droplet, where $\omega _x dA$ represents the flux of the streamwise vorticities through the cross-section. For $F_{L, z}^\omega$ , the combination of the volume and surface integrals ensures that the formulation remains conservative, making the specific choice of $\Omega _e$ less critical. In contrast, for $F_{L, x}^\omega$ , the selection of $A_{wake}$ downstream of the droplet is crucial. We have conducted a thorough study to ensure that the results are independent of the choice of $A_{wake}$ and the details of this study are presented in Appendix C. Note that the formulations in (6.1) are simplified for steady-flow conditions, so that scenarios involving vortex shedding are not addressed in this section. In addition, (6.1 b) is exactly the formulation proposed by Zenit & Magnaudet (Reference Zenit and Magnaudet2009) for lift force estimation on a zig-zagging bubble in experiments, while $\int \omega _x {\rm d} A$ denotes the strength of the circulation of the vortex threads.

The first set of numerical cases we investigate is $(Re, \chi ) = (100, 1.0)$ , with the viscosity ratio varying in the range of $0.01 \leqslant \mu ^* \leqslant 100$ , as examined in § 5.1 (see figure 9 b). We observe a continuous decrease in $C_L$ with increasing $\mu ^*$ , which eventually reverses at a viscosity ratio greater than $\mu ^* = 10$ . Notably, the counter-rotating streamwise vortices change signs before reaching $\mu ^* = 10$ (see figure 10). We compute the vorticity forces based on (6.1), with the results presented in figure 15(a). The lift coefficient $C_{L, num}$ , obtained from pressure and viscous integration (2.4), is shown by black squares. The components $C_{L, x}^\omega$ (blue circle) and $C_{L, z}^\omega$ (red circle) are computed using the vorticity moment method, with their combined sum, $(C_{L, x}^\omega + C_{L, z}^\omega )$ , indicated by the black circle. For $\mu ^* = 0.01$ , figure 15(b) shows the distribution of $-0.3 \leqslant \omega _z \leqslant 0.3$ on the $XOY$ plane around the droplet, where $r_e = 4R$ represents the outer radius of $\Omega _e$ used for integration. Figure 15(c) displays the iso-surfaces of $\omega _x = \pm 0.3$ in the droplet’s wake, with $l_A = 14R$ being the distance between $A_{wake}$ and the droplet centre, and $r_A = 4.4R$ as the corresponding radius for area integration. Appendix C confirms that the choices of $l_A$ and $r_A$ have already converged.

Figure 15(a) shows that, for any $\mu ^*$ , the sum of $C_{L, x}^\omega$ and $C_{L, z}^\omega$ , i.e. $(C_{L, x}^\omega + C_{L, z}^\omega )$ , aligns very closely with $C_{L, {num}}$ . This agreement validates the credibility of the force computation based on the vorticity moment in (6.1). It should still be noted that a slight discrepancy between the values of $C_{L,num}$ and $(C_{L, x}^{\omega} + C_{L, z}^{\omega} )$ is observed at high $\mu ^*$ , especially for $\mu ^* = 100$ . This discrepancy arises because the formulation in (6.1) is an approximation using the vorticity moment method, with further details discussed in Appendix C. Figure 15(a) also reveals that, for weakly to moderately viscous droplets ( $\mu ^* \leqslant 1$ ), $C_{L, x}^\omega$ is almost equal to $C_{L, {num}}$ , while $C_{L, z}^\omega$ is positive but negligible. This strongly suggests that, in this $\mu ^*$ regime, nearly all the lift originates from the streamwise vortices generated by the $L$ -mechanism. However, as $\mu ^*$ increases, the vorticity forces produced by $\omega _z$ grow sharply, corresponding to positive values, causing $C_{L, x}^\omega$ to deviate from $C_{L, {num}}$ . In particular, at $\mu ^* = 10$ , $C_{L, x}^\omega$ is a small negative value of $-5 \times 10^{-3}$ , compared with the positive lift of $C_{L, {num}} = 3 \times 10^{-2}$ . This deviation occurs because $C_{L, z}^\omega$ increases to $3.5 \times 10^{-2}$ , contributing more to the total lift. Consequently, the force decomposition explains why the streamwise vortices change sign at $\mu ^* = 10$ , yet the lift remains positive, addressing the question raised in the description of figure 10. This demonstrates that, in the moderate- to high- $Re$ regime for a moderately to highly viscous droplet, where the $S$ -mechanism becomes dominant, even if the $S$ -mechanism generates streamwise vortices that produce negative lift, the asymmetric distribution of $\omega _z$ in terms of the ‘extended Saffman mechanism’ may contribute more to the lift than $\omega _x$ , allowing the total lift to remain positive. As the viscosity ratio further increases to $\mu ^* = 50$ and $100$ , the streamwise vorticity-induced lift coefficient approaches $C_{L, x}^\omega = -5 \times 10^{-2}$ , while the total lift is $C_{L, {num}} = -5 \times 10^{-3}$ . This is because the positive lift generated by asymmetric $\omega _z$ also rises to $C_{L, z}^\omega = 4.5 \times 10^{-2}$ , offsetting most of the negative $C_{L, x}^\omega$ .

Furthermore, the contribution of $C_{L, z}^\omega$ becomes more pronounced at higher values of $\mu ^*$ , as the azimuthal vorticity generated around the droplet increases with $\mu ^*$ . This results in a more intense asymmetric convection of $\omega _z$ . This observation highlights a key difference from homogeneous flows past highly viscous droplets or rigid spheres, where the direction of the lateral lift force is predominantly determined by the counter-rotating streamwise vortices alone. In contrast, in shear flow, the streamwise vortices cannot solely dictate the magnitude or direction of the lift force because the azimuthal vortices also play a significant role. Additionally, figure 16 presents results for oblate droplets with $(Re, \chi ) = (100, 1.5)$ . A similar trend is observed in the variation of $C_L$ with respect to $\mu ^*$ . For weakly to moderately viscous droplets ( $\mu ^* \leqslant 1$ ), the lift coefficient $C_{L, num}$ is predominantly influenced by $C_{L, x}^\omega$ . However, for moderately to highly viscous droplets ( $\mu ^* \geqslant 1$ ), the influence of $C_{L, z}^\omega$ becomes significant. Figure 16(b) illustrates that the sign reversal of the streamwise vortices occurs at $\mu ^* = 6$ , leading to a negative $C_{L, x}^\omega$ . Despite this, the total lift coefficient $C_{L, num}$ remains positive at $\mu ^* = 6$ , as shown in figure 16(a), due to the positive contribution of $C_{L, z}^\omega$ .

Figure 16.

Decomposition of vorticity forces using (6.1), with droplets having fixed parameters of $(Re, \chi ) = (100, 1.5)$ . (a) Evolution of the vorticity forces with respect to the viscosity ratio $\mu ^*$ . (b) Iso-surfaces of $\omega _x = \pm 0.3$ for different viscosity ratios $\mu ^*$ . The results reveal that $C_{L, z}^\omega$ is less significant for weakly to moderately viscous droplets but becomes more important for moderately to highly viscous droplets. Additionally, focusing on the case of $\mu ^* = 6$ reveals that the sign reversal of the lift lags behind the sign reversal of $\omega _x$ .

Figure 17.

Evolution of the vorticity forces as a function of Reynolds number. (a) Results for spherical bubbles with parameters $(\mu ^*, \chi ) = (0.01, 1.0)$ . (b) Results for rigid spheres with parameters $(\mu ^, \chi ) = (100, 1.0)$ .

We further examine how $C_{L, x}^\omega$ and $C_{L, z}^\omega$ vary with the Reynolds number, keeping other parameters constant. This analysis includes both spherical bubbles and rigid spheres. For spherical bubbles, with parameters fixed at $(\mu ^*, \chi ) = (0.01, 1.0)$ , we vary the Reynolds number in the range of $50 \leqslant Re \leqslant 300$ and the results are presented in figure 17(a). It is evident that across all Reynolds numbers considered, $C_{L, x}^\omega$ aligns closely with $C_{L, num}$ , while $C_{L, z}^\omega$ remains positive and negligible, decreasing with $\mu ^*$ . This indicates that the streamwise vorticities induced by the $L\hbox{-}$ mechanism are the primary contributors to the lift force on the spherical bubble in this Reynolds number range. For rigid spheres, we fix the parameters at $(\mu ^*, \chi ) = (100, 1.0)$ and vary the Reynolds number in the range of $75 \leqslant Re \leqslant 200$ . The results, shown in figure 17(b), differ significantly from those of the inviscid bubbles. Specifically, $C_{L, x}^\omega$ remains negative due to the counter-rotating streamwise vortices originating from the $S\hbox{-}$ mechanism or the tilting eddies in the wake. However, $C_{L, num}$ is not necessarily negative, as $C_{L, z}^\omega$ consistently contributes positively, though it decreases with $\mu ^*$ . These findings suggest that, in the moderate- to high-Reynolds-number regime, the $L\hbox{-}$ mechanism induced streamwise vortices predominantly dictate the lift force on weakly to moderately viscous droplets. In contrast, for moderately to highly viscous droplets, the $S\hbox{-}$ mechanism does not solely determine the lift, but rather, the contribution from the extended Saffman mechanism remains significant.

Figure 18.

Evolution of the vorticity forces with respect to aspect ratio. (a) Results for bubbles with parameters $(Re, \mu ^*) = (300, 0.01)$ . (b) Results for rigid spheroids with parameters $(Re, \mu ^*) = (50, 100)$ .

We also examine the impact of the aspect ratio $\chi$ on the interplay between $C_{L, x}^\omega$ and $C_{L, z}^\omega$ in the moderate- to high-Reynolds-number regime. We begin by considering bubbles with parameters $(Re, \mu ^*) = (300, 0.01)$ , varying the aspect ratio in the range of $1.0 \leqslant \chi \leqslant 2.2$ . Figure 18(a) shows that $C_{L, x}^\omega$ closely matches $C_{L, num}$ for all aspect ratios studied. The contribution of $C_{L, z}^\omega$ becomes slightly noticeable only at larger aspect ratios ( $\chi \gt 2.0$ ), due to the azimuthal vorticities scaling with the aspect ratio as $\omega _{\phi , max } \propto \chi ^{8/3}$ . Thus, for oblate bubbles, the lift induced by streamwise vorticities, whether from the $L\hbox{-}$ or $S\hbox{-}$ mechanism, remains the primary source of the total lift. Conversely, for rigid spheroids with parameters $(Re, \mu ^*) = (50, 100)$ , figure 18(b) reveals that $C_{L, x}^\omega$ consistently deviates from $C_{L, num}$ . The positive values of $C_{L, z}^\omega$ are comparable to the negative values of $C_{L, x}^\omega$ , leading to a near-cancellation that results in small negative values for $C_{L, num}$ .

In summary, the decomposition of lift using the vorticity moment theory demonstrates that, in the moderate- to high-Reynolds-number regime where either the $L\hbox{-}$ or the $S\hbox{-}$ mechanism may dominate, counter-rotating streamwise vortices are not necessarily the sole or even the primary source of lift in the presence of shear. For droplets with moderate to high viscosity, the asymmetric convection of azimuthal vorticities, termed as the extended Saffman mechanism in the present study, plays a significant role in determining the lift, meaning that the reversal of lift is not synchronised with the reversal of $\omega _x$ . However, it remains challenging to quantify the different effects using the vorticity moment method. This is because, when using the momentum of the streamwise vorticity to obtain $F_{L, x}^\omega$ , we end up mixing two contributions, one from the $S\hbox{-}$ mechanism and other from the $L\hbox{-}$ mechanism. While we can evaluate which mechanism is dominating based on the sign of $F_{L, x}^\omega$ , we are unable to separate the contributions from the $S\hbox{-}$ and the $L\hbox{-}$ mechanisms individually. Similarly, the lift force induced by the azimuthal vorticities, which involves both the $S\hbox{-}$ and the extended Saffman mechnisms, cannot be accurately estimated using the current approach. Furthermore, in situations where axisymmetric standing eddies form in the wake of the droplet, the shear flow disrupts their axisymmetric structure, generating counter-rotating streamwise vortices. These tilted eddies also contribute to a negative $F_{L, x}^\omega$ , but this contribution cannot be clearly distinguished from the $S\hbox{-}$ or $L\hbox{-}$ mechanisms in the present analysis.

6.2. Influences of the shear rate on the lift

To this point, our analysis has been confined to a constant shear rate of $Sr = 0.2$ . However, previous studies on spherical bubbles (Legendre & Magnaudet Reference Legendre and Magnaudet1998) and rigid spheres (Kurose & Komori Reference Kurose and Komori1999) have demonstrated that the lift force experienced by a body is highly sensitive to variations in the shear rate. Consequently, we will extend our investigation to explore how the lift force evolves with different shear rates on various droplets. Additionally, we will continue to apply the force decomposition method (6.1), based on vorticity moments, to steady problems to elucidate the origin of the lift at these varying shear rates.

Figure 19.

Evolution of the vorticity forces with respect to shear rate. (a) Results for bubbles with parameters $(Re, \mu ^*, \chi ) = (300, 0.01, 1.0)$ . (b) Results for more viscous droplets with parameters $(Re, \mu ^*, \chi ) = (300, 10, 1.0)$ .

We first investigate the effects of shear rate on lift force by analysing bubbles with parameters $(Re, \mu ^*, \chi ) = (300, 0.01, 1.0)$ . The shear rate is varied from $0.01$ to $0.3$ , and the corresponding lift coefficient $C_L$ is depicted in figure 19(a). The data reveal a positive correlation between $C_{L, num}$ and the shear rate $Sr$ . Using (6.1), we decompose the lift into contributions from the streamwise vorticity $C_{L, x}^\omega$ and the azimuthal vorticity $C_{L, z}^\omega$ . It is evident that $C_{L, x}^\omega$ closely matches $C_{L, num}$ , still indicating that for bubbles, the enhancement in lift is primarily due to the $L\hbox{-}$ mechanism-induced streamwise vortices. For comparison, we consider a more viscous droplet with parameters $(Re, \mu ^*, \chi ) = (300, 10, 1.0)$ and examine how the lift coefficient changes with shear rate in figure 19(b). As $Sr$ increases, $C_{L, num}$ transitions from negative to positive values, with the transition threshold occurring at $Sr \approx 0.2$ . This shift suggests that, at higher shear rates, the contribution from azimuthal vorticities ( $C_{L, z}^\omega$ ) surpasses that from streamwise vorticities ( $C_{L, x}^\omega$ ). Although the negative contribution of $C_{L, x}^\omega$ increases with $Sr$ due to the intensified $S\hbox{-}$ mechanism and more pronounced tilted standing eddies, the positive contribution from $C_{L, z}^\omega$ increases even more sharply, leading to an overall positive lift force. This discrepancy between $C_{L, x}^\omega$ and $C_{L, num}$ underscores that for highly viscous droplets in shear flow, the counter-rotating streamwise vortices in the wake contribute only a portion of the total lift experienced by the droplet. In addition, we also investigate the influence of the shear rate on unsteady problems characterised by the vortex shedding in the wake of the droplet, and these results are presented and discussed in Appendix D.

In summary, increasing the shear rate ( $Sr$ ) consistently results in a larger positive lift force on the droplet. This dependency of the lift coefficient $C_L$ on $Sr$ helps explain the discrepancies observed in the Reynolds number ( $Re$ ) thresholds for lift reversal on a rigid sphere in different studies (Kurose & Komori Reference Kurose and Komori1999; Bagchi & Balachandar Reference Bagchi and Balachandar2002; Hölzer & Sommerfeld Reference Hölzer and Sommerfeld2009). The mechanisms underlying this lift enhancement vary between weakly to moderately viscous and moderately to highly viscous droplets. For weakly to moderately viscous droplets, a higher $Sr$ enhances the $L\hbox{-}$ mechanism, resulting in increased positive lift. Conversely, for moderately to highly viscous droplets, the lift enhancement is largely attributed to the intensified asymmetric convection of azimuthal vorticities, which is the extended Saffman mechanism that contributes significantly to the positive lift.