2.1. Shock tube platform and imaging techniques

In experiments, we employ the shock tube to generate a plane shock wave with a uniform airflow behind it. As the air stream impacts perpendicularly on the cylinder, the axisymmetric stagnation point flow can be constructed on the central position of the solid wall. Figure 1(a) shows the schematic diagram of the shock tube experiments. To ensure the symmetry of the wall-attached droplet, the wall surface and the droplet are positioned horizontally, and thus require the vertical motion of the shock wave. The shock wave with constant Mach number ( $Ma$ ) is generated through a vertical-positioned shock tube, followed by the post-wave air stream with a uniform velocity of $U_0$ . The radii of the droplet and the cylinder are denoted by $a$ and $R$ , respectively. The photograph of the experimental area inside the vertical shock tube is depicted in figure 1(b), showing the cylinder and the wall-attached droplet. The cylinder is placed inside the shock tube, which has a cross-sectional area of 100 mm $\times$ 28 mm and is equipped with glass windows on all sides to facilitate the observation of droplet evolutions. The radius of the cylinder ranges from $R = 1.5$ mm to 6 mm, which is very small relative to the cross-section of the shock tube. Therefore, the stagnation point flow can form close to the droplet as the air stream comes perpendicularly.

Figure 1.

(a) Schematic diagram of shock tube experiments, where the high-speed stagnation point flow is generated as the shock wave impacts vertically on the cylinder. The shock wave Mach number is denoted by $Ma$ , the velocity of the post-wave air stream is denoted by $U_0$ and the radii of the droplet and the cylinder are denoted by $a$ and $R$ . (b) Photograph of the experimental area inside the vertical-positioned shock tube.

In experiments, the shock wave and the post-wave airflow are generated by rupturing a membrane by compressed air inside the high-pressure section of the shock tube, which is located upstream of the experimental area. Through varying the thickness of membrane, we can control both the Mach number $Ma$ of the shock wave and the velocity of airflow $U_0$ . The relationship between $Ma$ and $U_0$ is

(2.1) \begin{equation} U_0=\frac {2c}{\gamma +1}(Ma-\frac {1}{Ma}), \end{equation}

where $\gamma$ is the gas adiabatic coefficient and can be chosen as $\gamma = 1.4$ for air, and $c$ represents the ambient speed of sound ahead of the shock wave. The shock wave Mach number in our experiments ranges from $Ma$ = 1.05 to 1.35, corresponding to the post-wave air stream velocity from $U_0$ = 27 m s⁻¹ to 174 m s⁻¹ and the post-wave air density from $\rho _g$ = 1.31 kg m $^-{^3}$ to 1.93 kg m $^-{^3}$ , respectively.

In our study, it is crucial to capture the high-speed surface evolution of millimetre-scaled droplets. To achieve high temporal and spatial image resolutions and also fine visualisation of droplet surface structures, the high-speed photography method is used to capture the aerodynamic deformation and breakup processes. The backlighting method is adopted to light up the experimental region, where the light source and high-speed camera are aligned on different sides of the droplet. This illuminating arrangement enables us to capture the droplet profile and the perturbation waves on the surface. The imaging set-up consists of a high-speed camera (Phantom V2012) with a macro lens (Canon MP-E 65 mm f/2.8 1–5X), realising a maximum temporal resolution of $1.52 \; \unicode {x03BC}\rm s$ and a maximum spatial resolution of $5.4\;\unicode {x03BC}\rm m$ per pixel.

As the high-speed photography cannot resolve the shock wave directly, high-speed schlieren imaging is also employed in our experiments to investigate the dynamics of the shock wave. Schlieren imaging enables the observation of density variations in air flow by detecting the displacement of light passing through fluids with different densities, which realises the visualisation of the shock wave and the direct calculation of shock wave Mach number.

To capture the flow characteristics of the post-wave air stream, particle image velocimetry (PIV) is also employed. This technique can quantitatively measure the flow field information by capturing the motion of tracing particles, and thus allows the validation of a stagnation point airflow close to the solid wall. In our experiments, tracing particles are dispersed into the experimental region of the shock tube from the top side. Once the particles fill the entire experimental region, compressed air is injected into the high-pressure section of the shock tube to rupture the membrane and generate the shock wave. Simultaneously, the motion of particles in the post-wave air flow is captured by high-speed photography. The tracing particles used in the experiment are a mixture of glycerin, triethylene glycol and water, with diameters of approximately 1–3 $ \;\unicode {x03BC} \rm m$ and a density of approximately 1.1 gmL⁻¹. The relaxation time of the particles approximates $t \sim 1 \;\unicode {x03BC} \rm s$ , which is much smaller than the characteristic time scale of the air flow across the experimental region. Therefore, these tracing particles can meet the testing requirements of high-speed air flows.

Using the PIV method, we first give a quantitative measurement of the axisymmetric stagnation point flow at the solid wall. The results are shown in figure 2(a), where the flow condition corresponds to $Ma$ = 1.07 and $U_0$ = 38 m s⁻¹. The background of the image exhibits the tracing particles, and the streamlines and velocity vectors obtained through cross-correlation arithmetic are given simultaneously. The motion of tracing particles is also given in supplementary movie 1 available at https://doi.org/10.1017/jfm.2025.085. Theoretically, for an axisymmetric stagnation flow field with a characteristic scale of $R$ (radius of the cylinder), the axial and radial dimensionless velocity components correspond to $\tilde {u}_z=-A\tilde {z}$ and $\tilde {u}_x=A\tilde {x}/2$ (Homann Reference Homann1936), where $\tilde {u}_z$ and $\tilde {u}_x$ are non-dimensionalised by $U_0$ , $\tilde {z}$ and $\tilde {x}$ are non-dimensionalised by $R$ , and $A$ is a constant. Figure 2(b) shows the experimental measurements of the velocity components along the $\tilde {z}$ and $\tilde {x}$ directions. It is observed that the values of velocity show a linear relationship with the coordinates, and the velocity in the $\tilde {z}$ direction is twice the velocity in the $\tilde {x}$ direction, indicating a good agreement between the theoretical predictions and the PIV measurements. Moreover, the coefficient $A$ can be decided by the slope of $\tilde {u}_z$ with $\tilde {z}$ in figure 2(b), which approximates to the constant value of 1.1.

Figure 2.

(a) PIV measurements of the axisymmetric stagnation point airflow field close to the solid wall under $Ma$ = 1.07 and $U_0$ = 38 m s⁻¹, where the velocity vectors and the streamlines are annotated simultaneously. (b) Variation of the dimensionless axial and radial velocities (denoted by $\tilde {u}_z$ and $\tilde {u}_x$ , respectively) with coordinates, which satisfies the linear relationship of $\tilde {u}_z=-A\tilde {z}$ and $\tilde {u}_x=A\tilde {x}/2$ . The coefficient $A$ is measured to be 1.1.

In our experiments, we employ water as the material of the droplet, with a density of $\rho _l = 997\,\rm kg\,m^-{^3}$ and a dynamic viscosity of $\mu _l = 0.001$ Pa $\cdot$ s. The surface tension coefficient between the droplet and the air is measured to be $\sigma = 72.1\,\rm mN\, m^-{^1}$ . All the droplets in experiments are in hemispherical shape with an apparent contact angle of 90 degrees on the solid wall. The aluminium alloy is chosen as the material of the cylinder, which ensures the desired contact angle after polishing and treatment with a low concentration of hydrophobic reagent. It is notable that as the contact line begins to spread on the solid wall under the continuous extrusion of vertical air flow, the apparent contact angle between the droplet and the wall must exceed the advancing contact angle (Sui et al. Reference Sui, Ding and Spelt2014), which is measured to be 120 degrees for the solid wall in our experiments. A pipette is used to control the droplet volume precisely, which varies from $2.5 \; \unicode {x03BC}\rm L$ to $19\;\unicode {x03BC}\rm L$ . Therefore, the radius of the droplets ranges from $a = 1.06$ mm to 2.09 mm. For the convenience of analysis, we define the radius ratio between the water droplet and the cylinder as $\eta = a/R$ . Since $R$ ranges from 1.5 mm to 6 mm, the value of $\eta$ falls within the range of 0.18–0.66 in our experiments.

2.2. Dimensionless parameters

During the droplet deformation and breakup in the air stream, the inertial force, the viscous force and the surface tension force act jointly. Choosing the droplet radius $a$ and the airflow velocity $U_0$ as the characteristic length and velocity, the Weber number ( $We$ ) and the Reynolds number ( $Re$ ) are defined to measure the relative importance between these forces, i.e.

(2.2) \begin{equation} We=\frac {\rho _gU_0^2a}{\sigma },\quad Re=\frac {\rho _gU_0a}{\mu _g}. \end{equation}

In our experiments, the Weber number ranges from $We=9$ to 1180. We also adopt a constant gas dynamic viscosity as the airflow velocity changes, which equals to $\mu _g=1.8 \times 10^{-5}$ Pa $\cdot$ s. Therefore, the Reynolds number ranges from $Re=1348$ to 39 500. Moreover, we define the Ohnesorge number ( $Oh$ ) to evaluate the viscous effect of the droplet, i.e.

(2.3) \begin{equation} Oh=\frac {\mu _l}{\sqrt {\rho _l\sigma a}}. \end{equation}

The value of $Oh$ ranges from 0.0026 to 0.0036 in experiments. Previous studies have indicated that the viscous force of a droplet is supposed to play an insignificant role below a critical value of $Oh$ ( $\approx 0.1$ ) (Guildenbecher et al. Reference Guildenbecher, López-Rivera and Sojka2009; Theofanous Reference Theofanous2011), and thus reveals that the viscous effect of the droplet can be ignored in our study as $Oh \sim O(10^{-3})$ .

It should be noted that the definition of $We$ and $Re$ cannot precisely represent the dynamics of the wall-attached droplet as the characteristic velocity of the local stagnation flow is different from $U_0$ . Since we have obtained the local flow fields of the droplet through PIV measurements, the dimensionless parameters that represent the aerodynamic deformation and breakup processes of the droplet can be obtained. In an axisymmetric stagnation point airflow with a characteristic size of $R$ , the characteristic gas velocity on the upper side of a droplet with a radius of $a$ could be represented by the local velocity, with the specific form of $U_{lo}=A \eta U_0$ . Consequently, the ratio of the local gas aerodynamic force to the droplet surface tension force, referred to as the local Weber number ( $We_{lo}$ ), can be expressed as

(2.4) \begin{equation} We_{lo}=\frac {\rho _gU_{lo}^{2}a}{\sigma }=A^2\eta ^2We, \end{equation}

where $A$ equals to the specific value of 1.1 according to the measurements in figure 2(b). Similarly, we can define the local Reynolds number, representing the ratio of local gas aerodynamic force to viscous force, i.e.

(2.5) \begin{equation} Re_{lo}=\frac {\rho _gU_{lo}a}{\mu _g}=A\eta Re. \end{equation}

In our experiments, $We_{lo}$ ranges from 0.35 to 622, indicating that both the surface tension force and the aerodynamic force can exert a predominant effect on droplet deformation according to the specific value of $We_{lo}$ . The value of $Re_{lo}$ ranges from 485 to 28 677 in experiments, indicating a weak effect played by the gaseous viscous force. Specifically, in axisymmetric stagnation point flow, the thickness $\delta$ of the local boundary layer on the solid wall can be calculated as $\delta =2.8aRe_{lo}^{-0.5}$ (Schlichting & Gersten Reference Schlichting and Gersten2016). Considering the magnitude of $Re_{lo}$ , the thickness of the wall boundary layer is at least one order smaller than the droplet radius. Therefore, the viscous effect of the airflow is ignored in our study. Based on the analysis above, $We_{lo}$ acts as the principal dimensionless parameter on droplet deformation, and the deformation modes and the droplet dynamics with the variation of $We_{lo}$ will be analysed in detail in later sections.

As the droplet deformation and breakup process are mainly characterised by the surface acceleration in a stagnation point airflow, the characteristic time scale is selected based on the acceleration of wall-attached droplets. Starting from the aerodynamic forces $F_d \sim {\rho _g}U_{lo}^{2}a^2$ and the droplet mass ${m_l} \sim {\rho _l}{a^3}$ , the characteristic acceleration of a droplet can be represented as $F_d/m_l$ . Therefore, the characteristic time scale of droplet deformation can be expressed as

(2.6) \begin{equation} T\sim \sqrt {\frac {a}{F_d/m_l}} \sim \dfrac {a}{U_{lo}}\sqrt {\dfrac {\rho _l}{\rho _g}}. \end{equation}

It should be noted that this characteristic time scale has been widely used in relevant studies considering the deformation of a free spherical droplet in uniform airflow (Theofanous Reference Theofanous2011; Xu et al. Reference Xu, Wang and Che2022), except that the uniform airflow velocity $U_0$ is replaced by the velocity $U_{lo}$ of the local stagnation airflow.

3.1. Interactions between shock wave and droplets

In the aerobreakup process of a wall-attached droplet induced by a shock wave, the droplet interacts with the shock wave first, and then deforms and breaks up in the post-wave air stream. In this section, we evaluate the effect of the shock wave on droplet evolution through the high-speed schlieren method, and the results are shown in figure 3(a). The droplet deformation process under approximate moments are also captured through high-speed photography, as shown in figure 3(b). In the experiments, the zero time instant ( $t=0 \;\unicode {x03BC}\rm s$ ) is chosen as the moment when the shock wave impacts the front edge of the droplet, as revealed by the high-speed schlieren image. The Mach number of the shock wave can be calculated by the displacement distance of the shock wave in two adjacent images, which equals to $Ma=1.2$ . As the incident shock wave (IS) interacts with the droplet, the reflected shock wave induced by the droplet upper surface (denoted by RS1) can be observed (see the panels at $t=1.52 \;\unicode {x03BC}\rm s$ and $3.04\;\unicode {x03BC}\rm s$ ). With the evolution of IS, the reflected shock wave induced by the solid wall (denoted by RS2) can also be observed ( $t=4.56 \;\unicode {x03BC}\rm s$ and $7.6\;\unicode {x03BC}\rm s$ ), resulting in the formation of complex wave structures close to the droplet and the solid wall. As the strength of RS2 approximates that of IS due to solid reflection, the airflow close to the droplet decelerates significantly and remains almost stagnant at the rear of RS2. However, since the cross-sectional area of the cylinder is much smaller than that of the shock tube, the divergent reflected shock waves dissipate rapidly ( $t=19.76 \;\unicode {x03BC}\rm s$ ), leading to the establishment of a post-shock air stream and stagnation point flow close to the wall-attached droplet. At $t=106.4\;\unicode {x03BC}\rm s$ , the droplet begins to exhibit visible deformation on the surface. Under the continuous extrusion of the air stream, the droplet deforms significantly and breaks up finally at $t=866.8\;\unicode {x03BC}\rm s$ . More details of droplet surface perturbations can be seen according to the images captured by high-speed photography (see figure 3 b).

Figure 3.

(a) Aerodynamic deformation and breakup process of a hemispherical droplet captured by high-speed schlieren, in which IS, RS1 and RS2 represent the incident shock wave, the reflected shock wave by the droplet upper surface and the reflected shock wave by the solid wall, respectively. (b) Aerodynamic deformation and breakup process of a hemispherical droplet captured by high-speed photography. The experimental parameters in panels (a) and (b) correspond to $Ma=1.2$ , $U_0=105\;\textrm m\,\rm s^-{^1}$ , $R= 4\;\textrm {mm}$ and $a=1.42\;\textrm {mm}$ .

According to the results in figure 3, the duration of shock–droplet interaction $(\sim 10 \;\unicode {x03BC} \rm s)$ is much shorter than the entire deformation time scale of the droplet $(\sim 1\;\rm ms)$ , indicating a secondary effect of the shock wave impacting compared with the extrusion effect of the post-wave air stream. Moreover, due to the large acoustic impedance of the liquid phase relative to the gas phase, the Richtmyer–Meshkov instability can be neglected. Consequently, our subsequent analyses mainly focus on the influence of high-speed airflow, but disregards the direct impact of shock wave, which is in accordance with the prior studies (Guildenbecher et al. Reference Guildenbecher, López-Rivera and Sojka2009; Theofanous Reference Theofanous2011).

3.2. Morphologies of droplets

Figure 4.

Side view of the aerodynamic deformation and breakup process of hemispherical droplets under (a) vibrating mode (VB), (b) compressing mode (CP), (c) sheet thinning mode (ST) and (d) shear-induced entrainment mode (SIE). The parameters are: (a) $We_{lo}=1.44$ $(\eta =0.21,We=27$ , corresponding to $U_0=38$ m s⁻¹, $R= 5$ mm and $a=1.05$ mm); (b) $We_{lo}=25.1$ $(\eta =0.36, We=160$ , corresponding to $U_0=75$ m s⁻¹, $R= 4$ mm and $a=1.42$ mm); (c) $We_{lo}=43.3$ $(\eta =0.36, We=276$ , corresponding to $U_0=105$ m s⁻¹, $R= 4$ mm and $a=1.42$ mm) and (d) $We_{lo}=274$ $(\eta =0.58, We=674$ , corresponding to $U_0=151$ m s⁻¹, $R= 2.5$ mm and $a=1.45$ mm).

In this section, the morphologies of hemispherical droplets during the deformation and aerobreakup process in a post-wave air stream are studied, in which high-speed photography is used to ensure the high spatial resolution of experimental images. It has been indicated in § 2.2 that the dynamics of droplet evolution is mainly decided by the local Weber number $We_{lo}$ . With the increase of $We_{lo}$ , four typical deformation and breakup modes of hemispherical droplets can be observed in the axisymmetric stagnation airflow, which include the vibrating (VB) mode, the compressing (CP) mode, the sheet thinning (ST) mode and the shear-induced entrainment (SIE) mode. The evolution characteristics of these modes from the side view are shown in figures 4(a)–4(d), respectively. The corresponding video of different modes are given as supplementary movie 2. For the convenience of analysis, we reselect the dimensionless zero time instant as the moment when the droplet surface presents visible deformation. It should be noted that according to (2.4), a specific value of $We_{lo}$ can be realised by modulating the air stream velocity $U_0$ (non-dimensionalised as $We$ ) and the relative size of droplet and cylinder (non-dimensionalised as $\eta =a/R$ ). Therefore, the corresponding values of $We$ , $\eta$ , $U_0$ , $a$ and $R$ are also given in the caption of figure 4.

In the VB mode, as shown in figure 4(a), the dynamical behaviour of the droplet is mainly characterised by the periodic vibration with finite surface deformations. At the early stage, the air stream compresses the upper surface of droplet, causing a small deformation of the droplet (see the panels from $t^*=$ 0–0.43). After reaching the maximal deformation, the upper surface begins to rebound due to surface tension, eventually returning to an approximately spherical shape at $t^* = 0.74$ . Subsequently, the droplet continues this vibrating behaviour periodically.

The dynamics of the CP mode (as shown in figure 4 b) exhibits significant differences from the VB mode, as the droplet undergoes substantial deformations in the entire process. The deformation features at the early stage mainly include the compression of the upper surface and the protrusion of the lateral surface, which results in the ‘ripple’ structure, as indicated by the sub-graph at $t^*=0.28$ . The ‘ripple’ evolves into a protruding ring-shaped structure (named as ‘rim’, see the sub-graph at $t^*=0.43$ ) between the droplet upper surface and lateral side. As the droplet continues to be compressed, the upper surface exhibits significant spreading characteristics in the horizontal direction. Due to the effect of contact angle hysteresis (Eral et al. Reference Eral, ’t Mannetje and Oh2013; Sui et al. Reference Sui, Ding and Spelt2014), the contact line between the droplet and the solid wall maintains static, which leads to the transformation of the trapezoidal droplet from a narrow top and wide bottom profile ( $t^*=0.43$ ) to a wide top and narrow bottom profile ( $t^*=0.74$ ). The droplet continues to be compressed until the contact angle between the droplet and the wall exceeds the advancing contact angle ( $t^*=0.94$ ). From there, the contact line spreads rapidly and the droplet eventually evolves into a liquid film on the solid wall.

The aerodynamic deformation and breakup processes of the droplet in the ST mode are shown in figure 4(c). Similar to the CP mode, the initial deformation features of the droplet mainly include the compression of the upper surface, and the formation of ‘ripple’ and ‘rim’ structures. However, different from the CP mode, the ‘rim’ structure in the ST mode does not maintain its initial profile but continues to extend outward and evolves into a liquid sheet (LS) structure (see the sub-graphs for $t^*=0.74$ and $t^*=0.94$ ). The LS then strips from the droplet main body and gradually thins, eventually detaching from the droplet main body and impacts the solid wall ( $t^*=1$ ). The contact line also begins to spread until the very later stage of droplet deformation.

In the SIE mode, as depicted in figure 4(d), the primary characteristics of the droplet are the vigorous peeling and breakup of liquid films and filaments on the droplet surface, resulting in the generation of a liquid mist. Although there still exists the ‘ripple’ structure on the upper surface in the early stages ( $ t^*=0.28$ ), the Kelvin–Helmholtz instability (KHI) develops rapidly at the droplet surface due to strong shear stress of the air flow ( $t^* = 0.43$ ). As the perturbation waves grow with time, the liquid mist peels from the droplet main body continuously, accompanied by the compression of the droplet upper surface ( $t^* = 0.74-1$ ).

Figure 5.

Oblique front view of the aerodynamic deformation and breakup process of hemispherical droplets under (a) vibrating mode (VB), (b) compressing mode (CP), (c) sheet thinning mode (ST) and (d) shear-induced entrainment mode (SIE). The parameters are: (a) $We_{lo}=1.44$ $(\eta =0.21,We=27$ , corresponding to $U_0=38$ m s⁻¹, $R= 5$ mm and $a=1.05$ mm); (b) $We_{lo}=7.8$ $(\eta =0.28,We=82$ , corresponding to $U_0=75$ m s⁻¹, $R= 5$ mm and $a=1.42$ mm); (c) $We_{lo}=36.9$ $(\eta =0.32,We=298$ , corresponding to $U_0=105$ m s⁻¹, $R= 5$ mm and $a=1.61$ mm) and (d) $We_{lo}=247$ $(\eta =0.37,We=1494$ , corresponding to $U_0=175$ m s⁻¹, $R= 5$ mm and $a=1.85$ mm).

To demonstrate the three-dimensional evolutive characteristics of wall-attached droplets, we also capture the aerodynamic deformation and breakup processes of droplets from the oblique front view, and the results are shown in figure 5. The corresponding video of different modes are given as supplementary movie 3. For the VB and CP modes (see figures 5 a and 5 b, respectively), it is observed that the droplet deformation remains axisymmetric during the entire evolution process. However, for the SIE mode shown in figure 5(d), the occurrence of KHI and the subsequent peeling of the secondary droplets can lead to significant three-dimensional evolution of the droplet surface. As for the ST mode shown in figure 5(c), although the stretching of the liquid sheet at a very later stage of deformation results in a slightly circumferential instability of the surface (see the panel at $t^*=1.16$ ), the droplet is able to maintain an overall axisymmetric profile during its deformation process. Therefore, for the CP and ST modes, we will build the axisymmetric theoretical models to predict the droplet deformation dynamics. The details will be give in a later section of § 4.3.

Through abundant experimental cases considering the variation of $We_{lo}$ , the parametric regions for VB, CP, ST and SIE modes can be determined, as shown in table 1. The transitions between VB/CP, CP/ST and ST/SIE modes are observed to occur close to the critical $We_{lo}$ of 2, 26 and 150, respectively. We emphasise that although a specific value of $We_{lo}$ can be realised by different combinations of $We$ and $\eta$ , as shown by (2.4), the individual values of $We$ or $\eta$ are found not to determine the boundary of modes transition. The detailed analysis on the variation of $We$ and $\eta$ can be seen in Appendix A, where the phase diagrams of different modes are given under the $We-1/\eta ^2$ , $We_{lo}-We$ and $We_{lo}-\eta$ planes. In our experiments, although the values of $We$ and $\eta$ vary within certain ranges, the critical values of $We_{lo}$ for modes transition hardly change, confirming that $We_{lo}$ is the key parameter on determining the droplet’s deformation and breakup mode.

Table 1.

Parametric regions of $We_{lo}$ for VB mode, CP mode, ST mode and SIE mode.

In the VB mode, when the air flow compresses the droplet upper surface, the surface energy increases due to droplet deformation. As the surface tension contributes to maintaining the approximately hemispherical profile of the droplet, the aerodynamic force is insufficient to continuously overcome the increasing surface tension force as the deformation continues. Therefore, the droplet exhibits a vibrating behaviour due to the competition between these two forces, corresponding to the typical range of $We_{lo} \sim O(1)$ . If the local aerodynamic force is large enough to dominate the surface tension, invariably, the upper surface of the droplet sustains the compressing behaviour, resulting in the CP mode. The criterion of modes transition between the CP mode and the ST mode approximates the critical value of $We_{lo} \sim O(10)$ , which indicates that the dominance of aerodynamic force relative to surface tension can result in different dynamic behaviours of the ‘rim’ structure developed on the droplet surface. When $We_{lo}$ increases to the values of $O(10^2)$ , the effect of surface tension is significantly weak compared with the aerodynamic force of the air stream, and the perturbation waves on the droplet surface tend to arise due to the strong shear effect of the airflow. Consequently, the Kelvin–Helmholtz instability dominates the breakup process and the droplet enters the SIE mode. It should be noted that in previous studies considering the free spherical droplets in a uniform airflow, the critical Weber number for the dominance of KHI also ranges of the order of $10^2$ (Hsiang & Faeth Reference Hsiang and Faeth1995; Guildenbecher et al. Reference Guildenbecher, López-Rivera and Sojka2009; Theofanous Reference Theofanous2011), which is consistent with the critical $We_{lo}$ obtained in our experiments.

Overall, due to the existence of stagnation point flow, the CP and ST modes exhibit distinctive deformation and breakup characteristics compared with the free spherical droplets in a uniform flow. This is mainly because the presence of the solid wall can suppress the development of Rayleigh–Taylor instability on the droplet upper surface and also prevent the formation of bubble structures. However, the VB mode and the SIE mode share similar evolutive characteristics with the free spherical droplets in an air stream. These similarities arise because the droplet deformation mechanisms under VB and SIE modes are not strongly related to the stagnation point flow close to the wall. Therefore, the critical values of $We_{lo}$ for the occurrence of VB and SIE modes almost coincide with those for free droplets in a uniform flow (Guildenbecher et al. Reference Guildenbecher, López-Rivera and Sojka2009; Theofanous & Li Reference Theofanous and Li2008). Considering the aforementioned reasons, we will analyse the aerodynamics of CP and ST modes as well as their transition mechanisms in detail in § 4, and then discuss dynamical behaviours of the VB mode and the SIE mode in § 5.

4.1. Temporal evolutions of droplets

At the early stage of evolution, the main deformation characteristics of the CP mode and the ST mode include the compression of the droplet upper surface and the occurrence of ‘ripple’ and ‘rim’ structures on the lateral side of the droplet. At the later stage of deformation, the droplets under these two modes both sustain a remarkable compression in the vertical direction and obvious spreading of the upper surface in the horizontal direction. The main difference between the CP mode and the ST mode is that a liquid sheet can be stretched from the droplet main body under the ST mode. These deformation characteristics have been shown clearly in figures 4 and 5. In this section, the temporal evolutions of droplets under these two modes are analysed quantitatively. Specific variables which characterise the aerodynamic deformation and breakup process of droplets are chosen, including the dimensionless height of the upper surface ( $h^*$ ), the dimensionless radius of the upper surface ( $l^*$ ), the radius of the protrusion measured from the centre of droplet to the rim ( $r^*_{max}$ ), and the angle between the protrusion position and the vertical direction ( $\theta$ ). The definitions of these variables can be seen in the insets of figures 6(a) and 7(a), and all the variables are non-dimensionalised by the length scale of droplet radius $a$ . It is notable that the values of $h^*$ , $l^*$ , $r^*_{max}$ and $\theta$ are interrelated when the upper surface has been compressed to a flat plane, with the quantitative relationship of $h^*=r^*_{max}\cdot$ sin $\theta$ and $l^*= r^*_{max}\cdot$ cos $\theta$ . Even so, we will analyse the temporal variations of these four characteristic parameters simultaneously as the value of $h^*$ represents the degree of compression along the vertical direction and the variations of $l^*$ represent the degree of droplet spreading at the horizontal direction. The evolutions of $r^*_{max}$ and $\theta$ are analysed to characterise the dynamics of ‘rim’ structures.

Figure 6.

Temporal evolutions of (a) the dimensionless height of droplet $h^*$ , (b) the dimensionless radius of upper surface $l^*$ , (c) the dimensionless radius of protrusion $r^*_{max}$ , and (d) the angle between the protrusion and the vertical direction $\theta$ under CP mode.

Figure 6 presents the temporal evolutions of the dimensionless quantities under the CP mode, in which the cases with different values of $We_{lo}$ are considered. The corresponding values of $We$ and $\eta$ for each case are also given. It is observed that the temporal evolutions of the droplet height $h^*$ and the radius of upper surface $l^*$ almost coincide under different values of $We_{lo}$ , as shown in figures 6(a) and 6(b), respectively. The temporal evolutions of the dimensionless radius of the protrusion $r_{max}^*$ are observed to vary slightly with $We_{lo}$ , as shown in figure 6(c). Specifically, the variations of $r_{max}^*$ present the accelerating process, with a larger $We_{lo}$ resulting in a faster protruding velocity, especially at the later stage of deformation. As for the temporal evolutions of the protrusion angle $\theta$ (see figure 6 d), it can be seen that before a critical instant $t^*=0.4$ , the values of $\theta$ remain almost constant and are slightly larger than $\pi /4$ , indicating that the protrusion on the droplet surface mainly goes through radial extension during this period. After the critical instant, the values of $\theta$ increase significantly, indicating that the protrusion mainly spreads along the horizontal direction as the droplet main body is compressed significantly. Therefore, the instant of $t^*=0.4$ can be chosen as the criterion to differentiate between the early stage and later stage of droplet deformation process.

Figure 7.

Temporal evolutions of (a) the dimensionless height of droplet $h^*$ , (b) the dimensionless radius of upper surface $l^*$ , (c) the dimensionless radius of protrusion $r^*_{max}$ , and (d) the angle between the protrusion and the vertical direction $\theta$ under the ST mode.

Figure 7 presents the temporal evolutions of the dimensionless quantities under the ST mode, considering the cases with different values of $We_{lo}$ (with the corresponding values of $We$ and $\eta$ given simultaneously). Clearly, the temporal evolutions of the dimensionless droplet height $h^*$ are similar to those of the CP mode, in which the results under different values of $We_{lo}$ almost coincide with each other, as shown in figure 7(a). Different from the CP mode, the evolution of the upper surface radius $l^*$ exhibits an obvious dependence with $We_{lo}$ , as shown in figure 7(b). Specifically, a faster spreading of the upper surface can be observed with the increase of $We_{lo}$ , especially beyond the early stage of droplet deformation (i.e. $t^* \geqslant 0.4$ ). The main reason lies in the formation of a liquid sheet structure. A larger $We_{lo}$ results in a stronger shear effect of the air flow close to the droplet, which stretches the droplet surface more quickly to form the liquid sheet, and thus results in a faster temporal increase of $l^*$ . It should be noted that different combinations of $We$ and $\eta$ values do not affect the droplet dynamics once the value of $We_{lo}$ is determined, as presented by the two cases with $We_{lo}=35$ . The formation of a liquid sheet at the later stage of deformation also leads to a significant increase of the radius of protrusion $r_{max}^*$ with the increase of $We_{lo}$ , as shown in figure 7(c). As for the angle $\theta$ of the protrusion on the lateral surface of a droplet, the evolution tendency under the ST mode resembles that observed under the CP mode, as illustrated in figure 7(d). Specifically, the values of $\theta$ remain nearly constant before the time instant of $t^* = 0.4$ and then increase rapidly until almost reaching $\pi /2$ , which signifies that the liquid sheet detaches from the droplet main body and impacts the wall eventually.

4.2. Transition mechanisms between CP and ST modes

As the presence of the LS brings in quantitative differences between CP and ST modes, the formation process of the LS is analysed in this section, which contributes to revealing the transition mechanisms between these two modes. According to the droplet morphologies in figures 4(b) and 4(c), whether the ‘rim’ structure continues to be stretched outward to form the LS structure is the main criterion for distinguishing the CP mode and the ST mode. Therefore, it is crucial to analyse the local gas flow field around the protrusion, which contributes to understanding how the aerodynamic force influences the formation of the LS structure. It is notable that the high velocity characteristics of the airflow and the small spatial scale make it challenging to use the PIV method to resolve the local airflow quantitatively. Therefore, we choose the schlieren method to analyse the local air flow field close to the ‘rim’ structure qualitatively. Due to the obvious local velocity gradient, a density gradient of air exists spatially, making the use of the schlieren method feasible.

Figure 8.

Droplet deformation processes captured by the high-speed schlieren and high-speed photography methods, showing the difference of flow separation between (a) the CP mode and (b) the ST mode. The images of high-speed photography are annotated with stream lines for sketch, and the separation points in schlieren images are marked with red arrows. The parameters are: (a) $We_{lo}=12.5$ $(\eta =0.22, We=259$ , corresponding to $U_0=105$ m s⁻¹, $R= 6$ mm and $a=1.33$ mm) and (b) $We_{lo}=38.9$ $(\eta =0.24, We=674$ , corresponding to $U_0=151$ m s⁻¹, $R= 6$ mm and $a=1.4$ mm).

Figures 8(a) and 8(b) show the droplet deformation processes under the CP mode and ST mode, respectively. The results captured by high-speed schlieren and high-speed photography are given in the right-half and left-half of each figure, respectively. According to the qualitative flow fields revealed by high-speed schlieren, we also annotate the sketch of streamlines on images obtained from high-speed photography. At the very early stage of deformation ( $t^* = 0.12$ ), the droplets in the CP mode and the ST mode exhibit nearly identical profiles, with only minor differences in the tendencies of the air flow field. Flow separation occurs near the lateral side of the droplet for both cases, as indicated by the schlieren images where the separation points are marked with arrows. The flow separation results in the formation of vortex regions downstream of the separation point, and the airflows upstream and downstream of the separation point converge. At $t^* = 0.28$ , although there is no significant difference between the droplet profiles under the CP and ST modes, obvious differences arise in the airflow field. In the CP mode, the flow separation point remains on the lateral side of the droplet, downstream of the ‘rim’ structure. In contrast, in the ST mode where the value of $We_{lo}$ is larger than that in the CP mode, the flow separation occurs directly at the ‘rim’ due to a higher local airflow velocity and a larger protrusion of the ‘rim’ structure. As the deformation continues (see panel at the instant of $t^* = 0.43$ ), if the flow separation point coincides with the location of the protrusion, the convergence of the airflow can promote the accumulation of liquid at the protrusion, increasing the protrusion’s amplitude and facilitating the formation of a liquid sheet (i.e. the ST mode), as shown by figure 8(b). However, the growth of a protrusion can induce more intense flow separation at the protrusion position, and thus creates a positive feedback loop for the stretching of the liquid sheet. Nevertheless, under the CP mode shown in figure 8(a), the liquid accumulation at the protrusion is relatively weak as the protrusion position locates upstream of the separation point. Therefore, a faster increase of $\theta$ in the ST mode can be observed compared with the CP mode, as has been shown in figures 6(d) and 7(d). At the later stage of deformation ( $t^* = 0.74$ ), the ‘rim’ under the ST mode is stretched significantly into an LS structure due to liquid accumulation. However, for the CP mode, although the continuous compression of upper surface leads to the coincidence of the ‘rim’ and the separation point, the extrusion is not sufficient enough to stretch the surface to form the liquid layer under a relatively small value of $We_{lo}$ . Therefore, the ‘rim’ remains unstretched, eventually being compressed together with the main body of the droplet by the airflow.

Based on the analysis above, it can be summarised that the greater aerodynamic force (represented by the larger $We_{lo}$ ) would lead to a larger protrusion at the droplet lateral side due to the superposition of the flow separation point and protrusion position, which contributes to the stretching of the liquid sheet and thus the transition from the CP mode to ST mode.

4.3. Theoretical analyses

To predict the dynamics of droplet evolution under the CP and ST modes quantitatively, we establish models to predict the droplet deformation process and propose theoretical analyses in this section. The CP mode and the ST mode are analysed by unified models, but the early stage and the later stage of droplet deformation are described separately. As has been indicated by figures 6(d) and 7(d), we choose the instant of $t^*=0.4$ to distinguish between the early stage and the later stage. For theoretical analyses, it is notable that the characteristic lengths chosen for the gas phase and the liquid droplet are different. For the airflow, its length scale is non-dimensionalised by the characteristic scale of the stagnation flow, i.e. $R$ , and the non-dimensional quantities are denoted with a superscript ( $\sim$ ). As for the droplet, the characteristic length is chosen as the droplet radius $a$ and the non-dimensional quantities are denoted with a superscript ( $*$ ).

4.3.1. Potential flow model at early stage of deformation

Up to now, the potential flow model has been successfully used to quantify the airflow around the droplet, enabling the theoretical predictions on the vibrating dynamics of a free droplet (Rimbert et al. Reference Rimbert, Escobar, Meignen, Hadj-Achour and Gradeck2020), and the deformation and breakup characteristics of a two-dimensional wall-attached hydrophilic droplet (Hooshanginejad et al. Reference Hooshanginejad, Dutcher, Shelley and Lee2020; Chen et al. Reference Chen, Hooshanginejad, Kumar and Lee2022). In this section, we also employ the potential flow model to investigate the axisymmetric stagnation airflow and the evolutive dynamics of the hemispherical droplet. At the early stage of deformation, the droplet goes through relatively small deformation on both the upper surface and the lateral side. For example, the height of upper surface is compressed by no more than 20 % as $t^* \leqslant 0.4$ (see figures 6 a and 7 a). Moreover, as the values of $Re_{lo}$ range from $O(10^3)$ to $O(10^4)$ under the parametric zone of the CP and ST mode, the viscous effect can be neglected at the upper surface even though the flow separation occurs close to the solid wall. Therefore, we employ the potential flow model within the framework of irrotational and incompressible flow to predict the airflow around the hemispherical droplet. Additionally, we assume that the external flow field of the air stream does not exhibit significant variations from the initial state at the early stage of droplet deformation, and thus only the steady potential flow solution considering the initial hemispherical droplet is employed to calculate the air flow field. While this assumption may introduce some deviations on the droplet lateral side due to the gradually pronounced droplet deformation and flow separation, it remains reasonably valid close to the upper surface, which can also be seen in a previous study by Hooshanginejad et al. (Reference Hooshanginejad, Dutcher, Shelley and Lee2020).

Figure 9.

Airflow fields obtained by PIV (left panel) and the potential flow model of (4.7) (right panel) under $We_{lo}=5.8$ $(\eta =0.27, We=69$ , corresponding to $U_0$ = 38 m s⁻¹, $R$ = 5 mm and $a$ = 1.33 mm). The spherical coordinates $\tilde {r}$ and $\tilde {\theta }$ accompanying the velocity components $\tilde u_r$ and $\tilde u_\theta$ are sketched, and the droplet surface is denoted by $\tilde {r}=\eta$ .

For solving the axisymmetric stagnation point airflow around the droplet, we establish a spherical coordinate system with the droplet centre as the origin, as shown in figure 9. On the incompressible framework of $\nabla \cdot \tilde u = 0$ , where $\tilde u$ denotes the flow velocity, there exists a stream function $\psi ( r,\theta )$ which satisfies

(4.1) \begin{equation} \frac {{\partial \psi }}{{\partial \tilde r}} = - \tilde r{\tilde u_\theta }\sin \theta , \quad \frac {{\partial \psi }}{{\partial \theta }} = {\tilde r^2}{\tilde u_r}\sin \theta , \end{equation}

where $\tilde u_\theta$ and $\tilde u_r$ (sketched in figure 9) represent the velocity components in the $\theta$ and $\tilde {r}$ directions, respectively. As the flow is also assumed to be irrotational, which can be written as $\nabla \times \tilde u = 0$ , the stream function $\psi ( r,\theta )$ satisfies

(4.2) \begin{equation} \frac {{{\partial ^2}\psi }}{{\partial {{\tilde r}^2}}} + \frac {1}{{{{\tilde r}^2}}}\frac {{{\partial ^2}\psi }}{{\partial {\theta ^2}}} + \frac {{\cot \theta }}{{{{\tilde r}^2}}}\frac {{\partial \psi }}{{\partial \theta }} = 0. \end{equation}

However, the stream function of an axisymmetric stagnation point flow without the droplet in spherical coordinates can be expressed as

(4.3) \begin{equation} {\psi _0} = - \frac {1}{2}A{\tilde r^3}\cos \theta {\sin ^2}\theta , \end{equation}

which has been verified under the cylindrical coordinates, as depicted in figure 2. When there exists a droplet, we introduce $\psi _1(r,\theta )$ as a correction of the stream function near the stagnation point, which can be separated into $R(\tilde r)$ and $\Theta (\theta )$ . In this situation, the stream function $\psi$ can be given as

(4.4) \begin{equation} \psi = {\psi _0} + {\psi _1} = - \frac {1}{2}A{\tilde r^3}\cos \theta {\sin ^2}\theta + R(\tilde r)\Theta (\theta ). \end{equation}

Combining (4.2), (4.3) and (4.4), and taking into account the boundary condition that the normal velocity equals zero at the droplet surface $\tilde r = \eta$ , we can obtain

(4.5) \begin{equation} \Theta (\theta ) = \cos \theta {\sin ^2}\theta , \quad R(\tilde r) = \frac {1}{2}A{\eta ^5}{\tilde r^{ - 2}}. \end{equation}

Therefore, we can finally obtain the forms of $\psi _1$ and $\psi$ as

(4.6) \begin{equation} \begin{aligned} {\psi _1}(\tilde r,\theta ) = R(\tilde r)\Theta (\theta ) = \frac {1}{2}A{\eta ^5}{\tilde r^{ - 2}}\cos \theta {\sin ^2}\theta ,\\ \psi = {\psi _0} + {\psi _1} = \frac {1}{2}A({\eta ^5}{\tilde r^{ - 2}} - {\tilde r^3})\cos \theta {\sin ^2}\theta . \end{aligned} \end{equation}

According to (4.1), we obtain the velocity components of the airflow, i.e.

(4.7) \begin{equation} \begin{aligned} {{\tilde u}_r} = \frac {1}{2}A({\eta ^5}{{\tilde r}^{ - 4}} - \tilde r)(1 - 3{\cos ^2}\theta ), \\ {{\tilde u}_\theta } = \frac {1}{2}A(2{\eta ^5}{{\tilde r}^{ - 4}} + 3\tilde r)\sin \theta \cos \theta . \end{aligned} \end{equation}

Particularly, at the surface of the droplet ( $\tilde r = \eta$ ), the radial and azimuthal velocities are

(4.8) \begin{equation} {\tilde u_r}(\eta ,\theta ) = 0,\quad {\tilde u_\theta }(\eta ,\theta ) = \frac {5}{4}A\eta \sin 2\theta . \end{equation}

To verify the potential flow solution, we compare the flow fields of PIV measurements with the theoretical results of (4.7), and the results are shown in figure 9. Qualitatively, the flow streamlines present similar profiles between the experimental results and the theoretical predictions. To give more quantitative validations, we also choose specific positions in the flow field and compare the velocity components of $\tilde u_r$ and $\tilde u_\theta$ between experimental and theoretical results, and the details can be found in Appendix B. According to (4.8), the minimum velocity on the droplet surface appears at the centre of the upper surface ( $\theta =0$ ), and the maximum airflow velocity occurs at the lateral side of $\theta =\pi /4$ . These tendencies can also be seen clearly on the velocity contour in figure 9.

Since we assume a steady potential flow model under small droplet deformation, according to the Bernoulli theorem, flow regions with higher velocities result in lower pressures, i.e.

(4.9) \begin{equation} {p_0}=p_1+ \frac {1}{2}{\rho _g}U_0^2 = {p_g}(r,\theta ) + \frac {1}{2}{\rho _g}(u_r^2 + u_\theta ^2), \end{equation}

where $p_0, p_1$ and $p_g$ represent the pressure at the stagnation point of the surface, at the far field and at arbitrary location within the flow field, respectively. Therefore, a minimum aerodynamic pressure is expected to occur at the position of $\theta =\pi /4$ on the droplet surface. This qualitatively explains why both CP and ST modes exhibit protrusions close to the angles of $\theta =\pi /4$ at the early stages, as shown in figures 6(d) and 7(d), respectively. It should be noted that the Bernoulli analysis has also been successfully applied in previous studies considering the deformation of a droplet in airflow (Villermaux & Bossa Reference Villermaux and Bossa2009; Rimbert et al. Reference Rimbert, Escobar, Meignen, Hadj-Achour and Gradeck2020) or the motion of a gas bubble in viscous liquid (Moore Reference Moore1959).

Similarly, the Navier–Stokes equations describing the evolution of a liquid droplet can also be simplified at the early stage of deformation. As the droplet exhibits relatively small deformation, the magnitude of the convective acceleration term in the Navier–Stokes equations is much smaller than that of the local acceleration term, and thus the convective acceleration term can be neglected. Based on this assumption, Yi et al. (Reference Yi, Zhu, Yang, Wang and Sun2017) simplified the inviscid Navier–Stokes equations and obtained the acceleration $a_l(r,\theta )$ at an arbitrary position inside the droplet with respect to pressure $p_l$ , with the form

(4.10) \begin{equation} a_l\left ( r,\theta \right ) = -\frac {1}{2\rho _l}\left ( \frac {1}{r}\frac {\partial ^2p_l}{\partial \theta ^2}+\frac {\cot \theta }{r}\frac {\partial p_l}{\partial \theta } \right ), \end{equation}

where $\rho _l$ is the density of the droplet and $p_l$ denotes the pressure inside the droplet. Through (4.10), the deformation of the droplet surface can be obtained by setting $r=a$ . By using the potential flow solution of (4.8) and the Bernoulli theorem of (4.9), we can obtain the aerodynamic pressure distribution on the droplet surface. Considering the Young–Laplace equation, we can obtain the pressure distribution inside the droplet surface as

(4.11) \begin{equation} p_l(a,\theta )-p_g(a,\theta )=\sigma (\frac {1}{R_1(\theta )}+\frac {1}{R_2(\theta )}), \end{equation}

where $R_1(\theta )$ and $R_2(\theta )$ represent the principal radii of curvature at the droplet surface. Finally, we can obtain the inviscid radial acceleration of the droplet surface after non-dimensionalisation, with the form

(4.12) \begin{equation} a^*_l(\theta )=-\frac {25}{16}\left ( \cos 4\theta +\cos ^2\theta \cos 2\theta \right )+\frac {1}{We_{lo}} \left (\frac {\partial ^2f(\theta )}{\partial \theta ^2}+\cot \theta \frac {\partial f(\theta )}{\partial \theta } \right ), \end{equation}

where $f(\theta )=a/R_1(\theta )+a/R_2(\theta )$ . Since the variation of the droplet curvature gradient is relatively small in the early deformation stage under the CP mode and the ST mode, the second term on the right-hand side of (4.12), representing the acceleration induced by surface tension force, remains much smaller than the first term in magnitude, which represents the acceleration induced by the aerodynamic force. For example, for the cases considered in figures 8(a) and 8(b), we can compare the magnitudes of the second term and the first term on the right-hand side of (4.12) through fitting the experimental surface profile of the droplet by a polynomial, which is able to provide the specific value of $f(\theta )$ . At the upper surface, it is found that the values of the second term remain at least two orders lower than those of the first term at different instants during the early stage, indicating that the effect of surface tension is relatively small. Therefore, our theoretical analysis neglects the contribution of the surface tension force for simplicity. Clearly, the extrema of droplet surface acceleration occur at azimuthal positions of $\theta =0$ and $\pi /4$ , corresponding to the diminution of droplet height and the extrusion of protrusion at the lateral side, respectively. The accelerations for $h^*$ and $r_{max}^*$ can be obtained as

(4.13) \begin{equation} a^*_l(0)=-\frac {25}{8},\quad a^*_l(\frac {\pi }{4})=\frac {25}{16}. \end{equation}

According to (4.13), the droplet height $h^*$ and the radius of protrusion $r^*_{max}$ both present uniform accelerating behaviours as the surface accelerations at $\theta =0$ and $\pi /4$ maintain constant values, and the corresponding relationship can be written as

(4.14) \begin{equation} h^*=1+a^*_l(0) \cdot t^{*2},\quad r^*_{max}=1+a^*_l(\frac {\pi }{4}) \cdot t^{*2}. \end{equation}

The quantitative comparison between experimental results and theoretical predictions will be given later in § 4.3.3.

4.3.2. Energy conservation model at later stage of deformation

Figure 10.

Sketch of the energy conservation model, where the surface profile is assumed as a trapezoid with a pinned contact line. As the time $t^*$ increases, the bottom edge of the trapezoid remains unchanged, while the top edge spreads horizontally and undergoes vertical compression under the action of aerodynamic force $F_d$ .

At the later stage of droplet evolution, as the droplet deforms significantly to an approximately trapezoidal profile, we employ an energy conservation approach to characterise the horizontal spreading and vertical compression of the upper surface. A sketch of the deformation process is given in figure 10. Upon the compression of the aerodynamic force of the airflow (denoted by $F_d$ ), the upper surface of the trapezoidal droplet evolves with the decrease of height $h$ and the increase of width $l$ . Furthermore, we assume a pinned contact line of the droplet for the convenience of analyses, which is consistent with the experimental observation that the contact line only exhibits minimal movement until the very later stage of deformation. For the energy conservation model of trapezoidal droplets, the aerodynamic work of airflow $W$ would be converted to the increase of droplet kinetic energy $E_k$ and surface energy $E_s$ , i.e.

(4.15) \begin{equation} W = E_k + E_s. \end{equation}

During a time interval from $t_0$ to $t_1$ , the work exerted by the air on the liquid droplet is given by

(4.16) \begin{equation} W=\int _{h(t_0)}^{h(t_1)}{F_d{\rm d}h}=\int _{t_0}^{t_1}{C_d \pi \rho _gU_h^2l^2h^{\prime}{\rm d}t}, \end{equation}

where the aerodynamic force is written as the form $F_d=C_d\pi \rho _gU_h^2l^2$ , $U_h$ ( $=AU_0h/R$ ) represents the characteristic velocity on droplet upper surface during the compression process and $C_d$ denotes the drag coefficient. As the droplet upper surface has been compressed to a horizontal plane with flow separation occuring at the lateral side when evolving to the later stage of deformation (see figure 8), the drag coefficient $C_d$ can be represented as $C_d=\pi /(\pi +4)$ , which has been derived analytically by Lamb (Reference Lamb1895) who considered the resistance generated by the perpendicular impacting of airflow on a flat plate through the free streamline theory. Considering the characteristic airflow velocity on the droplet surface $U_{lo}=A\eta U_0$ and the characteristic scale of droplet radius $a$ , we can derive the characteristic energy of droplet as

(4.17) \begin{equation} E \sim \rho _gU_{lo}^2a^3. \end{equation}

Therefore, the dimensionless form of (4.16) can be written as

(4.18) \begin{equation} W^*=\int _{t^*_0}^{t^*_1}{C_d \pi h^{*2}l^{*2}h^{\prime*}{\rm d}t^*}. \end{equation}

As the upper surface of droplet has been compressed to a plane in the later stage, the assumption of a trapezoidal profile of the droplet enables us to calculate the specific forms of the variations of the kinetic energy $E_k$ and the surface energy $E_s$ . The detailed deductions are given in Appendix C. After non-dimensionalisation, (4.15) finally becomes

(4.19) \begin{equation} \int _{t^*_0}^{t^*_1}{C_d h^{*2}l^{*2}h^{\prime*}{\rm d}t^*}=\frac {1}{20}h^*\left ( \frac {1}{6}+\frac {1}{2}l^*+l^{*2} \right ) \left ( 2h^{\prime*2}+l^{\prime*2} \right )+\frac { \left ( l^*+1 \right ) S^*+l^{*2}}{We_{lo}}, \end{equation}

where the left-hand side represents the work applied by the aerodynamic force on the droplet from initial moment $t_0^*$ to terminal moment $t_1^*$ , and the right-hand side represents the variations of the kinetic energy and the surface energy. Moreover, according to volume conservation, the droplet height $h$ and the upper surface radius $l$ satisfy the relationship

(4.20) \begin{equation} \frac {1}{3}\pi h\left ( l^2+a^2+la \right ) =\frac {2}{3}\pi a^3, \end{equation}

which corresponds to the dimensionless form

(4.21) \begin{equation} h^*\left ( l^{*2}+l^*+1 \right ) =2. \end{equation}

Starting from the instant $t^*=0.4$ , the temporal evolutions of $h^*$ and $l^*$ can be obtained through solving (4.19) in combination with (4.21) with the iteration of the time step.

4.3.3. Comparison between theoretical models and experiments

Figure 11.

Dynamic processes of droplets under (a) CP mode and (b) ST mode, where the experimental measurements (symbols) are compared with the potential flow model at the early stage (dashed lines) and the energy conservation model at the later stage (solid lines). The parameters are: (a) $We_{lo}=25.2$ $(\eta =0.36, We=160$ , corresponding to $U_0=75$ m s⁻¹, $R= 4$ mm and $a=1.42$ mm) and (b) $We_{lo}=59.3$ $(\eta =0.34, We=418$ , corresponding to $U_0=114$ m s⁻¹, $R= 4$ mm and $a=1.38$ mm).

Choosing the typical cases of the CP mode and the ST mode, we compare the theoretical models to the experimental results, as depicted in figures 11(a) and 11(b), respectively. According to the results shown in figures 6(d) and 7(d), the instant of $t^* = 0.4$ is chosen to divide the time intervals of the early stage and the later stage of deformation, in which different theoretical models are applied.

At the early stage of droplet deformation, the experiments are compared with the potential flow model described by (4.14), where the temporal evolutions of $h^*$ and $ r_{max}^*$ are given through the relationship of uniform acceleration. It is observed that the potential flow model can present favourable predictions on the evolutions of $h^*$ and $ r_{max}^*$ when $t^* \leqslant 0.4$ , both for the CP mode and the ST mode. It should also be emphasised that at time instants close to $t^*=0.4$ , the deviation between experiments and the theoretical model gradually appears, which is caused by the more obvious deformation of the droplet. Specifically, at $t^*=0.4$ , the deformations in $h^*$ and $r_{max}^*$ approximate 0.2 and 0.1, respectively (see the experiment results in figures 6 and 7). Therefore, the limitation of the potential flow model appears due to the assumption of small deformation.

At the later stage of droplet deformation, we start from the time instant $t_0^* = 0.4$ and use the energy conservation model (i.e. (4.19)) to predict the compressing and spreading dynamics of the droplet upper surface, as depicted in figures 11(a) and 11(b) for the CP mode and the ST mode, respectively. It is observed that before the critical instant around $t^* = 0.8$ , the energy conservation model provides relatively accurate predictions of the experimental results in both the CP mode (figure 11 a) and the ST mode (figure 11 b). This indicates that although the liquid sheet structure forms and stretches from the main droplet under the ST mode, it brings little influence on the prediction of droplet dynamical behaviours as the volume of the liquid sheet is relatively small compared with the main body of the droplet. At the very later stage of deformation (e.g. $t^* \gt 0.8$ ), the theoretical predictions under both the CP and ST modes diverge significantly from the experiments. Specifically, the values of $h^*$ predicted by theory decrease rapidly and become smaller than those measured in experiments, and the values of $l^*$ predicted by theory increase rapidly and remain larger than those in experiments. The reason for the deviations lies in the assumption of a pinned contact line in the theoretical model being not applicable at this very later stage, as revealed clearly by the experimental results in figure 4. The de-pinning of the contact line is related to the continuous compression of the droplet surface, which causes the apparent contact angle to exceed the advancing contact angle of the contact line (measured to be 120 degrees, as stated in § 2.1). When the contact line de-pins, the motion of the upper surface will also be affected by the bottom evolution of droplet due to volume conservation relationship, making the theoretical predictions depart from the experiments.

5.1. Spring–oscillator model of VB mode

As has been indicated before (see table 1, figures 4 and 5), the VB mode occurs under relatively low values of $We_{lo} (\sim O(1))$ , in which the competition between the aerodynamic force and the surface tension force leads to the periodic vibration of droplet surface. In this section, we use a spring–oscillator model to characterise the vibrating behaviour of droplets. The spring-oscillator model considers an object with mass $M$ connected to a spring with elastic coefficient $k$ . Under the action of an external force $F$ , the vibrating frequency $f$ of the oscillator maintains constant regardless of the magnitude of the external force (i.e. $f \sim \sqrt {k/M}$ ), while the vibrating amplitude $A$ is proportional to the magnitude of the external force (i.e. $A \sim F/k$ ). Likewise, during the vibration of droplets under the VB mode, the surface tension plays an analogous role to the spring while the aerodynamic force in the vertical direction acts as the external force, and the competition between the aerodynamic force and the surface tension is reflected by $We_{lo}$ . Therefore, we can infer that the droplet vibrating frequency is independent of $We_{lo}$ , while the vibrating amplitude is directly proportional to the value of $We_{lo}$ .

Figure 12.

(a) Dimensionless vibrating frequency and (b) dimensionless vibrating amplitude of droplets with the variation of $We_{lo}$ .

To verify the predictions of the spring–oscillator model, we measure the droplet vibrating frequency and amplitude at the central position of the upper surface, considering the variation of $We_{lo}$ . For the vibrating frequency, we non-dimensionalise the experimental results by the natural frequency of the droplet (denoted by $f_0$ ), which is calculated from Lamb’s correlation (Lamb Reference Lamb1895), i.e.

(5.1) \begin{equation} f_0=\sqrt {\frac {2\sigma }{\pi ^2\rho _la^3}}. \end{equation}

Figure 12(a) shows the dimensionless vibrating frequency of droplets $f^*$ with the variation of $We_{lo}$ , which indicates that the values of $f^*$ remain almost constant as $We_{lo}$ changes. The results are consistent with the spring–oscillator model, indicating that the VB mode can be recognised as a natural vibrating process of droplets. The variation of $We_{lo}$ can reflect the changes of droplet size or incoming airflow velocity, and thus the dimensionless vibrating frequency $f^*$ is independent of these factors. Moreover, the values of $f^*$ approximate 1.2 invariably, which are very close to the vibrating frequency of free spherical droplets studied by Shen et al. (Reference Shen, Li, Liu, Si, Liu, Tang and Fan2019), demonstrating the universality of the spring–oscillator model for predicting the droplet vibration. It should be noted that in a previous study by Rimbert et al. (Reference Rimbert, Escobar, Meignen, Hadj-Achour and Gradeck2020) where the vibration of a free droplet was considered, it has been observed that the dimensionless vibrating frequency would gradually decrease instead of maintain a constant value as the Weber number exceeds a critical value. However, this tendency has not been found for the wall-attached droplets in a stagnation flow as the transition from the VB mode to the CP mode would occur at a relatively small $We_{lo}$ close to 2. The vibrating amplitude of the droplet is non-dimensionalised by the droplet radius $a$ . For the spring–oscillator model, we substitute the external force $F$ and the elastic coefficient $k$ with the aerodynamic force $F_d$ $(=\rho _g U_{lo}^2 a^2)$ and the surface tension coefficient $\sigma$ . In this way, we can obtain the dimensionless vibrating amplitude $A^*$ of droplets as

(5.2) \begin{equation} A^*\sim \frac {F_d}{\sigma a}=We_{lo}. \end{equation}

Figure 12(b) shows the experimental measurements of $A^*$ under different values of $We_{lo}$ , which validate the proportional relation between $We_{lo}$ and $A^*$ , as predicted by the spring–oscillator model.

5.2. Kelvin–Helmholtz instability (KHI) of SIE mode

The aerodynamic deformation and breakup of a droplet present the SIE mode when the local Weber number is relatively large (e.g. $We_{lo}\sim O(10^2)$ , as indicated in table 1. The SIE mode is mainly characterised by the occurrence of KHI on the upper surface of the droplet, which further leads to the peeling and atomisation of secondary droplets. Figure 13 gives the entire process of the development of surface perturbation waves caused by KHI. It is observed that at the initial stage of deformation ( $t^*=0.269$ ), the upper surface forms the ‘ripple’ structure due to the uneven distribution of aerodynamic force, and the droplet behaviours are similar to those under the CP and ST modes. However, the KHI perturbation waves on the upper surface develop quickly under high- $We_{lo}$ conditions, and the amplitude of surface perturbation waves will increase to the approximate amplitude as the protrusion ( $t^*=0.318 \sim 0.367$ ). Similar to the protrusion which appears at the azimuthal position close to $\theta \approx \pi /4$ , the KHI wave also arises near this region due to the existence of maximum air flow velocity. Subsequently, the amplitude of the KHI perturbation waves increase and annihilate the protrusion under the shear effect of high-speed air flow ( $t^*=0.392$ ). The KHI waves gradually get peeled off from the main droplet, producing liquid filament structures ( $t^*=0.514$ ). Simultaneously, the high-speed air flow separates at the filament, and the local large air velocity accelerates the breakup of the filaments, eventually resulting in the formation of vast secondary droplets ( $t^*=0.612$ ).

Figure 13.

Deformation and breakup process of droplets on upper surface under SIE mode, where perturbation wavelength caused by KHI is defined as $\lambda _{KH}$ .

Up to now, a certain number of theoretical studies have been conducted to predict the perturbation wavelength of KHI, the most typical ones include the zero vorticity thickness model (Villermaux Reference Villermaux1998; Marmottant & Villermaux Reference Marmottant and Villermaux2004) and the thick vorticity thickness model (Marmottant & Villermaux Reference Marmottant and Villermaux2004). Schematic diagrams of these two models are given in figure 14. For the zero vorticity thickness model, as shown in figure 14(a), the velocity boundary layer and the shear force on the liquid–gas surface is totally ignored, and the piecewise uniform velocity profiles are assumed. For the thick vorticity thickness model (see figure 14 b), the surface shear force is taken into account by assuming a gas velocity boundary layer with a thickness of $\delta _g$ . Both models are established under the framework of linear instability analysis. In our study, we adopt these two models and compare the theoretical predictions with the experimental results. It is notable that during the entire process of perturbation development, obvious nonlinear behaviours can be observed as the filaments form and get peeled off from the main droplet. Therefore, we only measure the perturbation wavelength (denoted by $\lambda _{KH}$ ) at the initial stage of perturbation (e.g. $t^*=0.367 \sim 0.392$ in figure 13).

Figure 14.

Schematic diagram of (a) zero vorticity thickness model and (b) thick vorticity thickness model for KHI. (c) Comparison between the experimental results and the theoretical predictions of the dimensionless wavelength $\lambda _{KH}^*$ , in which the zero vorticity thickness model (5.6) and the thick vorticity thickness model (5.11) with $C$ = 0.1 are both considered.

In the zero vorticity thickness model, the analytical dimensional dispersion relationship of KHI can be written as (Villermaux Reference Villermaux1998)

(5.3) \begin{equation} \omega _{KH1}=\pm \frac {k_{KH1}}{\rho _l+\rho _g}\sqrt {\rho _l\rho _g\left ( U_l-U_g \right ) ^2-\left ( \rho _l+\rho _g \right ) \sigma _lk_{KH1}}, \end{equation}

where $\omega _{KH1}$ denotes the growth rate of perturbation, $k_{KH1}$ denotes the wavenumber which is inversely proportional to the wavelength $\lambda _{KH1}$ , $U_g$ denotes the local velocity of airflow and $U_l$ denotes the local velocity of the liquid. According to (4.8) obtained from the potential flow model, we choose the local gas velocity at azimuthal position $\theta =\pi /4$ where the unstable waves arise, which gives $U_g=5/4 U_{lo}$ . As $\rho _l \gg \rho _g$ and $U_g \gg U_l$ , the dispersion relationship can be non-dimensionalised by the characteristic time and length scales, i.e.

(5.4) \begin{equation} \omega _{KH1}^*=\dfrac {R}{AU_0}\sqrt {\dfrac {\rho _l}{\rho _g}}*\omega _{KH1},\quad k_{KH1}^*=ak_{KH1}. \end{equation}

Therefore, we can obtain the dimensionless form of the dispersion relationship as

(5.5) \begin{equation} \omega _{KH1}^*=\pm k_{KH1}^* \sqrt {\frac {25}{16}-\frac {k_{KH1}^*}{We_{lo}}}. \end{equation}

The dimensionless wavelength corresponding to the maximum perturbation growth rate can be written as

(5.6) \begin{equation} \lambda _{KH1}^{*}=\frac {2\pi }{k_{KH1\max }^{*}}=\frac {5\pi }{2We_{lo}}, \end{equation}

indicating the power law relationship of $\lambda _{KH1}^* \sim We_{lo}^{-1}$ . It is notable that the negative first power law between the wavelength and the Weber number was also adopted in previous studies considering the KHI of a liquid–liquid interface (Abe et al. Reference Abe, Kizu, Arai, Nariai, Chitose and Koyama2004; Hadj-Achour et al. Reference Hadj-Achour, Rimbert, Gradeck and Meignen2021).

In the thick vorticity thickness model, the dispersion relationship of KHI can be obtained as (Marmottant & Villermaux Reference Marmottant and Villermaux2004)

(5.7) \begin{equation} e^{-2\kappa }=\left [ 1-\kappa \left ( 2\varOmega +1 \right ) \right ] \frac {1+0.5\kappa \left ( \rho _l/\rho _g+1 \right ) \left ( 2\varOmega -1 \right )}{1+0.5\kappa \left ( \rho _l/\rho _g-1 \right ) \left ( 2\varOmega -1 \right )}. \end{equation}

Here, $\kappa =k_{KH2}\delta _g$ and $\varOmega =\omega _{KH2}/\kappa ( U_g-U_l )$ , where $\delta _g$ denotes the thickness of the gas boundary layer with the specific form of

(5.8) \begin{equation} \delta _g=C\frac {a}{\sqrt {Re_{lo}}}. \end{equation}

Here, $C$ refers to a constant coefficient and $Re_{lo}$ represents the local Reynolds number at the upper surface of azimuthal position $\theta =\pi /4$ . Therefore, we also choose $U_g=5/4 U_{lo}$ as the local airflow velocity. For the convenience of further analysis, $Re_{lo}$ can be defined in terms of $We_{lo}$ , i.e.

(5.9) \begin{equation} Re_{lo}=\frac {5\rho _gU_{lo}a}{4\mu _g}=\frac {5\sqrt {We_{lo}}}{4Oh}\frac {\mu _l}{\mu _g}\left ( \frac {\rho _g}{\rho _l} \right ) ^{\frac {1}{2}}, \end{equation}

where $Oh$ $(=\mu _l/\sqrt {\rho _l \sigma a})$ is the Ohnesorge number defined by the liquid phase. Considering $\rho _l \gg \rho _g$ and $U_g \gg U_l$ , we non-dimensionalise the dispersion relationship according to the characteristic scales in (5.4), and the perturbation wavenumber with the largest growth rate can be obtained as

(5.10) \begin{equation} k_{KH2\max }^{*}=1.5\sqrt {\frac {\rho _g}{\rho _l}}\frac {1}{\delta _g}. \end{equation}

Thus, the dimensionless wavelength of the most unstable KHI wave is

(5.11) \begin{equation} \lambda _{KH2}^{*}=\frac {2\pi }{k_{KH2\max }^{*}}=\frac {4\pi }{3}C\sqrt {\frac {\rho _l}{Re_{lo}\rho _g}}=\frac {4\pi }{3}C\sqrt {\frac {4Oh}{5\sqrt {We_{lo}}}\frac {\mu _g}{\mu _l}\left ( \frac {\rho _l}{\rho _g} \right ) ^{\frac {3}{2}}}, \end{equation}

indicating the power law relationship of $\lambda _{KH2}^* \sim We_{lo}^{-1/4}$ .

The comparison between the perturbation wavelengths in experiments and those predicted by theoretical models are given in figure 14(c), where the variation of $We_{lo}$ is considered. Clearly, the dimensionless wavelength $\lambda _{KH}^*$ decreases with the increase of $We_{lo}$ , which is qualitatively consistent with the theoretical models. However, there are significant discrepancies in the quantitative comparison. Specifically, the wavelengths measured in experiments are much larger than those of the theoretical predictions of the zero vorticity thickness model (see (5.6)), indicating that it is inappropriate to totally ignore the gas velocity boundary layer close to the droplet surface. As for the thick vorticity layer model, a good agreement between the theoretical predictions and the experimental measurements can be reached when choosing the constant coefficient $C= 0.1$ for (5.11). According to (5.8), the boundary layer thickness at $C=0.1$ ranges at $\delta _g \sim 10^{-3} a$ , which is significantly small compared with the size of the droplet. It is notable that the specific value of $C=0.1$ is also consistent with other studies when approximating the thickness of boundary layer (e.g. $C=0.1$ for Guo et al. (Reference Guo, Kang, Mu, Li and Si2023) and $C=0.12$ for Qaddah et al. (Reference Qaddah, Chapelle, Bellot, Jourdan, Kewalramani, Deborde, Hammes and Rimbert2024)), indicating the universality of the thick vorticity layer model on predicting the wavelength of KHI.