2.1. High-speed shadowgraphy and deformation measurement

Shadowgraph images of the droplets were captured using a high-speed camera (Photron Fastcam SA-5 1000K-M3) with backlighting provided by a constant-source floodlight (Dedolight DLHM4-300U). The settings at which the camera were operated were varied based on the requirements for each set of flow conditions and are tabulated in table 1. The lenses used provided spatial resolutions of approximately 0.0431 and 0.0390 mm pixel$^{-1}$ for the ‘side-view’ and ‘end-on’ sets, respectively. This was determined by comparing measurements of both the droplet and air jet needle outer diameters with measurements of each from the images where both were found to give the same value within the error of detecting the edges in the image (${\approx }1$ pixel per edge).

Table 1. Camera settings.

Two independent sets of high-speed shadowgraphy video were captured: the side view conventionally shown in similar experiments, and an end-on view. The end-on view has been shown previously for low viewing angles (${<}30^{\circ }$ to the side view) and primarily for the ST and catastrophic regimes by Theofanous et al. (Reference Theofanous, Mitkin, Ng, Chang, Deng and Sushchikh2012) and at ${<}45^{\circ }$ for the BS regime by Zhao et al. (Reference Zhao, Liu, Xu, Li and Lin2013), but has not been published before to the authors’ knowledge for high viewing angles or the other morphologies of breakup. The side-view set facilitated the measurement of the cross-stream (CS) and streamwise (SW) deformation of the droplets in time as has been done by others, while the end-on images allowed for an additional visualization of the deformation process. For the side-view experiments, the droplets were held ${<}1$ mm from the air-needle tip, with the line-of-sight of the camera at a $90^{\circ }$ angle to the jet centreline. For the end-on-view experiments, the line-of-sight of the camera was adjusted to approximately $30^{\circ }$ from the jet centreline such that the jet was directed towards the camera. In this case, the droplets were suspended approximately 10 mm from the tip of the air needle to allow backlighting of the droplet as the air needle would otherwise interfere. Examples of the images taken from the side and end-on views are shown in figure 3.

Figure 3. Example shadowgraphy images showing side view (a) and end-on view (b). The top row shows the initial undeformed droplet, while the bottom row shows the deformed droplets at the initiation time, $T_i$. (c) Illustrates the camera's alignment to the air flow for the side-view (top) and end-on (bottom) cases.

To find the CS and SW deformation of the droplets in time, the extents of the leeward, windward, top and bottom of the drop were identified for each frame of the video from the side-on experimental set. Image processing was carried out using an in-house Python code. The shadowgraph contour of the droplet was found using Canny edge detection in Python's scikit-image package. For the leeward, top and bottom extents of the droplet, the image was scanned from the outer edge of the image towards the droplet contour until the $x$-$y$ location of each extent was found. The irregularity of the windward side of the drop in the late stages of the deformation makes the definition of the windward extent difficult. To make the measurement more consistent within each test, the windward extent is found by tracing straight across the droplet from the leeward extent. Example results of this method are shown in figure 4. The SW dimension is then found as the difference in the $x$ locations of the leeward and windward extents, while the CS dimension is found as the difference in the $y$ locations of the top and bottom extents. Figure 5 shows an example of the CS and SW measurement for $We=13$. This method has been found to identify the edges of the droplets within $\pm$1 pixel ($\pm$0.04 mm).

Figure 4. Example of the deformation measurement for $We \approx 13$.

Figure 5. Droplet deformation in time for SW and CS directions for $We = 13$. Insets show measurement axes (top) and images of the droplet at the approximate deformation times. The initiation time corresponds to the point at which the droplet has reached its minimal streamwise deformation. Every fourth data point is plotted for clarity.

2.2. Flow measurement

Since the instantaneous air-flow rate is difficult to measure when the solenoid is first opened, the steady-state air flow was measured and used to characterize the flow. The air-flow rate was measured using a gas-flow rotameter (Matheson FM 1050 E700 for flow rates 13.75–30 l.p.m., Cole Parmer 03217-34 for flow rates 30–68 l.p.m.) and a pressure transducer (WIKA A-10). The measured standard air-flow rate was corrected to real conditions following the manufacturer's recommended method. The volumetric flow rate of the jet issuing from the air needle at atmospheric temperature and pressure, $Q_{jet}$, was then found by mass conservation and used to solve for the average velocity of the jet issuing from the air needle by $U_{ave} = Q_{jet}/A_{g}$ where $A_g$ is the outlet flow area of the gas jet.

The centreline (CL) velocity of the jet, $U_{CL}$, which is used as the characteristic gas speed, $U$, was found by assuming that the velocity profile of the flow in the needle follows the power law velocity profile, $U_{CL} = U_{ave}(n+1)(2n+1)/2n^2$, with $n = -1.7 + 1.8 \log Re_{CL}$ (Pritchard & Leylegian Reference Pritchard and Leylegian2011). The value of the power law exponent $n$ was found iteratively, using $Re_{ave}$ as an initial guess for $Re_{CL}$. It should be noted that this equation is valid for $Re_{CL} > 20\,000$ ($n>6$) and that in the present study $Re_{CL} = 5200\text {--}25\,000$ ($5 < n < 6.2$), however, good agreement in $We$ transitions and in various deformation measurements between previous works and the present experiments suggest that the approximation holds (see appendix B).

For the side-view cases, the pendant drop was positioned close enough to the air needle that the air jet in the vicinity of the deforming drop is in the near-field region, thus the centreline velocity of the jet does not change downstream (Schlichting & Gersten Reference Schlichting and Gersten2017). As in previous experiments, it was assumed that the upstream air-flow field is unaffected by the presence of the drop (for example, Flock et al. Reference Flock, Guildenbecher, Chen, Sojka and Bauer2012). Therefore, $U_{CL}$ is the characteristic air-flow velocity near the drop, $U$. This was also assumed for the end-on view experiments where the drop is farther from the nozzle tip, however, it is likely that the predicted velocities are high as the average velocity in the vicinity of the drop will be lower than when the drop is closer to the nozzle exit where the velocity profile near the centreline is flatter. For this reason, the end-on view experiments were used primarily for qualitative description of the process. The few quantitative measurements that were taken from the end-on views are indicated separately from the side-view experiments.

4.1. Internal flow

The modelling of the internal droplet flow first carried out by Villermaux & Bossa (Reference Villermaux and Bossa2009) is replicated in this section for completeness of the present analysis. However, no assumptions about the shape of the droplet will be made that constrain the surface curvatures in the present section as to allow for greater generality of the resulting equation. These assumptions will be left to later sections. The coordinate system used in the present analysis is illustrated in figure 16(a).

Figure 16. Illustrations of (a) the coordinate system of the present analysis, (b) the actual droplet shape during the balancing period and (c) the approximated ellipsoidal geometry during the balancing period.

The liquid flow field inside the drop is found using the Navier–Stokes and conservation of mass equations in cylindrical coordinates for incompressible, axisymmetric radial flow,

(4.1)\begin{gather} \rho_{l}\left(\frac{\partial u_{r}}{\partial t}+u_{r} \frac{\partial u_{r}}{\partial r}\right)={-}\frac{\partial p_{l}}{\partial r}+\mu_l\left[\frac{1}{r} \frac{\partial}{\partial r}\left(r \frac{\partial u_{r}}{\partial r}\right)-\frac{u_{r}}{r^{2}}\right], \end{gather}

(4.2)\begin{gather}r \frac{\partial h}{\partial t}+\frac{\partial\left(r u_{r} h\right)}{\partial r}=0, \end{gather}

where $p_l(r)$ is the radial pressure distribution inside the drop and $h = h(t)$ is the SW thickness of the drop as it flattens.

Defining the droplet's CS radius as $R=R(t)$ and assuming that its volume is proportional to $R^2h$ (true for either ellipsoidal or cylindrical geometries), (4.2) gives the radial velocity profile inside the drop as $u_r(r)=r\dot {R}/R$. Upon substitution into (4.1), the viscous terms cancel. Integrating over $R$ leads to (4.3).

(4.3)\begin{equation} R\ddot{R} = \frac{2}{\rho_l}\left(p_l(0) - p_l(R)\right). \end{equation}

The liquid pressure inside the drop is related to the surrounding air pressure by the Laplace pressure at the interface and is given by (4.4) where $\sigma$ is the interfacial surface tension and $\kappa$ is the local surface curvature.

(4.4)\begin{equation} p_l(r)=p_g(r)+\sigma \kappa. \end{equation}

The air-flow field is approximated as an inviscid stagnation-point flow against a plane at the windward pole of the drop with the axial component of the air velocity given by $U_z=-a z U/d_0$, where $U$ is the free-field air-flow speed, $z$ is the axial direction measured from the windward pole of the droplet and $a$ is the stretching rate of the air flow (i.e. the spacing between the streamlines of the flow). Here, $a$ essentially relaxes the planar stagnation-point flow assumption to approximate the windward droplet face topology. The limiting shapes of the deformation process are the initial sphere and a final disk, which correspond to $a=6$ and ${\rm \pi} /4$, respectively (Villermaux & Bossa Reference Villermaux and Bossa2009). The pressure field in $r$ and $z$ of the surrounding air flow is found by conservation of mass and momentum assuming inviscid, incompressible, quasi-steady and axisymmetric flow, which is solved at $z=0$ to find the pressure at the windward face of the drop given by (4.5).

(4.5)\begin{equation} p_{g}(r)=p_{g}(0)-\rho_{g} \frac{a^{2} U^{2}}{8 d_{0}^{2}} r^{2}. \end{equation}

Substituting (4.4) and (4.5) into (4.3) and non-dimensionalizing gives

(4.6)\begin{equation} \frac{\ddot{R}}{R}=\frac{1}{\tau^2}\left(\frac{a}{2}\right)^2 \left[1-\frac{8}{a^2We}\frac{(\kappa_p - \kappa_w)d_0^3}{R^2}\right], \end{equation}

where $\kappa _p$ and $\kappa _w$ are the curvatures at the periphery and windward pole of the drop, respectively. It is clear from (4.6) that the deformation rate of the droplet depends on its shape as it deforms as the aerodynamic forces depend on the stretching rate while the surface tension forces depend upon the curvature at the windward face and periphery of the drop, all of which change as the droplet deforms. In order to avoid the cyclical argument that to solve the deformation of the droplet one must know the deformation of the droplet, assumptions must be made to describe the droplet's shape during the deformation to determine the values of $\kappa _p$, $\kappa _w$ and $a$.

Equation (4.6) is the most general form of the model for the internal liquid flow of a droplet in a high-speed air stream. Similar approaches have been taken by Villermaux & Bossa (Reference Villermaux and Bossa2009) and Kulkarni & Sojka (Reference Kulkarni and Sojka2014), however, these attempts did not predict the constant radial growth rate that is found here. These analyses are discussed further in appendix C.

4.2. Droplet deformation

Using the equations for the internal flow of the droplet, the problem now turns to determining the constant radial growth rate that is observed during the initiation period of the breakup.

Upon inspection of (4.6), a constant radial growth rate would occur when the term in square brackets becomes zero, i.e. the aerodynamic and surface tension forces balance. This is expected as the deformation causes the droplet's windward face to flatten, thus $a$ decreases from $6$ towards ${\rm \pi} /4$ while $\kappa _w$ tends towards zero as the windward face flattens and $\kappa _p$ increases as the periphery thins, as illustrated in figure 16(b). This process is described as ‘flow balancing’ as the forces acting on the droplet change through the deformation until they balance resulting in the constant growth rate ($\ddot {R}=0$). This description is consistent with the division of the deformation stage into two periods; the initial windward pole flattening period and the constant radial expansion period as described in § 3.

The constant radial expansion rate is achieved once the flows have balanced, which occurs at the flow balancing time, $T_{bal}$. Since the radial flow near the pole at early $T$ does not greatly affect the flow at the periphery of the drop where the CS deformation measurement is taken, $T_{bal}$ is defined by the point at which the SW deformation becomes approximately constant. From the present experiments, $T_{bal}\approx 1/8$ appears to be a good estimate as shown in figure 17 where the constant SW deformation rate begins at $T\approx 1/8$. In other words, $T_{bal}$ is taken as the time at which the disk first begins to form on the windward face of the drop.

Figure 17. (a) CS and SW droplet deformation $L/d_0$ vs. time $T$ for $We = 9.6$ (Bag) and $We = 22.5$ (BS) over the initiation period. (b) Magnification of the region given by the box in (a). The constant SW deformation rate occurs at $T_{bal} \approx 1/8$ for both cases, indicated by the vertical dashed line.

Initially, the droplet is spherical, therefore $\kappa _w = \kappa _p$ and the surface tension terms cancel. However, as the windward face flattens, $\kappa _w$ rapidly diminishes and becomes negligible. The curvature of the periphery is composed of two curvatures: the thickness curvature, $\kappa _{h}=2/h$, and the azimuthal curvature, $\kappa _{\theta }=1/R$. In previous analyses, $\kappa _{\theta }$ was neglected as it was assumed to be much smaller than $\kappa _{h}$ (see appendix C). However, this is not the case for the levels of deformation that occur up to $T_{i}$ where $2R/d_0 < 1.6$, as shown in figure 18 where the ratio between the curvatures is less than an order of magnitude for deformations less than $2R/d_0 \approx 2$. This is likely a significant reason that 2-D numerical simulations, which neglect this curvature, perform poorly compared to their 3-D equivalents (Strotos et al. Reference Strotos, Malgarinos, Nikolopoulos and Gavaises2016). In figure 18, $h/d_0$ is found by conservation of mass for both ellipsoidal and cylindrical geometries ($h/d_0 = (d_0/2R)^2$ and $h/d_0 = 2/3(d_0/2R)^2$, respectively). The azimuthal curvature therefore cannot be neglected in the analysis of the droplet deformation for $T<T_{i}$, including over $T_{bal}$.

Figure 18. Ratio of the thickness curvature to the azimuthal curvature vs. extent of deformation for ellipsoid and disk geometries ($h/d_0 = (d_0/2R)^2$ and $h/d_0 = 2/3(d_0/2R)^2$, respectively).

Equation (4.6) is then,

(4.7)\begin{equation} \frac{\ddot{R}}{R}=\frac{1}{\tau^2}\left(\frac{a}{2}\right)^2\left[1-\frac{64}{a^2We} \left( \left(\frac{d_0}{h}\right) \left(\frac{d_0}{2R}\right)^2 +\left(\frac{d_0}{2R}\right)^3 \right)\right]. \end{equation}

To finalize the definition of the curvatures, an assumption must be made about the geometry of the droplet during the flow balancing such that $h$ and $R$ can be related. It is assumed that during the flow balancing period the droplet can be approximated by an ellipsoid as illustrated in figure 16(c) such that, by conservation of mass with the initially spherical droplet, $h/d_0 = (d_0/2R)^2$ (this assumption is discussed further in appendix C). This simplifies (4.7) to 4.8, where the only remaining unknown is $a$

(4.8)\begin{equation} \frac{\ddot{R}}{R}=\frac{1}{\tau^2}\left(\frac{a}{2}\right)^2\left[1-\frac{64}{a^2We} \left( 1 +\left(\frac{d_0}{2R}\right)^3 \right)\right], \end{equation}

where $a$ relates to the approximation of the flow around the deforming object as a planar stagnation-point flow. While its limits are defined as $a=6$ for the initial sphere and $a={\rm \pi} /4$ for a disk, the relationship between $a$ and the extent of deformation (i.e. the path of $a$ from $6$ to ${\rm \pi} /4$) is unknown. While it is possible that a relationship between $a$ and the shape of the droplet could be determined from numerical simulation or by highly detailed particle image velocimetry (PIV) measurements, it will be assumed that the dominant value of $a$ over the deformation period is its initial value as this results in the highest $\ddot {R}$ during the balancing period. Additionally, it is observed that the shape of the droplet as a whole is more spherical than disk like during the balancing period. Thus, over the balancing period, it is assumed that $a$ takes a constant value of $a=6$. The constant radial expansion rate is then found by integrating (4.8) over the flow balancing time. As described earlier in this section, the CS dimension of the drop remains roughly constant during the flow balancing, thus $R(T<T_{bal}) \approx d_0/2$. This integration leads to (4.9)

(4.9)\begin{equation} \frac{\dot{R}}{d_0/2} = \frac{1}{\tau}\left(\frac{a}{2}\right)^2\left(1-\frac{128}{a^2We}\right) T_{bal}. \end{equation}

The prediction of (4.9) is compared to the experiments in figure 19, where the constant radial growth rate is measured by performing a linear fit to the cross-sectional diameter measurements over the period $0.25<T<T_{i}$ for all cases. Good agreement between the theory and experiments is found for Bag and BS breakup morphologies (figure 19b). MB and ST morphologies exhibit greater deviation due to the early breakup as described in § 3, in particular the breakup of the liquid fragment at the top of the drop which results from the suspending needle which becomes the windward most point of the drop which affects the measurement as in § 2. This has a greater effect on the fitting procedure as there are fewer data points during the deformation phase which are undisturbed owing to the faster progression of the deformation relative to the frame rate used.

Figure 19. Comparison of predicted expansion rate from (4.9) to all experiments (a) and to NB, Bag and BS morphologies only (b). The box in (a) shows the region that is magnified in (b).

Figure 20 compares the present deformation prediction as well as the previous deformation models (the TAB model of O'Rourke & Amsden (Reference O'Rourke and Amsden1987) and the models of Villermaux & Bossa Reference Villermaux and Bossa2009; Kulkarni & Sojka Reference Kulkarni and Sojka2014) to the present experiments for a range of $We$ covering Bag and BS morphologies. Note that the present model does not rely on exactly predicting the temporal deformation and instead focuses only on predicting the constant deformation rate. As such, for figure 20, the present model is plotted assuming that the CS expansion first passes the periphery of the drop at $T \approx 0.25$, with the slope given by (4.9). Excellent prediction of the present model is shown, while the previous models show significant deviation and indeed a misrepresentation of the constant trend.

Figure 20. Comparison of TAB (O'Rourke & Amsden Reference O'Rourke and Amsden1987), Villermaux & Bossa (Reference Villermaux and Bossa2009) and Kulkarni & Sojka (Reference Kulkarni and Sojka2014) models and (4.9) to the CS measurement of the present experiments during initial deformation.

The present prediction also agrees with the numerical results of Jain et al. (Reference Jain, Tyagi, Prakash, Ravikrishna and Tomar2018), who studied the effect of density ratio, $\rho ^{\ast }=\rho _l/\rho _g$, on the breakup. In their work, $\dot {R}(T=1)=4.8$, $2.0$ and $0.8$ for $\rho ^{\ast } = 10$, $150$, and $1000$, respectively, at $We=20$. The corresponding values using the present analysis are $\dot {R}(T=1)=5.9$, $1.5$, and $0.6$. The prediction of the present analysis is higher than that of Jain et al. (Reference Jain, Tyagi, Prakash, Ravikrishna and Tomar2018) at $\rho ^{\ast }=10$ as the resistance due to the gas phase density on the expanding periphery is not accounted for. The present prediction is lower for the other cases as the measurements of Jain et al. (Reference Jain, Tyagi, Prakash, Ravikrishna and Tomar2018) are taken after $T_i$ ($T_i<1$), when the flow dynamics begins to change, as will be described in § 4.4. Chou & Faeth (Reference Chou and Faeth1998) empirically found the growth rate to be $2 \dot {R} / d_0 = 1.06$ and $1.22$ for their water experiments at $We=13$ and $20$, respectively ($2\dot {R} / d_0\tau = 0.5$, with $\tau =0.472$ and $0.41$ ms for $We=13$ and $20$, respectively). This result is low compared to the present predictions of $1.8$ and $2.26$; however, several high-viscosity cases were included in the results of Chou & Faeth (Reference Chou and Faeth1998). Adjusting their analysis for their water experiments only gives $2\dot {R}/d_0=1.31$ and $1.95$. These results agree relatively well with those of the present study considering the much lower time resolution of Chou & Faeth (Reference Chou and Faeth1998), which may cause a higher uncertainty in the linear fit.

4.3. Rim formation and prediction of initiation geometry

The majority of the child droplets that form from the breakup of the drops come from the rim geometries described in § 3. Furthermore, the rim size can be used to find the initiation geometry in order to determine the breakup morphology. It is therefore necessary to predict the rim size when it forms, which can be found by solving for the internal flow inside the windward disk.

As the droplet's windward face expands radially, the disk that forms has thickness $h$. At the instant of initiation, the periphery of the disk ‘rolls’ as the surface tension pinches the sides of the disk periphery to form a thick rim of major and minor diameters $2R_i/d_0$ and $d_{i} \approx h_i$, respectively. This geometry is illustrated for all breakup morphologies in figure 21. The pressure balance across the peripheral interface (4.4) can be used to find $h_i$, giving (4.10)

(4.10)\begin{equation} d_{i} = h_{i} =\left[\frac{1}{4 \sigma}\left(\rho_{l} \dot{R}^{2}+\rho_{a} \left(\frac{a}{2}\right)^{2} U^2\left(\frac{R_i}{d_{0}}\right)^{2}\right)-\frac{1}{2 R_i}\right]^{{-}1}. \end{equation}

Figure 21. Illustration of the (a) initial droplet and the deformed droplet geometry at $T_i$ for (b) Bag, (c) BS, (d) MB and ST morphologies, indicating measurement locations of $h_i$ and $R_i$. The rim is marked by cross-hatching while the undeformed core is marked by grey fill.

While $2R_i/d_0$ and $h_i$ both serve to balance the pressure across the interface, one must be known in order to uniquely determine the other. At this time, the mechanism which determines the unique values is unknown, therefore, an assumption will be made for $2R_i/d_0$ based on the present experiments so that $h_i$ may be determined. Figure 22(a) shows the measured CS diameter of the droplets at $T_{i}$, which for all breakup modes is determined from the average of the Bag and BS morphologies as $2R_{i}/d_0 \approx 1.6$. This value agrees with the breakup criterion for the TAB model, $2R/d_0 =1.5$, as well as the experimental findings of Zhao et al. (Reference Zhao, Liu, Cao, Li and Xu2011a), $2R/d_0 =1.8$, who used the older definition of $T_i$ that is slightly later than the definition presently used, resulting in a slightly higher value as discussed in § 3.

Figure 22. Measurements of droplet periphery thickness (a) and diameter (b) at $T_{i}$. Periphery thickness measurement is compared to the present model (4.11). Droplet CS diameter reaches a value of $2R_{i}/d_0 \approx 1.6$ for all breakup modes.

In this late stage of deformation, $a$ will satisfy the balancing of (4.3) for the initiated geometry and is found to be $a \approx \sqrt {80/We}$. This suggests a shift in dominance of the aerodynamic forces over the inertia of the expanding rim when $We > 80$, which coincides with the transition to the ST morphology as described by Guildenbecher et al. (Reference Guildenbecher, López-Rivera and Sojka2009). This then gives that $We_{c,ST}=80$, which agrees with the present experiments as well as with Guildenbecher et al. (Reference Guildenbecher, López-Rivera and Sojka2009) and Zhao et al. (Reference Zhao, Liu, Li and Xu2010). Using these assumptions, (4.10) simplifies to

(4.11)\begin{equation} \frac{d_{i}}{d_0} = \frac{h_{i}}{d_0} \approx \left[\frac{1}{4} \left( \frac{\rho_l \dot{R}^2 d_0}{\sigma} \right) + 2.6 \right]^{{-}1} = \frac{4}{\left(We_{rim}+10.4\right)}, \end{equation}

where the term $\rho _l \dot {R}^2 d_0 / \sigma$ is grouped as a ‘rim Weber number’, $We_{rim}$, which describes the competition between the radial momentum induced at the droplet periphery and the restoring surface tension of the stable droplet. For $We_{rim} < 1$, the radial flow momentum is not strong enough to significantly overcome surface tension and the drop will oscillate with NB occurring. However, for $We_{rim}>1$, the radial momentum of the drop's periphery is great enough to overcome this restoring force, allowing the periphery to continuously extend, facilitating the thinning and blowout of the interior. This criticality defines the onset of the Bag morphology, where $We_{rim}=1$ and (4.9) lead to $We_{c,Bag}\approx 8.8$, which agrees well with the present experiments, as shown in figure 23(a), and to the previous results discussed in § 1. Since the derivation of (4.11) does not account for any slowing in $\dot {R}$, which results from the thinning periphery, the present prediction is slightly high near $We_{c,Bag}$ where $We_{rim}$ is low.

Figure 23. Measurements of $We_{rim}$ vs. $We$ (a) and periphery thickness vs. $We_{rim}$ (b) compared to the prediction of (4.11). The criticality $We_{rim}=1$ agrees with the transition to the Bag morphology and defines $We_{c,Bag}\approx 8.8$ (a), and periphery thickness vs. $We_{rim}$ compared to the prediction of (4.11) (b).

The prediction of (4.11) is compared to the experiments with $We$ and $We_{rim}$ in figures 22(b) and 23(b), respectively. The measurements shown in these figures are taken at the periphery of the droplet in order to measure the thickness of the rim, as measurements taken across the centre of the drop result in larger values due to the undeformed core. Examples of this measurement are shown in figure 24. Although the trend agrees well with the experiments, the prediction of (4.11) consistently over-predicts the measurements by $\approx 0.05$ (note that a conservative measurement error of $\pm$1 pixel results in ${\pm }0.025$). A possible cause of this overprediction is that the present analysis neglects flow separation around the periphery of the windward disk. Flow separation would cause a lower pressure at the periphery than is presently modelled, especially near $T_i$ when the windward disk is at its thinnest. This affects the balance of the surface tension and pressure terms that govern the rim's formation and would cause it to be thinner than predicted.

Figure 24. Example of $h_i$ measurement from droplet periphery for low (a) and high (b) $We$.

With a good prediction of the dimensions of the initiation geometry, the core volume can be determined and used to predict the initiation geometry. The windward disk is modelled as a flat disk with rounded edges having a radius of $R_i$ and a thickness of $h_i$ (i.e. the intersection of a short cylinder and a torus). The leading coefficient of the $h_i^3$ term in the volume equation for this geometry is very small compared to the other coefficients and can therefore be neglected. The volume of the disk, $V_d$, is thus given by (4.12). The undeformed core is then the volume that remains from the initially spherical droplet ($V_0 = 4/3 {\rm \pi}(d_0/2)^3$), given by $V_{c}=V_{0}-V_{d}$ (4.13)

(4.12)\begin{gather} \frac{V_{d}}{V_0} = \frac{3}{2} \left[ \left(\frac{2R_i}{d_0}\right)^2 \left(\frac{h_i}{d_0}\right) - 2 \left(1-\frac{\rm \pi}{4}\right) \left(\frac{2R_i}{d_0}\right) \left(\frac{h_i}{d_0}\right)^2 \right] \end{gather}

(4.13)\begin{gather}\frac{V_{c}}{V_0} = 1 - \frac{3}{2} \left[ \left(\frac{2R_i}{d_0}\right)^2 \left(\frac{h_i}{d_0}\right) - 2 \left(1-\frac{\rm \pi}{4}\right) \left(\frac{2R_i}{d_0}\right) \left(\frac{h_i}{d_0}\right)^2 \right]. \end{gather}

Figure 25 plots $V_d$ (a) and $V_c$ (b) for the present experiments with the predictions. As with the prediction of $d_{i}$ (4.11), the trend matches well; however, $V_d$ is slightly overpredicted resulting in the underprediction of $V_c$, delivering nonsensical results such as $V_d > V_0$ for low $We$. The over/underprediction of these volumes is the consequence of the overprediction of $d_{i}$, therefore, to improve the model for the volume distribution at initiation, the prediction of $d_{i}$ must be further improved. However, the experimental results can still be used to define the criterion for the transition between the various breakup morphologies. These trends also agree reasonably well with the experimental results of Dai & Faeth (Reference Dai and Faeth2001), considering the potential errors of their measurement method (i.e. measuring the volume of the drops after the breakup, as well as relatively low spatial resolution).

Figure 25. Volume of initiation disk (a) and undeformed core (b) vs. $We$. Dashed and dash-dotted lines indicate where all and half of the droplet volume is contained in the disk, respectively. Solid lines give the predictions of (4.12) (a) and 4.13 (b).

For Bag breakup, it is expected that there is no undeformed core, therefore this morphology will occur when $V_d/V_0 \approx 1$ ($V_c/V_0 \approx 0$), which matches the experiments as indicated by the dashed line in figure 25.

As described in § 3, BS, MB and ST morphologies occur when $V_c/V_0 > 0$ ($V_d/V_0 < 1$), however, the criterion for $V_c/V_0$ required for the transition from BS to MB is not as obvious as the transition from Bag to BS. From figure 25 it can be seen that this transition occurs at $V_c/V_0 \approx 0.5$ (indicated by the dash-dotted line), i.e. exactly half of the drop is undeformed. The mass per unit frontal surface area of the undeformed core is maximal at this point (assuming it takes roughly the shape of a spherical cap). This gives the core the greatest resistance to being pushed downstream by the air flow and allows the sheet to be blown out before the core deforms significantly. This explanation agrees with the physical description of the MB breakup process in § 3.

The transition to ST from MB is governed by the aerodynamic forces and not $V_c$, therefore there is no expected insight into the ST transition from this volume-based analysis.

With this rough prediction of the breakup morphology and of the rim dimensions at initiation, a prediction for the child drop sizes which result from the breakup of the rim can now be formulated.

4.4. Rim expansion and bag blowout

The breakup sizes from the rim will be proportional to its diameter at breakup. It is therefore necessary to find the rim diameter after it has grown to the breakup size. From conservation of mass on the rim, which formed at $T_{i}$ with an initial minor and major diameter ($d_{i}$, and $2R_i$, respectively) and the final rim minor and major diameter ($d_{f}$ and $2R_f$), $d_{f}$ is found by (4.14)

(4.14)\begin{equation} d_{f} = d_{i} \sqrt{\frac{R_i}{R_f}} \end{equation}

Since the initiation geometry, $R_i$ and $d_{i}$, have already been defined in § 4.3, the remaining factor to be determined is $R_f$, i.e. the extent of the rim's expansion. The formation of the rim coincides with the formation of the sheet between the rim and the undeformed core which becomes the bag structures after $T_{i}$. The initial stage of the rim's growth continues at the constant $\dot {R}$ found by (4.9). During this stage, the bag size, $\beta$, is smaller in the SW direction than in the CS. Once the bag grows to and past the size of the rim, it pulls on the rim and accelerates its growth until the bag breaks, after which the rim continues to grow at a constant rate. This is seen in figure 5, and the geometries of interest are illustrated in figure 26.

Figure 26. Illustration of droplet geometry at (a) $T_0$, (b) $T_i$, (c) $T_i^+$ and (d) $T_f$ for bag breakup with the relevant dimensions for the rim expansion and bag growth indicated. The rim is marked by cross-hatching, while the remainder of the windward disk is marked by the white area.

When the bag eventually ruptures, it quickly recedes around the bag (as described by Lhuissier & Villermaux Reference Lhuissier and Villermaux2012) and collides with the droplet's rim, instigating its breakup. Therefore, it is assumed that the final size of the rim is proportional to the final size of the bag when it breaks. In the cases of BS and MB breakup, the bag is confined further by the undeformed core. By determining the bag size at burst, $\beta _f$, one can approximate $R_f$ and thus determine the child drop size, $d_c$. Vledouts et al. (Reference Vledouts, Quinard, Vandenberghe and Villermaux2016) modelled the breakup of accelerating liquid shells by supposing that the interface is susceptible to the Rayleigh–Taylor instability with a growth rate, $\omega$, given by (4.15),

(4.15)\begin{equation} \omega(t^{{\ast}})=\sqrt{\frac{\rho_l \ddot{\beta}^{2}(t^{{\ast}}) h(t^{{\ast}})}{2 \sigma}}, \end{equation}

where $\beta (t^{\ast })$ is the SW dimension of the bag and $h(t^{\ast })$ is the membrane thickness at time $t^{\ast }=t-t_i$, as illustrated in figure 26(c). The breakup time was estimated by (4.16) where $C$ is a constant of $O(1)$

(4.16)\begin{equation} \int_{0}^{t^{{\ast}}_{b}} \omega(t^{{\ast}}) \,\mathrm{d} t^{{\ast}} = C. \end{equation}

Performing a force balance on the tip of the bag as in Villermaux & Bossa (Reference Villermaux and Bossa2009), the bag tip acceleration is found as

(4.17)\begin{equation} \ddot{\beta}(t^{{\ast}})=\frac{p_0}{\rho_l h(t^{{\ast}})}. \end{equation}

Although Kulkarni & Sojka (Reference Kulkarni and Sojka2014) included surface tension in the force balance, its effect is negligible. At its largest value (when the bag size is the size of the rim), it effectively lowers the pressure by $\sigma /2 << p_0$ for all cases.

The value of $h(t^{\ast })$ is determined by conservation of mass of the expanding disk, approximating it as a cylinder of radius $R-h_i$ and thickness $h$, giving $h(t^{\ast })\approx V_b / {\rm \pi}(R-h_{i})^{2}$ where $V_b$ is the volume of the bag. The volume of the bag is the volume of the windward disk that is not in the rim ($V_b = V_d - V_r$). $V_d$ was given earlier by (4.12). The rim is modelled as a torus of major and minor diameters $2R_i$ and $h_i$, with $V_r$ given by (4.18). It is worthwhile to note that for the Bag morphology, (4.18) estimates that the rim contains 50 %–80 % of the droplet volume which agrees with the findings of previous works (see appendix B), and is thus a good approximation of the droplet geometry. Since the disk expands at a constant rate, the radius can be estimated by $R=R_i + \dot {R}t^{\ast }$

(4.18)\begin{gather} \frac{V_{r}}{V_{0}}=\frac{3 {\rm \pi}}{2}\left[\left(\frac{2 R_i}{d_{0}}\right) \left(\frac{h_i}{d_{0}}\right)^{2}-\left(\frac{h_i}{d_{0}}\right)^{3}\right] \end{gather}

(4.19)\begin{gather}\frac{V_{b}}{V_{0}}=\frac{3}{2}\left[\left(\frac{2 R_i}{d_{0}}\right)^{2} \left(\frac{h_i}{d_{0}}\right)-\left(\frac{2 R_i}{d_{0}}\right)\left(\frac{h_i}{d_{0}}\right)^{2}\left(2+\frac{\rm \pi}{2}\right)+ {\rm \pi}\left(\frac{h_i}{d_{0}}\right)^{3}\right] \end{gather}

The breakup time, $t^{\ast }_{b}$, is then found by substituting $h$ into (4.15) and performing the integral in (4.16), leading to (4.20) where the only term that remains to be determined in the constant $C$

(4.20)\begin{equation} t^{{\ast}}_{b}=\frac{\left[\left(\dfrac{2 R_{i}}{d_{0}}\right)- 2\left(\dfrac{h_{i}}{d_{0}}\right)\right]}{\left(\dfrac{2 \dot{R}}{d_{0}}\right)} \left[{-}1+\sqrt{1+C \frac{8 \tau}{\sqrt{3 W e}} \sqrt{\frac{V_{b}}{V_{0}}} \frac{\left(\dfrac{2 \dot{R}}{d_{0}}\right)}{\left[\left(\dfrac{2 R_{i}}{d_{0}}\right)-2 \left(\dfrac{h_{i}}{d_{0}}\right)\right]}}\right]. \end{equation}

A much more detailed analysis of the instability growth in the bag is required in order to determine $C$ analytically, however, the present data are not well suited for this as the instability cannot be resolved on the bag. To determine $C$, (4.20) is fit to the breakup time data of the present experiments and is found to be $C \approx 9.0$, as shown in figure 27(b). In Vledouts et al. (Reference Vledouts, Quinard, Vandenberghe and Villermaux2016), (4.16) was applied directly to an already thin shell where the instability grows for the entire duration of $t^{\ast }_b$, thus the breakup time is relatively short and $C=O(1)$. In the present analysis, the definition of $t^{\ast }_b$ to which (4.20) is fit includes the initial formation of the bag after $T_i$ when $\ddot {\beta }$ is small and $h$ is relatively large. In this case, the instability does not develop appreciably for much of $t^{\ast }_b$, and $C>O(1)$. The bag size at time $t^{\ast }$, $\beta (t^{\ast })$, is found by integrating (4.17) with the substitution of $h$, leading to (4.21)

(4.21)\begin{align} \frac{\beta\left(t^{*}\right)}{d_{0}} &= \frac{3}{4}\left(\frac{V_{0}}{V_{b}}\right) \frac{1}{\tau^{2}}\left\{\left[\left(\frac{2 R_{i}}{d_{0}}\right)-2 \left(\frac{h_{i}}{d_{0}}\right)\right]^{2} \frac{t^{* 2}}{2} \right. \nonumber\\ &\quad \left.+\left(\frac{2 \dot{R}}{d_{0}}\right) \left[\left(\frac{2 R_{i}}{d_{0}}\right)-2\left(\frac{h_{i}}{d_{0}}\right)\right] \frac{t^{* 3}}{3}+\left(\frac{2 \dot{R}}{d_{0}}\right)^{2} \frac{t^{* 4}}{12}\right\}. \end{align}

Figure 27. Breakup time (a) and bag size at burst (b) vs. $We$. The solid lines are (4.20) which is fit to the data giving $C \approx 9$ (a), and (4.21) evaluated at $t^{\ast }_b$ (b).

Solution of (4.21) at time $t^{\ast }_{b}$ (4.20) gives the bag size at the time it ruptures, $\beta _f$, which is compared to the experiments in figure 27(b). Although the prediction of $\beta$ and $\beta _f$ match the experiments within the errors typical of these experiments (see appendix B), it is clear that the model does not capture all of the dynamics of the bag. In particular, the model underpredicts $\beta _f$ at the higher $We$ range of the BS morphology ($20 < We < 30$). It is likely that this is the result of the shape of the bag around the stamen. This shape causes the bag to be thicker between the stamen and the tip of the bag, thus the thinnest point of the bag is between the bag tip and the rim. This would result in a slightly longer breakup time, as the most susceptible point on the bag is farther away from the tip where the acceleration is expected to be the greatest.

Equation (4.21) evaluated over $0<t^{\ast }<t^{\ast }_b$ is compared to experiments in figure 28, showing excellent agreement over a range of $We$ for each morphology, with the exception of the upper $We$ limit for MB (f). In this case, the SW measurement appears to grow faster than the prediction. At the high $We$ limit of the MB morphology, the ST mechanics which draw the periphery downstream begin to become significant. As a result, the rim is drawn somewhat downstream during this phase and dominates the SW measurement in the early stages of the bag growth, causing apparent disagreement between the present model and experiments. Even though this partial ST occurs, the breakup ultimately persists as in MB.

Figure 28. Experimental CS and SW droplet deformation $L/d_0$ vs. time $T$ for various $We$ and morphologies (annotated in figure). Equation (4.21) is given by the solid black lines, where the solid black dots indicated the breakup point as calculated from (4.20). Every fourth data point is plotted for clarity.

The relationship between $\beta _f$ and $2R_f$ depends on the morphology of the breakup. The breakup geometries are illustrated in figure 29, showing the assumed relationships between $\beta _f$ and $2R_f$. For Bag breakup, the bag is only constrained by the rim, and since in the late stage of the bag growth the bag size matches the rim diameter, the relationship $2R_f \approx \beta _f$ is suitable (a) (see figure 7). For BS, the stamen causes a central fold of the bag, thus $2R_f \approx 2\beta _f$ (b) (see figure 8). This is similar for the MB morphology, however, the undeformed core plays a more significant role, and $2R_f \approx d_0 + 2\beta _f$ ($c$) (see figures 10 and 11). Finally, for ST, the rim of the droplet is pulled downstream right after it forms at diameter $2R_i/d_0$ (d) (see figure 13). The ratio of $2R_f$ to $\beta _f$ is plotted for all morphologies in figure 30(a). Figure 30 generally agrees with the above assumptions, showing that for Bag breakup $2R_f/\beta _f \approx 1$ and for MB $2R_f/\beta _f > 2$, however, it is clear that there is a continuous transition between the Bag and BS morphologies which is not captured here. As with the prediction of $\beta _f$, this is expected to be the result of the thickness modulations which occur due to the Bundt cake shape of the bag in the BS morphology. The prediction of $2R_f$ is shown in figure 30(b). While the prediction generally agrees with the average of each morphology, it is again clear that a derivation which takes into account the effect of the continuity in the change of the core size as a result of the growing undeformed core on $2R_f/\beta _f$ will provide better agreement across the Bag and BS morphologies.

Figure 29. Illustration of the breakup geometry for (a) Bag, (b) BS, (c) MB and (d) ST morphologies. Measurements of and relationships between $R_f$ and $\beta _f$ are indicated. The rim is marked by cross-hatching while the undeformed core is marked by grey fill.

Figure 30. (a) Rim vs. bag size at burst ($2R_f/\beta _f$) for all breakup morphologies. The solid black line gives $2R_f/\beta _f = 1$ and the dashed black line gives $2R_f/\beta _f=2$, indicating the assumed value of $2R_f/\beta _f$ for Bag and BS morphologies, respectively. (b) Rim major diameter at burst ($2R_f/d_0$) with prediction of (4.14).

The prediction of the rim expansion ratio, $R_i/R_f$, can now be substituted into (4.14), which is compared with the experiments in figure 31(a), showing good agreement with the results despite not fully capturing the continuous transitions as described previously.

Figure 31. Measured rim minor diameter (a) and rim-thinning ratio (b) vs. $We$. Solid and open markers are for side and end-on experiments, respectively. Lines give the prediction of (4.14) (a) and the median rim-thinning ratio (b).

In the interest of formulating a continuous prediction, the median rim thinning ratio, $d_r/d_i \approx 0.64$, is estimated from the experiments (figure 31b) to be used as a comparison in the next section. Note that the spread of the rim-thinning ratio in figure 31 is fairly large, thus this term will only be used as a rough estimation in order to assess the effects of the prior assumptions on the breakup prediction without the complexity of the morphology variation.

Chou & Faeth (Reference Chou and Faeth1998) performed a similar analysis using several simplifications. In their analysis, the growth of the rim was determined directly by performing the momentum balance between the rim and bag assuming the rim and bag thicknesses to be constant and the rim volume to be $V_r/V_0=0.65$. This was then fit to their results using a coefficient analogous to $C_d$, which accounts for variations in the relationship of the stagnation pressure and rim expansion. The expansion was then related to the rim thinning by conservation of mass as was done in the present work, although the terminal rim size was determined empirically. The rim was also assumed to continue growing after the breakup of the bag, with bag and rim breakup times determined empirically. The rim breakup was modelled by Rayleigh–Plateau instability, as in § 4.5. Due to these simplifying assumptions, as well as a relatively small data set with which to fit, their work was not able to predict beyond the narrow range of $15<We<20$. Essentially, the work of the present section accounts for the variations in these assumptions, allowing it to achieve good prediction across a range of $We$ and morphologies.

4.5. Rim breakup

Confounding factors such as the number and location of the rupture sites and collisions of the rupture fronts with each other and with the main rim of the droplet are numerous (see, for example, figure 32). While many of these phenomena have been studied previously (Bremond & Villermaux Reference Bremond and Villermaux2005; Lhuissier & Villermaux Reference Lhuissier and Villermaux2012; Vledouts et al. Reference Vledouts, Quinard, Vandenberghe and Villermaux2016; Poulain, Villermaux & Bourouiba Reference Poulain, Villermaux and Bourouiba2018), predictions of when each mechanism will be preferred in the breakup are still illusive. The present experiments cannot be compared directly to these works due to the limitations of the high-speed camera. The bag breakup period occurs on a very short time scale ($T < 0.2$) and at a very small spatial scale ($d_b$ of the order of $\mathrm {\mu }$m). Since it is not possible to achieve the desired time and spatial resolution while maintaining the wide field of view that is necessary to capture the entire breakup process, the present experiments are unable to provide insights into the breakup of the bags formed in droplet breakup.

Figure 32. Time-series images of the breakup of the bag showing (a) a single rupture and (b) multiple ruptures. Here, $T$ gives the time of the first frame, while ${\rm \Delta} T$ gives the incremental time between each frame; $T$ increases from left to right and continues in the bottom row.

As mentioned previously, the rims that form in droplet breakup contain much more of the droplet's mass than the bags. Additionally, much of the volume of the bag merges with the rim of the droplet once the bag recedes. The result is that most of the droplet's volume ends up being part of the rim breakup process and the droplets formed from the bag are far less important in terms of volume than those formed by the rim. It is therefore of higher importance to model the breakup of the rim.

Once the bag has completed its breakup, the rim of major and minor diameters $2R_f$ and $d_{r,f}$ corrugates and breaks into many child droplets of size $d_c$, as shown and illustrated in figure 33. The rim is given a significant initial perturbation when the receding bag collides with it, thus the instability grows and fractures the rim very quickly and the change in the diameters of the rim during its breakup can be neglected. The rim, approximated by a slender column of diameter $d_{r,f}$, is unstable by the capillary Rayleigh–Plateau instability which results in a child droplet size given by (4.22) (Ashgriz Reference Ashgriz2011)

(4.22)\begin{equation} \frac{d_c}{d_0} = 1.89 \frac{d_{f}}{d_0}. \end{equation}

Figure 33. Illustration (a) and end-on images (b) of rim breakup.

The validity of this mechanism is established by comparing its prediction of the breakup size given the experimental rim size just before its corrugation and the experimental final breakup sizes, shown in figure 34. The child droplet sizes were measured from regions of the images where droplets from the breakup of the rims could be clearly identified using the ImageJ particle size measurement tool. The reported sizes are the geometric mean of the measured child droplets for each experimental condition. Droplets from the large nodes such as those seen in figure 33 were excluded as their size does not directly relate to the final size of the rim. Figure 34 shows good agreement between the theory and the experiments given the potential measurement error for these small geometries. However, it can be seen that the prediction is somewhat higher than the measurements. This is likely the result of neglecting the formation of larger nodes on the rim during the rim and bag growth phase which would result in less mass in the ligament portions of the rim between the nodes, and thus a smaller breakup size.

Figure 34. Rayleigh–Plateau breakup prediction (4.22) based on experimental rim minor diameter vs. experimental rim child drop size. Solid and open markers are for side and end-on experiments, respectively.

The overall prediction of the model is compared to experiments in figure 35(a). As with the prediction of $d_f/d_0$ in figures 30 and 31, the breakup prediction generally agrees with the average behaviour for each breakup morphology, however, the breakup size for the Bag morphology is overpredicted. This is because the model does not account for the beading that occurs on the rim which causes large nodes to accumulate in equal intervals around the rim, which can be seen in figure 33. Zhao et al. (Reference Zhao, Liu, Li and Xu2010) suggest that these nodes are the result of a Rayleigh–Taylor instability acting on the liquid rim. The instability results in a preferential distribution of the droplet's volume near the peak of the waves and as the rim expands the segments in between the nodes are drawn thin while the nodes remain largely intact. Since some of the rim volume is stored in these nodes, the ligament portions of the rim grow thinner than the prediction as the mass conservation which led to (4.14) no longer holds. As the beading effect is much more prominent for the thicker rims, it has a much greater effect on the Bag morphology than the others.

Figure 35. Measured rim child droplet size vs. $We$. Solid and open markers are for side and end-on experiments, respectively. Lines give the prediction of (4.22).

Other causes of the overprediction of $d_c/d_0$ in figure 35(a) are the overprediction of $d_i/d_0$ from (4.10) (see figure 22), and the slight overprediction of the Rayleigh–Plateau breakup mechanism ((4.22) and figure 34). The effects of these overpredictions are illustrated in figure 35(b), where the prediction using the median rim-thinning ratio (see figure 31b) is compared to the same using a shift of $-0.05$ in $d_i/d_0$ which accounts for the overprediction of (4.10), as discussed in § 4.3. The shifted result better predicts the experiments, however, there is still a slight overprediction which is similar to that of the Rayleigh–Plateau mechanism seen in figure 34.