1. Introduction
The hydrodynamic deformation of a droplet driven by a laser pulse is a process of high interest for the nanolithography industry, in which the state-of-the-art technique for producing extreme ultraviolet (EUV) light involves deforming a tin droplet into a thin sheet that is subsequently ablated into a light-emitting plasma (Bakshi Reference Bakshi2009; Fomenkov et al. Reference Fomenkov2017; Versolato Reference Versolato2019; Sizyuk & Hassanein Reference Sizyuk and Hassanein2020). Similar deformation processes can also be found in scenarios without a laser, in which the deformation is triggered by the droplet impact onto a large solid surface (Scheller & Bousfield Reference Scheller and Bousfield1995; Bergeron et al. Reference Bergeron, Bonn, Martin and Vovelle2000; Josserand & Thoroddsen Reference Josserand and Thoroddsen2016) or a solid pillar (Rozhkov, Prunet-Foch & Vignes-Adler Reference Rozhkov, Prunet-Foch and Vignes-Adler2002; Villermaux & Bossa Reference Villermaux and Bossa2011; Wang & Bourouiba Reference Wang and Bourouiba2017) for example. In these cases, the droplet is also expanded into a thin sheet, which eventually can fragment or retract due to the effects of surface tension. These processes are also industrially relevant in a wide range of applications, such as spray coating and printing (Derby Reference Derby2010; Lohse Reference Lohse2022). Hence, understanding their fluid mechanics is essential for further optimising current systems.
In a typical nanolithography scenario, a spherical tin droplet is irradiated by a nanosecond laser pulse in vacuum. This interaction ablates part of the droplet surface, which becomes a tin plasma that rapidly expands also on the laser pulse time scale (
$\tau _p \sim 10$
ns). As this plasma expands, it transfers momentum to the liquid droplet, which propels with a velocity
$U_z \sim 100\,{\textrm{m s}}^{-1}$
and expands radially with velocity
$\dot {R}_0 \sim 100\,{\textrm{m s}}^{-1}$
. Therefore, this expansion happens hydrodynamically on the inertial time scale (
$\tau _i = R_0 / \dot {R}_0 \sim 100$
ns), with
$R_0 \sim 10\,\unicode{x03BC} {\textrm{m}}$
being the droplet radius. In the third stage, when the droplet is already significantly deformed into a sheet, surface tension slows down its expansion, eventually forming a rim, fragmentation and sheet retraction. This final stage happens on the capillary time scale
$\tau _c = \sqrt {\rho R_0^3/\gamma } \sim 10\,\unicode{x03BC}$
s, where
$\rho = 6900\,{\textrm{kg m}}^{-3}$
is the liquid tin density and
$\gamma = 0.55\,$
N m–1 its surface tension (Liu et al. Reference Liu, Kurilovich, Gelderblom and Versolato2020).
Previous works by Klein et al. (Reference Klein, Bouwhuis, Visser, Lhuissier, Sun, Snoeijer, Villermaux, Lohse and Gelderblom2015) and Gelderblom et al. (Reference Gelderblom, Lhuissier, Klein, Bouwhuis, Lohse, Villermaux and Snoeijer2016) have numerically simulated laser-induced droplet deformation by assuming that the laser (and the generated plasma) interact with the droplet by applying a near-instantaneous pressure on its surface, setting it in incompressible motion. As illustrated in figure 1, this pressure is commonly approximated by a function
$f(\theta )$
defined over the parametrised surface of the droplet, where
$\theta = 0$
and
$\theta =\pi$
correspond to the side of the droplet hit by the laser and the opposite side, respectively. Different profiles have been suggested for this pressure function, most commonly a Gaussian profile is used, in which the standard deviation
$\sigma$
is tuned to mimic variations in the focus size of the laser beam. A truncated cosine has also been suggested by Gelderblom et al. (Reference Gelderblom, Lhuissier, Klein, Bouwhuis, Lohse, Villermaux and Snoeijer2016), and the authors show that this is the only profile capable of generating a perfectly symmetric sheet.
Figure 1. Initialisation strategy to obtain a velocity field from a pressure profile given on the droplet surface. (a) Sketch of the droplet and laser. (b) Examples of possible initial pressure profile functions on the droplet surface. (c) Radial and axial velocity components on three sheets obtained from different pressure profile types. The values of
$u_z$
represent the droplet velocity relative to its centre-of-mass translational velocity.
Kurilovich et al. (Reference Kurilovich, Basko, Kim, Torretti, Schupp, Visschers, Scheers, Hoekstra, Ubachs and Versolato2018) and Hernandez-Rueda et al. (Reference Hernandez-Rueda, Liu, Hemminga, Mostafa, Meijer, Kurilovich, Basko, Gelderblom, Sheil and Versolato2022) utilise the full radiation-hydrodynamic code RALEF-2D to simulate the interaction between a laser and the tin, leading to the plasma formation that transfers momentum to the droplet. While their code could not simulate late-time droplet deformation, the early time interaction between the laser and tin was reported, and the authors observed pressure profiles that can match Gaussian functions in some aspects, such as propulsion-to-deformation energy ratios. However, no attention is given to the actual morphology of the deformed sheet coming from these pressure profiles.
While Gaussians and the truncated cosine profiles have been successfully used for the purpose of these previous studies, they were not capable of reproducing one striking feature of some experiments: a sheet that curves in the opposite direction to the incident laser beam, as also illustrated in figure 1. This sheet morphology, which we refer to as a ‘positive curvature’, is intriguing as it is characteristic, for example, of a droplet being pushed by a constant wind tunnel, which is not the case in a near-instantaneous focused laser interaction. Depending on the experiment’s parameters, a ‘negative’ curvature can also be observed. This naturally leads to the following question: Is it possible to obtain a single functional form for the pressure profile to explain positive and negative curvatures, as we see both in experiments? In this work we go beyond the early time dynamics and show that this is possible by suggesting a pressure profile based on a raised cosine function. As will be demonstrated later, this profile was chosen due to its capability of mimicking laser beam focus as the Gaussian profile and providing a forward curvature that is more comparable with experiments. Correctly predicting the morphology and curvature of the sheet provides important information for EUV-light generation applications. The shape and thickness of the sheet will directly influence how much energy is necessary to ablate the liquid and generate light, such good predictions of the sheet morphology will allow for better tuning for the optimal experimental parameters. It has also been shown by Engels et al. (Reference Engels, Meijer, Schubert, van der Zande, Ubachs and Versolato2023) that the vaporisation of a thin sheet creates a vapour with a shape that follows the local curvature of the sheet, indicating that controlling the morphology of the sheet is also essential for shaping the vaporisation (and plasma creation) process.
The raised cosine profile is illustrated in figure 1, along with the more classical Gaussian and cosine profiles. We also show examples of typical sheet shapes that can be numerically obtained using the three pressure profiles discussed in this paper. In the Gaussian example, a sheet that curves towards the laser is shown (negative curvature). In the raised cosine example, the sheet curves away from the laser origin (positive curvature). The centre image in figure 1(c) shows the truncated cosine case that presents a symmetric sheet and, according to our definition, presents zero curvature. The expressions for the three pressure profiles are shown in table 1. Similarly to the Gaussian function, the proposed raised cosine also has a tuneable parameter
$W$
, which controls the width of the profile.
Table 1. Different initial pressure profiles used throughout this work. In the cosine and raised cosine expressions,
$H$
represents the standard Heaviside step function.
While in this work we obtain positive curvatures by using raised cosine pressure profiles, we note that this is not the only function capable of giving this sheet morphology. Other functions can also be chosen or carefully constructed to obtain similar deformation dynamics. During an exploratory study of different pressure profiles, we observed that, in order to obtain a realistic positive curvature, the following characteristics should be present in a chosen function: first, the function should be peaked at
$\theta = 0$
; second, the peak should be wide enough so that not all the pressure is focusing at
$\theta = 0$
, which would give negative curvature; and third, while the peak needs to be wide, the function cannot fully wrap the droplet to
$\theta = \pi$
, which would stop expansion. This requires the function to have a wide peak, but short tail. The second condition can be easily achieved with the classical Gaussian profile by increasing the
$\sigma$
parameter, which widens the peak. However, Gaussian profiles have long tails, such that widening the peak would automatically break the third condition, as the profile would wrap the whole droplet. The ‘tail size’ of a distribution is often characterised by excess kurtosis, which is the distribution’s kurtosis minus the reference Gaussian kurtosis. A Gaussian function has, therefore, an excess kurtosis of exactly 0 and is referred to as a mesokurtic function. Leptokurtic functions are those with a positive excess kurtosis, and these are not desired in this study since they lead to fat-tailed profiles. Platokurtic functions have negative excess kurtosis, which leads to a short tail. Therefore, platokurtic functions are optimal for our study in order to guarantee the third condition above. Out of many known statistical distributions, the raised cosine is one of the few smooth functions that provide a constant and negative excess kurtosis (
$ \textit{Kurt} = -0.59$
) that does not depend on the parameter
$W$
. Due to its simplicity, while providing the necessary kurtosis condition, we chose to use the raised cosine in this study. Throughout this work we show that the raised cosine can, indeed, be a suitable choice to simulate situations in which the deformed sheet presents a positive curvature.
This paper is organised as follows. Section 2 introduces the experimental set-up used. In § 3 the problem is mathematically stated, and the numerical approach to simulate the droplet deformation is described. Section 4 discusses the numerical and experimental results obtained. Finally, § 5 concludes the results and presents future perspectives. Additional numerical and experimental details can be found in the appendices.
2. Experiment
Our experimental set-up has previously been described in detail (Kurilovich et al. Reference Kurilovich, Klein, Torretti, Lassise, Hoekstra, Ubachs, Gelderblom and Versolato2016; Kurilovich Reference Kurilovich2019; Liu et al. Reference Liu, Meijer, Hernandez-Rueda, Kurilovich, Mazzotta, Witte and Versolato2021). In short, in these experiments, a kilohertz sequence of liquid tin droplets (270
$^{\circ }\textrm{C}$
) is vertically dispensed in a vacuum environment (
$10^{-7}$
mbar) by a droplet generator. The pressure and frequency of the generator are adjusted to obtain droplets with diameters ranging from
${27}\,{\unicode{x03BC} \textrm{m}}$
to
${59}\,{\unicode{x03BC} \textrm{m}}$
. The droplets fall at approximately
${10}{\,{\textrm{m}} \, {\textrm{s}}^{-1}}$
through a horizontal light sheet generated by a continuous-wave HeNe laser. The scattered light is then detected by a photomultiplier tube and downsampled to
${10}\,{\textrm{Hz}}$
, serving as the start trigger for the experiment.
Firstly, a droplet is hit by the Gaussian intensity pre-pulse (
$\lambda = {1064}\,{\textrm{nm}}$
, a full width at half-maximum value (FWHM) of 10 ns, circularly polarised). This laser interacts with the droplet by ablating part of its surface and generating a plasma, exerting pressure on the remaining liquid tin. This pressure rapidly propels and expands the tin droplet on the order of 100
${\textrm{m}} \, {\textrm{s}}^{-1}$
to a thin axisymmetric sheet (Kurilovich et al. Reference Kurilovich, Klein, Torretti, Lassise, Hoekstra, Ubachs, Gelderblom and Versolato2016, Reference Kurilovich, Basko, Kim, Torretti, Schupp, Visschers, Scheers, Hoekstra, Ubachs and Versolato2018). The propulsion velocity
$U_z$
aligns with the propagation direction of the laser, while the expansion occurs radially with an initial velocity
$\dot {R}_{0}$
. The radial expansion speed subsequently decreases due to the surface tension that is exerted on the edge of the sheet (Villermaux & Bossa Reference Villermaux and Bossa2011; Gelderblom et al. Reference Gelderblom, Lhuissier, Klein, Bouwhuis, Lohse, Villermaux and Snoeijer2016; Kurilovich et al. Reference Kurilovich, Klein, Torretti, Lassise, Hoekstra, Ubachs, Gelderblom and Versolato2016).
The time scales governing propulsion and expansion accelerations are similar to the duration of the laser pulse (ns) and are much shorter than the time scale of the subsequent fluid dynamics deformation (Gelderblom et al. Reference Gelderblom, Lhuissier, Klein, Bouwhuis, Lohse, Villermaux and Snoeijer2016; Kurilovich et al. Reference Kurilovich, Klein, Torretti, Lassise, Hoekstra, Ubachs, Gelderblom and Versolato2016). Figure 2(a) illustrates the typical experimental response of droplets to the pre-pulse impact, characterised by a beam width of approximately 100
$\unicode{x03BC} {\textrm{m}}$
FWHM and showcasing a sheet that curves away from the laser beam. Conversely, figure 2(b) depicts the response to a pre-pulse with a beam diameter of approximately 20
$\unicode{x03BC} {\textrm{m}}$
, exhibiting negative curvature.
In order to observe the evolution of the liquid tin, we obtain shadowgraphs as depicted in figure 2. The shadowgraphy imaging set-up is described in detail by Kurilovich (Reference Kurilovich2019) and briefly summarised here. It consists of a dye-based illumination source and a charge-coupled camera coupled to a long-distance microscope, granting a spatial resolution of approximately 5
$\unicode{x03BC} {\textrm{m}}$
. The illumination source emits pulses with a duration of 5 ns (FWHM) and a spectral bandwidth of 12 nm (FWHM) centred at 560 nm. We utilise these shadowgraphy pulses (SP) for backlight illumination of the side-view acquisitions (at 90
$^{\circ }$
concerning the laser axis) to capture the curvature dynamics of the liquid tin sheet. We select a discrete number of time steps in the experiment and record 20 frames in a stroboscopic manner for each step, each frame representing a distinct laser--droplet interaction event. This methodology allows us to apply post-filtering techniques, such as selecting optimally aligned laser-to-droplet events.
Figure 2. Laser pulse schemes with specific irradiation geometries lead to either (a) positive or (b) negative curvature, visualised by shadowgraphs.
The non-dimensional parameters for each individual experiment can be seen in tables 3 and 4 in Appendix A.
3. Problem formulation and numerical approach
We perform fluid simulations to numerically simulate the laser-induced deformation of a tin droplet into a sheet. The governing equations for the isothermal incompressible bi-phase (droplet and ambient) flow are the continuity and momentum conservation, given by
\begin{align} \rho \left ( \frac {\partial \boldsymbol {u}}{\partial t} + \boldsymbol{\nabla }\boldsymbol{\cdot }( \boldsymbol {u}\boldsymbol {u} ) \right ) & = - \boldsymbol{\nabla }p + \boldsymbol{\nabla }\boldsymbol{\cdot }\left ( 2 \mu \boldsymbol {D} \right ) + \boldsymbol{f}_{\gamma }, \end{align}
\begin{align} \boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol {u} & = 0, \\[9pt] \nonumber \end{align}
where
$\boldsymbol {u}$
and
$p$
are the velocity and pressure fields,
$\boldsymbol {D} = [\boldsymbol{\nabla }\boldsymbol {u} + (\boldsymbol{\nabla }\boldsymbol {u})^T ]/2$
is the deformation rate tensor,
$\rho$
and
$\mu$
are the density and viscosity of the fluids, respectively. We note that, in this one-fluid formulation,
$\rho$
and
$\mu$
are functions with values that change across the droplet--ambient interface. The expression for these functions is given later in this paragraph. In the numerical method used here, the surface tension force is defined as a body force
$\boldsymbol{f}_{\gamma } = \gamma \kappa \delta _s \boldsymbol {n}$
, where
$\kappa$
is the local curvature of the interface,
$\gamma$
the constant surface tension coefficient,
$\boldsymbol {n}$
is the unit vector normal to the interface and
$\delta _s$
is the Dirac delta function centred on the interface (Popinet Reference Popinet2009; Tryggvason, Scardovelli & Zaleski Reference Tryggvason, Scardovelli and Zaleski2011). The droplet interface is tracked using a volume of fluid (VOF) scheme (Hirt & Nichols Reference Hirt and Nichols1981), in which a scalar colour function
$c(\boldsymbol {x}, t)$
indicates the fraction of droplet fluid contained in each numerical cell. The local density and viscosity are obtained by linearly interpolating using the value of
$c$
. So
$\rho$
and
$\mu$
from (3.1) are defined as
\begin{align} \rho (c) = & \ c \ \rho _d + (1 - c)\rho _a, \end{align}
\begin{align} \mu (c) = & \ c \ \mu _d + (1 - c)\mu _a, \end{align}
where the indices
$d$
and
$a$
refer to the properties of the droplet and ambient fluids, respectively. While in experiments the droplet is contained within a vacuum chamber, due to numerical limitations, we keep the ambient fluid properties set to
$\rho _a = 10^{-4} \rho _d$
and
$\mu _a = 10^{-4} \mu _d$
. We note that, in the experimental vacuum chamber, the actual density ratio is
$\rho _a \approx 10^{-13} \rho _d$
. In Appendix B we provide a small sweep of this ratio parameter to show that our solution is not significantly affected by it with our chosen value of
$10^{-4}$
.
Equations (3.1) and (3.2) can be non-dimensionalised by rescaling variables with the following choices:
\begin{equation} \boldsymbol {x} = R_0\bar {\boldsymbol {x}}, \hspace {10pt} t = \frac {R_0}{U_z}\bar {t}, \hspace {10pt} \boldsymbol {u} = U_z\bar {\boldsymbol {u}}, \hspace {10pt} p = \rho _d \,U_z^2\bar {p}, \hspace {10pt} \kappa = \frac {1}{R_0}\bar {\kappa }, \hspace {10pt} \delta _s = \frac {1}{R_0}\bar {\delta _s}.\end{equation}
Here
$R_0$
is the initial radius of the droplet and
$U_z$
is the droplet propulsion velocity obtained after the laser hit.
Substituting (3.5) into (3.1) and (3.2) , we obtain the non-dimensional version of the governing equations, given by
\begin{align} \bar {\rho }\left ( \frac {\partial \bar {\boldsymbol {u}}}{\partial t} + \boldsymbol{\nabla }\boldsymbol{\cdot }( \bar {\boldsymbol {u}}\bar {\boldsymbol {u}} ) \right ) & = - \boldsymbol{\nabla }\bar {p} + \frac {1}{\textit{Re}}\boldsymbol{\nabla }\boldsymbol{\cdot }\left ( 2 \bar {\mu } \bar {\boldsymbol{D}} \right ) + \frac {1}{\textit{We}} \bar {\kappa } \bar {\delta _s} \boldsymbol {n}, \end{align}
\begin{align} \boldsymbol{\nabla }\boldsymbol{\cdot }\bar {\boldsymbol {u}} & = 0, \end{align}
where
\begin{equation} \textit{Re} = \frac {\rho _d \, U_z \, R_0}{\mu _d}, \quad {\textrm{and}} \quad \textit{We} = \frac {\rho _d \, U_z^2 \, R_0}{\gamma } \end{equation}
are the Reynolds and Weber numbers, respectively. Later in the text we also use a version of these two numbers based on the droplet initial expansion velocity
$\dot {R}_0$
as a velocity scale, that is,
$\textit{Re}_{\textit{def}} = ({\rho _d \dot {R}_0 R_0}/{\mu _d})$
and
$\textit{We}_{\textit{def}} = ({\rho _d \dot {R}_0^2 R_0}/{\gamma })$
. The non-dimensional density and viscosity functions are obtained from (3.3) and (3.4) and are given by
\begin{align} \bar {\rho }(c) = & \ c + (1 - c)\rho _a/\rho _d, \end{align}
\begin{align} \bar {\mu }(c) = & \ c + (1 - c)\mu _a/\mu _d. \\[9pt] \nonumber \end{align}
Equations (3.6) and (3.7) require an appropriate velocity field as an initial condition to be solved. To obtain an initial condition for the velocity, we adopt the same strategy proposed by Gelderblom et al. (Reference Gelderblom, Lhuissier, Klein, Bouwhuis, Lohse, Villermaux and Snoeijer2016), which we briefly describe here. We know that the process is nearly inviscid and that the time scale of this laser-plasma kick
$(\tau _p \approx 10^{-8}\,{\textrm{s}})$
is much smaller than the inertial
$(\tau _i\approx 10^{-7}\,{\textrm{s}})$
and capillary
$(\tau _c\approx 10^{-5}\,{\textrm{s}})$
time scales, that is,
$\tau _{p} \ll \tau _{i} \ll \tau _{c}$
. While the droplet surface is pressurised by the plasma, shock waves travel across the droplet within the acoustic time scale
$\tau _a = R_0/c_s \approx 8\times 10^{-9}\,{\textrm{s}}$
, where
$c_s$
is the speed of sound in the liquid tin. It was shown by Grigoryev et al. (Reference Grigoryev2018), Meijer et al. (Reference Meijer, Kurilovich, Eikema, Versolato and Witte2022) that for short laser pulses, the focusing of these waves in the droplet’s centre can generate a cavitating bubble with dynamics that can affect the droplet morphology (also see Jalaal et al. Reference Jalaal, Li, Klein Schaarsberg, Qin and Lohse2019a
). However, the laser parameters chosen here fall within a range where these effects are not pronounced and, in fact, no effects of early time cavitation are observed within our droplets. Therefore, since the acoustic and laser-plasma time scales are much smaller than that of inertia, we assume that the tin plasma applies a near-instantaneous pressure kick on the droplet surface, which, in turn, instantaneously develops a pressure field within the whole droplet volume. This pressure kick will be provided as a function
$f(\theta )$
over the parametrised surface of the droplet as in figure 1. With these assumptions, the Navier–Stokes equations can be reduced into the following Laplace equation for the pressure:
\begin{equation} {\nabla} ^2 \bar {p} = 0. \end{equation}
Equation (3.11) is solved in spherical coordinates assuming symmetry in the azimuthal angle
$\varphi$
. Therefore, it is solved in the domain
$(r, \theta ) \in [0, 1] \times [0, \pi ]$
with boundary condition
$\bar {p}(1, \theta ) = f(\theta )$
. A solution to (3.11) can be obtained by decomposing the pressure field in Legendre polynomials
$P_n$
, which will give
\begin{equation} \bar {p}(r, \theta ) = \sum _{n=0}^{\infty } A_n r^n P_n( \cos \theta ), \end{equation}
with coefficients
\begin{equation} A_n = \frac {2n + 1}{2}\int _0^\pi f(\theta ) P_n( \cos \theta ) \sin \theta \ {\textrm{d}}\theta . \end{equation}
The initial velocity field within the droplet can then also be obtained from the simplified Navier–Stokes equations by
\begin{equation} \bar {\boldsymbol {u}}_{0} = - \alpha \boldsymbol{\nabla }\bar {p} = -\alpha \left [ \sum _{n=0}^{\infty } n A_n r^{n-1} P_n( \cos \theta ), \ \ - \frac {1}{r}\sum _{n=0}^{\infty } A_n r^n P^{\prime}_n( \cos \theta ) \sin \theta \right ], \end{equation}
where
$\alpha$
can be numerically tuned to change the magnitude of the velocity field. Since we opted to use the propulsion speed as our velocity scale in (3.5), we always take
$\alpha$
such that the non-dimensional propulsion velocity is
\begin{equation} \bar {U}_z = \frac {3}{2} \int _0^\pi f(\theta ) \cos \theta \sin \theta \ \text{d}\theta = 1. \end{equation}
The velocity field
$\bar {\boldsymbol {u}}_{0}$
from (3.14) is then used as an initial condition for (3.6) and (3.7).
Equations (3.6) and (3.7) are solved numerically using the open-source free-software language Basilisk C (Popinet & Collaborators Reference Popinet2013–2021). The droplet is created at the centre of a square domain
$[-10R_0, \ 10R_0] \times [-10R_0, \ 10R_0]$
that is fully discretised with a non-uniform quadtree grid (Popinet Reference Popinet2003, Reference Popinet2009). To accurately resolve the flow structure inside the droplet and its shape, we apply increased refinement levels for the liquid phase and also at the interface. The maximum quadtree level of refinement used is 14, resulting in grid cells with a minimum size of
$\varDelta = 20R_0/(2^{14}) = 0.0012R_0$
.
The numerical code then solves the governing equations using a projection method and a multilevel Poisson solver. We refer to the works of Popinet (Reference Popinet2009, Reference Popinet2015) for more details of the VOF implementation and Appendix C for validation tests of our code. We note that the software language Basilisk C used here has also been previously validated in many other works involving deformable surfaces, such as the works by Sanjay, Lohse & Jalaal (Reference Sanjay, Lohse and Jalaal2021) and by Popinet (Reference Popinet2009).
4. Results and discussion
4.1. Predicting curvature based on the initial velocity field
We begin our study by looking at the initial velocity calculated at time
$t = 0$
from (3.14) and attempting to predict how the sheet will curve after deformation.
In order to do this, we first subtract the centre-of-mass velocity from
$\bar {\boldsymbol {u}}_{0}$
such that we obtain only the velocity field responsible for the droplet deformation. From (3.14) and (3.15), this field is given by
\begin{equation} \bar {\boldsymbol {u}}_{\textit{def}}(r, \theta ) = \bar {U}_z\left [\cos \theta , \sin \theta \right ] + \! \left [\! - \sum _{n=0}^{\infty } n A_n r^{n-1} P_n( \cos \theta ), \frac {1}{r}\sum _{n=0}^{\infty } A_n r^n P^{\prime}_n( \cos \theta ) \sin (\theta )\! \right ]. \end{equation}
An indicative of how the droplet will deform is the radial component of
$\bar {\boldsymbol {u}}_{\textit{def}}$
at the droplet surface, that is,
\begin{equation} \bar {u}_{r_{\textit{def}}}(1, \theta ) = \bar {U}_z\cos \theta - \sum _{n=0}^{\infty } nA_nP_n( \cos \theta ). \end{equation}
Figures 3(a) and 3(b) show the value of the radial velocity (4.2) as a function of
$\theta$
for three Gaussian and three raised cosine profiles. Intuitively, we see negative velocities for low and high
$\theta$
, with positive values at intermediate
$\theta$
. This is expected since the droplet is being ‘squeezed’ in the laser-hit direction while it expands upwards. A striking feature of the Gaussian curves is that the high values of
$\bar {u}_{r_{\textit{def}}}$
always happen before
$\pi / 2$
, regardless of the choice for the parameter
$\sigma$
. This is also shown in figure 3(d–f), where the arrows represent
$\bar {u}_{r_{\textit{def}}}$
on the droplet surface. This observation indicates that the maximum radial expansion for a Gaussian profile will always happen behind the droplet centre of mass, such that a negative curvature is to be expected. We note that this had also previously been observed by Gelderblom et al. (Reference Gelderblom, Lhuissier, Klein, Bouwhuis, Lohse, Villermaux and Snoeijer2016), where the authors show a convergence of the direction of maximum expansion to an angle close to
$\pi /2$
.
For a raised cosine profile, on the other hand, we can see that the maximum value of
$\bar {u}_{r_{\textit{def}}}$
can be positioned beyond
$\pi /2$
if
$W$
is chosen to be large enough. This can be seen in figures 3(b), 3(g–i). Since the maximum radial expansion happens after the droplet centre of mass, we expect the deformed sheet to have a positive curvature. It is worth noting that the negative curvature can also be expected from raised cosine profiles, as long as
$W$
is chosen small enough. This highlights the versatility of this pressure profile choice, as both outcomes can be achieved with the correct tuning of a single parameter.
Figure 3. Initial droplet radial velocity
$\bar {u}_{r_{\textit{def}}}(r, \theta )$
in spherical coordinates at time
$t = 0$
for different profiles and profile parameters. (a–b) Velocity at the droplet surface over
$\theta$
for three Gaussian and three raised cosine profiles, respectively. The insets show the normalised pressure profiles over the droplet surface. (c) Value of
$\theta _{\textit{max}}$
as a function of profile parameter. Panels (d–i) show the
$\bar {u}_{r_{\textit{def}}}$
field within the droplet for the six cases shown in panels (a) and (b). The arrows represent the magnitude and sign of
$\bar {u}_{r_{\textit{def}}}(1, \theta )$
and the dashed line indicates the point of maximum
$\bar {u}_{r_{\textit{def}}}$
for each
$r$
.
To quantitatively predict if a profile can provide positive or negative curvature, we can calculate the angle where
$\bar {u}_{r_{\textit{def}}}(1, \theta )$
is maximised, let us call that angle
$\theta _{\textit{max}}$
. We can find this local maximum by derivating (4.2) and solving the equation
\begin{equation} {\bar{u}}^{\prime}_{r_{\textit{def}}}(1, \theta ) = - U_z\sin \theta + \sum _{n=0}^{\infty } n \sin \theta A_n P^{\prime}_n( \cos \theta ) = 0. \end{equation}
The solution to this equation is shown in figure 3(c) for the Gaussian and raised cosine profiles as a function of their respective parameters,
$\sigma$
and
$W$
. Once again, we see more clearly that the Gaussian solution cannot cross the
$\theta _{\textit{max}} = \pi /2$
limit, while the raised cosine does. We note that the Gaussian results are only shown up to
$\sigma \approx 1.2$
, since for higher values, the Gaussian pressure profile wraps the entire droplet such that no significant expansion happens.
We can also show that a Gaussian pressure profile can never provide
$\theta _{\textit{max}} \gt \pi /2$
. If we calculate the derivative (4.3) at the point
$\theta = \pi /2$
, we have
\begin{equation} {\bar{u}}^{\prime}_{r_{\textit{def}}}\left (1, \frac {\pi }{2}\right ) = - U_z + \sum _{n=0}^{\infty } n A_n P^{\prime}_n(0). \end{equation}
If the value of (4.4) is positive, we know the function is still increasing at
$\pi /2$
, such that the maximum velocity will happen at
$\theta _{\textit{max}}\gt \pi /2$
. Analogously, if the value is negative, we will have
$\theta _{\textit{max}}\lt \pi /2$
. If we substitute the Gaussian profile into
$A_n$
and
$\bar {U}_z$
in (4.4), we obtain
\begin{align} {\bar{u}}^{\prime}_{r_{\textit{def}}}\left (1, \frac {\pi }{2}\right ) &= - \frac {3}{2}\int _0^\pi \cos \theta \sin \theta e^{ -\theta ^2/(2\sigma ^2) }\text{d}\theta \nonumber\\ & \quad + \sum _{n=0}^{\infty } n \frac {2n + 1}{2} P^{\prime}_n(0) \int _0^\pi P_n( \cos \theta ) \sin \theta e^{ -\theta ^2/(2\sigma ^2) } \text{d}\theta , \end{align}
which is always negative.
The radial deformation velocity from (4.2) can also be defined for other fixed radial positions (instead of only at the droplet surface). By doing that for multiple values of
$r$
and storing the angle at which the maximum velocity is attained, we obtain the dashed lines shown in figure 3(d–i). This line indicates the direction of maximum radial expansion not only at the droplet surface but also within the fluid. This shows that the expansion direction smoothly varies from the droplet centre to the surface. An almost-flat curve is obtained for a Gaussian with
$\sigma = 1.2$
and
$W\approx 2$
, indicating that these droplets will expand almost perfectly vertically for all
$\theta$
cross-sections.
4.2. Curvature measurements after droplet deformation
The parameter
$\theta _{\textit{max}}$
was introduced as an initial condition indicator of the sheet curvature characteristics in the previous section. In this section we attempt to demonstrate that this parameter can indeed be correlated with the actual later-time sheet curvature. Since we have demonstrated in § 4.1 that Gaussian profiles can never provide positive curvature, we now focus our study only on raised cosine profiles.
Multiple simulations were performed varying the width parameter in the range
$W \in [1.0, 3.0]$
. Values of
$W$
lower than these result in numerical difficulties in the convergence of (3.14), while for higher values, no significant expansion can be observed, as the pressure fully wraps the droplet as also noted by Gelderblom et al. (Reference Gelderblom, Lhuissier, Klein, Bouwhuis, Lohse, Villermaux and Snoeijer2016) for wide Gaussian profiles. The deformation Reynolds and Weber numbers were always kept high enough (
$\textit{Re}_{\textit{def}}\gt 1000$
and
$\textit{We}_{\textit{def}} \gt 1000$
) such that no significant viscous effects are present, and significant capillary effects are only observed at later times in the simulations (particularly the formation of a rim at the edge of the sheet). None of these effects are expected to affect the bulk curvature formation of the sheet.
To quantify the curvature of the deformed sheet, we extract the sheet cross-section from the simulations and perform circular fits. This process is illustrated in figure 4. The two sides of the sheet are separately extracted, and individual circles are fitted to each side, resulting in two curvatures
$\kappa _{1}$
and
$\kappa _{2}$
. Each value will be defined as negative if the corresponding circle centre is in the same direction as the laser or positive otherwise. We then define the simulated sheet curvature as
\begin{equation} \kappa _{\textit{avg}} = \frac {\kappa _{1} + \kappa _{2}}{2}. \end{equation}
With this definition, we obtain
$\kappa _{\textit{avg}} = 0$
if the sheet is perfectly symmetric, as in the case of the cosine profile shown previously in figure 1. It is positive when the sheet curves away from the laser origin and negative when it curves towards the laser origin. These three possibilities are illustrated in figure 4, showing the fitted circles used to estimate the curvatures.
Figure 4. Method used to define and extract the sheet curvature of experiments and simulations based on fitted circles. The black lines represent the sheet interface, as extracted from a simulation or an experiment. The orange and red lines show the circles that were fitted to the laser side and back side of the sheet, respectively.
We begin by looking at how the numerical curvatures
$\kappa _1$
,
$\kappa _2$
and
$\kappa _{\textit{avg}}$
develop over time. Figure 5(a) presents
$\kappa _{\textit{avg}}$
over non-dimensional time for different raised cosine profiles. The full evolution of all simulations in this figure can also be seen in supplementary movie 1 available at https://doi.org/10.1017/jfm.2025.10665. Initially, all curves start from the origin since
$\kappa _{\textit{avg}} = 0$
in the case of a perfect sphere. Over time, each simulation reaches a peak curvature as the droplet experiences its initial deformation, which then tends again to zero as the sheet expands and becomes flatter. While some curves show positive values and others negative, none present a change of sign over time. This means that
$\kappa _{\textit{avg}}$
is a good theoretical quantity to classify if a simulation results in positive or negative curvature since this classification will not depend on the chosen measurement time.
Figure 5. Sheet curvature over time for different values of the raised cosine parameter
$W$
. (a) Average curvature
$\kappa _{\textit{avg}}$
, (b)
$\kappa _1$
and (c)
$\kappa _2$
. The inset shows the peak average curvature
$\kappa ^{\textit{max}}_{\textit{avg}}$
for each simulation as a function of
$W$
. The transparent bands show the uncertainty in the curvature measurement for simulations as defined in Appendix D. (d,e) Droplet/sheet interface at selected time stamps for two simulations with
$W = 1.25$
and
$W = 2.75$
, respectively. The components of the velocity field in cylindrical coordinates
$u_x$
and
$u_z$
are also shown in the snapshots. The exact parameters for these simulations are presented in table 2 of Appendix A.
For low values of
$W$
, as illustrated in figure 5(d), the focused profile quickly pierces the laser side of the droplet, changing it from its initial convex shape to a concave shape, which results in negative values for
$\kappa _{\textit{avg}}$
. Due to the focused profile, we can see that most of the velocity field is concentrated at the laser side of the droplet, which deforms very fast and creates the typical negative curvature shape. This initial curvature formation happens in approximately half the expansion (inertial) time scale, as can be seen in figure 5(a) by the moment when a peak in
$\kappa _{\textit{avg}}$
is obtained. After this peak is reached, most of the deformation observed will be simple sheet expansion, such that the curvature will again reduce in magnitude as the sheet becomes more extended. This initial curvature deformation is seen very clearly by looking at the initial field
$u_z$
, also shown in figure 5(d) for
$W = 1.25$
. We observe that
$u_z$
has only a small concentrated area of strong positive values at the laser side of the droplet. The opposite side of the droplet maintains mostly its spherical shape during this stage. Over time, due to mass conservation, the
$u_z$
component is converted into
$u_x$
as the droplet is progressively ‘squeezed’ in the
$z$
direction and expands in
$x$
. This causes the
$u_z$
field to become closer to uniform within the droplet, as seen in figure 5(d). By observing the velocity component field
$u_x$
for
$W =1.25$
we can also see a good indication that negative curvature will be obtained, since the higher values of
$u_x$
are strongly concentrated on the laser side of the droplet at initial time. Due to this, the droplet expands towards the laser, and the curvature becomes increasingly negative until it eventually reaches its peak.
On the other hand, the unfocused pressure profile wraps the droplet so that the laser side of its surface does not experience strong deformation and maintains its approximate initial spherical shape. This can also be seen clearly through the velocity field
$u_z$
in figure 5(e) for
$W = 2.75$
. Most of the early deformation happens on the opposite side of the droplet, where a strong negative velocity field pushes the surface inwards, creating a positive curvature. The velocity field
$u_x$
also indicates positive curvature, since the highest values of
$u_x$
are concentrated on the droplet side not illuminated by the laser. Similarly to the previous case, we see also that the curvature experiences a peak. Once again, this is set by the time it takes for the initial velocity field to establish the shape of the sheet, after which only expansion with a fixed shape will be present.
The inset in figure 5(a) presents the value of this peak (
$\kappa ^{\textit{max}}_{\textit{avg}}$
) as a function of the profile parameter
$W$
. Interestingly, we see a linear growth of the peak curvature within this range of
$W$
. Another interesting observation is that the time scale of the peak formation is not heavily dependent on
$W$
and approximately half the inertial time scale. Due to the similar time scale magnitude, the processes of sheet curving and sheet expansion happen simultaneously and cannot be easily separated over time.
Figure 5(b,c) also shows the individual curvatures
$\kappa _1$
and
$\kappa _2$
for the laser side and opposite side of the droplet, respectively. Unlike
$\kappa _{\textit{avg}}$
, these individual curvatures can change sign over time. The curvature
$\kappa _1$
will always present a sign change for simulations with small
$W$
. This happens as initially
$\kappa _1 = 1$
due to the spherical shape and eventually becomes negative as the laser side of the droplet is strongly pushed by a focused pulse. However,
$\kappa _2$
will not change sign for simulations of small
$W$
since the opposite side of the droplet does not experience strong pressure to deform, and its initial curvature does not change direction. On the other hand, for simulations with high
$W$
, the exact opposite happens, that is,
$\kappa _2$
changes sign over time while
$\kappa _1$
does not.
As suggested by figure 5, a direct correlation between
$\kappa _{\textit{avg}}$
and
$W$
will vary with time; however, the sign change will not be time dependent. Since our main interest is to know when a sheet will be positively or negatively curved, we now study this correlation between
$\kappa _{\textit{avg}}$
and
$W$
at a fixed time.
In figure 6 we present the measured value of
$\kappa _{\textit{avg}}$
as a function of the profile width
$W$
at a specific time of the simulations
$t \boldsymbol{\cdot }\dot {R}_0/R_0 = 1$
, where
$\dot {R}_0$
is the initial expansion velocity experienced by the sheet. The time scale
$R_0 / \dot{R}_0$
is then an inertial time scale related to how fast the sheet expands. We chose to match this specific non-dimensional time for all simulations so that all sheets have the same amount of expansion at the moment of measurement. As expected, we see that the curvature starts at negative values for focused pressure profiles (low
$W$
) and then monotonically increases with
$W$
. In the same figure we visualise the values of
$\theta _{\textit{max}}$
introduced in figure 3. We can observe that both
$\theta _{\textit{max}}$
and
$\kappa _{\textit{avg}}$
follow qualitatively a very similar trend. Moreover, the
$\theta _{\textit{max}}$
curves indicate that the curvature will flip at
$W = 1.9$
, and the
$\kappa _{\textit{avg}}$
curve indicates the flip around
$W \approx 2$
, which is a reasonable agreement.
To visualise the correlation between the predicted and measured curvatures, we plot in figure 6(b) the value
$\theta _{\textit{max}}$
vs
$\kappa _{\textit{avg}}$
. We observe that the curve passes not too far from the point
$(\theta _{\textit{max}}, \kappa ) = (\pi /2, \ 0)$
, which is the point that indicates a perfect prediction of the curvature flip. For illustration, in the bottom row of figure 6, we show the sheet cross-section for eight of the simulations used in the other panels. One can easily see from these snapshots that the curvature sign flips around
$W = 2$
. Other interesting characteristics can also be seen in these snapshots outside the curvature measurements. We observe, for example, that the thickness profile along the sheet is closer to constant for high values of
$W$
, while the sheet is very thick at its centre and thin at the edges for low
$W$
. These thin edges will eventually numerically disconnect from the sheet, forming rings, as our grid is not fine enough to resolve it properly. We note that, in experiments, we do indeed observe small droplets fragmenting near the edge of the sheet. However, this is not what is observed here since our simulation assumes axisymmetry and asymmetric fragments cannot be reproduced. The near-constant thickness profile obtained from high
$W$
is beneficial for the experimental application of EUV-light generation since the tin over the whole sheet can be evenly ablated into plasma by a following laser pulse.
Figure 6. Simulated average curvature and
$\theta _{\textit{max}}$
for different raised cosine parameters
$W$
. The curvature is always obtained at time
$t \boldsymbol{\cdot }\dot {R}_0/R_0 = 1$
. Panel (a) shows
$\theta _{\textit{max}}$
and
$\kappa _{\textit{avg}}$
over
$W$
. (b) Correlation between
$\theta _{\textit{max}}$
and
$\kappa _{\textit{avg}}$
. (c) Snapshots for data points shown in panels (a) and (b). The coloured lines represent the two circle fits in figure 4. Error bars estimate the uncertainty from fitting circles using different sections of the sheet; see Appendix D for details. The exact parameters for these simulations are presented in table 2 of Appendix A.
4.3. Comparison with experiments
We have seen so far that the raised cosine profile can provide deformed sheets with a curvature that goes from negative to positive as
$W$
is changed. We now attempt to visualise how well this can be correlated to the sheet shapes obtained experimentally when the ratio between laser focus and droplet size is changed.
To quantify the sheet curvature experimentally, we use the side-view shadowgraphs obtained following the experimental method described in § 2. For a given shadowgraph, we apply an edge-extraction algorithm to obtain the interface of the sheet. We then perform two circle fits as described previously in figure 4 for simulations. However, unlike in the simulations, we note that the shadowgraphs do not show a cross-section of the sheet, they only show a projection of the three-dimensional sheet into the camera plane. Consequently, if the sheet has a curvature, the concave side of the interface will look flat, and the actual curvature of that side will be hidden in the projection. With this limitation in mind, in this section we opted to define our curvature differently from the simulation case in (4.6). We consider only the side of the sheet that displays the largest curvature in absolute value, while the opposite side is ignored since it is likely hiding an internal curvature that cannot be seen. Therefore, the curvature in this section will be defined as
\begin{equation} \kappa _{\textit{max}} = \begin{cases} \kappa _1, &\quad |\kappa _1| \geq |\kappa _2|, \\ \kappa _2, &\quad \text{otherwise}, \end{cases} \end{equation}
where
$\kappa _1$
and
$\kappa _2$
refer to the curvature of each side of the sheet as illustrated previously in figure 4.
We now sweep over the ratio between the diameters of the laser beam and the droplet
$(d/D_0)$
. This is done by keeping fixed the beam diameter
$d = 20\,\unicode{x03BC} \text{m}$
and using droplets with a diameter in the range
$D_0 \in [27, \ 59]\,\unicode{x03BC} \text{m}$
. For each experiment, side-view snapshots were taken at the times
$[0, \ 500, \ 1000, \ 2000, \ 2500]\,\text{ns}$
. The first two snapshots were used to estimate the expansion velocity of the sheet
$\dot {R_0}$
, which is then used to obtain the non-dimensional expansion time
$t_{\textit{exp}} = {t \boldsymbol{\cdot }\dot {R_0}/R_0}$
. In order to measure the curvature at similar expansion times, we then select the snapshot of each experiment that has a time closest to
$t_{\textit{exp}} = 5$
. The actual selected time will vary between experiments since the experimental snapshots are only available at limited time steps.
Figure 7 shows the experimental images at the selected time. As expected, we observe that a focused laser pulse (with respect to the droplet size) pierces the centre of the droplet, resulting in a sheet with negative curvature. On the other hand, a wide beam results in a sheet that curves positively. Since the droplet is small in comparison to the beam size, we hypothesise that the laser generates a plasma that wraps around most of the droplet, such that the tips also experience a pressure that can push them forward. A detailed description of this plasma-wrapping phase is given by Hemminga et al. (Reference Hemminga, Poirier, Basko, Hoekstra, Ubachs, Versolato and Sheil2021), where this process is numerically investigated and the authors verify that the plasma expands across the whole droplet surface. The value of
$\kappa _{\textit{max}}$
is plotted in panel (a) as a function of
$d/D_0$
. The triangular points indicate the cases in which
$\kappa _1$
was selected as
$\kappa _{\textit{max}}$
, while for the circle points,
$\kappa _2$
was selected. Similarly to the simulation results, the curvature flips from negative to positive as the ratio
$d/D_0$
increases, eventually approaching a plateau. This flip from negative to positive curvature happens around
$d / D_0 \approx 0.575$
, indicated as a vertical line in figure 7. For comparison, we also plot the simulated curvature as a function of
$W$
in the same panel using a separate horizontal axis for
$W$
. The linear fit used to correlate the
$W$
axis and the
$d/D_0$
axis is given by
$W = 4.37 \ (d/D_0) - 0.48$
. This fit, however, is specific for the range of ratios shown here
$d/D_0 \approx 0.5$
and would not realistically hold for more focused laser pulses, where a nonlinear scaling would be necessary. We note that to fairly compare numerical and experimental results, the numerical curvature shown here is also
$\kappa _{\textit{max}}$
, which is different from the curvature
$\kappa _{\textit{avg}}$
used in the previous sections. For simulations, we note that the flip between
$\kappa _1$
or
$\kappa _2$
being chosen for
$\kappa _{\textit{max}}$
happens at
$W \approx 2$
, which is also indicated by the vertical line in figure 7 that separates the triangular and circular simulation points. Overall, a good agreement between the experimental and simulated curvatures can be seen, which motivates the usage of the raised cosine parameter
$W$
as a numerical analogous of the experimental ratio
$d/D_0$
.
Figure 7. Experimentally measured sheet curvature for different ratios
$d/D_0$
. For all cases shown, we use
$E_{pp}=0.8$
mJ,
$d=20\,\unicode{x03BC}$
m and
$\tau _p=10$
ns, while
$D_0$
is varied from
$59$
to
$27\,\unicode{x03BC}$
m. (a) the experimental curvature is plotted as a function of
$d/D_0$
. Grey markers are the individual repetitions of the experiment and black markers/bars show their mean. Error bars estimate the uncertainty from using circle fits; see Appendix D for details. As a qualitative comparison, the numerical curvature is shown as a function of the raised cosine parameter
$W$
(red). (b) a shadowgraph for each
$d/D_0$
is shown, going from
$59$
(top left) to
$27\,\unicode{x03BC}$
m (bottom right). The orange and red curves indicate the circles that were fitted to each side of the sheet, as described in figure 4. The bright spots are generated by plasma light that overexposes the camera chip. In some frames, this plasma is offset relative to the droplet centre as a result of vertical motion of the droplet over time. Our images are flipped vertically, so the droplet is displaced upward in relation to the plasma location. The exact parameters for these simulations and experiments are presented, respectively, in tables 2 and 3 of Appendix A.
In order to visually compare the sheet deformation between experiments and simulations, in figure 8 we show a time lapse of three different scenarios involving positive and negative curvatures. Each frame contains a shadowgraph of an experiment overlayed with the corresponding simulated sheet cross-section. The inset of each frame shows the projection of the three-dimensional simulated sheet, which would be the equivalent comparison to what is seen in the experimental images. The curvature
$\kappa _{\textit{exp}}$
is measured over time for all three cases and plotted at the top panel of the same figure.
Figure 8(b) contains shadowgraphs at five time stamps of a tin droplet experiment that exhibits negative curvature. To match this specific experiment, we perform a simulation with a focused raised cosine profile of
$W = 1$
and an expansion Weber number of
$\textit{We}_{\textit{def}} = 2022$
. This overlayed comparison between simulation and experiment accentuates the main advantage of the numerical simulations: the actual sheet thickness can be seen in the simulation, while it is completely hidden in the projected two-dimensional experimental images. Very good agreement is observed in the bulk area of the sheet between experiment and simulation, while some discrepancy is observed only as we get closer to the edges. Experimentally, strong asymmetric sheet fragmentation is observed at the edge of the sheets as a result of a violent expansion. As previously discussed, this behaviour cannot be captured by our simulations, which are axisymmetric and only present the expected formation of a rim. At the edges of the experimental images, we also notice an area in which the sheet deviates from the bulk curvature, becoming nearly flat abruptly. This behaviour is also not captured by our initialisation approach, which creates a smooth curvature change along the droplet and sheet.
Figure 8(c) shows the same experimental--simulation visual comparison now for a case with positive curvature. We chose time stamps that result in expansion times similar to those in the previous case. At the earlier time of
$t_{\textit{exp}} = 1.36$
we observe a very good agreement between the simulated and experimental sheets. Note that no onset of experimental or numerical fragmentation has yet been observed at the sheet edge since the sheet volume is better distributed with a near-constant thickness profile in this case. At later times, while the positive curvature is maintained for both the simulated and experimental sheets, a discrepancy can be observed: the laser side of the sheet remains uniformly curved in the simulation, while in the experiments, it flattens out and only the edges curve forward. This is a striking feature of the experimental sheets that is still numerically unreproducible with our current choice of pressure profile. Numerically, we also see the formation of a small rim at later times that curves slightly back. The rim is formed due to capillary effects acting at a late time, and we believe this rim is slightly pulled back due to drag in the simulations (the outer medium is not a vacuum in simulations due to numerical limitations).
Figure 8(d) showcases the laser-induced deformation of a larger water droplet, reproduced from the work of Klein et al. (Reference Klein, Bouwhuis, Visser, Lhuissier, Sun, Snoeijer, Villermaux, Lohse and Gelderblom2015). In this experiment, a large water droplet with
$D_0 = 0.9\,\text{mm}$
is used, and due to its size, the velocity received by the droplet is small, such that only a small Weber number (
$\textit{We} = 16$
or
$\textit{We}_{\textit{def}} = 141$
) is attained. We note that this is a striking difference from our experiments, in which a micrometre-sized droplet is used and deformed with a much higher Weber number
$(\textit{We}_{\textit{def}} \gt 1000)$
. Not only do the size, velocity and time scales differ from our experiments, but also some aspects of the mechanism that induces the initial pressure. In this particular water droplet, only a vapour cloud is created, as opposed to the tin plasma that quickly expands and wraps the droplet in our experiments. We note, however, that with appropriate laser parameter choices, plasma can also be generated from water droplets as reported, for example, by Klein et al. (Reference Klein, Bouwhuis, Visser, Lhuissier, Sun, Snoeijer, Villermaux, Lohse and Gelderblom2015) and Marston et al. (Reference Marston, Mansoor, Thoroddsen and Truscott2016), but these are not the conditions relevant to the comparison below. Besides all these differences, we show that the raised cosine approach can also be successfully used to describe this flow. In this case, a focused profile with
$W = 1.25$
was chosen to obtain the negative curvature observed in the experiment. The good visual agreement is observed in the frame-by-frame comparison, which is also confirmed in the
$\kappa _{\textit{max}}$
measurements over time shown in the top plot. This confirms that this approach can be used in a wide range of droplet impact scenarios, with different time, length and velocity scales and different mechanisms of propulsion.
Figure 8. Frame-by-frame comparison between sheets obtained experimentally and numerically for three different sets of parameters. The background of each plot shows the experimental side-view shadowgraph, while the red/green/yellow outlines show the simulated sheet cross-section. The insets show the side-view projection of the simulated sheet. (a) Numerical and experimental curvature over time for each case in the following panels. Circles and error bars represent the experimental results, and dashed lines represent the numerical results. The transparent bands show the uncertainty in the curvature measurement for simulations as defined in Appendix D. The uncertainty for experimental points is defined in the same manner. (b) Case from a focused laser beam resulting in negative curvature. (c) Unfocused beam resulting in positive curvature. (d) Experiment with a large water droplet reproduced from Klein et al. (Reference Klein, Bouwhuis, Visser, Lhuissier, Sun, Snoeijer, Villermaux, Lohse and Gelderblom2015). The exact parameters for these simulations and experiments are presented in table 4 of Appendix A.
5. Discussion, conclusion and perspective
We numerically and experimentally investigated the morphology of a liquid sheet obtained by laser-induced droplet deformation. Direct numerical simulations are performed to understand how the initial pressure profile exerted on the droplet surface can lead to different late-time sheet morphologies. We propose a pressure profile based on a raised cosine function and demonstrate that, by tuning its parameter
$W$
, we can obtain better morphological agreement with experiments compared with previous functions proposed in the literature.
Our results show that using an instantaneous pressure impulse we can still obtain both negative and positive sheet curvatures as long as the pressure profile function is correctly chosen. The previously proposed profiles in the literature (Gaussian and cosine functions) could not provide both sheet curvature types, which illustrates again the importance of carefully choosing or determining the pressure profile. The actual determination of this pressure profile from experiments is a difficult task since it comes from the complex interaction between a laser-generated plasma and a droplet. As this plasma expands and wraps the droplet, the determination of the exerted pressure is not trivial. Therefore, while we show in this work that we can obtain positive curvatures as long as a suitable pressure impulse is given, we note that more work on determining the real shape of this impulse from an experimental point of view is still required. Computationally, determining the ‘correct’ pressure impulse would require complex plasma physics modelling of the interaction between tin plasma and the droplet, which requires different numerical techniques that are out of the scope of the current work.
We have analysed the corresponding initial velocity field for a given pressure profile and proposed a method to predict whether a given profile will produce a positive or negatively curved sheet. This method shows that a positively curved sheet can be obtained with raised cosine functions of
$W\gt 1.9$
, which can never be achieved with the traditional Gaussian pressure. We observe that in order to obtain positive curvature, a pressure profile needs to present a wide peak and a short tail, which the raised cosine can achieve due to its negative excess kurtosis but not the Gaussian (zero excess kurtosis). The initial prediction is then compared with the curvature obtained after the simulation of the droplet deformation. Good agreement is obtained as we see that the sheet curves forward for approximately
$W\gt 2$
. Predicting sheet morphology from the initial condition at time zero can be very beneficial for future works on probing different pressure profiles and their corresponding deformation dynamics. Without having to perform expensive dynamical fluid simulations, it is possible to test different profiles and tune them accordingly to obtain the desired morphology. This allows for efficient and careful design of pressure profiles to obtain different sheet morphologies.
Focusing on the raised cosine profile, we have shown that the numerical parameter
$W$
can be correlated to the experimental beam-to-droplet size ratio
$d/D_0$
. A negative sheet curvature is obtained for low values of
$W$
and
$d/D_0$
, while a positive curvature is obtained for high values. This effect is associated with the pressure profile being either focused at the centre of the droplet or wrapping around most of the droplet surface, which can experimentally happen due to the fast plasma expansion.
While the raised cosine profile is capable of answering the riddle of curvature inversion observed in experiments, some other experimental characteristics are still not fully matched, notably the expansion-to-propulsion velocity ratio highlighted in Appendix E. Therefore, additional research on determining the optimised pressure profile for specific experiments could still be performed. We believe some insights on this optimised profile could be obtained numerically using full radiation-hydrodynamics solvers, such as by Hernandez-Rueda et al. (Reference Hernandez-Rueda, Liu, Hemminga, Mostafa, Meijer, Kurilovich, Basko, Gelderblom, Sheil and Versolato2022), experimentally by using retarding field analysers to measure energy of flying ions at different angles as by Poirier et al. (Reference Poirier, Hemminga, Lassise, Assink, Hoekstra, Sheil and Versolato2022), or analytically with strongly simplified plasma expansion models as by Murakami et al. (Reference Murakami, Kang, Nishihara, Fujioka and Nishimura2005), Bayerle et al. (Reference Bayerle2018). All of these approaches rely heavily on studying the dynamics of plasma expansion, and are beyond the scope of the present work.
The approach used in this work can be used more generally than only in laser--droplet applications for EUV-light generation. We have shown that the raised cosine function also correctly captures the laser-induced expansion dynamics of a large millimetre-sized water droplet performed by Klein et al. (Reference Klein, Bouwhuis, Visser, Lhuissier, Sun, Snoeijer, Villermaux, Lohse and Gelderblom2015). These deformation dynamics have also already been shown to be similar to the case of droplets impacting narrow solid pillars, such that this numerical simulation approach could also be used for droplet impact problems to design or optimise such systems. Moreover, the combination of the present simulation methods and the counterpart experimental insights may be leveraged to shed further light on a wide range of other laser--material interactions, such as laser-induced forward transfer (Jalaal et al. Reference Jalaal, Klein Schaarsberg, Visser and Lohse2019b ; Serra & Piqué Reference Serra and Piqué2019), laser-induced breakdown spectroscopy (Fortes et al. Reference Fortes, Moros, Lucena, Cabalín and Laserna2013) and pulsed laser deposition (Shepelin et al. 2023), paving the way for more predictive models and achieving optimised process control.
Supplementary movie
Supplementary movie is available at https://doi.org/10.1017/jfm.2025.10665.
Funding
This work was conducted at the Advanced Research Center for Nanolithography (ARCNL), a public–private partnership between the University of Amsterdam (UvA), Vrije Universiteit Amsterdam (VU), Rijksuniversiteit Groningen (UG), the Dutch Research Council (NWO) and the semiconductor equipment manufacturer ASML, and was partly financed by ‘Toeslag voor Topconsortia voor Kennis en Innovatie (TKI)’ from the Dutch Ministry of Economic Affairs and Climate Policy.
Declaration of interests
The authors report no conflict of interest.
Appendix A. Parameters from simulations and experiments
In total, eight simulations are used in figures 5, 6 and 7, with parameters shown in table 2. The input values for
$\textit{Re}$
and
$\textit{We}$
change significantly between simulations as we keep the values of
$\textit{Re}_{\textit{def}}$
and
$\textit{We}_{\textit{def}}$
close to 2000 while changing
$W$
. As mentioned previously, with the large values of
$\textit{Re}_{\textit{def}}$
and
$\textit{We}_{\textit{def}}$
, no significant effects of viscosity and surface tension are visible within the inertial-expansion time scale.
Table 2. Non-dimensional groups for all simulations shown in figures 5, 6 and 7. We note that the propulsion-based numbers
$\textit{We}$
and
$\textit{Re}$
are the actual inputs provided to our solver, while
$\textit{We}_{\textit{def}}$
and
$\textit{We}_{\textit{def}}$
are measured during post-process since they are the actual relevant numbers to study droplet expansion.
Table 3. Non-dimensional groups for the experiments shown in figure 7. For a fixed laser beam energy and size, the droplet size is varied as a parameter in our experimental set-up. The Reynolds and Weber numbers (both for propulsion and expansion) are measured as an output.
In figure 7, eleven experiments are shown as we perform a sweep in the droplet-to-beam size ratio
$d/D_0$
by changing
$D_0$
. The parameters for each of these experiments is shown in table 3. Note that, similarly to what was done in the simulations, the propulsion Weber and Reynolds numbers change drastically between experiments, but
$\textit{We}_{\textit{def}}$
and
$\textit{We}_{\textit{def}}$
always stay large.
Finally, the parameters for the three simulations and three experiments used in figure 8 are presented in table 4. We note that the three experiments in this figure span a wide range of parameters and set-ups, including also large (milimetric) water droplets.
Appendix B. Effect of other parameters in the droplet expansion
From (3.6)–(3.10), we can identify that, for a given pressure profile, our system is determined by four parameters: the Reynolds number (
$\textit{Re}$
), Weber number (
$\textit{We}$
), viscosity ratio (
$\mu _a/\mu _d$
) and density ratio (
$\rho _a/\rho _d$
). In the main body of this paper we focus on the parameter ranges relevant to obtaining liquid sheets in a setting relevant to application: high Reynolds and Weber numbers, density and viscosity ratios close to zero. In this appendix we move away from these limits and provide a brief overview of how these parameters can affect the droplet expansion and the sheet morphology over time. For the simulations in this appendix, we assume the base values
$\textit{Re} = 1000$
,
$\textit{We} = 1000$
,
$\mu _a/\mu _d = 10^{-4}$
,
$\rho _a/\rho _d = 10^{-4}$
. In each parameter sweep, one of these values will be changed while the other three are kept fixed as the base value.
Table 4. Non-dimensional groups for the simulations and experiments shown in figure 8.
We begin by studying the effect of the Reynolds number in our system by running simulations with
$\textit{Re} = \{10, \ 100, \ 1000\}$
. In figure 9 we show simulations for these three values of
$\textit{Re}$
while also sweeping over the pressure profile parameter
$W$
. Similarly to what was previously done in figure 6, we measure the sheet curvatures at the expansion time scale
$t = R_0/\dot {R_0}$
and plot the two individual curvatures in panels (d) and (e). We observe that, for the two higher values of
$\textit{Re}$
, the influence in the sheet morphology is minimal at the time of measurement, since the time scale of viscous dissipation is much larger than the relevant inertial time scale of sheet expansion in which the measurement was taken. Some effect is observed only in the case of
$W = 2.75$
, due to the small expansion velocity of this pressure profile that leads to small effective inertia. Significant effects are only observed for the simulations of
$\textit{Re} = 10$
. In this scenario, the viscous dissipation time scale is comparable to the time of expansion, particularly in the cases of high
$W$
(small inertia), in which strong dissipation is quickly observed that significantly slows down the droplet expansion to the point that it never becomes a sheet. This reduction in expansion velocity is also made quantitatively clear in panel (f).
Figure 9. Effect of the Reynolds number in the sheet curvature at time
$t \boldsymbol{\cdot }\dot {R_0}/R_0 = 1$
. The values
$\textit{We} = 1000$
,
$\mu _a/\mu _d = 10^{-4}$
,
$\rho _a/\rho _d = 10^{-4}$
are kept fixed, and a sweep is performed in the Reynolds number and pressure profile parameter
$W$
. Minimal viscous effects are observed for
$\textit{Re}=1000$
and
$\textit{Re}=100$
, with significant dissipation being present only for
$\textit{Re} = 10$
.
Figure 10. Effect of the Weber number in the sheet curvature at time
$t \boldsymbol{\cdot }\dot {R_0}/R_0 = 1$
. The values
$\textit{Re} = 1000$
,
$\mu _a/\mu _d = 10^{-4}$
,
$\rho _a/\rho _d = 10^{-4}$
are kept fixed, and a sweep is performed in the Weber number and pressure profile parameter
$W$
. Almost no capillary effects are observed for
$\textit{We}=1000$
. For
$\textit{We}=100$
and, particularly,
$\textit{We}=10$
, strong capillary-driven droplet retraction and oscillations are present for the high values of
$W$
.
We follow by looking into the influence of the Weber number through simulations with
$\textit{We} = \{10, 100, 1000\}$
. This parameter sweep is shown in figure 10. Differently to the results in the Reynolds sweep, we already observe very significant capillarity effects for
$\textit{We} = 100$
. This is due to the square root dependency between the inertial and capillary time scale
$\tau _{c} / \tau _{i} = \sqrt {\textit{We}}$
, which keeps the two time scales closer to each other. Therefore, strong surface tension driven deformation and retraction can already be seen for the higher values of
$W$
at
$\textit{We} = 100$
. For the more extreme case of
$\textit{We} = 10$
, we see very clear rim formation for the low values of
$W$
and, once again, a very fast suppression of droplet expansion at high
$W$
, which leads to capillary-driven droplet oscillations. These oscillations are underdamped due to the choice of
$\textit{Re} = 1000$
made here. We note that the actual measurement of a ‘sheet curvature’ makes little sense in the cases shown for
$\textit{We} = 10$
and high
$W$
, since no sheet is actually formed, but we still show the measurements for consistency.
In figure 11 we analyse the effects of the viscosity ratio by sweeping over four orders of magnitude
$\mu _a/\mu _d = \{10^{-4}, 10^{-3}, 10^{-2}, 10^{-1}\}$
. This sweep corresponds to changing the amount of viscous dissipation only within the ambient fluid. With the very large Reynolds number used here (which also affects the ambient fluid, according to our (3.6)), no effects of viscous dissipation in the ambient phase is observed for this sweep. Therefore, for viscous dissipation in the outer fluid to be relevant, either the other three parameters need to be different or the viscosity ratio must be significantly bigger, such that the outer fluid is actually more viscous than the droplet itself.
Figure 11. Effect of the viscosity ratio in the sheet curvature at time
$t = R_0/\dot {R_0}$
. The values
$\textit{Re} = 1000$
,
$\textit{We} = 1000$
,
$\rho _a/\rho _d = 10^{-4}$
are kept fixed, and a sweep is performed in the viscosity ratio and pressure profile parameter
$W$
. Due to the very small density and high Reynolds number, no effects are observed here even though the viscosity ratio is within four orders of magnitude.
Figure 12. Effect of the density ratio in the sheet curvature at time
$t = R_0/\dot {R_0}$
. The values
$\textit{Re} = 1000$
,
$\textit{We} = 1000$
,
$\mu _a/\mu _d = 10^{-4}$
are kept fixed, and a sweep is performed in the density ratio and pressure profile parameter
$W$
. Almost no effects are observed for the two smaller values of density ratio. For the higher values, significant drag is present, particularly for the cases with small
$W$
, resulting in changes in the measured curvature and reducing expansion velocity.
We finalise by investigating the effect of the density ratio in this system by sweeping through the values
$\rho _a / \rho _d = \{ 10^{-4}, \ 10^{-3}, \ 10^{-2}, \ 10^{-1} \}$
, as shown in figure 12. This parameter sweep is particularly important since it can also demonstrate how close we are from the experimental ‘vacuum’ limit that we are interested in in the main paper. The simulations for the first two density ratios (
$10^{-4}$
and
$10^{-3}$
) look practically identical. This verifies that the density ratio chosen in the main body of this paper (
$10^{-4}$
) is sufficiently small to be used as an approximation to the real experimental scenario (vacuum). For the third density ratio
$(10^{-2})$
, some effects of drag are observed in the sheet morphology for the cases of small
$W$
. Since drag is proportional to the velocity squared, it is understandable that the cases of small
$W$
present the most disturbances, since the expansion velocity is high for small
$W$
. For the extreme case in which the ratio is
$10^{-1}$
, effects of drag can be observed for all values of
$W$
. For high
$W$
, the tips of the sheet get dragged inwards, which slightly changes the measured curvature. For small
$W$
, the surface instabilities become extremely large with the high expansion velocity leading to cases in which the surface is so unstable that no actual curvature measurements can be performed
$(W = 1)$
. We note that, for consistency with the main body of the paper, in all of these simulations the centre-of-mass propulsion velocity
$U_z$
is subtracted from the droplet at the initial condition, such that only expansion is present. In a real scenario, not only drag coming from the expansion would be present, but also from the centre-of-mass horizontal propulsion. To clarify how important the propulsion drag in the droplet morphology is, we show in figure 13 the same density ratio sweep without the removal of propulsion velocity at time zero. We initially observe that, for density ratio
$10^{-4}$
, these simulations are identical to those without propulsion, validating that the removal of the propulsion velocity is an acceptable strategy for the situations considered in this paper. As the density ratio is increased, the expansion drag is observed but additional drag also comes from the horizontal translation, which curves the tips of the sheet towards the laser (left side). In the extreme case of ratio
$10^{-1}$
, the horizontal drag is strong enough that a negative curvature is observed for all values of
$W$
.
Figure 13. Effect of the density ratio in the sheet curvature at time
$t = R_0/\dot {R_0}$
. This case differs from figure 12 as the propulsion velocity is not removed at the initial condition, so drag will also be observed due to the
$z$
-direction velocity. The values
$\textit{Re} = 1000$
,
$\textit{We} = 1000$
,
$\mu _a/\mu _d = 10^{-4}$
are kept fixed, and a sweep is performed in the density ratio and pressure profile parameter
$W$
. Almost no effects are observed for the two smaller values of density ratio. For the higher values, significant drag is present, particularly for the cases with small
$W$
, resulting in changes in the measured curvature and reducing expansion velocity. For the highest density ratio, the horizontal drag is strong enough that the sheet always displays a negative curvature, regardless of chosen
$W$
.
Appendix C. Numerical implementation and validation
Figure 14. Radius and thickness of the droplet/sheet over time for four different simulations. All simulations are performed with Gaussian pressure profiles, and the Gaussian width is varied in the set
$\sigma \in \{ \pi /8, \ \pi /6, \ \pi /4, \ \pi /3 \}$
. Results from this work (solid lines) are compared against the results from the boundary integral simulations by Gelderblom et al. (Reference Gelderblom, Lhuissier, Klein, Bouwhuis, Lohse, Villermaux and Snoeijer2016) (circles) for numerical validation of our implementation.
Basilisk (Popinet & Collaborators Reference Popinet2013-2021) and its Navier–Stokes and VOF solvers have been widely used and validated in other published works, including many problems with droplet deformation dynamics, e.g. in recent studies of França et al. (Reference França, Jalaal and Oishi2024), Li & Cheng (Reference Li and Cheng2023) and Sanjay et al. (Reference Sanjay, Zhang, Lv and Lohse2025). In our work, we use the standard Navier–Stokes solver within Basilisk with no custom modifications. The only custom implementation is the calculation of an initial condition for the velocity field as given in (3.14). In summary, this initialisation is performed according to the following steps.
-
(i) One of the pressure profile functions
$f(\theta )$
from table 1 is chosen by the user. -
(ii) The coefficients of the first 30 Legendre polynomials
$P_n$
are loaded from a spreadsheet. -
(iii) The first 30 coefficients
$A_n$
are calculated according to (3.13). The definite integral in this equation is numerically evaluated using the composite Simpson’s 1/3 rule. -
(iv) The propulsion velocity
$U_z$
is obtained according to (3.15). At this stage, this propulsion velocity is not yet unity, since the function
$f(\theta )$
used is not necessarily normalised for this. -
(v) For each grid cell within the droplet, the following steps are employed.
-
(a) The components of the velocity field in spherical coordinates (
$u_r$
and
$u_\theta$
) are set according to (3.14). We simply use
$\alpha = 1$
in this equation, which means the resulting field will not yet be normalised correctly. -
(b) A conversion is made to obtain the components of the velocity field in cylindrical coordinates (
$u_x$
and
$u_z$
). We also divide the velocity field by the propulsion velocity calculated in step (iv), which now guarantees the correct normalisation. In most cases (except figure 13), we also subtract the normalised propulsion velocity from the velocity field, which means the droplet’s centre of mass will not move in the simulations. Therefore, the components
$u_x$
and
$u_z$
are calculated through(C1)
\begin{align} u_x(x, z) & = \frac {u_r\sin {\theta } + u_\theta \cos {\theta }}{U_z}, \end{align}
(C2)
\begin{align} u_z(x, z) & = - \frac {u_r\cos {\theta } - u_\theta \sin {\theta }}{U_z} - 1. \\[2pt] \nonumber \end{align}
-
After the steps above, the equations are evolved in time with the standard Basilisk solver (centered.h header file). Surface tension is imposed with the standard tension.h header file and the interface is tracked and evolved using the two-phase.h header. In the first few time steps a very small time-step size is used
$(10^{-6}$
–
$10^{-7})$
while the solver slightly corrects the solution close to the droplet boundaries in order to guarantee full agreement with the equations in that area. After this initial adaptation, a larger time-step size is used. We set the tolerance for the Basilisk’s internal Poisson–Helmholtz solver to
$5 \boldsymbol{\cdot }10^{-5}$
and reduce the coefficient for the Courant--Friedrichs--Lewy (CFL) time-step restriction to
${\textrm{CFL}} = 0.5$
.
To validate our implementation of this approach with Basilisk C, we compare our results against the boundary integral simulations performed by Gelderblom et al. (Reference Gelderblom, Lhuissier, Klein, Bouwhuis, Lohse, Villermaux and Snoeijer2016). Figure 14 shows the radius and central thickness of the sheet over time for four simulations with different Gaussian pressure profiles. All four simulations were performed with fixed propulsion Reynolds and Weber numbers
$\textit{Re} = 5000$
and
$\textit{We} = 800$
, respectively. The pressure profile was always kept as a Gaussian, and the Gaussian width was varied for each simulation within the set
$\sigma \in \{ \pi /8, \ \pi /6, \ \pi /4, \ \pi /3 \}$
. We note that the curves for both the radius and the thickness agree well with the reference data from Gelderblom et al. (Reference Gelderblom, Lhuissier, Klein, Bouwhuis, Lohse, Villermaux and Snoeijer2016), serving as validation of our solutions.
Appendix D. Uncertainty estimation in the curvature measurements
In all figures of this work in which the curvature is estimated, we provide error bars for each measured point from both simulations and experiments. These bars estimate the uncertainty in the curvature measurement due to choices made in the selection of points to be used in the circle fit. We note in simulations and experiments that the sheet does not always presents the same curvature within its bulk area and its tips. Therefore, the measured curvature could vary if one fits a circle using the entire sheet or by using only points of the sheet that are close to the bulk. In order to make this choice, we define a parameter
$\phi \in [0, 1]$
, such that
$\phi = 1$
indicates the entire sheet will be used in the circle fit and
$\phi \rightarrow 0$
indicates only points close to the sheet centre will be used. We proceed to perform various circle fits by sweeping values of
$\phi$
, which will produce a distribution of curvatures. We note that, if the sheet was indeed a perfect circle, this distribution would contain a single value. As the sheet deviates from a circle, the distribution would widen and have a deviation around its mean. We select the mean of this distribution as the measured curvature (and corresponding circle) that represents our sheet. The standard deviation of the distribution is plotted as an error bar. The circles visualised in the snapshots of figures 6 and 7 are those obtained in this manner. This strategy is applied to create the error bars in figures 6 and 7 as well as the transparent bands in figures 5 and 8.
For simulation results, we note that all error bars are relatively small, indicating that the entire sheet is close to a circle and that fitting choices do not significantly impact our results. Slightly larger uncertainties are obtained only for simulations with small
$W$
due to the curved edges of the sheet as can be seen in figure 6.
For the experimental results, the error bars are significantly larger as can be seen in figure 7. This reflects in part the fact that the sheets in the experiment are, compared with the simulation, not equally well described by a circle curve. Additional research on determining the optimised pressure profile for specific experiments could still be performed to bring simulations and experiment in more complete agreement. In all figures with experimental results, we also show multiple semi-transparent points showing the individual repetitions of the same experiments. The opaque point in the foreground is the average of all repetitions, and we note that the standard deviation from this distribution is smaller than the one coming from the circle fitting strategy.
We note that, after choosing a circle with the procedure described above, we have also calculated a root-mean-square (RMS) error for all experiments and simulations, comparing the actual points of the extracted sheet to the chosen circle. This RMS leads to error bars that are significantly smaller than those coming from the standard deviation of the distribution of the circle fitting results. With this in mind, we chose to use only the standard deviation as our error bars in all figures.
Appendix E. Expansion-to-propulsion velocity ratio
While the raised cosine profile can provide a more realistic direction for curvature, it still presents a discrepancy from experimental observations. This discrepancy comes in the form of the expansion-to-propulsion velocity ratio
$\dot {R_0} / U_z$
. This quantity indicates how fast the droplet expands compared with how fast it propels forward.
In figure 15(a) the blue circles show
$\dot {R}_0/U_z$
as a function of the raised cosine width
$W$
. A very sharp diverging increase is observed for
$W \rightarrow 0$
, which indicates that almost all the energy is being directed to deforming the droplet and not to propel it. For higher values of
$W$
, the ratio decreases and eventually vanishes into a plateau where the droplet only propels and never expands. To validate if the raised cosine exhibits the expected velocity ratio for a given
$W$
, we compare our results with those previously reported by Hernandez-Rueda et al. (Reference Hernandez-Rueda, Liu, Hemminga, Mostafa, Meijer, Kurilovich, Basko, Gelderblom, Sheil and Versolato2022) using the full radiation-hydrodynamic code RALEF-2D. In that software, the pressure profile is not given as an input. Instead, the actual laser parameters are provided, and the full laser interaction with the tin droplet is simulated, such that the pressure profile and
$\dot {R}_0/U_z$
are both outputs of the simulation. We fitted the pressure profiles provided by Hernandez-Rueda et al. (Reference Hernandez-Rueda, Liu, Hemminga, Mostafa, Meijer, Kurilovich, Basko, Gelderblom, Sheil and Versolato2022) using raised cosines and plotted the velocity ratio provided by their code as a function of the fitted
$W$
. We can see that the velocity ratios obtained are a good match with the initialisation approach used in our simulations, which validates the initialisation method.
As we have seen previously in figure 6, in order to obtain significant positive curvature with a raised cosine, the width parameter needs to be
$W \gt 2$
. In this range of widths, we see in figure 15 that the velocity ratio is always below 1, such that we cannot obtain fast-expanding sheets that exhibit a positive curvature. This limitation is directly visible in figure 15(b), where we plot the predicted sheet curvature as a function of
$\dot {R}_0/U_z$
, and we can clearly see that we are only able to obtain simulations with positive curvature along with velocity ratios that are below 1. This is still a discrepancy from experiments that needs to be further studied in the future, as somehow, the experimental sheets exhibit both positive curvature and fast expansion simultaneously.
Figure 15. Measurements of the velocity ratio
$\dot {R}_0 / U_z$
for raised cosine profiles. (a) Velocity ratio as a function of
$W$
for the initialisation approach in this work and as predicted with the code RALEF by Hernandez-Rueda et al. (Reference Hernandez-Rueda, Liu, Hemminga, Mostafa, Meijer, Kurilovich, Basko, Gelderblom, Sheil and Versolato2022). (b) Numerical average curvature
$\kappa _{\textit{avg}}$
as a function of the velocity ratio in our simulations.
Open access
