Non-Galilean Taylor–Culick law governs sheet dynamics in unsteady fragmentation

Abstract

We present the results of a combined experimental and theoretical investigation of sheet evolution, expansion and retraction, under unsteady fragmentation upon drop impact on a surface of comparable size to that of the drop. We quantify and model the effect of the continuous time-varying – unsteady – shedding of droplets from the sheet via its bounding rim. We present and validate especially developed advanced image processing algorithms that quantify, with high accuracy, the key quantities involved in such unsteady fragmentation, from sheet, to rim, to ligaments, to droplet properties. With these high precision measurements, we show the important effect of continuous unsteady droplet shedding on the sheet dynamics. We combine experiments and theory to derive and validate governing equations of the sheet that incorporate such continuous shedding – associated with continuous loss of momentum and mass – from unsteady fragmentation. Combining this theory with the universal unsteady rim dynamics discovered in Wang et al. (Phys. Rev. Lett., vol. 120, 2018, 204503), we show that the governing equation of the sheet can be reduced to a continuous-shedding, non-Galilean Taylor–Culick law, from which we deduce new analytical expressions for the time evolution of the sheet radius. We show the robustness of the predictions to changes of fluid properties, including surface tension and moderate fluid viscosity and elasticity, including use of physiological mucosalivary fluid. We also reconcile prior literature's inconsistent experimental results on the sheet dynamics upon drop impact.

1. Introduction and goals

Sheet dynamics, expansion and retraction, coupled with its fragmentation into secondary droplets arises in industrial, environmental and physiologic processes, from printing, spray cooling of surfaces, to delivery of nutrients or pesticides, decontamination, or disease transmission (Rein 1993; Yarin 2006; Traverso et al. 2013; Bourouiba, Dehandschoewercker & Bush 2014; Gilet & Bourouiba 2014; Laan et al. 2015; Comiskey et al. 2016; Josserand & Thoroddsen 2016; Wang & Bourouiba 2018a). A large family of such fragmentation processes are unsteady, where droplets can be generated at the very early time depending on substrate properties (Stow & Stainer 1977; Ting & Keller 1990; Xu, Zhang & Nagel 2005; Xu, Barcos & Nagel 2006) and also continuously shed with properties that vary with time from impacts of spray drops or intermittent jets on surfaces or films, for example. These include crown splash upon drop impact on dry or wet surfaces (Josserand & Thoroddsen 2016) or crescent-moon fragmentation upon drop-on-drop interaction on surfaces (Gilet & Bourouiba 2015; Wang & Bourouiba 2018a).

All these processes involve the transition of the impacting fluid bulk to a sheet, expanding and then retracting, which is typically surrounded by a rim, itself destabilizing into droplets (figure 1). Recent work (Wang & Bourouiba 2018b) showed that both size and speed distributions of the droplets ejected in unsteady fragmentation from a sheet, upon drop impact on a surface of comparable size to the drop, are shaped by the unsteadiness of the sheet dynamics. Moreover, this prior study showed that, contrary to prior discussions of the problem, droplets are shed continuously during the sheet expansion and retraction, with most shed, in fact, prior to maximum expansion of the sheet (figure 2). The fundamental questions we address in this paper are the following:

1. (i) Is the volume shed via continuous droplet shedding from the rim during sheet expansion comparable to the rim thickness surrounding the sheet?

   

   Figure 1. $(a)$ Sequence of unsteady sheet expansion from drop impact on a target of comparable size, $d_r = 6.3\ \textrm {mm}$, to that of the drop. The top image shows the expanding sheet viewed from the side at the time of maximum radial expansion. The bottom panels show the top view; the time interval is $\delta t = 2.5\ {\rm ms}$. ($b$) Definition of physical quantities.

   

   Figure 2. ($a$) Snapshot of unsteady sheet fragmentation upon drop impact on a target of comparable size, for $We = 963$, at time $t=0.2\tau _{cap}$ prior to maximum extension, reached at $t_m=0.43\tau _{cap}$. Secondary droplets are continuously shed from the expanding sheet prior to $t_m$. ($b$) Diameter distribution of the droplets shed throughout unsteady sheet fragmentation upon drop impact, for $We = 963$. The different colours represent the partition of different groups of data. The total area of the two colours is equal to unity. The colours do not overlap. This shows that most droplets are shed before the time of maximum sheet extension.

   We especially developed advanced algorithms allowing for high precision in quantification of all parts of such fragmenting system (§ 3) and show that the answer to this question is positive (§ 4).

2. (ii) How does such continuous shedding affect the mass and momentum of the unsteady sheet?

   We find that such an effect cannot be neglected when modelling the sheet dynamics (§§ 4–6).

3. (iii) How should the sheet dynamics be modelled to accurately account, in a tractable manner, for the essential physics elucidated in the prior two questions, namely the continuous unsteady fragmentation?

   Combining our recent insights on the unsteady rim self-adjustment selecting a thickness equal to the time-varying local capillary length determined by the sheet acceleration (Wang et al. 2018), with the global sheet dynamics derived herein, we derive and validate the governing equation for the sheet which incorporate unsteady continuous fluid shedding. We show that these equations can be reduced to a compact non-Galilean Taylor–Culick law (§ 5), with analytical solutions that capture the essential physics of the sheet dynamics under continuous unsteady fragmentation well (§ 7).

4. (iv) Is the prediction of the sheet dynamics under unsteady continuous fragmentation robust to changes of fluid properties? And does our prediction allow to reconcile prior literature's inconsistent results on the sheet dynamics in the air, upon drop impact?

   We show that the answers to these two last questions are also positive in § 8.

To address the above four questions, we focus on the canonical unsteady fragmentation process from drop impact on a surface of comparable size, $d_r$, to that of the impacting drop $d_0$. Upon impact on such surface, the impacting bulk fluid falling at velocity $u_0$ is transformed into a sheet radially expanding in the air (figure 1$a$) with velocity and thickness profiles $u(r,t)$ and $h(r,t)$, respectively, and with sheet radius $r_s(t)$, where $t$ is the time since impact and $r$ is the radial position of the sheet (figure 1$b$). We start with a review of prior studies on such sheet dynamics and discuss the knowledge gaps.

2. Background

2.1. Discrepancies on sheet evolution upon drop impact

Prior experimental studies aimed at rationalizing the canonical physics of drop impacts on solid surfaces (Rein 1993; Yarin 2006; Josserand & Thoroddsen 2016). The focus was on predicting the maximum radius of the sheet, $r_m$, as a function of the Weber number $We=\rho u_0^2 d_0/\sigma$ and Reynolds number $Re=u_0 d_0/\nu$, where $\rho$ is the density of the fluid, $\sigma$ its surface tension and $\nu$ its kinematic viscosity. However, discrepancies between the scaling laws proposed (Madejski 1976; Scheller & Bousfield 1995; Clanet et al. 2004; Laan et al. 2015; Lee et al. 2016; Wildeman et al. 2016) persist. A number of factors can be at the origin of such discrepancies. Among them are the inconsistent surface properties used in the different studies coupled with the difficulty in accurately capturing the associated distinct viscous stresses. To clarify this potential source of inconsistency and derive the essential physics involved, our focus is on the inviscid limit. In this inviscid limit, first, Yarin & Weiss (1995), Roisman, Berberovi & Tropea (2009) and Eggers et al. (2010) studied theoretically and numerically drop impacts on solid surfaces. Second, Lastakowski et al. (2014) and Riboux & Gordillo (2016) studied experimentally and theoretically drop impacts on superheated surfaces. The vapour layer, from the Leidenfrost effect, effectively reducing surface stresses (Wachters & Westerling 1966), mimicking an inviscid sheet expansion. Third, drop impacts on targets of comparable size to the impacting drop were used to create a sheet expanding in the air (figure 1). This last approach was used to study stationary sheets from liquid jet impact on small targets by Savart (1833) and Taylor (1959a). Rozhkov, Prunet-Foch & Vignes-Adler (2002) first investigated such drop impacts and examined the dynamics of sheet expansion in the air theoretically and experimentally. The experimental and theoretical study of Villermaux & Bossa (2011) led to the derivation of an analytic prediction of the entire time evolution of the radius $R(t)$ of the expanding sheet in the air (figure 1$b$).

Without viscous effects, Rozhkov et al. (2002), Villermaux & Bossa (2011) and Lastakowski et al. (2014) developed their respective theoretical models based on control volume analysis using mass and momentum balance at the rim. In each study, the numerical or theoretical predictions are in good agreement with the study's own experimental results. Thus, it would first appear that, in the inviscid regime, the physics of sheet expansion from inviscid drop impact is well understood. However, figure 3 compares the experimental data of the sheet radius evolution over the dimensionless time $T = t/\tau _{cap}$, which reveals large discrepancies between the above studies. Here $\tau _{cap}=\sqrt {\rho \Omega _0/{\rm \pi} \sigma } = \sqrt {\rho d_0^3/6\sigma }$ is the capillary time scale characteristic of the sheet expansion with $\Omega _0 = {\rm \pi}\rho d_0^3/6$ the volume of the impacting drop. In particular, note the non-monotonic increase of maximum radius $r_m$ with $We$. For the same values of the Weber number reported, the non-dimensionalized values of maximum sheet radius vary by up to a factor of 2 and non-dimensionalized times of maximal extension can vary by more than a factor of 1.5. Discrepancies of such magnitude are important as they lead to inconsistent predictions of coating ranges and poor control of $r_m$, of the order of $0.5\sqrt {We}d_0$ to $1.5\sqrt {We}d_0$, or secondary droplet production relevant for a wide range of important spray, electronics, printing or decontamination applications. The discrepancies between these studies reflect that ingredients are still missing from the current physical picture of sheet expansion and unsteady fragmentation upon drop impact. We discuss three ingredients next.

Figure 3. Comparison of the time evolution of the dimensionless sheet radius $R_s = r_s/d_0$ upon drop impact on surfaces from prior studies. The sheet radius, $r_s$, is non-dimensionalized by the diameter of the impacting drop $d_0$. The Weber number is $We =\rho u_0^2 d_0/\sigma$, with $u_0$ the impacting velocity of the drop, $\sigma$ the surface tension and $\rho$ the drop fluid density. Time is non-dimensionalized with the capillary time scale, $\tau _{cap}=\sqrt {\rho \Omega _0/{\rm \pi} \sigma } = \sqrt {\rho d_0^3/6\sigma }$, characteristic of the sheet expansion, with $\Omega _0 = {\rm \pi}d_0^3/6$ the volume of the impacting drop. The inset shows the non-dimensional time at which the maximum sheet radius is reached, $T_m = t_m/\tau _{cap}$, as a function of $We$ in different studies. Experimental conditions for each study are given in the table, including the target-to-drop size ratio $\eta = d_r/d_0$ (Wang & Bourouiba 2017), with $d_r$ the size of the small surface, i.e. the target.

2.2. Sheet thickness

The radius, $r_s$, of the sheet is governed by both the velocity and thickness profiles of the expanding sheet. The velocity profile was derived as $u(r,t) = u/t$ (Yarin & Weiss 1995; Roisman et al. 2009; Eggers et al. 2010; Villermaux & Bossa 2011) and verified experimentally (Lagubeau et al. 2012; Lastakowski et al. 2014; Wang & Bourouiba 2017). However, the thickness profile of the expanding sheet proposed for drop impacts on solid surfaces and on small targets differ. For a drop impact on a small target, Rozhkov, Prunet-Foch & Vignes-Adler (2002, 2004) proposed a sheet thickness profile $h(r,t) \sim t/r^3$ assuming that the impacting drop has a decreasing volume flux $q_s \sim 1/r$ during expansion, while Villermaux & Bossa (2011) proposed a sheet thickness profile analogous to that of the sheet from a steady liquid jet impacting a surface, but with time variation such that $h(r,t) \sim f(t)/r$, leading to $h(r,t) \sim 1/rt$. Wang & Bourouiba (2017) proposed and validated a unified self-similar thickness profile reconciling these two profiles

(2.1a)\begin{equation} H(R,T) = \frac{1}{T^2}F \left(\frac{R}{T} \right), \end{equation}

with

(2.1b)\begin{equation} F(X) = \frac{1}{a_1X + a_2X^2 + a_3X^3}, \end{equation}

where $a_1=24$, $a_2=-38$ and $a_3=-34$ are constants independent of $We$ derived from the initial and boundary conditions of the impact. Here, $R = r/d_0$ is the radial position in the sheet non-dimensionalized by the diameter $d_0$ of the impacting drop and $T = t/\tau _{imp}$ is the time non-dimensionalized by the impact time $\tau _{imp}=d_0/u_0$, rather than the capillary time $\tau _{cap}$ as used in figure 3. The physics underlying the difference between these two time scales is discussed later (§ 2.3). $H(R,T) = h(r,t)/d_0$ is the sheet thickness, non-dimensionalized by $d_0$, at position $R$ and time $T$. The unified profile (2.1) was shown to be robust and to capture well the experimental data of independent studies (Vernay, Ramos & Ligoure 2015; Wang & Bourouiba 2017). For drop impacts on inviscid solid surface, the thickness profile also follows a similarity profile $H = (F/{T^2})({R}/{T})$ but with a different expression of $F(X)$ from that of (2.1) (Roisman et al. 2009; Eggers et al. 2010). Such distinction of sheet thickness profile between drop impact on small targets or infinite inviscid solid surfaces, with respect to the impacting drop size, could explain the discrepancy between the results of Lastakowski et al. (2014) and those of others shown in figure 3.

2.3. Importance of multiple time scales

Physically, the self-similarity profile of the sheet thickness (2.1) is governed by the impact time scale $\tau _{imp} = d_0/u_0$, since the formation of the thickness profile of the expanding sheet is only governed by the inertia of the fluid emanating from the impacting drop close to onset of impact. However, the time evolution of the sheet radius, as discussed in detail in § 4.3, is governed by the balance of inertia and surface-tension forces acting on the rim over a longer time scale, the capillary time scale, resulting from

(2.2)\begin{equation} \rho d_0^3 \frac{d_0}{\tau_{cap}^2} \sim \sigma d_0 \quad \Longrightarrow \quad \tau_{cap} \sim \sqrt{\frac{\rho d_0^3}{ \sigma}}. \end{equation}

Here, we choose a conventional definition of the capillary time as used in prior studies (Villermaux & Bossa 2011; Vernay et al. 2015; Wang & Bourouiba 2018b; Wang et al. 2018) with

(2.3)\begin{equation} \tau_{cap} = \sqrt{\frac{\rho \Omega_0}{{\rm \pi} \sigma}} = \sqrt{\frac{\rho d_0^3}{6\sigma}} = \sqrt{\frac{We}{6}}\tau_{imp}. \end{equation}

Choosing $d_0$ and $\tau _{cap}$ (2.3), the non-dimensional variables are then

(2.4a–e)\begin{equation} R = \frac{r}{d_0}, \quad T = \frac{t}{\tau_{cap}}, \quad R_s = \frac{r_s}{d_0}, \quad H = \frac{h}{d_0} \quad \text{and} \quad U = \frac{u}{d_0/\tau_{cap}}. \end{equation}

The sheet thickness profile (2.1) expressed in terms of the dimensionless time non-dimensionalized by the capillary time scale $\tau _{cap}$ reads

(2.5)\begin{equation} H(R,T) = \frac{1}{T^2}\bar{F}\left(\frac{R}{T}\right), \end{equation}

with

(2.6)\begin{equation} \bar{F}(X) = \frac{6}{{We}}F\left(\sqrt{\frac{6}{{We}}}X\right) = \frac{\sqrt{{6We}}}{6a_3X^3 +a_2X^2\sqrt{{6We}} + a_1X{We}}, \end{equation}

where $a_1 = 24$, $a_2 = -38$ and $a_3 = 34$ are the same coefficients as used in (2.1), and $We = \rho u_0^2 d_0/\sigma$ is the impacting $We$.

2.4. Rod-to-drop size ratio, $\eta = d_r/d_0$

The sheet radius $r_s$ upon drop impact on a target of diameter $d_r$ is affected by the ratio of the rod-to-drop size, $\eta = d_r/d_0$. Wang & Bourouiba (2017) showed that when $\eta$ is too small, the sheet cannot expand horizontally. When $\eta$ is too large, surface stresses induce a sufficiently high dissipation to alter the velocity profile of the expanding sheet in the air. The authors showed that a range of $1.4 \leq \eta \leq 1.9$ is needed to obtain a horizontal sheet while also being able to neglect surface stresses. Several studies shown in figure 3 are in fact outside this range, which can explain the large discrepancies among them. Indeed, $1.4 \leq \eta \leq 1.9$ is satisfied for the experiments in Rozhkov et al. (2002), Vernay et al. (2015) and the present study, and for all three, the maximum sheet radius $r_{s}$ monotonically increases with Weber number, with a maximum extension reached at the same dimensionless time around $T_m = t_m/\tau _{cap} = 0.43$ (figure 3-inset), recalling that $\tau _{cap}=\sqrt {\rho \Omega _0/{\rm \pi} \sigma } = \sqrt {\rho d_0^3/6\sigma }$ is the capillary time scale characteristic of the sheet expansion, with $\Omega _0 = {\rm \pi}\rho d_0^3/6$ the volume of the impacting drop.

2.5. Continuous droplet shedding

Prior theoretical models predicting the radius of the inviscid expanding sheet (Eggers et al. 2010; Villermaux & Bossa 2011; Lastakowski et al. 2014; Riboux & Gordillo 2016) were developed considering a control volume on the rim and applying rim mass conservation and momentum balance. To simplify the problem, the fluid shed from the rim during expansion prior to the maximum radius was systematically neglected. The underlying assumption being that the sheet fragmentation and secondary droplet shedding occurs primarily after the maximum expansion of the sheet is reached. However, using specifically developed advanced image processing algorithms (§ 3.2), Wang & Bourouiba (2018b) showed that droplets are shed continuously throughout the sheet expansion with most shedding occurring prior to the sheet's maximum expansion, and indeed even continuing during its retraction (figure 2). Moreover, Wang et al. (2018) showed that the rim destabilizes and sheds fluid, in the form of ligaments and droplets, due to the local coupled interfacial and inertial instability driven by the continuous deceleration of the rim (Roisman 2010; Agbaglah & Deegan 2014), imposing an instantaneous self-adjustment of the rim thickness, $b$, to satisfy a robust local and instantaneous Bond number $Bo = \rho b^2 (-\ddot {r_s})/\sigma = 1$ condition. In turn, this condition has fundamental implications for the continuous volume shed from the rim over the course of the full sheet dynamics. These implications are the focus of the present paper.

3. Experimental approach and advanced image processing (AIP) algorithms

3.1. Experimental set-up

Figure 4($a$) shows the schematic diagram of our experimental set-up, with two high-speed cameras used to record the experiments from side and top views, simultaneously. The frame rate of the top-view and side-view cameras are 20 000 and 8000 frames per second, respectively. The pixel resolution of videos recorded from top and side views are ${\simeq }50 \ \mathrm {\mu }\textrm {m}\,\textrm {pixel}^{-1}$ and ${\simeq }30\ \mathrm {\mu }\textrm {m}\,\textrm {pixel}^{-1}$, respectively. An impacting drop is released from a needle from different heights selecting a range of different impact velocities $u_0$, that are measured directly using high-speed imaging. For most of this paper, we present results using de-ionized water drops with Nigrosin dye of concentration $1.2\ \textrm {g}\,\textrm {l}^{-1}$, with density $\rho =1.0\times 10^3\ \textrm {kg}\,\textrm {m}^{-3}$, surface tension $\sigma =72\times 10^{-3}\ \textrm {N}\,\textrm {m}^{-1}$, and kinematic viscosity $\nu =1.0\times 10^{-6}\,\textrm {m}^2\,\textrm {s}^{-1}$ (table 1). In the second part of this paper, other fluids are used to assess robustness of the results to change of fluid viscosity, elasticity and surface tension (table 2). The room temperature where the experiments are conducted is $21\pm 2\,^{\circ }\textrm {C}$.

Figure 4. ($a$) Schematic diagram of the experimental set-up. $(b)$ Side-view image taken prior to drop impact. ($c$) Time evolution of the diameter of an impacting drop measured by various approaches to estimate its diameter. For all approaches, the measurement error of the drop diameter is ${\approx }0.03\ \textrm {mm}$. Only a volume-based measurement is free of influences induced by the drop deformation, while the other approaches induce measurement errors of $10\,\%\unicode{x2013}25\,\%$ for the impacting drop volume; $t = 0$ is the time at which the droplet makes contact with the surface.

Table 1. Initial conditions used for the experiments of water drop impacts, including the impact drop volume-based diameter $d_0$ (figure 4), drop aspect ratio $d_x/d_y$ at the time of impact, impacting velocity $u_0$, Weber number, $We = \rho u_0^2 d_0/\sigma$, and Reynolds number, $Re = u_0d_0/\nu$, where $\rho = 1.0\times 10^3\ \textrm {kg}\,\textrm {m}^{-3}$, $\nu = 1.0\times 10^{-6}\ \textrm {m}^2\,\textrm {s}^{-1}$ and $\sigma = 72\ \textrm {mN}\,\textrm {m}^{-1}$, are the density, kinematic viscosity and surface tension of the drop, respectively. $N_{exp}$ is the number of experiments for each condition. For each drop diameter, the rod size is chosen to maintain the rod-to-drop size ratio, $\eta$, that is optimal to ensure an inviscid two-dimensional sheet horizontal formation (Wang & Bourouiba 2017).

Table 2. Properties of the fluids used in this study, with corresponding initial conditions and dimensionless Weber, $We$, Reynolds, $Re$, and Deborah, $De$, numbers. Glycerol–water mixtures were used to vary the fluid viscosity, with percentage indicating the mass fraction used. Mixtures of poly(ethylene oxide) (PEO, $1\times 10^6\ \textrm {g}\,\textrm {mol}^{-1}$) and poly(ethylene glycol) (PEG, $1\times 10^4\ \textrm {g}\,\textrm {mol}^{-1}$) were used to examine the role of fluid elasticity, quantified by the relaxation time $\tau _E$ in $De = \tau _E/\sqrt {\rho d_0^3/8\sigma }$, consistent with the definition in Wang (2018). The percentages shown indicate the mass concentration of each compound in the mixture. A naturally elastic biofluid was also used: mucosalivary fluid (the rheology of which is characterized in Bourouiba 2021). Dodecane and Dimethyl sulfoxide (DMSO) were used to vary the fluid surface tension. Six experiments per impact condition and fluid are conducted. Impact target diameter used for the first five fluids in the table is $d_r = 6.3\ \textrm {mm}$, and for the last two fluids is $d_r = 6\ \textrm {mm}$, to ensure a drop-to-target size ratio, $\eta$, between 1.5 and 1.7.

To generate the sheet, the drops impact stainless-steel cylindrical rods of different diameters $d_r$ chosen to ensure formation of a horizontal expanding sheet (Wang & Bourouiba 2017). The rod surface is polished by a 1000-Grit sand paper, associated with a root-mean-square roughness $R_{rms} < 0.25\ \mathrm {\mu }\textrm {m}$.

Each group of experiments is repeated at least six times (table 1). We ensure that the impacts are at the centre of the rod, enabling an axisymmetric sheet expansion (figure 1). The difficulty of studying fragmentation lies in its multi-scale nature and the high precision of experimental measurements required at various scales to decipher the underlying deterministic processes from fluctuations and thus gain fundamental insights on the underlying physics. A limited number of physical quantities, such as the diameter, $d_0$, of the impacting drops and the radius, $r_s$, of the expanding sheet can be measured manually accurately, but with obvious low efficiency. Thus, we next present specifically developed algorithms that allow us to measure key multi-scale critical physical quantities systematically and, to our knowledge, with unprecedented precision.

3.2. Advanced AIP algorithms and validation

3.2.1. Volume of the impacting drop

From the side view, the diameter $d_0$ of the impacting drop can be measured directly in one direction or based on the visible drop area, $A$, with $d_0 = \sqrt {4A/{\rm \pi} }$ (figure 4$b$). For spherical drops, both methods give the same diameter. However, drops naturally oscillate. Figure 4($c$) shows the time evolution of drop volume $\Omega _0= {\rm \pi}d_0^3/6$ calculated based on the drop diameter measured in one direction directly from the side view versus that computed from the measured area of the drop. The drop volume should be invariant in time, yet neither one of the two methods satisfy this invariance. To address this issue, we proceed to measure the drop diameter $d_0$ based on the measurement of its volume, $\Omega _0$, instead: $d_0=(6\Omega _0/{\rm \pi} )^{{1}/{3}}$. Indeed, despite deformation, the shape of the drop remains typically axisymmetric aligned with its trajectory and its cross-section in this direction remains spheroidal. Thus, the drop volume, $\Omega _0$, is

(3.1)\begin{equation} \Omega_0 = \int_{0}^{d_2} \frac{\rm \pi}{4}d_c^2(y) \, \mathrm{d}y, \end{equation}

where $d_c(y)$ is the diameter of the cross-section at the position $y$ and $d_2$ is the length of the drop in the vertical direction (figure 4$b$). The measurements obtained satisfy volume conservation (figure 4$c$). We use the same volume-of-revolution method to measure the volume of fluid remaining on the rod, $\Omega _c$, during the impact, and volume of the secondary droplets, $\Omega _d$, from continuous sheet fragmentation (figure 5$a$).

Figure 5. ($a$) Steps of our AIP algorithms analysing unsteady sheet fragmentation. The algorithms first separate the fluid into several parts. ($b$) Time evolution of the fluid volume in each part of the system is measured. ($c$) The sum of the volume of fluid in all parts recovers the volume of the initial impacting drop, showing the high accuracy of our AIP algorithms. The error bar on the total volume curve shows that the cumulative error of measurement for each part of the fluid is ${<}10\,\%$. The data shown are for $We = 963$.

Drop oscillation can lead to variation in aspect ratio of the drop at the time of impact. Prior studies (Durst 1996; Saylor & Grizzard 2004; Yun & Lim 2014) reported an effect of this aspect ratio on the dynamics of impact. We maintain the aspect ratio of the drop at the time of impact on the surface for all groups of experiments within the range $d_x/d_y = 1 \pm 0.05$ (table 1). The convergence of the evolution of the diameter for different measurement methods in figure 4($b$) indicates that the impacting drop becomes spherical at the time of impact. Occasionally, the aspect ratio $d_x/d_y$ deviates from unity by $15\,\%$, for which the variation of the maximum sheet radius is around $3\,\%\unicode{x2013}5\,\%$.

3.2.2. Tracking and total volume conservation

We derived and developed multi-step AIP algorithms to capture all key physical quantities of the sheet, rim, ligaments and secondary droplets listed in figure 1($b$) despite their change in morphology and topology with time. For a given frame, one algorithm first removes the background of the image to isolate the fluid fragmenting (figure 5$a$i) and then detects the inner and outer contours of the rim and ligaments, and the expanding sheet (figure 5$a$ii). Between the inner and outer contours, the rim-ligament connection is detected (figure 5$a$iii). Outside the outer contour, the ejected secondary droplets are also detected (figure 5($a$iv) and § 3.2.1).

For the expanding sheet, the radius of the sheet $r_s$ is captured by measuring the total area $A_s$ within the inner contour and calculating $r_s =\sqrt {A_s/{\rm \pi} }$ (figures 5($a$ii) and 5($a$vi)). The thickness profile $h(r,t)$ of the sheet is measured by light-absorption as shown in Wang & Bourouiba (2017). With the radius $r_s$ and thickness profile $h(r,t)$, the volume of fluid in the sheet is

(3.2)\begin{equation} \Omega_s(t) = \int_{d_r/2}^{r_s(t)} 2{\rm \pi} h(r,t)r\, \mathrm{d}r. \end{equation}

For the rim–ligament connection, the extruded parts are the ligaments, which can be systematically and precisely separated from the rim by our algorithms (figure 5$a$vii). Upon separation, the thickness of the rim, $b$, is measured locally at each angular position and averaged along the entire rim. The length of each ligament $\ell$ is measured as the distance from its root to its tip. The width of each ligament $w$ is measured locally at each cross-section and averaged along the centreline of the ligament (figure 5$a$vii). Using cylindrical coordinates, the volume of fluid in the rim, $\Omega _r$, and total ligaments, $\Omega _{\ell }$, can be evaluated as

(3.3a,b)\begin{equation} \Omega_r(t) = \int_{0}^{2{\rm \pi}} \frac{\rm \pi}{4} b^2(\theta,t) R(t)\, \mathrm{d}\theta \quad \text{and} \quad \Omega_{\ell}(t) = \sum_{n=1}^{N(t)} \int_{0}^{\ell_n(t)} \frac{\rm \pi}{4} w_n^2(\ell,t) \, \mathrm{d}\ell, \end{equation}

where $b(\theta ,t)$ is the local thickness of the rim, $N(t)$ is the number of ligaments, along the rim, $\ell _n(t)$ is the major-axis length of the $n$th ligament and $w_n(\ell ,t)$ is the local width of the $n$th ligament along the major axis at time $t$. The major axis of each ligament is computed from the ligament's contour detected by our AIP algorithms (figure 5$a$vii). The tip of each ligament detected at each frame is linked to the location of the tip in another frame with the ligament-tracking algorithm described in Wang & Bourouiba (2018b), enabling us to track the evolution of each ligament throughout the sheet expansion.

The position of the secondary droplets is also linked at each frame, by our AIP algorithms (Wang et al. 2018), to construct their trajectory (figure 5$a$viii). Hence, the number of trajectories represents accurately the number of droplets ejected throughout the sheet dynamics undergoing continuous fragmentation. The travelling speed of a droplet is measured by taking the difference in its position along its own trajectory. The droplet ejection speed is defined as the initial speed along its trajectory, measured at the time of shedding at which an associated change in the length of the ligament of origin of the droplet is detected. The volume of each droplet, $\Omega _{d}$, is measured using (3.1) along its trajectory, with the cumulative droplet volume shed at time $t$, $\Omega _{d}$, being

(3.4)\begin{equation} \Omega_d(t) = \sum_{m=1}^{M(t)}\Omega_{dm}, \end{equation}

with $\Omega _{dm}$ the volume of the $m$th ejected droplet and $M(t)$ the number of droplets ejected by time $t$. Here, we note that using (3.1) to measure the secondary droplet volume is based on the fact that secondary droplets oscillate along their trajectory and their cross-section in axisymmetric in that direction. Wang & Bourouiba (2018b) verified experimentally that ${>}10\,\%$ of droplets shed from the sheet could have other degrees of motion due to the complexity of the necking. Thus they have a small contribution to the estimation of total volume ejected. In addition, Thoroddsen, Takehara & Etoh (2012) observed very tiny and fast-moving micro-splashing droplets of size $5\unicode{x2013}20\ \mathrm {\mu }\textrm {m}$ at the very early stage of impact when the impacting Reynolds number $Re > 14\,000$. Such tiny micro-splashing droplets are smaller than our pixel resolution (§ 3.1). We do not capture these in this study. Their total volume is negligible with respect to the total volume of secondary droplets. Their speed is, however, very high, but their deceleration induced by the air drag on their motion, limiting their range, is also very high.

Finally, figure 5($c$) shows the time evolution of the fluid volume in each of the different parts of the system discussed above, clearly showing that our AIP algorithms ensure both cumulative and instantaneous volume conservation. Indeed, all parts that are detected and tracked sum up, at all times, to the initial volume of the impacting drop $\Omega _0$. This conservation shows the high accuracy of our AIP algorithms in tracking and extracting key features of the system despite their time-varying complex morphology, such as merging and bending of ligaments. With this high precision in data collection, we can now tackle the analysis of the quantitative results obtained.

4. Shedding throughout the sheet expansion

4.1. Importance of continuous shedding

Using the AIP algorithms discussed and validated in § 3.2, figure 6($b$) shows the experimental measurement of the volume of the rim $\Omega _r$ compared to the volume shed from the rim as a function of time. The cumulative volume shed from the rim at time $t$ includes the volume of the ligaments $\Omega _{\ell }$ at that time and the cumulative volume of droplets $\Omega _d$ ejected up to that time (figure 6). Clearly, the volume shed from the rim is of the same order of magnitude as that of the rim throughout the entire sheet dynamics. Thus, the volume shed from the rim cannot be neglected when modelling the sheet evolution. This is true for all $We$ (figure 6$b$–$d$). In this section, we discuss further the implications of this robust finding for the dynamics of the sheet expansion. All critical quantities involved in the physical models hereafter are examined both experimentally and theoretically.

Figure 6. $(a)$ Sheet fragmentation with (red) the detected fluid shed from the rim isolated from the fluid in the rim (blue). The fluid shed from the rim is either transferred into ligaments or secondary droplets. Scale bar is 4 mm. ($b\unicode{x2013}d$) Time evolution of the volume in the rim compared to the cumulative volume shed from the rim, including the volume in all ligaments and the cumulative volume of all ejected droplets, for three different $We$. All quantities are non-dimensionalized. For all different $We$, both volumes have the same order of magnitude throughout the sheet dynamics, indicating that the fluid shed by the rim cannot be neglected when modelling the sheet dynamics.

4.2. Validity of the model of sheet dynamics without shedding from the rim

Upon impact on the surface, the drop is transformed into an expanding sheet bounded by a rim (figure 1). The fluid of the drop is first transferred into the sheet and then accumulated into the rim. Mass conservation in the sheet, without the rim, reads

(4.1a)\begin{equation} \int_{r_0}^{r_s(t)} 2{\rm \pi} rh(r,t) \, \text{d}r=\int_{t_0}^{t} 2{\rm \pi} q_s(r_0,t) \, \text{d}t -\int_{t_0}^t 2{\rm \pi} q_{in}(t)\, \text{d}t, \end{equation}

with

(4.1b,c)\begin{equation} q_s(r,t) = u(r,t)h(r,t)r \quad \text{and} \quad q_{in}(t) = h(r_s,t)r_s(t)(u(r_s,t)-\dot{r_s}), \end{equation}

where $u(r,t)$ and $h(r,t)$ are the velocity and thickness profiles of the sheet, respectively (figure 1$b$); $r_s(t)$ is the radius of the sheet and $\dot {r_s}$ is its time derivative, which is the sheet expanding speed, as well as the rim velocity. The exact, and complex, initial condition of sheet formation and expansion at early time is inaccessible; hence, here $t_0$ does not refer to the time of first contact between the drop and the surface but to the time after which a sheet forms. A good choice is $t_0 = \tau _{imp}$, where $\tau _{imp} = d_0/u_0$ is the drop impact time scale, with $r_0 = r_s(t_0)$. $q_s(r,t)$ is the volume per radian and unit time entering the sheet at radial position $r$ and time $t$. $q_{in}(t)$ is the volume per radian and unit of time entering the rim from the sheet.

Note that with the velocity profile $u(r,t) = r/t$ and the self-similar thickness profile (2.1) of the sheet, (4.1) holds identically. Indeed, the derivation of the sheet profiles (Wang & Bourouiba 2017) already accounted for mass conservation (4.1). Considering mass conservation at the rim and neglecting shedding from the rim, rim mass conservation per radian reads

(4.2a)\begin{equation} \frac{\textrm{d}}{\textrm{d}t}(ar_s) = q_{in}, \end{equation}

with

(4.2b)\begin{equation} a = \frac{{\rm \pi} b^2}{4}, \end{equation}

where $a$ is the cross-section area of the rim and the volume rate entering the rim, $q_{in}(t)$, is given by (4.1). We also assume the cross-section of the rim to be circular, and the thickness (diameter) of the rim, $b$, to be much smaller than the diameter of the sheet $2r_s$ (figure 1$b$). Momentum balance on the rim gives the total momentum change in the rim equal to the net momentum influx into the rim in addition to the forces exerted on the rim. The only force acting on the rim in the direction of its radial motion is the restoring surface tension force exerted by the sheet. Momentum balance, per radian, on the control volume following the rim then reads

(4.3)\begin{equation} \frac{\textrm{d}}{\textrm{d}t}\left(\rho ar_s\cdot\dot{r_s}\right) = \rho q_{in}(t)u(r_s,t)- 2\sigma r_s, \end{equation}

where the dot $\langle \ \cdot \ \rangle$ indicates the time derivative of the variable. Combining (4.1)–(4.3) leads to the governing equation

(4.4)\begin{equation} \underbrace{\rho(ar_s)\ddot{r_s}}_{\textit{inertia}} - \underbrace{\rho h(r_s,t)r_s(u(r_s,t) - \dot{r_s})^2}_{\textit{momentum influx}} \quad +\underbrace{2\sigma r_s}_{\textit{surface tension}} = 0, \end{equation}

with

(4.5)\begin{equation} ar_s = \int_{t_0}^t q_{in}(t)\, \text{d}t + a_0r_0 = \int_{t_0}^{t} q_s(r_0,t) \, \text{d}t - \int_{r_0}^{r_s} rh(r,t) \, \text{d}r + a_0r_0, \end{equation}

where $a_0$ is the initial cross-section area of the rim at initial time $t_0$. The above integro-differential equations are difficult to solve analytically. To simplify, the common assumption made in prior studies (Roisman, Rioboo & Tropea 2002; Villermaux & Bossa 2011) was that the inertial effects associated with the acceleration are small and negligible compared to the momentum input from the sheet into the rim. This assumption leads to the elimination of the highest-order derivative term in (4.4), thus, simplifying it to

(4.6)\begin{equation} - \rho h(r_s,t)(u(r_s,t) - \dot{r_s})^2 + 2\sigma = 0 \quad \Longrightarrow \quad u(r_s,t) - \dot{r_s} = \sqrt{\frac{2\sigma}{\rho h(r_s,t)}}, \end{equation}

which is the Taylor–Culick speed (Taylor 1959b; Culick 1960). This would mean that the difference between the velocity of fluid entering the rim $u(r_s,t)$ and the velocity of the rim $\dot {r_s}$ is fully determined by the sheet thickness at the rim $h(r_s,t)$. Similar to the absence-of-shedding assumption, this assumption of negligible inertial effects in the rim has not been verified experimentally. With our algorithms discussed in § 3.2, we measure the quantities involved experimentally.

Figure 7($a$) shows the measured time evolution of the absolute value of the inertial term associated with the acceleration $\rho (ar_s)\ddot {r_s}$ compared to the momentum-influx term $\rho h(r_s,t)r_s(u(r_s,t)-\dot {r_s})^2$ and the surface-tension term $2\sigma r_s$. Although the inertial term is smaller, the three terms remain of the same order of magnitude throughout the entire sheet dynamics. Figure 7($b$) shows that the ratio between the absolute value of the inertia term and the surface-tension force term is approximately 0.4. Hence, the inertia of the rim associated with the rim acceleration is non-negligible. Given that the acceleration is changing with time, the rim dynamics would thus have to be considered in a non-Galilean frame of reference. Moreover, figure 7($b$) shows that the ratio between the inertia and surface-tension terms remains constant throughout the entire sheet evolution, and this is robust to the change of $We$. We explain such constant ratio theoretically in subsequent sections (§ 5.1). It is important to note that the measured inertia term $(\rho ar_s)\ddot {r_s}$ from our experiments accounts for the fluid shed from the rim. In other words, if no shedding was occurring from the rim, more fluid would accumulate in it and the inertial effect would be even greater. Thus, our measurements do not support the assumption of absence of shedding prior to the maximum sheet expansion. They also do not support neglecting the inertia of the rim, the assumption of which leads to a Taylor–Culick dynamics, (4.6).

Figure 7. ($a$) Measured time evolution of the three terms of the rim momentum balance (4.4), for $We = 679$. Each term is non-dimensionalized by $\rho d_0^4/\tau _{cap}^2 = 6\sigma d_0$. ($b$) Time evolution of the ratio of surface-tension force to inertia term in (4.4), for a range of $We$ values. This term remains constant throughout the sheet dynamics, and is robust for the change of $We$. This result is further discussed in § 5.1, with the experimental value of the ratio being in good agreement with the theoretical prediction ${\rm \pi} /8$ (5.2) derived from the local rim constraint $Bo = 1$ from (4.14) (Wang et al. 2018).

4.3. Rim dynamics incorporating continuous fluid shedding and inertial effects

We proceed to develop and validate a rim theoretical model that enables us to incorporate both the non-negligible shedding (figure 6) and non-negligible inertial effects on the rim (figure 7), while ensuring mass conservation for the full system, starting with

(4.7)\begin{equation} \frac{\textrm{d}}{\textrm{d}t}(ar_s) = q_{in}-q_{out}, \end{equation}

where $q_{out}$ is the volume per radian and unit time leaving the rim. The shedding from the rim also induces momentum leaving the rim per unit of time. Momentum balance with a control volume following the rim reads

(4.8)\begin{equation} \frac{\textrm{d}}{\textrm{d}t}\left(\rho ar_s\cdot\dot{r_s}\right)= \rho q_{in}u(r_s,t)- \rho q_{out}v_{out}-2\sigma r_s, \end{equation}

where $v_{out}$ is the speed of fluid leaving the rim. Equations (4.8) and (4.3) differ by the term $\rho q_{out}v_{out}$ representing the momentum leaving the rim. Combining (4.7) and (4.8), the modified governing equation for the sheet expansion becomes

(4.9)\begin{equation} (\rho ar_s)\ddot{r_s} - \rho q_{in}(u(r_s,t) - \dot{r_s}) + \rho q_{out}v_\ell + 2\sigma r_s = 0, \end{equation}

where $v_\ell = v_{out} - \dot {r_s}$ is the speed of fluid leaving the rim in the reference frame of the moving rim. The combination of (4.1), (4.7) and (4.9) describes the sheet expansion dynamics, in which two quantities remain unknown: the rate of volume loss from the rim, $q_{out}$, and its associated fluid outward speed $v_{\ell }$.

Clearly the time evolution of the expanding sheet (4.9) results from the balance of the inertia and surface tension acting on the rim, confirming that the time evolution of the sheet radius is indeed governed by the capillary time scale $\tau _{cap}$ (2.3), rather than the impact time scale $\tau _{imp} = d_0/u_0$. For the remainder of this paper, unless specified otherwise, we use lower case letters for dimensional variables and capital letters for variables non-dimensionalized using $d_0$ and $\tau _{cap}$.

4.4. Volume rate leaving the rim: $q_{out}$

We now determine the rate of fluid volume leaving the rim, $q_{out}$. The rim mass conservation enables us to link $q_{out}(t)$ to the rate of fluid volume entering the rim from the sheet, $q_{in}(t)$. Using (4.1) for $q_{in}(t)$, the cumulative volume $\Omega _{in}(t)$ entering the rim is

(4.10)\begin{equation} \Omega_{in}(t) = \Omega_{in0} + \int_{t_0}^{t} 2{\rm \pi} h(r_s,t)r_s(t)(u(r_s,t)-\dot{r_s})\,\text{d}t, \end{equation}

where $\Omega _{in0}$ is the volume of fluid that enters the rim prior to $t=t_0$. Non-dimensionalizing by the initial volume of the impacting drop $\Omega _0 = {\rm \pi}d_0^3/6$, leads to $\mathcal {V}_{in} = \Omega _{in}/\Omega _0$ which reads

(4.11)\begin{equation} \mathcal{V}_{in}(T) = \mathcal{V}_{in0} + 12\int^T_{T_0} H(R_s,T)R_s(T)(U(R_s,T)-\dot{R_s}(T)) \, \text{d} T, \end{equation}

where $\mathcal {V}_{in0} = \Omega _{in0}/\Omega _0$ is the dimensionless volume entering the rim prior to $T_0$, and the other non-dimensional variables are defined in (2.4a--e). Using the sheet velocity profile $U(R_s,T)=R_s(T)/T$, and the self-similar sheet thickness profile (2.5), we can re-write the expression of the cumulative volume $\mathcal {V}_{in}$ as

(4.12)\begin{equation} \mathcal{V}_{in}(T) = \mathcal{V}_{in0} + 12\int^{T}_{T_0} - \frac{R_s(T)}{T}\bar{F}\left(\frac{R_s(T)}{T}\right) \, \text{d} \left(\frac{R_s(T)}{T}\right). \end{equation}

By changing the integral variable $U(T) = R_s(T)/T$ and using the expression of $\bar {F}(X)$ (2.5), we can rewrite the expression of the cumulative volume $\mathcal {V}_{in}$ entering the rim as

(4.13)\begin{equation} \mathcal{V}_{in}(T) = \mathcal{V}_{in0} + 12\int^{{R_0}/{T_0}}_{{R_s(T)}/{T}} \frac{\sqrt{{6We}}}{6a_3U^2 +a_2U\sqrt{{6We}} + a_1{We}} \, \mathrm{d}U. \end{equation}

Experimentally, the time evolution of $\Omega _{in}(t)$ can be measured by summing the volume in the rim $\Omega _r(t)$, the ligaments $\Omega _{\ell }(t)$ and the cumulative volume in the ejected droplets $\Omega _d(t)$ using (3.3a--e) and (3.4). Figure 8($b$) shows the experimental measurement of the cumulative volume entering the rim $\Omega _{in}$ as a function of time. The estimation of $\Omega _{in}$ using (4.13) based on the experimental data of the sheet radius $r_s$ matches the experimental data very well for all $We$ (figure 8$b$, inset). Such a match further supports the accuracy of the sheet velocity and thickness profiles (Wang & Bourouiba 2017) and the theoretical expression of the volume influx (4.13). Thus, we can use these validated predictions to pursue our derivation of the sheet evolution.

Figure 8. ($a$) Schematic diagram of the sheet. Experimental measurements of ($b$) the cumulative fluid volume $\Omega _{in}$ entering the rim from the sheet compared with the prediction (4.13), of ($c$) the rim thickness $b$ compared with the prediction (4.14) and of ($d$) the cumulative volume $\Omega _{out}$ shed from the rim compared with the prediction (4.18). The solid lines are the theoretical predictions and the main graphs are shown for $We = 679$. The insets in (b) and (d) show the ratio of the experimental values with the theoretical predictions of the corresponding quantity for all $We$. All ratios remain equal to one, indicating good agreement between predictions and measurements, and robustness to change of $We$ values.

As described in § 2.5, the shedding of fluid from the rim is due to rim destabilization, governed by the local instantaneous Bond number of the rim remaining equal to one (Wang et al. 2018), inducing a self-adjustment that maintains the rim thickness $b(t)$ equal to the local and instantaneous capillary length $l_c(t)$, with

(4.14)\begin{equation} Bo = \frac{\rho b^2 (-\ddot{r_s})}{\sigma}=1 \quad \Longrightarrow \quad b = l_c = \sqrt{\frac{\sigma}{\rho (-\ddot{r_s})}}. \end{equation}

Figure 8($c$) shows a very good match, and robustness with respect to change of $We$, of the time evolution of the measured rim thickness with the local instantaneous capillary length $l_c$ using (4.14) based on the experimental data of the rim acceleration, as expected from (Wang et al. 2018). Using (2.4a--e), the rim thickness in dimensionless form reads

(4.15)\begin{equation} B = \frac{b}{d_0} = \sqrt{\frac{\sigma}{\rho d_0/\tau_{cap}^2 (-\ddot{R_s})}} = (-6\ddot{R_s})^{-1/2}. \end{equation}

Taking a cylindrical coordinate system, the volume of the rim can be expressed as

(4.16)\begin{equation} \Omega_{r}(t) = 2{\rm \pi} r_s \left(\frac{\rm \pi}{4}b^2\right) = 2{\rm \pi} d_0^3 \frac{{\rm \pi} R_s}{24 (-\ddot{R_s})}. \end{equation}

Similarly, non-dimensionalized by the initial volume of the impacting drop, $\Omega _0 = {\rm \pi}d_0^3/6$, the dimensionless rim volume $\mathcal {V}_r$ reads

(4.17)\begin{equation} \mathcal{V}_{r}(T) = \frac{\Omega_{r}(T)}{\Omega_0} = \frac{{\rm \pi} R_s}{2 (-\ddot{R_s})}. \end{equation}

The rim mass conservation (4.7) imposes that the cumulative volume of fluid shed from the rim is $\Omega _{out}(t) = \Omega _{in}(t) - \Omega _r(t)$. Using (4.13) and (4.17), the dimensionless cumulative volume $\mathcal {V}_{out}$ shed from the rim reads

(4.18)\begin{gather} \mathcal{V}_{out}(T) = \frac{\Omega_{out}}{\Omega_0} = \mathcal{V}_{in0} + 12\int^{{R_0}/{T_0}}_{{R_s(T)}/{T}} \frac{\sqrt{{6We}}}{6a_3U^2 +a_2U\sqrt{{6We}} + a_1{We}} \, \mathrm{d}U - \frac{{\rm \pi} R_s}{2 (-\ddot{R_s})}. \end{gather}

Similar to $\Omega _{in}$, the time evolution of $\Omega _{out}(t)$ can be measured experimentally by our AIP algorithms, by taking the sum of the volume in the ligaments $\Omega _{\ell }(t)$ and the cumulative volume in the ejected droplets $\Omega _{out}(t)$. Figure 8($d$) shows the experimental measurement of the cumulative volume $\Omega _{out}(t)$ shed from the rim compared to its estimation from (4.18). The estimation matches the experimental data very well, and is robust to changes of $We$, indicating that the fluid shed from the rim is captured entirely by the criterion $Bo = 1$ (4.14). In other words, using (4.1) and (4.7), the volume per radian and unit of time shed by the rim is

(4.19)\begin{equation} q_{out}(t) = h(r_s,t)r_s(t)(u(r_s,t)-\dot{r_s}) - \frac{\textrm{d}}{\textrm{d}t}\left(\frac{\rm \pi}{4}b^2r_s \right). \end{equation}

Using the length scale $d_0$ and time scale $\tau _{cap}$, we choose the volume rate to scale as $d_0^3/\tau _{cap}$ leading to the dimensionless volume rate per radian

(4.20)\begin{equation} Q_{out}(T) = \frac{q_{out}}{d_0^3/\tau_{cap}} = \frac{R_s}{T^2}\bar{F}\left(\frac{R_s}{T}\right)\left(\frac{R_s}{T}-\dot{R_s}\right) - \frac{\textrm{d}}{\textrm{d}T}\left(\frac{\rm \pi}{4}B^2R_s \right), \end{equation}

expressed as a function of the variables (2.4a--e). Note that the same expression can also be obtained by taking the derivative of (4.18) and recalling that $\mathcal {V}_{out}$ is non-dimensionalized by the initial volume of the drop $\Omega _0 = {\rm \pi}d_0^3/6$, while the $Q_{out}(T)$ is the volume outward per radian and time, non-dimensionalized by $d_0^3/\tau _{cap}$. Thus,

(4.21)\begin{equation} Q_{out}(T) = \frac{{\rm \pi}/6}{2{\rm \pi}} \frac{\textrm{d} \mathcal{V}_{out}}{\textrm{d}T} = \frac{1}{12}\frac{\textrm{d} \mathcal{V}_{out}}{\textrm{d}T}. \end{equation}

With the volume leaving the rim per unit of time determined, we can now turn to determine the associated fluid outward speed $v_{\ell }$.

4.5. Speed of fluid leaving the rim: $v_{\ell }$

The mechanism that governs the fluid emanating from the rim into ligaments is discussed in detail in a separate study (Wang & Bourouiba 2021). In the non-inertial reference frame of the rim (figure 9$a$), which decelerates, a fluid parcel is subject to the fictitious force $\rho \ddot {r_s}$, inducing a fluid flow from the junction – the region at the foot of the ligament – into the ligament. This outward speed $v_\ell$ is estimated using a local analysis involving Bernoulli applied between point (1), at the centre of the rim, and (2), at the centre of the ligament (figure 9$a$). Taking the centreline of the rim as the reference level for this body force, the fluid in the rim has a zero velocity at point (1), the rim-curvature-induced pressure is thus $2\sigma /b$, with the potential of the fictitious force being null. At point (2), the outward speed is $v_\ell$, the ligament-curvature-induced pressure is $2\sigma /w$, while the potential of the fictitious force is $\rho (-\ddot {r_s})(r_c + b/2)$, where $r_c$ is the effective radius of curvature at the rim–ligament junction. The Bernoulli equation applied between points (1) and (2) then reads

(4.22)\begin{equation} 0 + \frac{2\sigma}{b} = \frac{1}{2}\rho v_\ell^2 + \frac{2\sigma}{w} - \rho\left(-\ddot{r_s}\right)\left(r_c+\frac{b}{2}\right). \end{equation}

In a separate investigation (Wang & Bourouiba 2021) we conducted a detailed and systematic study of isolated ligaments along the rim for different impact $We$ and found that the effective radius of curvature at the rim-ligament junction, $r_c$, is in fact close to the width of the ligament $w$. The outward speed is thus

(4.23)\begin{equation} v_\ell = \sqrt{(2w+b)(-\ddot{r_s}) + \frac{4\sigma}{\rho}\left(\frac{1}{b}-\frac{1}{w}\right)}. \end{equation}

The experimental measurement of $v_\ell$, at each time, is indirect using $v_\ell (t) = q_{e\ell }(t)/a_\ell (t)$, where $q_{el}(t)$ is the population average volume leaving the rim to enter each ligament, per unit of time. This is obtained from the derivative of the cumulative volume shed from the rim into each ligament measured by our AIP algorithms (§ 3). Here, $a_\ell (t)={\rm \pi} w^2(t)/4$ is the cross-sectional area of each ligament averaged along the ligament (figure 9$b$). Figure 9($c$) shows the experimental measurement of outward speed $v_\ell$ as a function of time for $We = 679$, compared with the estimation (4.23), which are in good agreement, and robust to change of impact $We$ (figure 9$c$, inset).

Figure 9. ($a$) Schematic of the junction between the rim and a ligament as shown in the inset of ($c$). ($b$) Quantities used to compute the experimental value of the outward velocity $v_\ell$ from the rim to the ligaments. ($c$) Experimental data of the outward velocity $v_{\ell }$ as a function of time for $We = 963$. The average experimental outward velocity of fluid shed from the rim is measured from the average volume rate entering each ligament, $q_{e\ell }$, divided by the average ligament cross-sectional area $a_{\ell }= {\rm \pi}w^2/4$. The cumulative measurement error of $v_{\ell }$ is around 15 %. The solid line is the prediction of the outward velocity $v_\ell$ from (4.23), which is in good agreement with experimental measurements. The upper-right inset shows the ratios between the theoretical prediction and the data of $v_{\ell }$ for different $We$, all remaining equal to unity.

Here, we note that the effect of aerodynamic forces on the sheet expansion can be estimated. The pressure force $F \sim \rho _{a} (\dot {r_s})^2 b$, where $\rho _a$ is the air density, acting on the rim-ligament system during the sheet evolution is ${<}4\,\%$ of the capillary force, and is comparable to our measurement error. Thus, we neglect aerodynamic forces in our analysis.

5. Unified closed-form theory of sheet dynamics under continuous unsteady fragmentation

5.1. The $Bo = 1$ criterion: balance of inertial and surface-tension forces and first system reduction

With (4.1) and (4.9), and the expression of quantities $q_{out}$ (4.19) and $v_\ell$ (4.23), we now have a closed system governing the sheet evolution, incorporating inertial and shedding effects. Using the local rim condition $Bo = 1$ (4.14), the first term in (4.9) is re-written as

(5.1)\begin{equation} (\rho ar_s) \ddot{r_s} = \left(\rho \frac{\rm \pi}{4}b^2 r_s \right) \ddot{r_s} = -\frac{\rm \pi}{4} \sigma r_s, \end{equation}

which is analogous to a surface-tension force. Thus, the ratio between the absolute value of the inertia and the surface-tension force terms in (4.9) is constant throughout the sheet dynamics and independent of the impact $We$, with

(5.2)\begin{equation} \frac{|(\rho ar_s ) \ddot{r_s}|}{|2\sigma r_s|} = \frac{\rm \pi}{8} \approx 0.39. \end{equation}

This constant is in excellent agreement with the experimental measurements of this ratio (figure 7$b$), confirming that the inertia associated with the rim acceleration is, indeed, non-negligible.

More importantly, the $Bo = 1$ criterion means that the local instantaneous inertial force acting on the rim is always equivalent, in magnitude, to the instantaneous local capillary force. Remarkably, this balance enables us to reduce the complexity of the inertial term, a second-order derivative of $r_s$ and the cross-sectional area of the rim $a$ which would require taking the integral of (4.7). Instead, the equivalent surface-tension term is linear in $r_s$, thus dramatically, reducing the complexity of (4.9) while ensuring that both inertial and continuous-shedding effects remain captured. The reduced system reads

(5.3)\begin{equation} \underbrace{-\,\rho q_{in}\left(\underbrace{u(r_s,t) - \dot{r_s}}_{v_{in}}\right)}_{p_{in}}+ \underbrace{\rho q_{out}v_\ell}_{p_{out}} + \left(2-\frac{\rm \pi}{4}\right)\sigma r_s = 0. \end{equation}

Recall that

(5.4a,b)\begin{equation} q_{in}=\rho h(r_s,t)r_s\left(\frac{r_s}{t} -\dot{r_s}\right) \quad \text{and} \quad u(r_s,t) = \frac{r_s}{t}, \end{equation}

which are functions of sheet radius $r_s$ and time $t$ only. Thus, the only remaining term to close the reduced system governing the sheet evolution is the momentum leaving the rim per unit of time $\rho q_{out}v_\ell$, tackled next.

5.2. Closed system of equations governing the sheet evolution, $r_s(t)$

The expressions of $q_{out}$ and $v_\ell$ are given by (4.21) and (4.23), respectively, as

(5.5a,b)\begin{equation} q_{out} = q_{in} - \frac{\rm \pi}{4}\frac{\textrm{d}}{\textrm{d}t}\left(b^2 r_s\right) \quad \text{and} \quad v_\ell = \sqrt{(2w+b)(-\ddot{r_s}) + \frac{4\sigma}{\rho}\left(\frac{1}{b}-\frac{1}{w}\right)}, \end{equation}

where $q_{in}$ is the volume entering the rim per unit of radian and time (4.1) and $b =\sqrt {{\sigma }/\rho (-\ddot {r_s})}$ is the rim thickness determined by $Bo = 1$ (4.14). Thus, the only unknown quantity in (5.3) is the ligament width, $w$.

Figure 10($a$) shows the time evolution of ligament width, $w$, compared to rim thickness, $b$, with inset showing that the ratio $w/b$ is approximately 1–1.4 throughout the sheet evolution, and for all different $We$. Thus $w$ and $b$ are considered approximately equal, at first order. Using $w\approx b$ and $Bo = 1$, $v_\ell$ simplifies to

(5.6)\begin{equation} v_\ell = \sqrt{3(-\ddot{r_s})b} = \sqrt{\frac{3\sigma}{\rho b}}. \end{equation}

We quantify the error induced by the approximation $w\approx b$ on the sheet dynamics, by comparing the magnitudes of the rates of momentum inward and outward in (5.3). Figure 10($b$) shows the time evolution of the measured fluid speed from the rim to the ligaments, $v_\ell$, compared with the fluid speed from the sheet into the rim, $v_{in} = u(r_s,t)-\dot {r_s}$, both shown in the rim reference frame: $v_\ell$ is systematically smaller than $v_{in}$ during the entire sheet dynamics, with $v_{in}/v_{out} \approx 4$ and constant. Figure 10($c$) shows the comparison between the volume shed from the rim per unit of time and radian, $q_{in}$, and the volume entering the rim per unit of time and radian, $q_{out}$: $q_{in}/q_{out}$ ranges from 1 to 3 (figure 10$c$, inset). Thus, the rim momentum loss, $p_{out} = \rho q_{out}v_\ell$, is systematically smaller than the momentum gain, $p_{in} = \rho q_{in}(u(r_s,t) - \dot {r_s})$, in (5.3) (figure 10$d$), with an instantaneous ratio $p_{in}/p_{out} = 3-10$, and a time-averaged ratio $\approx 5$ over the course of the sheet's life (figure 10$d$, inset). Thus, recalling that $w/b = 1\unicode{x2013}1.4$ (figure 10$a$, inset), the approximation $w \approx b$ introduces a time-averaged error of approximately $10\,\%$ on $v_\ell$ (4.23). With the time-averaged ratio $p_{in}/p_{out} \approx 5$ in (5.3), the error on the prediction of $r_s(t)$ introduced by the approximation $w\approx b$ and plugging the expression of $v_\ell$ (5.5a,b) into $p_{out}$ in (5.3) is $0.1\times (1/5) = 0.02$, namely a negligible error of $2\,\%$.

Figure 10. ($a$) Comparison of the measured time evolution of the width of a ligament, $w$, and rim thickness, $b$, for $We = 679$, non-dimensionalized by $d_0$. The ratio $w/b$ is mostly invariant to change of $We$ (inset). ($b$) Comparison of measured time evolution of fluid speed $v_{in} = u(r_s,t)-\dot {r_s}$ entering the rim and outward speed $v_\ell$ leaving the rim, in the rim reference frame and for $We = 679$. Both are non-dimensionalized by $d_0/\tau _{cap}$. The ratio $v_{in}/v_{out}$ remains close to 4 throughout the sheet evolution (inset). ($c$) Comparison of measured time evolution of the volume rate $q_{in}$ entering the rim and the volume rate $q_{out}$ leaving the rim in the rim reference frame, for $We = 679$, both non-dimensionalized by $d_0^3/\tau _{cap}$. The ratio $q_{in}/q_{out}$ is mostly invariant to change of $We$ (inset). ($d$) Comparison of measured time evolution of the rate of momentum from the sheet to the rim, $p_{in} = \rho q_{in} v_{in}$, and the rate of momentum from the rim into the ligaments/droplets, $p_{out} = \rho q_{out}v_\ell$, in the rim reference frame, for $We = 679$. Both quantities are non-dimensionalized by $\rho d_0^4/\tau _{cap}^2 = 6\sigma d_0$. The ratio of $p_{in}/p_{out}$ ranges from 10 to 3, showing that $p_{out}$ remains systematically smaller than $p_{in}$.

We emphasize that, physically, the width of a ligament, $w$, and the rim thickness, $b$, should not be considered to be equal for all purposes. They are in fact selected by different underlying mechanisms. The rim thickness is determined by the $Bo = 1$ criterion and a self-adjustment of the rim thickness selecting its width to be the local capillary length based on the sheet instantaneous acceleration; while the width of a ligament is determined by the fluid mass and momentum shed by the rim combined with the ligament population growth dynamics. In a separate study (Wang & Bourouiba 2021), we elucidate the details of the ligament evolution and prediction of the time evolution of the ligament width, $w$, important for prediction of the properties of the droplets shed over time. Here, for the purpose of elucidating and capturing the sheet evolution, at first order, $w\approx b$ is an acceptable approximation (introducing $2\,\%$ error on $r_s$) that enables us to decouple the sheet from the ligament dynamics, thus close the governing system, with the rim momentum loss now reading

(5.7)\begin{equation} p_{out} = \rho q_{out}v_\ell = \rho \left[ q_{in} - \frac{\rm \pi}{4}\frac{\textrm{d}}{\textrm{d}t}\left(b^2 r_s\right) \right] \sqrt{\frac{3\sigma}{\rho b}}. \end{equation}

Figure 11($a$) confirms the validity of (5.7), showing that the time evolution of the measured $p_{out} = \rho q_{out}v_\ell$ matches well (5.7) computed using the experimental measurements of $r_s$ and $b$. The closed system governing the sheet radius $r_s$ is thus

(5.8a)\begin{equation} \displaystyle -\rho q_{in}\underbrace{\left(\frac{r_s}{t} - \dot{r_s}\right)}_{v_{in}} + \underbrace{\rho\left[ q_{in} - \frac{\rm \pi}{4}\frac{\textrm{d}}{\textrm{d}t}\left(b^2 r_s\right) \right] \sqrt{\frac{3\sigma}{\rho b}}}_{\rho q_{out} v_\ell} + \left(2-\frac{\rm \pi}{4}\right)\sigma r_s = 0, \end{equation}

with

(5.8b,c)\begin{equation} q_{in}= h(r_s,t)r_s\underbrace{\left(\frac{r_s}{t} -\dot{r_s}\right)}_{v_{in}} \quad \text{and} \quad b =\sqrt{\frac{\sigma}{\rho(-\ddot{r_s})}}. \end{equation}

This closed system is a nonlinear differential equation involving only $r_s(t)$ and its derivatives, while incorporating the inertial effects and continuous fluid shedding from the rim. With $p_{out} = \rho q_{out} v_\ell$ in (5.8) being a second-order derivative of the sheet radius $r_s$ via its dependence on $b$; (5.8) is a nonlinear third-order ordinary differential equation on $r_s$. This system is self-contained and can be directly solved numerically to predict $r_s$, using known initial conditions.

Figure 11. ($a$) Experimental measurement of the momentum leaving the rim per unit of radian $p_{out} = \rho q_{out}v_\ell$ as a function of time. The experimental $p_{out}$ is indirectly measured by multiplying the measured volume rate outward $q_{out}$ and the measured outward speed $v_\ell$ shown in figure 9($b$). The experimental measure of $p_{out}$ is compared with the estimation (5.7) based on the experimental data of $r_s$ and $b$, and also compared with the surface-tension force per radian exerted by the sheet $2\sigma r_s$. Each term is non-dimensionalized by $\sigma d_0$. The inset of ($a$) shows the ratio of the measured $p_{out}$ with $2\sigma r_s$, multiplied by $We^{3/8}$. The dash shows the theoretical value of the $\text {ratio} = 2.1$ (5.13), in good agreement with the experimental data. ($b$) Time evolution of the measured critical length $\lambda (t)$ between neighbouring growing ligaments along the rim, relative to the instantaneous rim thickness $b$. The inset of ($b$) shows the ratio $\lambda /b$ normalized by $We^{3/8}$ for different $We$ collapses on a single curve. The dash line shows the theoretical value of the normalized ratio equal to 0.56 (5.11), in good agreement with the experimental data.

To gain further physical insights, an analytical closed-form solution is also desirable. We discuss next how a key result about the ligament dynamics and momentum outward allows us to further reduce (5.8) to a first-order equation on $r_s(t)$ more suitable to derive an analytical solution for the sheet evolution for different $We$ impact regimes, while continuing to capture continuous-shedding and inertial effects.

5.3. Reduced sheet evolution system, with continuous shedding: non-Galilean Taylor–Culick law

The fluid shed from the rim enables the growth of the ligaments along it. Wang & Bourouiba (2021) show that there is a critical minimum distance $\lambda$ between two adjacent ligaments required to maintain the growth of a fully formed ligament. The minimum distance between ligaments is determined by the constraint that the rate of volume leaving the rim within a distance $\lambda$ must equal the rate of volume entering one ligament at that time, $q_{e\ell }(t)$, namely,

(5.9a)\begin{equation} \frac{q_{out}(t)}{r_s(t)}\lambda(t) = q_{e\ell}(t) = a_\ell(t) v_\ell(t), \end{equation}

with

(5.9b)\begin{equation} a_\ell(t) = \frac{\rm \pi}{4}w(t)^2, \end{equation}

where, $q_{out}$ is the volume rate leaving the rim per radian, $q_{out}(t)/r_s(t)$ is the volume rate leaving the rime per unit arclength, $a_\ell (t)$ is the average cross-sectional area of the ligaments and $v_\ell (t)$ is the average outward fluid speed (4.23). Using (5.9) and the approximation $w\approx b$, the momentum loss term (5.7) is re-expressed as

(5.10)\begin{equation} p_{out} = \rho q_{out}v_\ell = \rho \frac{r_s}{\lambda}a_\ell v_\ell^2 = \frac{\rho r_s}{\lambda} \frac{\rm \pi}{4}b^2 \left(\frac{3\sigma}{\rho b}\right) = \frac{3{\rm \pi}}{4}\left( \frac{b}{\lambda} \right)\sigma r_s. \end{equation}

Figure 11($b$) shows that $\lambda (t)/b(t)$ remains approximately constant throughout the sheet life. Wang & Bourouiba (2021) derived that, in fact, the ratio of the minimum distance $\lambda$ with the rim thickness $b$ is a constant in time, with dependence on $We$ such that

(5.11)\begin{equation} \frac{\lambda}{b} \approx 0.56 We^{3/8}, \end{equation}

consistent with the experimental data in figure 11($b$). Substituting (5.11) into the momentum loss term (5.10) gives

(5.12)\begin{equation} p_{out} = \rho q_{out}v_\ell = \frac{1.34{\rm \pi}}{We^{3/8}}\sigma r_s. \end{equation}

Thus, similarly to the inertia term (5.1), the momentum exiting the rim, $p_{out}$, introduces an equivalent surface-tension force. Hence, the ratio between the magnitude of $p_{out}$ and the surface-tension force term in (5.3) is

(5.13)\begin{equation} \frac{|\rho q_{out}v_\ell|}{|2\sigma r_s|} = \frac{0.67{\rm \pi}}{We^{3/8}} \approx \frac{2.1}{We^{3/8}}. \end{equation}

Figure 11($a$) shows the measured time evolution of the momentum exiting the rim $p_{out} = \rho q_{out}v_\ell$, compared with the surface-tension force $2\sigma r_s$. A systematic gap between the two curves is observed with the ratio between them matching the theoretical value given by (5.13), a prediction that remains robust to the change of $We$ (figure 11$a$, inset), confirming that the approximate expression (5.10) captures, at first order, $p_{out}$ well. Using (5.10), the governing equation of sheet radius $r_s$, (5.3), can be further reduced to

(5.14a)\begin{equation} -\rho h(r_s,t)\left(\frac{r_s}{t} - \dot{r_s}\right)^2 + \left(2-\beta\right)\sigma = 0, \end{equation}

with

(5.14b)\begin{equation} \beta = \frac{\rm \pi}{4}-\frac{1.34{\rm \pi}}{We^{3/8}}, \end{equation}

which is analogous to a Taylor–Culick law, with an additional term $\beta \sigma$. We emphasize that the reduction of the governing equation of sheet radius to an analogue Taylor–Culick law is not because the inertia of the rim is negligible. Instead, it is due to the continuous shedding of the fluid from the rim, with the instantaneous volume rate leaving the rim, $q_{out}(t)$, governed by the nonlinear rim destabilization, which enables the rim to self-adjust its thickness to remain governed by a universal instantaneous $Bo = 1$ criterion, with an instantaneous speed of fluid leaving the rim, $v_\ell (t)$, of shedding determined by a local rim–ligament dynamics (Bernoulli).

It is the combination of these two local mechanisms underlying the rim destabilization and the ligament growth that enables us to convert the rim inertia ($m\ddot {r_s}$) and $p_{out} = \rho q_{out}v_\ell$ of the shedding to effective additional surface-tension forces acting on the rim, the additional term $-\beta \sigma$. The original Taylor–Culick law (4.6) is derived for steady flow in a Galilean frame (non-accelerating), while (5.14) is derived for an unsteady sheet dynamics, where the rim constantly decelerates and sheds fluid, i.e. moving in a non-Galilean frame. Thus, we name (5.14) the non-Galilean – or non-inertial – Taylor–Culick law for rims undergoing unsteady fragmentation with continuous shedding. The great advantage of the non-Galilean Taylor–Culick law is that it reduces the complex and high-order rim governing equation (4.9) to a first-order differential equation in addition to enabling accurate prediction of the droplet shedding occurring continuously throughout the sheet expansion and retraction. Using the length scale $d_0$ and capillary time scale $\tau _{cap}$ with the dimensionless variables introduced by (2.4a--e), the final governing equation in dimensionless form reads

(5.15a)\begin{equation} -6H(R_s,T)\left(\frac{R_s}{T}-\dot{R_s}\right)^2 + \left(2-\beta\right) = 0, \end{equation}

with

(5.15b)\begin{equation} \beta = \frac{\rm \pi}{4}-\frac{1.34{\rm \pi}}{We^{3/8}} \end{equation}

with

(5.15c)\begin{equation} H(R_s,T) = \frac{T\sqrt{{6We}}}{6a_3R_s^3 +a_2R_s^2T\sqrt{{6We}} + a_1R_sT^2{We}}, \end{equation}

where $H(R_s,T)$ is the thickness profile of the sheet at the rim derived from (2.5) and Wang & Bourouiba (2017).

6. Solutions of the non-Galilean Taylor–Culick law

Before deriving the exact solution of (5.15), we first review the general behaviour of its solutions, requiring a choice of initial conditions. Both theoretical and experimental studies of early-time dynamics of the sheet expansion (Thoroddsen et al. 2012; Gordillo, Lhuissier & Villermaux 2014; Philippi, Lagrée & Antkowiak 2016) showed that, at very early times, the rim velocity $\dot {r_s}(0)$ is proportional to the drop impact speed $u_0$. This is reasonable as, at very early times prior to the formation of the sheet, the surface tension force is not dominant, with an initial rim horizontal momentum proportional to the incoming vertical momentum of the impacting drop, thus, $\dot {r_s}(0) \sim u_0$. However, the exact factors remain debated, from $\dot {r_s}(0) = u_0$ (Villermaux & Bossa 2011) to $\dot {r_s}(0) \sim 1/\sqrt {T}$ with singularity at $T = 0$ (Thoroddsen et al. 2012; Riboux & Gordillo 2014). Without loss of generality, and for the purpose of illustration of the properties of the solution of (5.15) governing the sheet reaching the air after expansion on finite surface, we use $\dot {r_s}(0) \approx 2u_0$ as the initial condition.

Non-dimensionalizing the initial velocity $\dot {r_s}(0) = 2u_0$ using the length scale $d_0$ and time scale $\tau _{cap}$ (2.3), the dimensionless initial condition reads

(6.1)\begin{equation} \dot{R_s}(0) = \frac{\dot{r_s}(0)}{d_0/{\tau_{cap}}} = \frac{2u_0}{d_0}\sqrt{\frac{\rho d_0^3}{6\sigma}} = \sqrt{\frac{2}{3}We}. \end{equation}

Figure 12 shows the solution of (5.15) for impact $We$ ranging from 500 to 500 000, with the sheet radius $R_s$ normalized by $\sqrt {We}$. Note that we chose to explore the full range of $We$ values to illustrate the regimes of solutions (5.15). However, phenomena such as prompt splash and sheet disintegration with hole formation would occur in experiments for such high values of Weber number (Mundo, Sommerfeld & Tropea 1995; Yarin 2006; Josserand & Thoroddsen 2016) under certain conditions of air pressure, surface properties, fluid properties and purity (Poulain, Villermaux & Bourouiba 2018). Assuming that these effects are suppressed, for such wide range of $We$, the maximum radius $R_m$, normalized by $\sqrt {We}$, and the time of maximum radius $T_m$ of the solutions of (5.15) vary slowly, $10\,\%$ change, over a wide range of impact $We$ spanning three orders of magnitude. We show next that such slow variations of $R_m/\sqrt {We}$ and $T_m$ are clearly embedded in the structure of (5.15).

Figure 12. ($a$) Full solution of the non-Galilean Taylor–Culick law (5.15) for impact $We$ ranging from 50 to 500 000. The radius of the expanding sheet is normalized by $\sqrt {{We}}$. ($b$) The equivalent maximum radius $Y_m$ and ($c$) the time of maximum radius $T_m$ as a function of $We$, normalized by the maximum radius $Y_m^{\infty }$ and the time of maximum radius $T_m^{\infty }$ as $\textit{We}\to \infty$. The inset of ($a$) shows the snapshot of an expanding sheet at the time of maximum radius for $We = 250$, for which the rim does not fragment.

We introduce a new variable $Y_m = R_m/\sqrt {We}$. Substituting the new variable into the non-Galilean Taylor–Culick equation (5.15) gives

(6.2a)\begin{equation} -6H(Y_s,T)\left(\frac{Y_s}{T}-\dot{Y_s}\right)^2 + \left(2-\beta\right) = 0, \end{equation}

with

(6.2b)\begin{equation} \beta = \frac{\rm \pi}{4}-\frac{1.34{\rm \pi}}{We^{3/8}}, \end{equation}

with

(6.2c)\begin{equation} H(Y_s,T) = \frac{T\sqrt{\mathrm{6}}}{6a_3Y_s^3 +a_2Y_s^2T\sqrt{{6}} + a_1Y_sT^2}, \end{equation}

from which, it is clear that the only explicit dependence on $We$ is $1.34{\rm \pi} /We^{3/8}$, which corresponds to the momentum outward discussed in § 5.3. Given that $\dot {R_s}(0) \sim \sqrt {We}$, the structure of the equation on $Y_m$ (6.2) has initial conditions free of explicit $We$ dependence. Moreover, $\lim _{{We} \to \infty }$ using the $1.34{\rm \pi} /We^{3/8} = 0$. In other words, as $We$ becomes large (6.2) becomes independent of $We$, thus a loss of explicit dependence of $Y_m = R_m/\sqrt {We}$ and $T_m$ on $We$ as $We$ increases. This is clearly confirmed by the full solution (figure 12$a$).

Even when $We$ is not very high, such as $We \sim O(10^2)\unicode{x2013}O(10^3)$, the variations of $Y_m$ and $T_m$ with $We$ remain limited. To better quantify the slow variation of $Y_m$ and $T_m$ for the full range of $We$, we denote $Y_m^{\infty }$ and $T_m^{\infty }$ the normalized maximum radius and time of maximum radius, respectively, as $We \to \infty$. Figures 12($b$) and 12($c$) show the ratio of normalized maximum radii $Y_m/Y_m^{\infty }$ and the ratio of times of maximum radius $T_m/T_m^{\infty }$, respectively, as a function of $We$. The full range of Weber values can be separated into three different regimes. In the high-$We$ regime of $We >10^4$, both $Y_m$ and $T_m$ can reach $95\,\%$ of the upper limit value, thus, in this range, $1.34{\rm \pi} /We^{3/8}$ in (5.15) is negligible, with a contribution to the solution of the order of $5\,\%$ only. In this high-$We$ regime, (5.15) simplifies to

(6.3a)\begin{equation} -6H(R_s,T)\left(\frac{R_s}{T}-\dot{R_s}\right)^2 + \left(2-\frac{\rm \pi}{4}\right) = 0, \end{equation}

with

(6.3b)\begin{equation} H(R_s,T) = \frac{T\sqrt{{6We}}}{6a_3R_s^3 +a_2R_s^2T\sqrt{{6We}} + a_1R_sT^2{We}}. \end{equation}

In the low-$We$ regime of $We<250$, the Weber-dependent term $1.34{\rm \pi} /We^{3/8}$ in $Y_m$ and $T_m$ is non-negligible (figure 12$b$,$c$). However, in this low-$We$ regime (figure 12$a$, inset), no fragmentation occurs typically as the time of rim destabilization into corrugations and their growth into ligaments is of the same order of magnitude as that of the sheet expansion itself (Clanet et al. 2004; Wang et al. 2018). In this regime, prior models not incorporating droplet shedding throughout the sheet dynamics, such as (4.4), would be valid. However, given the lack of fluid fragmentation inducing ejection of droplets beyond the impact area this regime is not the focus of this study.

In the last intermediate-$We$ regime of $250< We <10\,000$, the Weber-dependent term $1.34{\rm \pi} /We^{3/8}$ continues to affect $Y_m$ and $T_m$ within $10\,\%$, with a variation of both $Y_m/Y_m^{\infty }$ and $T_m/T_m^{\infty }$ from 0.85 to 0.95, (figure 12$b$,$c$). In this regime, (5.15) can be simplified by maintaining the weak dependence on $We$ as an approximate constant term of mean value 0.9. Namely, $1.34{\rm \pi} /We^{3/8} \approx 3{\rm \pi} /28$ in this $We$ range used in (5.15) to read as

(6.4a)\begin{equation} -6H(R_s,T)\left(\frac{R_s}{T}-\dot{R_s}\right)^2 + \left(2-\frac{\rm \pi}{7}\right) = 0, \end{equation}

with

(6.4b)\begin{equation} H(R_s,T) = \frac{T\sqrt{{6We}}}{6a_3R_s^3 +a_2R_s^2T\sqrt{{6We}} + a_1R_sT^2{We}}. \end{equation}

In this formulation, over the range of $We$ of interest, the $Y_m$ and $T_m$ are effectively independent of Weber, at first order. The error introduced by this approximation in the prediction of $r_s$ is limited to $5\,\%$ (figure 12), of the same order of magnitude as our measurement errors, thus negligible.

6.1. Application-specific selection of relevant $We$ regime of impact

Common unsteady fragmentation processes upon drop impact in nature are in the intermediate-$We$ regime (blue in figure 12$b$,$c$). Indeed, the size of a drop is restricted by the capillary length, the definition of which can be generalized from the gravitational acceleration to any resulting acceleration induced by other forces. Taking $a_c$ as drop acceleration, the capillary length is $d_c \sim \sqrt {{\sigma }/{\rho a_c}}.$ The terminal velocity, obtained by balancing the drag and the acceleration, is

(6.5)\begin{equation} \rho d_c^3 a_c \sim \frac{1}{2}C_d\rho_a u_c^2 d_c^2 \quad \Longrightarrow \quad u_c \sim \sqrt{\frac{\rho}{\rho_a}\cdot\frac{d_c a_c}{C_d}}, \end{equation}

where $\rho _a$ is the air density and $C_d$ is the $Re$-dependent drag coefficient (figure 13$a$). The resulting maximum $We$ is

(6.6)\begin{equation} We_{max} = \frac{\rho u_c^2 d_c}{\sigma} \sim \frac{\rho}{\rho_a} \cdot \frac{\rho d_c^2 a_c}{\sigma C_d} = \frac{\rho}{\rho_a} \frac{1}{C_d} \approx 10 \frac{\rho}{\rho_a}, \end{equation}

where the minimum value of $C_d = 0.1$ was used (figure 13$a$), showing that the density ratio between liquid and gas phases sets the maximum value of the impact $We$ of drops naturally or artificially generated. Figure 13($b$) shows that the density ratio of most common fluids, including alcohol, oil, water-based solution, acid and alkaline, with air would be approximately $O(10^3)$, thus an associated $We_{max} \sim O(10^4)$, same as the boundary between the intermediate-$We$ and the high-$We$ regimes (figure 12). In extreme conditions, such as a mercury fluid at high temperature, falling in air or another fluid falling in lower pressure surrounding gas, $We_{max}$ could reach $O(10^5)$. Note that for impacts of drops on high-speed moving objects such as airplanes, the $We_{max}$ could, however, be higher, thus producing sheets governed by (6.3). Hereafter, we focus on the sheet evolution in the intermediate-$We$, with $We\sim 250-10^4$, most relevant for environmental and health applications, and with (6.4) as governing equation.

Figure 13. ($a$) Drag coefficient, $C_d$, for a sphere as a function of Reynolds number $Re_{air} = ud/\nu _a$, where $\nu _a$ is the air kinematic viscosity. ($b$) Density values for common liquids at room temperature.

6.2. Independence of the temporal dynamics from the impact energy

Focusing on the sheet dynamics in the intermediate-$We$ regime, we introduce the new variable $Y(T) = R_s(T)/\sqrt {{We}}$ into the reduced non-Galilean Taylor–Culick equation (6.4) leading to

(6.7)\begin{equation} \frac{6\sqrt{6}T}{6a_3Y^3 +\sqrt{6}a_2Y^2T + a_1YT^2}\left(\frac{Y}{T}-\dot{Y}\right)^2 = \left(2-\frac{\rm \pi}{7}\right). \end{equation}

Since the Weber-dependent term corresponding to the momentum outward term, in this regime, is approximately constant, (6.7) is independent of $We$. Moreover, recalling that $\dot {R_s}(0) \sim \sqrt {{We}}$ (6.1), the initial condition for $Y(T)$ is also independent of $We$. With both the governing equation on $Y(T)$ and its initial condition being independent of $We$, $Y(T)$ is also independent of $We$. More precisely, the temporal evolution of the sheet radius $R_s(T)/\sqrt {{We}}$ is governed by a characteristic time scale, the capillary time $\tau _{cap}$, that is independent of the initial impact energy or $We$. This is a property typical of linear dynamical systems, such as harmonic oscillators, with their temporal evolution, or period, independent of their amplitude (energy). However, such property is not generally true for nonlinear systems, even when they are as simple as a weakly nonlinear pendulum for which the time to peak depends on its amplitude or energy. We would thus certainly not expect, a priori, this property to be true for a complex nonlinear system such as (6.7). In a separate study, we show that (6.7) indeed satisfies this peculiar property (Wang & Bourouiba 2022), implying that two drops with the same $\tau _{cap}$ impacting a surface at the same time, but with different energies, or $We$, would in fact evolve in synchrony, including reaching their maximum sheet expansion and final retraction at the same non-dimensional time. This property is clear in our data (figures 14 and 15$a$) and is also notable in the data of prior studies, which show the time of maximum radius $T_m$ to be independent of $We$ (Villermaux & Bossa 2011; Lastakowski et al. 2014), while the maximum radius $R_s$ is proportional to $\sqrt {{We}}$ (Eggers et al. 2010; Villermaux & Bossa 2011; Riboux & Gordillo 2016).

Figure 14. ($a$) Temporal evolution of the non-dimensional radius, $R_s$, of the expanding sheet, for different $We$. The vertical dashed line shows that independently of the $We$, the sheet reaches its maximum extension at the same dimensionless time $T_m = t_m/\tau _{cap}$. ($b$) The temporal evolution of the sheet radius is shown for different impacting drop diameters $d_0$, confirming that the magnitude of the sheet radius is governed by $We$ and not $d_0$; while the time of maximum radius is reached at the same non-dimensional time independent of $We$ (inset), i.e. $T_m = const$.

Figure 15. ($a$) Measured temporal evolution, with time normalized by $\tau _{cap}$, of the dimensionless sheet radius, normalized by $\sqrt {{We}}$, collapses on a single curve. The solid line is the full solution of (6.4), and the dashed line is its approximate analytic expression (7.3), both matching the data very well. The inset shows that the measured, and predicted maximum radius (7.2), of the expanding sheet scales as $\sqrt {{We}}$, robust to the change of impacting drop diameter. ($b$) Measured temporal evolution of the sheet radius with different viscosity and elasticity values within the range of validity of the rim thickness $Bo = 1$ criterion, which is well captured by the prediction (7.3).

7. Full and approximate expression of the sheet radius $R_s(T)$

7.1. Initial condition and full solution

The non-Galilean Taylor–Culick law imposes that the sheet radius reaches its maximum at a fixed dimensionless time regardless of initial $We$. Our focus being on the unsteady sheet dynamics incorporating the shedding, not just the very early contact dynamics, without loss of generality, we use the time of maximum radius $T_m$ as the condition to fully determine (6.4). Experimentally, $T_m = 0.43$ (figure 14$b$ inset), consistent with the prior literature (figure 3 inset). Using $T_m = 0.43$, the full numerical solution of (6.4) shown in figure 15($a$), captures the data well. Recalling the properties of invariance of time evolution of the sheet radius § 6, the full solution predicts and confirms that the sheet radius $R_s(T)$, normalized by $\sqrt {{We}}$, follows a universal time evolution independent of impact $We$. Thus, confirming that the maximum radius $R_m$ is proportional to $\sqrt {We}$ (figure 15$a$, inset).

7.2. Approximate explicit solution

Having verified that the full solution captures the data well, we now turn to the derivation of the approximate explicit solution of (6.4). Our focus is the dynamics of shedding throughout the sheet dynamics, not just its very early time, thus we leverage what we already learned from the full solution of (6.4) confirmed by data (figure 14), namely that the maximum radius $R_s(T_m)$ is reached at a constant dimensionless time $T_m = 0.43$, independent of impact $We$ (figure 14$b$ inset). This peculiar, but robust condition is sufficient to derive an explicit approximate expression for the solution of (6.4), which we here expand around time $T = T_m$, associated with $\dot {R_s}(T=T_m) = 0$, the maximum radius $R_m = R_s(T_m)$ determined by

(7.1)\begin{equation} a_3T_m\left(\sqrt{\frac{6}{{We}}}R_m\right)^2 + (a_2T_m^2-\alpha)\left(\sqrt{\frac{6}{{We}}}R_m\right) + a_1T_m^3 = 0, \end{equation}

where $\alpha = 6/(2-{\rm \pi} /7)$. $a_1$, $a_2$ and $a_3$ are constant coefficients of the sheet thickness profile, independent of $We$, determined in Wang & Bourouiba (2017) and $T_m = 0.43$, thus (7.1) is algebraic on $R_m$, with

(7.2)\begin{equation} R_m = \frac{1}{2a_3T_m}\sqrt{\frac{6}{{We}}}\left(\alpha - a_2T_m^2 - \sqrt{(a_2T_m^2-\alpha)^2 - a_1a_3T_m^4}\right) \approx 0.12\sqrt{{We}}, \end{equation}

which is close to the full solution of (6.4) and matching the experiments (figure 14$b$ inset). However, different from the maximum radius $R_m$, the derivation of the explicit analytical solution for the time evolution of sheet radius $R_s(T)$ is not trivial, due to the nonlinearity of the thickness profile, $H(R_s,T)$. Leveraging the property of independence of the temporal evolution of the sheet from the initial energy, or $We$, (§ 6.2), and using the condition of maximum radius reached at $T_m = 0.43$, we derive an approximate analytical solution using a power series expansion (appendix A) which ultimately lead to

(7.3a)\begin{equation} R_s(T) = \sqrt{{We}}\left[0.15(T-T_m)^3-0.4(T-T_m)^2 + \mathcal{R}_m\right] \end{equation}

with

(7.3b,c)\begin{equation} T_m = 0.43 \quad \text{and} \quad \mathcal{R}_m = R_m/\sqrt{{We}} = 0.12 \end{equation}

for $250<We<10\,000$, which also captures the data and full solution very well (figure 15$a$). Note that the two coefficients in the approximate analytical solution (7.3) are not free nor fitted parameters. They were derived from the power series expansion of the non-Galilean Taylor–Culick system (appendix A), showing that the two coefficients are universal, independent of impact conditions in the intermediate-$We$ regime: $250 < We < 10\,000$.

8. Robustness of the prediction of the sheet evolution, $r_s(t)$

8.1. Comparison with experiments using fluids of various properties

To further verify the robustness of our model (5.15) and (7.3) for the prediction of the sheet radius, we conducted experiments using fluids of various properties (table 2). The non-Galilean Taylor–Culick equation (5.15) that governs the sheet radius is derived for the inviscid limit, in which the sheet evolution is independent of the fluid viscosity, as well as the fluid elasticity, when the impact $We =\rho u_0^2 d_0/\sigma$ and $Re =u_0 d_0/\nu$ are much larger than 1. However, Wang et al. (2018) shows that as the fluid viscosity and elasticity increase, the $Bo = 1$ criterion for the local dynamics of the rim is violated when the local rim $\hat {Re} = \sqrt {\sigma b/(\rho \nu ^2)} < 6\sqrt {2}$ and local rim $De = \tau _E/\sqrt {\rho b^3/\sigma } > 16$, where $\tau _E$ is the relaxation time of the elastic fluid. When the $Bo = 1$ criterion breaks, the simplification of the inertial term associated with the rim acceleration in the momentum conservation (5.1) is violated. In addition, the expressions of the volume outward $q_{out}$, (4.19), determined based on the $Bo = 1$ criterion also becomes invalid. Thus, the prediction (7.3) of the time evolution of the sheet radius would be violated when the fluid viscosity and elasticity are high enough to push the system outside the region of validity of the rim $Bo = 1$ condition. On the other end, the prediction (7.3) of the sheet radius should be robust within the range of the fluid viscosity and elasticity where the rim $Bo = 1$ criterion holds, which we verify next.

We use the glycerol–water mixtures to increase the fluid viscosity and the PEO–PEG solution to increase the fluid elasticity (table 2). Figure 15($b$) shows the time evolution of the sheet radius of $50\,\%\unicode{x2013}50\,\%$ glycerol–water mixture and PEO0.1 %–PEG1 % water solution, compared with that of water. It shows that within the regime of validity of the rim thickness’ $Bo = 1$ constraint, the radius of the sheet of the fluids with different viscosity and elasticity continue to be captured by (7.3) well. Finally, we use dodecane and Dimethyl sulfoxide (DMSO) to vary the fluid surface tension. The surface tension of dodecane $\sigma = 25\pm 2$ ${\rm mN}\ {\rm m}^{-1}$ and DMSO $ \sigma = 42 \pm 2\ {\rm mN}\ {\rm m}^{-1}$ are measured by tensiometer, leading to a time evolution of the sheet radius that remains well captured by our prediction, (7.3).

8.2. Comparison of the prediction with the prior literature results

Having verified the robustness of our prediction of the sheet radius evolution (7.3) for fluids of various properties, we compare our prediction with experimental data and predictions from prior studies. Figure 16($a$) shows the experimental data of the time evolution of sheet radius from Vernay et al. (2015), which used the same impact condition as this study, compared with our data and prediction (7.3). It shows that their data match well our experimental data and both of them match the prediction (7.3). Figure 16($a$) also shows the comparison of our experimental data with theoretical and numerical predictions in the literature, for various solid surfaces.

Figure 16. (a) Comparison of our sheet radius prediction, (7.3), with other models in the literature. ($b$) Comparison of our prediction (7.3) with the experimental data from prior studies. ($c$ and $d$) Experimental data of the maximum sheet radius, $R_{m}$, and the time of maximum radius, $T_{m}$, as functions of $We$ from different studies. The values of $T_{m}$ remain remarkably fixed at ($\approx 0.43$), for most studies with our prediction, (7.2), capturing also $R_{m}$ well. Abbreviations for different studies are defined in figure 3.

Figure 16($c$) shows the experimental data of the time evolution of sheet radius from other studies, for inviscid sheet expansion upon drop impact on solid surfaces, but with slightly different impact conditions, which can explain the discrepancies as discussed in § 2. However, after normalizing by $\sqrt {{We}}$, the curves from most studies collapse. Figure 16($b$) shows the dimensionless time $T_{m}$, non-dimensionalized by $\tau _{cap}$, of the maximum sheet radius as a function of Weber number for all these studies, showing that for all but one study $T_{m} = 0.43$, consistent with our data. Figure 16($d$) shows that the maximum radius of the sheet also follows our prediction $R_m = 0.12\sqrt {{We}}$ well, further confirming the robustness of (7.3) for the prediction of the time evolution of the radius of sheet dynamics in the air.

The deviation of the experimental data of Villermaux & Bossa (2011) from the others may be caused by differences in experimental conditions as a small drop-to-rod size ratio is used, $\eta = 1$ (figure 3), which is out of the range (Wang & Bourouiba 2017) that guarantees a complete momentum transfer from the vertical to horizontal direction upon drop impact. An insufficient transfer can reduce the extension of the sheet, consistent with the lower $R_m$ in Villermaux & Bossa (2011). In addition, fluid loss on the edges of the target used to control the sheet ejection angle – a brass jacket fitted to the impacted cylinder – combined with initial potential prompt splash for a fluid such as ethanol may also explain the smaller time of maximum radius, $T_m$, in Villermaux & Bossa (2011). An effective diameter $d_{0e} = 0.8\ d_0$ and an effective impact Weber number $We_e = 0.4\ We $. We would reconcile the data with the rest of the studies shown in figure 16. These effective values of $d_{0e}$ and $We_e$ would correspond to an impact speed $u_{0e}$ reduced to $70\,\%$ of the measured impacting drop $u_0$. This is similar to the required adjustment in initial condition, in particular initial energy used, discussed in Klein et al. (2020).

The experimental condition from Lastakowski et al. (2014), which involves a drop impact on a super-heated infinite solid surface (figure 3), is more distinct than that of our experimental set-up. However, the deviation of their data from ours remains small. As mentioned in § 2.2, we attribute the discrepancy to the expected difference between the sheet thickness profile of a sheet on an infinite surface compared to that of a sheet in the air. Other potential factors, such as the effect of the vapor layer or the different volume rate of fluid shed from the rim, could also explain the deviation and are to be determined in future work.

9. Conclusion

We presented the results of a combined experimental and theoretical investigation of sheet dynamics – expansion and retraction in the air – under unsteady fragmentation upon drop impact on a surface of comparable size to that of the drop. We used especially developed AIP algorithms to quantify, with high accuracy, key quantities involved in unsteady fragmentation, from sheet, to rim, to ligaments, to droplet properties. First, the AIP method enabled us to show the importance of continuous shedding from the rim bounding the sheet: such shedding not only removes a non-negligible amount of mass from the rim (figure 6), but also a non-negligible amount of momentum (appendix B). Quantifying and modelling the continuous shedding occurring during both expansion and retraction of the sheet is not only important for the prediction of the resulting spray properties, but also for the prediction of the sheet dynamics itself as it is influenced by the continuous loss of the rim's mass and momentum, from shedding.

The sheet formed from drop impact expands and retracts, continuously decelerating throughout the process, while it is bounded by a rim from which the ligaments and drops are shed. The inherent unsteadiness of this dynamics associated with continuous deceleration introduces a fictitious force, in the non-Galilean reference frame of the rim, acting outward during both phases of sheet extension and retraction in the air. This fictitious force induces a self-adjustment of the rim thickness (Wang et al. 2018), with its thickness continuously adjusting to be equal to the local instantaneous capillary length defined with the rim's acceleration: a robust $Bo = 1$ condition on the rim (4.14). In this paper we derived the governing equation capturing the intricate interplay between such a non-Galilean dynamics of the sheet – and its bounding rim in particular – coupled with the continuous shedding of fluid from the rim into ligaments and drops governed by the $Bo = 1$ condition. At each step of the derivations, we conducted detailed comparisons between theory and experimental measurements to validate our approach and intermediate assumptions. The shedding induces sinks of mass and momentum on the rim (4.8).

In the range of impact Weber numbers relevant for applications involving spray production, $We> 250$, we show that the resulting rim's nonlinear equation has a remarkable property of time to peak invariance to initial impact $We$ – or energy (§ 6.2 and figure 12). This is a property typical of linear systems, such as simple harmonic oscillators, that reach their maximum amplitude at a time independent of the initial energy injected into the system at rest (§ 6.2). Leveraging this robust property of the nonlinear equation of motion of the rim, we demonstrate that the full governing equation of the sheet radius, (4.9), can be reduced to a continuous-shedding non-Galilean Taylor–Culick law (5.15) and more specifically takes the form (6.4) for common ranges of applications with $250 < We < 10^4$ (§ 6.1).

We showed that the continuous-shedding non-Galilean Taylor–Culick law (6.4), governs the sheet dynamics while accounting for continuous ligament formation and droplet shedding. We showed that the predictions from (6.4), both in the form of a full solution and approximate analytical solution (7.3), capture very well the sheet dynamics of our experiments and those of the prior literature on sheet expansion in the air (figure 15$a$). Finally, we showed the robustness of our predictions to changes of fluid properties, including moderate fluid viscosity and elasticity, including biological mucosalivary fluid (as characterized in Bourouiba 2021), and reconciled prior literature's inconsistent experimental results on the sheet dynamics upon drop impact (figure 15$b$).