1 Introduction
Observing marks of a brief rain shower on the beach or on mud is a common experience to many of us. Some of these rain imprints may even have been preserved over geological time scales and are discovered as sedimentary structures that were first reported back in 1839 and are known as raindrop impressions (Cunningham Reference Cunningham1839). These impressions were suggested as a means to determine the air density of remote eras (Lyell Reference Lyell1851), an idea which has been explored recently (Som et al. Reference Som, Catling, Harnmeijer, Polivka and Buick2012). As their name suggests, raindrop impressions are expected to be concave pits and, therefore, the presence of convex sedimentary features with underneath cavities has led to speculations on alternative causes, e.g. they could be air bubbles that dried in mud without bursting (Lyell Reference Lyell1851). These speculations give rise to questions regarding the actual origin of raindrop impressions (Moussa Reference Moussa1974; Metz Reference Metz1981). In this paper, we will demonstrate that under some conditions a cavity is formed naturally underneath a droplet impacting on dry sand.
The formation of three morphologies; (a) truffle, (b) pie and (c) pancake. Each row represents five consecutive snapshots of an impact. The leftmost snapshot is at the moment of reaching the maximum expansion of the droplet,
$t=t^{\ast }$
(namely,
$t^{\ast }=5.2~\text{ms}$
, 5.6 ms and 5 ms, respectively from top to bottom). The rightmost one is the last frame,
$t=t^{\,f}$
, showing the final steady state. The snapshots in between correspond to
$t^{\ast }+5~\text{ms}$
,
$t^{\ast }+10~\text{ms}$
,
$t^{\ast }+15~\text{ms}$
. In each row snapshots are taken from the same viewpoint. Scale indicators of 1 mm by 1 mm are given at the left bottom in the leftmost image of each row.
From a physics perspective the morphology of the marks left by droplet impact on granular targets depends on the impact dynamics. If the grains of the target are wettable to the impacting liquid, the final morphology of the impact crater consists of two parts: the deformed target surface and the liquid–grain residual. While the deformed target surface is always concave, various morphologies of the liquid–grain residual have been reported (Katsuragi Reference Katsuragi2010, Reference Katsuragi2011; Delon et al. Reference Delon, Terwagne, Dorbolo, Vandewalle and Caps2011; Zhao et al. Reference Zhao, Zhang, Tjugito and Cheng2015a ; Zhao, de Jong & van der Meer Reference Zhao, de Jong and van der Meer2015b ). Loosely speaking, the final morphology is a consequence of the dynamics after the droplet reaches its maximum spreading diameter. If the liquid–grain mixture is able to retract, the final residual is strongly curved and concentrates in the centre of the crater (cf. figure 1 a), the so-called truffle and doughnut morphology (Zhao et al. Reference Zhao, de Jong and van der Meer2015b ). Otherwise the mixture is ‘frozen’ in its expansion state and covers most of the concave crater (cf. figure 1 c). However, even without retraction, the liquid–grain mixture could still undergo vertical deformation while its horizontal dimension is unchanged – the centre is lifted up (cf. figure 1 b), and in some cases even protrudes from the sandy surface. Furthermore, we discovered that this lift creates an unexpected feature associated with the last morphology, a cavity underneath, which is only observable with X-ray tomography (cf. figure 2 h). Immediate questions that arise are where this phenomenon is located in the parameter space and what is the underlying physical mechanism that creates the cavity. To answer these questions, we conduct a series of experiments varying the impacting liquid, the impact speed, the grain size and the wettability of individual grains.
2 Experimental methods
In our experiments the impacting droplet is composed of either water or ethanol mixed with food dye (mass fraction
${<}2\,\%$
) for visualization purposes. The radius of the droplet
$R_{d}$
is fixed to 0.9 mm for ethanol. The radius of the water droplets is fixed to 1.4 mm for most experiments and to 1.75 mm occasionally. The impacting droplet is released from a nozzle above the substrate. Small oscillations may be induced during the pinch-off of the droplet. Nevertheless the effect of these oscillations may safely be neglected in the inertia-driven spreading regime we focus on in our experiments. The impact speed,
$U$
, reaches from
$1.5~\text{m}~\text{s}^{-1}$
to
$5.5~\text{m}~\text{s}^{-1}$
by altering the falling height. The granular target consists of ceramic beads with the specific density
$\unicode[STIX]{x1D70C}_{g}=6000~\text{kg}~\text{m}^{-3}$
. Three different grain sizes are used, and the hydrophilicity of the grains is enhanced after cleaned with a piranha solution (see table 1) (Zhao, de Jong & van der Meer Reference Zhao, de Jong and van der Meer2017). The packing density of the bed
$\unicode[STIX]{x1D719}$
is prepared in the range of 0.59–0.63 by air fluidization and subsequent mechanical taps. While the droplet deformation is visualized with a high-speed camera, at the same instance, the deformation of the target surface is measured by an in-house built high-speed laser profilometer (Zhao et al.
Reference Zhao, de Jong and van der Meer2015b
). The impacting droplet appears to be a dark ellipse (due to perspective) in the video. When detecting the spreading diameter cross-correlations are computed between dark ellipse templates and the (background subtracted) image, and the template size corresponding to the maximum correlation gives the spreading diameter. The evolution of the crater profile with time,
$t$
, is quantified as a function
$z(r,t)$
. It is good to stress that here, as well as in all of the analyses in this study, axisymmetry of the impact crater is assumed. Here,
$z$
represents the vertical coordinate, and
$r$
is the in-plane radial coordinate with the impact centre at
$r=0$
. When the location of impact is far from the laser lines, the part of
$z(r,t)$
close to the impact centre is unresolved. We discard the experiments where the unresolved region is beyond
$r=1~\text{mm}$
, and extrapolate the profile to
$r=0$
with a parabolic fit otherwise, where a parabola centred at
$r=0$
is fitted to the last 20 datapoints closest to the centre. The parabolic fitting function is given in appendix C. This fit is implemented on the available data in the range of
$r<2~\text{mm}$
. To estimate the potential error the same method is applied to a fully resolved experimental profile, and the difference between the fitted crater depth and the measured one is very small (0.0043 mm), and is representative of potential errors made by using the fit. Note that the corresponding error in the crater area is even smaller, since the fitted region of
$r<1~\text{mm}$
, the crater centre, is typically smaller than the region of interest (see below).
Summary of the contact angle of water and ethanol,
$\unicode[STIX]{x1D703}_{w}$
and
$\unicode[STIX]{x1D703}_{e}$
, and the grain size,
$d_{g}$
.
The leftmost column (a–c) illustrates the crater profile at the moment of maximal droplet expansion
$t=t^{\ast }$
(dashed lines) and the moment the final static state of the crater is reached,
$t=t^{\,f}$
(solid lines). For (a,b), the differences of crater profiles at these two moments result from the deformation of the liquid–grain mixtures, whereas the discrepancy at the crater rim in (c) is caused by grain splashing. The middle column (d–f) shows the time evolution of crater area,
$\unicode[STIX]{x1D6F4}$
(see text for definition). The moment
$t^{\ast }$
of maximum droplet expansion is indicated by vertical dashed lines in (d–f). In (e) the dash dotted line denotes
$\unicode[STIX]{x1D6F4}(t^{\,f})$
. Both columns correspond to the impacts in figure 1. The rightmost column (g–i) displays vertical cross-sections from X-ray tomograms through the impact centre where the residual is highlighted with yellow colour. These tomograms are scanned for droplet impacts on glass beads, different experiments from that in (a–f). Details can be found in appendix A.
Figure 1 exemplifies three residual morphologies: truffle, pie and pancake. Those residual shapes result from different post-impact processes, as can be seen in the experimental movies in the supplementary materials available online at https://doi.org/10.1017/jfm.2019.742. With the acquired crater profile,
$z(r,t)$
, differences between those processes can be quantified. We first compare the profiles at two moments in time. One moment is at
$t=t^{\ast }$
, where the impacting droplet reaches its maximum expansion radius
$R_{d}^{\ast }$
. The other moment is the last frame of the experiment,
$t=t^{\,f}$
, where the crater has obtained its final static state. Minutes or hours afterwards the experiment, the liquid–grain residual is dried. Its shape may only change slightly from that at
$t^{\,f}$
. The crater profiles of experiments shown in figure 1 are presented in figure 2(a–c). While the profiles at
$t^{\ast }$
and
$t^{\,f}$
are largely identical for the pancake shape (cf. figure 2
c), both truffle and pie shapes display prominent evolution (cf. figure 2
a,b). The difference between the truffle and pie shapes is clear as well. For the truffle shape the evolution occurs through the whole crater, from rim to the centre, and the final profile is highly curved. In contrast, for the pie shape, the part of the profile far from the impact centre remains pinned, and the final profile is relatively flat, with the centre lifted upwards.
To further quantify the crater evolution we calculate the reduced crater area by integrating the crater profile between
$r=0$
and
$r=R_{0}$
,
$\unicode[STIX]{x1D6F4}(t)=\unicode[STIX]{x03C0}\int _{0}^{R_{0}}\sqrt{1+z^{\prime 2}}\,\text{d}r^{2}-\unicode[STIX]{x03C0}R_{0}^{2}$
. Here,
$z^{\prime }=\unicode[STIX]{x2202}z(r,t)/\unicode[STIX]{x2202}r$
is the derivative of the profile with respect to
$r$
, and
$R_{0}$
is the radial distance where the crater profile intersects with
$z=0$
at
$t=t^{\ast }$
(cf. figure 2
b). The results corresponding to figure 2(a–c) are given in (d–f) in the same figure. As the flat area,
$A_{0}=\unicode[STIX]{x03C0}R_{0}^{2}$
, is subtracted,
$\unicode[STIX]{x1D6F4}(t)$
always increases from 0 upon impact to a local maximum at
$t=t^{\ast }$
for all three cases (cf. figure 2
d–f). Afterwards the time evolution
$\unicode[STIX]{x1D6F4}(t)$
exhibits interestingly distinct signatures corresponding to the three typical crater profile evolutions described above. The profile of pancake shape does not change between
$t^{\ast }$
and
$t^{\,f}$
, therefore,
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}=\unicode[STIX]{x1D6F4}(t^{\,f})-\unicode[STIX]{x1D6F4}(t^{\ast })=0$
(cf. figure 2
f). In contrast, for the truffle shape, the liquid–grain mixture contracts into a highly curved residual, more curved than the crater profile at
$t^{\ast }$
, which results in a positive area difference,
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}>0$
(cf. figure 2
d). For the pie shape, as the centre of the crater is lifted, the profile becomes flatter, and
$\unicode[STIX]{x1D6F4}$
decreases. Therefore, the area difference is negative, i.e.
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}<0$
(cf. figure 2
e). The lift of the crater centre creates a cavity revealed by the X-ray tomography scan after
$t^{\,f}$
in figure 2(h) (refer to appendix A for details of the X-ray tomogram). Hence, the sign of
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}$
can be used to indicate morphology:
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}=0$
for the pancake shape,
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}>0$
for the truffle shape and
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}<0$
for the pie shape. In experiments the oscillation of the remaining liquid surface introduces fluctuations into
$\unicode[STIX]{x1D6F4}(t)$
after
$t>t^{\ast }$
. When computing
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}$
we average the data of
$\unicode[STIX]{x1D6F4}(t)$
in the last 2 ms to obtain the value
$\unicode[STIX]{x1D6F4}(t^{\,f})$
at
$t^{\,f}$
.
Figures 1 and 2 indicate that the final morphology is the outcome of the post-impact dynamics of the liquid–grain mixture. What leads to the different post-impact dynamics? From the formation of truffle morphology it is clear that in this case the dynamics is dominated by the Laplace pressure, the product of the surface tension and local curvature of the liquid–air interface, which are developed during the impact. This also give us clues to the formation of the other two morphologies.
3 Mixing ratio
Similar to impacts on a solid surface, upon impact the droplet starts to spread, and the circumference curvature is developed along its edges. Meanwhile, a central curvature also arises, which is unique to impacts on deformable substrates, like the granular target studied here (cf. figure 3 inset). These two curvatures determine the post-impact dynamics. According to the corresponding crater dimensions (Zhao et al.
Reference Zhao, de Jong and van der Meer2015b
, Reference Zhao, de Jong and van der Meer2017; de Jong, Zhao & van der Meer Reference de Jong, Zhao and van der Meer2017), the circumference curvature is always larger than the central one. Nevertheless, their relative significance can be altered by the degree of mixing between liquid and grains during the impact. When there is little mixing, the circumference curvature is dominant and leads to retraction of the liquid–grain mixture, which finally results in a ball-like shape, a truffle. In the other extreme of complete mixing, where the whole droplet penetrates into the granular substrate during the impact, both curvatures are lost, and no further dynamics occurs after the impact. The final morphology preserves the maximum expanded shape, the so-called pancake shape. However, in the intermediate range of mixing, where more than half of the droplet is mixed with grains, the circumference curvature is suppressed, which prevents the retraction. The radial extension of the residual is then frozen in the expanded state. Nevertheless, the central curvature of the top liquid–air interface may still lift the liquid–grain mixture and create the underneath cavity displayed in figure 2(h), the pie morphology. In this dual-curvature model the morphology appears to be a function of mixing ratio at
$t^{\ast }$
, namely the ratio of the liquid volume that is mixed with grains,
${\mathcal{V}}_{m}$
, to that of the whole droplet,
$V_{d}=4\unicode[STIX]{x03C0}R_{d}^{3}/3$
.
Main plot: crater area difference,
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}$
, plotted versus the estimated mixing ratio
${\mathcal{V}}_{m}/V_{d}$
. The sign of
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}$
, here normalized by the corresponding flat area
$A_{0}$
, indicates the final morphology. The morphology develops with the estimated mixing ratio
${\mathcal{V}}_{m}/V_{d}$
. Three regimes can be distinguished which are in accordance with the residue morphology shown in figure 1: a, truffle; b, pie; c, pancake. See text for details and for the definition of
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}$
and
$A_{0}$
. Inset: a sketch of the two curvatures occurring in a spread droplet. In practice, the droplet (blue) is mixed with grains which are not drawn here.
To estimate
${\mathcal{V}}_{m}$
, we first decompose it into the contact area,
${\mathcal{A}}_{c}$
, and the penetration depth of the liquid into the substrate,
$L$
. Although in reality the contact area evolves with time in a complicated manner, for the purpose of this analysis it is approximated by a constant, namely the geometric mean of the initial droplet radius
$R_{d}$
and its maximum spreading radius
$R_{d}^{\ast }$
, i.e.
${\mathcal{A}}_{c}=\unicode[STIX]{x03C0}R_{d}R_{d}^{\ast }$
. The development of
$L$
is modelled by Darcy’s law (refer to appendix B), and its solution is
$$\begin{eqnarray}L(t)=\sqrt{\frac{2\unicode[STIX]{x1D705}P}{\unicode[STIX]{x1D707}_{l}(1-\unicode[STIX]{x1D719})}\,t}.\end{eqnarray}$$
In the above equation the permeability of the substrate,
$\unicode[STIX]{x1D705}=(1-\unicode[STIX]{x1D719})^{3}d_{g}^{2}/(180\unicode[STIX]{x1D719}^{2})$
, is defined by the Carman–Kozeny relation (Carman Reference Carman1956), and
$\unicode[STIX]{x1D707}_{l}$
is the dynamic viscosity of the liquid. The driving pressure
$P$
consists of two parts: the inertia pressure and capillary pressure, where the effect of wettability of grains on the liquid penetration is included in the latter. The complete expression of
$P$
can be found in appendix B and Zhao et al. (Reference Zhao, de Jong and van der Meer2017). In the model the spreading time is estimated as half of the intrinsic capillary oscillation time of the droplet (Richard, Clanet & Quere Reference Richard, Clanet and Quere2002; Okumura et al.
Reference Okumura, Chevy, Richard, Quéré and Clanet2003; Delon et al.
Reference Delon, Terwagne, Dorbolo, Vandewalle and Caps2011),
$\unicode[STIX]{x1D70F}_{c}={\textstyle \frac{1}{2}}\sqrt{(\unicode[STIX]{x03C0}/6)(\unicode[STIX]{x1D70C}_{l}D_{0}^{3}/\unicode[STIX]{x1D70E})}$
, where
$\unicode[STIX]{x1D70E}$
is the surface tension of the liquid. With all the quantities in (3.1) defined, we could finally evaluate
${\mathcal{V}}_{m}=(1-\unicode[STIX]{x1D719}){\mathcal{A}}_{c}L(\unicode[STIX]{x1D70F}_{c})$
. Note that the factor
$(1-\unicode[STIX]{x1D719})$
accounts for the presence of grains in the mixture. The resultant
${\mathcal{V}}_{m}$
increases with impact speed
$U$
, grain size
$d_{g}$
and wettability,
$\cos \unicode[STIX]{x1D703}_{c}$
.
The normalized crater area development
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}/A_{0}$
is plotted against mixing ratio
${\mathcal{V}}_{m}/V_{d}$
in figure 3. Three regimes corresponding to individual morphologies can be defined according to the sign of
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}$
: a.
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}>0$
, the truffle regime; b.
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}<0$
, the pie regime; and c.
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}\approx 0$
, the pancake regime. As morphology develops with the mixing ratio, the first cross-over, between regime a and b, is found at
${\mathcal{V}}_{m}/V_{d}\approx 0.55$
, where approximately half of the droplet is mixed with grains, and the circumference curvature is lost. As a result, surface tension is not able to retract the droplet anymore, as discussed in Zhao et al. (Reference Zhao, de Jong and van der Meer2015b
). The second cross-over, between b and c, is at
${\mathcal{V}}_{m}/V_{d}\approx 1$
, where the whole droplet penetrates into the granular substrate, and the liquid–grain mixture is completely frozen at its maximum expansion. Because of mass conservation,
${\mathcal{V}}_{m}/V_{d}>1$
is physically clearly not expected, and it simply indicates that the whole droplet penetrates into the granular target before reaching
$t=\unicode[STIX]{x1D70F}_{c}$
. In fact, this is consistent with a transition of the droplet spreading dynamics. It has been shown that the viscous dissipation of liquid penetrating into the granular target becomes the largest energy sink in this regime (Zhao et al.
Reference Zhao, de Jong and van der Meer2017). Hence, regime c may also be interpreted from an energy perspective: no further dynamics would happen after the impact energy is largely dissipated by the viscosity of the liquid moving in the pores between the grains.
4 Lifting criterion
Both the wettability and the size of the grains play a role in the degree of mixing, and thus in the formation of final morphologies. In general, increasing the grain size or the wettability enhances the mixing ratio. For instance, impacts on
$d_{g}=257~\unicode[STIX]{x03BC}\text{m}$
result in the pancake shape, while impacts on
$d_{g}=98~\unicode[STIX]{x03BC}\text{m}$
do not enter that regime. Although the dual-curvature model introduced above captures this general tendency of the morphology development, in regime b (
$0.55<{\mathcal{V}}_{m}/V_{0}<1$
) there are impacts of relatively large
$d_{g}$
and low
$U$
that actually end in the pancake shape (
$\unicode[STIX]{x0394}\unicode[STIX]{x1D6F4}\approx 0$
) rather than the expected pie morphology. Therefore, the mixing ratio alone is not sufficient to explain the formation of pie morphology. The missing factor is the force balance. If the crater is too shallow, the Laplace pressure offered by the central curvature cannot overcome the gravitational pressure, then no lifting would occur.
The lifting of the liquid–grain mixture is a process where surface energy,
$E_{s}$
is converted into gravitational energy
$E_{g}$
. The condition under which lifting can be initiated is given by
$\text{d}(E_{s}+E_{g})/\text{d}t|_{t=t^{\ast }}<0$
. This condition is path dependent. To obtain a criterion in terms of accessible experimental quantities we introduce a virtual lifting process, where the concave liquid–grain mixture at
$t^{\ast }$
is raised to be flat. Following this deformation the centre depth of the crater
$Z_{c}$
is decreased from
$Z_{c}^{\ast }$
to 0. Note that here depth is used for the downwards direction. Using the chain rule and
$\text{d}Z_{c}/\text{d}t<0$
the lifting criterion can be rewritten now as
$$\begin{eqnarray}\left.\frac{\text{d}(E_{s}+E_{g})}{\text{d}Z_{c}}\right|_{Z_{c}=Z_{c}^{\ast }}>0.\end{eqnarray}$$
If the crater profile is parabolic,
$z(r,Z_{c})=Z_{c}(r^{2}/R_{0}^{2}-1)$
, one can readily obtain the proportionalities
$E_{g}\propto V(Z_{c})\propto Z_{c}$
and
$E_{s}\propto \unicode[STIX]{x1D6F4}(Z_{c})\propto Z_{c}^{2}$
(details in appendix C). Furthermore, the excess surface energy
$E_{s}$
decreases from
$\unicode[STIX]{x1D70E}\unicode[STIX]{x1D6F4}^{\ast }$
to 0, and the gravitational energy
$E_{g}$
increases from 0 to
$\unicode[STIX]{x1D70C}_{m}\bar{L}gV^{\ast }$
(cf. figure 4). Here,
$g=9.8~\text{m}~\text{s}^{-2}$
is the gravitational acceleration,
$\unicode[STIX]{x1D70C}_{m}=\unicode[STIX]{x1D70C}_{l}(1-\unicode[STIX]{x1D719})+\unicode[STIX]{x1D70C}_{g}\unicode[STIX]{x1D719}$
is estimated as the density of the mixture and
$\bar{L}={\mathcal{V}}_{m}/[\unicode[STIX]{x03C0}{R_{d}^{\ast }}^{2}(1-\unicode[STIX]{x1D719})]=LR_{d}/R_{d}^{\ast }$
is the average mixing layer thickness. Also,
$\unicode[STIX]{x1D6F4}^{\ast }=\unicode[STIX]{x1D6F4}(t^{\ast })$
and
$V^{\ast }=\unicode[STIX]{x03C0}\int _{0}^{R_{0}}z(r,t^{\ast })\,\text{d}r^{2}$
are the crater area and the crater volume at
$t^{\ast }$
or when
$Z_{c}=Z_{c}^{\ast }$
. Let
$\unicode[STIX]{x0394}E_{s}=\unicode[STIX]{x1D70E}\unicode[STIX]{x1D6F4}^{\ast }$
and
$\unicode[STIX]{x0394}E_{g}=\unicode[STIX]{x1D70C}_{m}\bar{L}gV^{\ast }$
denote the total virtual change of
$E_{s}$
and
$E_{g}$
respectively. Both energies can then be written explicitly as a function of
$Z_{c}$
$E_{g}=\unicode[STIX]{x0394}E_{g}(1-Z_{c}/Z_{c}^{\ast })$
and
$E_{s}=\unicode[STIX]{x0394}E_{s}(Z_{c}/Z_{c}^{\ast })^{2}$
. Substituting this into (4.1) results in the lifting criterion,
$$\begin{eqnarray}\frac{\unicode[STIX]{x0394}E_{s}}{\unicode[STIX]{x0394}E_{g}}>0.5.\end{eqnarray}$$
Note that
$\unicode[STIX]{x0394}E_{s}$
and
$\unicode[STIX]{x0394}E_{g}$
can be computed from the measured crater profile at
$t^{\ast }$
.
A sketch of the virtual lifting process. The lifting process is denoted by the decrease of
$Z_{c}$
. The lifting may occur if the total energy
$E_{s}+E_{g}$
decreases at the beginning of the lifting, i.e. equation (4.1).
The data in regimes b and c in figure 3 are plotted versus the ratio of the virtual change of the surface energy,
$\unicode[STIX]{x0394}E_{s}$
, and the gravitational energy,
$\unicode[STIX]{x0394}E_{g}$
. See text for their definitions. The same symbols and colours as that in figure 3 are used here. The dashed line indicates the predicted lifting criterion,
$\unicode[STIX]{x0394}E_{s}/\unicode[STIX]{x0394}E_{g}=0.5$
, which separates the pie and pancake shapes in the overlapping range of mixing ratio.
The measured energy ratio of the data in regime b and regime c in figure 3 is plotted against the mixing ratio in figure 5. It can be seen that the data points corresponding to the pancake shape in regime b (
$0.55<{\mathcal{V}}_{m}/V_{d}<1$
) are at or below the proposed lifting criterion
$\unicode[STIX]{x0394}E_{s}/\unicode[STIX]{x0394}E_{g}=0.5$
, whereas the pie data all lie well above. For a given mixing ratio the energy ratio is proportional to the product of two impact parameters: the area-to-volume ratio,
$\unicode[STIX]{x1D6F4}^{\ast }/V^{\ast }$
, and squared the maximum expansion radius of the droplet,
${R_{d}^{\ast }}^{2}$
, i.e.
$\unicode[STIX]{x0394}E_{s}/\unicode[STIX]{x0394}E_{g}\propto \unicode[STIX]{x1D6F4}^{\ast }{R_{d}^{\ast }}^{2}/V^{\ast }$
. Both
$\unicode[STIX]{x1D6F4}^{\ast }/V^{\ast }$
and
$R_{d}^{\ast }$
increase with impact speed
$U$
(refer to appendix C). Therefore, the impacts with large
$d_{g}$
and/or higher wettability but small
$U$
in regime b fail to meet the lifting criterion.
With the established criterion we could discuss its possible application to the raindrop impression mentioned in the beginning. In nature, raindrops fall at their terminal velocity. For a droplet radius of 1.4 mm, the terminal velocity is approximately
$7~\text{m}~\text{s}^{-1}$
. First, we evaluate the mixing ratio as a function of grain size
$d_{g}$
using (3.1). Parameters used in the evaluation are the packing density of the target
$\unicode[STIX]{x1D719}=0.59$
, the contact angle of the grains
$\cos \unicode[STIX]{x1D703}_{c}=0.4$
and the maximum spreading radius of the droplet
$R_{d}^{\ast }/R_{d}=4$
. It is possible to reach regime b of the mixing ratio for the grain size range of
$35{-}75~\unicode[STIX]{x03BC}\text{m}$
. Silt falls in such a size range (Boggs Reference Boggs1995), and its mineral origin is quartz and feldspar, lighter than ceramic beads used here. Therefore, the density of the mixture
$\unicode[STIX]{x1D70C}_{m}$
would be smaller. Furthermore, the larger impact speed results in a larger
$R_{d}^{\ast }$
, and the area volume ratio,
$\unicode[STIX]{x1D6F4}^{\ast }/V^{\ast }$
would be comparable to the impacts shown here (appendix C). It is therefore very likely that the resultant
$\unicode[STIX]{x0394}E_{s}/\unicode[STIX]{x0394}E_{g}$
would have satisfied the lift criterion. In conclusion, it is possible for raindrop impact on a substrate like silt to result in the pie morphology with a underneath cavity. It is however worth pointing out that the curvature of the final pie shape is unlikely to exceed that of the liquid–grain mixture at
$t^{\ast }$
(although the sign is opposite). The radius of curvature of the final morphology is larger than both the horizontal and vertical dimensions,
$R_{d}^{\ast }$
and
$Z_{c}^{\ast }$
. Therefore, the underneath cavity created by droplet impacts should be well distinguishable from that resulting from a nearly spherical bubble in the mud.
5 Conclusion
In this work, we first discuss the discovery of a counterintuitive cavity between the impacting droplet and the deformed granular target using X-ray tomography in figure 2. Its presence in the parameter space is further characterized by the distinct crater profile evolution measured by laser profilometry. It turns out that the morphologies of the liquid–grain residual develop depending on the mixing ratio between the liquid and grains. The cavity creation phenomenon is located in the regime of parameter space where the mixing ratio is intermediate (cf. figure 3). In this regime the circumference curvature of the deformed droplet is suppressed by liquid–grain mixing, while the central curvature of the liquid layer on the top leads to the lift of the liquid–grain mixture and the creation of the cavity underneath. The dual-curvature model gives the necessary condition for the lifting, but turns out not to be sufficient, i.e. there are experiments in this regime that do not show any lifting. The missing condition is the force balance. We derive a simple lifting criterion in terms of experimentally accessible quantities in (4.2), which successfully separates lifting from the no-lifting morphologies in this regime (cf. figure 5). Towards the end of this exploration we would like to point out one more detail that needs further attention. In the derivation of the lifting criterion in (4.2) we assume a parabolic crater profile. A parabola has the maximum curvature at its centre which implies that the lifting process is led by the centre of the mixture. In practice, the profile may deviate from a parabola and may have the maximum curvature away from centre (refer to appendix D). In spite of the difference between the model and experiments, the validity of the resultant criterion appears to be robust.
Acknowledgements
The authors thank Dr M. Schröter and Professor S. Herminghaus for offering the opportunity to perform X-ray tomography measurements at the Max-Planck Institute of Dynamics and Self-Organization. This work is financed by the Netherlands Organization for Scientific Research (NWO) through VIDI grant no. 68047512.
Supplementary movies
Supplementary movies are available at https://doi.org/10.1017/jfm.2019.742.
Appendix A. X-ray tomography
The X-ray tomography scanning was performed with the Nanotom® at the Max-Planck Institute of Dynamics and Self-Organization in Göttingen, Germany. The droplet impacts shown in figure 2(g–i) were conducted on substrates consisting of soda-lime beads with various grain sizes, droplet diameters and impact speeds. Figures 2(g) and 2(h) are of the same grain size
$d_{g}\sim 90~\unicode[STIX]{x03BC}\text{m}$
and the same droplet diameter 2.73 mm, but with impact speeds
$U=2.6~\text{m}~\text{s}^{-1}$
and
$U=4.2~\text{m}~\text{s}^{-1}$
. Figure 2(i) is a scan of an impact with
$U=4.2~\text{m}~\text{s}^{-1}$
,
$d_{g}\sim 250~\unicode[STIX]{x03BC}\text{m}$
and droplet diameter 3.2 mm. To distinguish the liquid–grain mixture and the rest of the substrate the impacting droplet is mixed with a suspension of caesium nanoparticles which appears bright in the tomogram. The liquid–grain residual is highlighted in figure 2(g–i) accordingly. Note that the initial packing density may influence the morphology. For the tomograms shown in figure 2(g–i) the initial packing density is 0.55, 0.57 and 0.58 respectively. In figure 6 a cross-section of a tomogram is shown for the same impact parameters as figure 2(h), but with a reduced initial packing density of 0.55. Similar to figure 2(h) a cavity can be seen. The details of the morphology display differences, for instance, the centre is concave.
A vertical cross-section of a tomogram scan for water droplet impact on a packing of soda-lime beads. The impact speed, the grain size and the droplet size are the same as figure 2(h). However, the packing density is lower, namely equal to 0.55.
Appendix B. The evaluation of mixing ratio
To evaluate the volume of the liquid mixed with grains,
${\mathcal{V}}_{m}$
, we use Darcy’s law to model the penetration of the droplet into the granular substrate,
$$\begin{eqnarray}\boldsymbol{Q}=\frac{\unicode[STIX]{x1D705}{\mathcal{A}}_{c}}{\unicode[STIX]{x1D707}_{l}}\unicode[STIX]{x1D735}P.\end{eqnarray}$$
In the above equation the permeability of the substrate,
$\unicode[STIX]{x1D705}=(1-\unicode[STIX]{x1D719})^{3}d_{g}^{2}/(180\unicode[STIX]{x1D719}^{2})$
, is defined by the Carman–Kozeny relation (Carman Reference Carman1956), where
$\unicode[STIX]{x1D735}P$
is the pressure gradient,
${\mathcal{A}}_{c}$
is the contact area between the droplet and the substrate and
$\unicode[STIX]{x1D707}_{l}$
is the dynamic viscosity of the liquid. Since the pressure gradient is mainly in the vertical direction, equation (B 1) can be reduced to a scalar equation. Denote the penetration depth of the liquid as
$L$
, and the pressure gradient can be estimated as
$P/L$
. On the other hand, volume conservation gives
$(1-\unicode[STIX]{x1D719}){\mathcal{A}}_{c}\,\text{d}L/\text{d}t=\text{d}{\mathcal{V}}_{m}/\text{d}t=|\boldsymbol{Q}|$
, when taking the presence of the grains in the mixing layer into account. All together (B 1) becomes an ordinary differential equation for the mixing layer thickness
$L$
with respect to time
$t$
,
$$\begin{eqnarray}\frac{\text{d}L}{\text{d}t}=\frac{\unicode[STIX]{x1D705}P}{\unicode[STIX]{x1D707}_{l}(1-\unicode[STIX]{x1D719})L}.\end{eqnarray}$$
Its solution is (3.1) in the main text.
The driving pressure
$P$
consists of two components: the capillary pressure
$P_{c}$
and the inertial pressure
$P_{i}$
. The capillary pressure is independent of the details of the impact dynamics as soon as the liquid contacts the sand and can be obtained from experimental parameters,
$P_{c}=4\unicode[STIX]{x1D70E}\cos \unicode[STIX]{x1D703}_{c}/d_{c}$
, where
$d_{c}=(2(1-\unicode[STIX]{x1D719})/3\unicode[STIX]{x1D719})d_{g}$
is the average diameter of capillaries between grains derived from the Carman–Kozeny relation. On the other hand, the inertial pressure needs to be corrected with the deformation of the substrate,
$P_{i}=\unicode[STIX]{x1D70C}_{l}U^{2}R_{d}/(R_{d}+Z_{c}^{\ast })$
, as the deceleration of the droplet is reduced in comparison with the impact on a solid surface. Here,
$\unicode[STIX]{x1D70C}_{l}$
is the liquid density, and
$Z_{c}^{\ast }$
is the maximum vertical deformation of the substrate measured by the crater profile. This correction relates to the effective energy transfer during impact (Zhao et al.
Reference Zhao, de Jong and van der Meer2015b
) and has been validated by the analysis of
$R_{d}^{\ast }$
(Zhao et al.
Reference Zhao, de Jong and van der Meer2015b
, Reference Zhao, de Jong and van der Meer2017) and crater formations (de Jong et al.
Reference de Jong, Zhao and van der Meer2017). Another difference between
$P_{i}$
and
$P_{c}$
is their action period;
$P_{i}$
acts during an inertial time scale,
$\unicode[STIX]{x1D70F}_{i}=2(R_{d}+Z_{c}^{\ast })/U$
. In contrast,
$P_{c}$
lasts as long as the contact between liquid and grains exists. As we are interested in the liquid–grain mixing until the maximum droplet spreading is reached, the spreading time is estimated as half of the intrinsic oscillation time of the droplet,
$\unicode[STIX]{x1D70F}_{c}={\textstyle \frac{1}{2}}\sqrt{(\unicode[STIX]{x03C0}/6)(\unicode[STIX]{x1D70C}_{l}D_{0}^{3}/\unicode[STIX]{x1D70E})}$
, where
$\unicode[STIX]{x1D70E}$
is the surface tension of the liquid. These two time scales provide relative weights for
$P_{i}$
and
$P_{c}$
in the spreading phase of the droplet, and the average effect of the total pressure is evaluated as
$P=(\unicode[STIX]{x1D70F}_{i}/\unicode[STIX]{x1D70F}_{c})P_{i}+P_{c}$
(Zhao et al.
Reference Zhao, de Jong and van der Meer2017). The two components of
$P$
also correspond to two length scales;
$R_{d}$
is the length scale of the action of the inertial pressure. Due to droplet spreading, the capillary pressure acts over a larger length scale of
$R_{d}^{\ast }$
. The contact area
${\mathcal{A}}_{c}$
used in the main text accounts for this difference.
Appendix C. Virtual lifting process
In the main text we consider a virtual lifting process, where the crater profile is a parabola,
$$\begin{eqnarray}z(r,Z_{c})=Z_{c}\left(\frac{r^{2}}{R_{0}^{2}}-1\right).\end{eqnarray}$$
The centre depth
$Z_{c}$
is decreased during the lifting, while the position where the profile intersects with 0,
$r=R_{0}$
, is fixed. This process is sketched in figure 7.
The profile represents the top surface of the liquid–grain mixture. The area of the top surface can be written as an explicit function of
$Z_{c}$
,
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6F4}(Z_{c}) & = & \displaystyle \unicode[STIX]{x03C0}\int _{0}^{R_{0}}\sqrt{1+z^{\prime 2}}\,\text{d}r^{2}-\unicode[STIX]{x03C0}R_{0}^{2}\nonumber\\ \displaystyle & = & \displaystyle \frac{\unicode[STIX]{x03C0}R_{0}^{4}}{6Z_{c}^{2}}\left[\left(1+4\frac{Z_{c}^{2}}{R_{0}^{2}}\right)^{3/2}-1\right]-\unicode[STIX]{x03C0}R_{0}^{2}\nonumber\\ \displaystyle & {\approx} & \displaystyle 2\unicode[STIX]{x03C0}Z_{c}^{2}.\end{eqnarray}$$
The last step is obtained by the expansion to the second order of
$Z_{c}^{2}/R_{0}^{2}$
, when the crater is shallow, i.e.
$4Z_{c}^{2}\ll R_{0}^{2}$
. The corresponding surface energy
$E_{s}=\unicode[STIX]{x1D70E}\unicode[STIX]{x1D6F4}(Z_{c})$
is thus quadratic in
$Z_{c}$
.
On the other hand, the lift of the mixture increases the gravitational energy,
$E_{g}(Z_{c})$
. Denote the density of the liquid–grain mixture as
$\unicode[STIX]{x1D70C}_{m}$
, and assume the average thickness of the mixture,
$\bar{L}$
, is constant during lifting. Given
$E_{g}(Z_{c}^{\ast })=0$
, the gravitational energy can be computed by the integral of volume elements lifted from
$z(r,Z_{c}^{\ast })$
to
$z(r,Z_{c})$
$$\begin{eqnarray}\displaystyle E_{g}(Z_{c}) & = & \displaystyle \unicode[STIX]{x1D70C}_{m}g\bar{L}\unicode[STIX]{x03C0}\int _{r=0}^{r=R_{0}}[z(r,Z_{c})-z(r,Z_{c}^{\ast })]\,\text{d}r^{2}\nonumber\\ \displaystyle & = & \displaystyle \unicode[STIX]{x1D70C}_{m}g\bar{L}[V(Z_{c}^{\ast })-V(Z_{c})].\end{eqnarray}$$
Substituting the parabolic profile in (C 1) into the above definition gives
$$\begin{eqnarray}E_{g}(Z_{c})=\unicode[STIX]{x1D70C}_{m}g\bar{L}\frac{\unicode[STIX]{x03C0}}{2}R_{0}^{2}(Z_{c}^{\ast }-Z_{c}).\end{eqnarray}$$
The gravitational energy in (C 4) is seen to be linear with
$Z_{c}$
.
A sketch of the virtual lifting process.
Applying (C 2) and (C 4) to the lifting criterion (see the main text),
$$\begin{eqnarray}\left.\frac{\text{d}E_{s}+\text{d}E_{g}}{\text{d}Z_{c}}\right|_{Z_{c}=Z_{c}^{\ast }}>0,\end{eqnarray}$$
one obtains the lift criterion in terms of total virtual change of
$E_{s}$
and
$E_{g}$
,
$\unicode[STIX]{x0394}E_{s}/$
$\unicode[STIX]{x0394}E_{g}>0.5$
. One may easily check that the obtained criterion is equivalent to the force condition where the surface tension force exceeds that of gravity.
This energy ratio can be further written explicitly,
$$\begin{eqnarray}\displaystyle \frac{\unicode[STIX]{x0394}E_{s}}{\unicode[STIX]{x0394}E_{g}} & = & \displaystyle \frac{\unicode[STIX]{x1D70E}}{\unicode[STIX]{x1D70C}_{m}g}\frac{\unicode[STIX]{x1D6F4}(Z_{c}^{\ast })}{V(Z_{c}^{\ast })}\frac{1}{\bar{L}}\nonumber\\ \displaystyle & \propto & \displaystyle \frac{\unicode[STIX]{x1D6F4}(Z_{c}^{\ast })}{V(Z_{c}^{\ast })}\frac{{R_{d}^{\ast }}^{2}}{{\mathcal{V}}_{m}}=\frac{\unicode[STIX]{x1D6F4}^{\ast }}{V^{\ast }}\frac{{R_{d}^{\ast }}^{2}}{{\mathcal{V}}_{m}}.\end{eqnarray}$$
Both
$\unicode[STIX]{x1D6F4}^{\ast }$
and
$V^{\ast }$
are measurable in experiments. The average thickness of the liquid–grain mixture is evaluated,
$\bar{L}={\mathcal{V}}_{m}/{R_{d}^{\ast }}^{2}$
. For a given mixing ratio
${\mathcal{V}}_{m}$
, the larger the spreading of the droplet,
$R_{d}^{\ast }$
, and/or the larger area volume ratio,
$\unicode[STIX]{x1D6F4}^{\ast }/V^{\ast }$
, the larger the energy ratio is. The value of
$R_{d}^{\ast }$
is a monotonic increasing function of the impact speed
$U$
. The area volume ratio measured in experiments increases with
$U$
for small impact speed but largely saturates for high speed (cf. figure 8). Therefore, roughly speaking, for a given mixing ratio,
$\unicode[STIX]{x0394}E_{s}/\unicode[STIX]{x0394}E_{g}$
increases with
$U$
.
The analysis above is performed on a parabolic crater shape and its virtual lifting process. For axisymmetry shapes of higher order the lifting path is expected to be more complex. Nevertheless, a similar procedure could still be applied to obtain the lifting criterion. The resultant lifting criterion in that case may however not be expressed as the ratio of the total virtual change of energies.
The ratio of the area and volume
$\unicode[STIX]{x1D6F4}^{\ast }/V^{\ast }$
is measured from the experimental profile at
$t^{\ast }$
. The scattering for a given
$U$
reflects the influence of the varying packing fraction of the granular target.
Appendix D. Experimental lifting process
(a) Snapshots of the crater height profile at the moment reaching maximum depth and 5 ms and 10 ms afterwards. (b) Time evolution of the curvature in the
$(r,z)$
-plane in the centre and at a distance equal to the original droplet radius,
$r=R_{d}$
. The data are from an experiment in regime b for a water droplet impact with
$U=4.1~\text{m}~\text{s}^{-1}$
,
$R_{d}=1.4~\text{mm}$
,
$d_{g}=98~\unicode[STIX]{x03BC}\text{m}$
and
$\cos \unicode[STIX]{x1D703}=0.6$
. The crater profiles shown here are fully resolved by the laser profilometry and therefore do not contain any fitting.
In practice the lifting process may be different from that sketched in figure 7, and the profile may deviate from a parabola. An example can be seen in figure 9, where the curvature around
$r=R_{d}$
is larger than that at the centre. As a consequence, lifting is initiated by decreasing the curvature at
$R_{d}$
, and the central curvature is enhanced to further lift the residual. This lifting process is different from the virtual one depicted in figure 7. However, this difference is not significant enough to invalidate the lifting criterion introduced here.
Open access
