1. Introduction
Evaporating droplet phenomena are ubiquitous in a wide range of problems in nature, ordinary life and industry. In past decades, following the pioneering work of Deegan et al. (Reference Deegan, Bakajin, Dupont, Huber, Nagel and Witten1997), significant progress has been made in understanding the evaporation of single-component droplets, as explained by reviews by Cazabat & Guena (Reference Cazabat and Guena2010), Erbil (Reference Erbil2012) and Gelderblom, Diddens & Marin (Reference Gelderblom, Diddens and Marin2022). Not only single-component droplets but also multi-component droplets are quite important for many industrial applications, such as inkjet printing (Kuang, Wang & Song Reference Kuang, Wang and Song2014; Lohse Reference Lohse2022), spray cooling (Kim Reference Kim2007) and biological assays (Garcia-Cordero & Fan Reference Garcia-Cordero and Fan2017). On the contrary to the single-component case, for a droplet containing two or more components with different evaporation rates, the concentration of components near the liquid–gas interface can be changed during the evaporation. Then the concentration gradient along the interface can cause the surface tension gradient that results in the formation of solutal Marangoni instabilities (Sternling & Scriven Reference Sternling and Scriven1959; Lohse & Zhang Reference Lohse and Zhang2020; Wang et al. Reference Wang, Orejon, Takata and Sefiane2022), creating complex features such as chaotic flows (Christy, Hamamoto & Sefiane Reference Christy, Hamamoto and Sefiane2011), uniform pattern formation (Kim et al. Reference Kim, Boulogne, Um, Jacobi, Button and Stone2016) and phase segregation (Kim & Stone Reference Kim and Stone2018; Li et al. Reference Li, Lv, Diddens, Tan, Wijshoff, Versluis and Lohse2018). In addition to evaporation, a surface tension gradient can be caused by the external vapour, generating the vapour-driven solutal Marangoni effect. Using a vapour source, Hegde et al. (Reference Hegde, Chakraborty, Kabi and Basu2018) and Park et al. (Reference Park, Ryu, Sung and Kim2020) controlled the flow pattern inside the droplet, Cira, Benusiglio & Prakash (Reference Cira, Benusiglio and Prakash2015) and Malinowski, Parkin & Volpe (Reference Malinowski, Parkin and Volpe2020) achieved a contactless droplet movement, and Malinowski et al. (Reference Malinowski, Volpe, Parkin and Volpe2018) manipulated a particle deposition pattern after droplet evaporation.
Typically, in the evaporation of a multi-component droplet, the surface tension gradient is determined by the evaporation flux of more volatile components. Due to the high evaporation flux at the contact line of a sessile droplet with contact angle less than 90
$^\circ$
(Deegan et al. Reference Deegan, Bakajin, Dupont, Huber, Nagel and Witten1997; Hu & Larson Reference Hu and Larson2002), as the droplet evaporates, it fails to maintain its initial concentration, leading to significant concentration gradients along the droplet interface. This results in surface tension differences along the interface, generating Marangoni flow. However, in the case of the external vapour source, the surface tension gradient is defined with the vapour distribution surrounding the droplet interface, often without considering the adsorption or absorption flux (Hegde et al. Reference Hegde, Chakraborty, Kabi and Basu2018; Malinowski et al. Reference Malinowski, Volpe, Parkin and Volpe2018; Kabi, Pal & Basu Reference Kabi, Pal and Basu2020; Park et al. Reference Park, Ryu, Sung and Kim2020). This approach assumes the rapid establishment of vapour–liquid equilibrium between vapour in the air and the adsorbed molecule in the droplet by considering the fast adsorption dynamics of water-soluble vapour (Wilson & Pohorille Reference Wilson and Pohorille1997), indicating that the vapour concentration at the droplet interface is proportional to the vapour pressure above the interface according to Henry’s law (Henry Reference Henry1803). However, under the uniform distribution of external vapour, Marangoni flow occurs due to the high vapour absorption flux near the contact line (Majumder et al. Reference Majumder2012; Wang et al. Reference Wang, Karapetsas, Valluri, Sefiane, Williams and Takata2021; Kim Reference Kim2022; Yang et al. Reference Yang, Pahlavan, Stone and Bain2023). Additionally, while vapour molecules may adsorb on the droplet interface within tens of picoseconds (Wilson & Pohorille Reference Wilson and Pohorille1997), the Marangoni flow induced by external vapour concentration gradients typically takes several minutes to fully develop, as reported by Park et al. (Reference Park, Ryu, Sung and Kim2020). This indicates that vapour–liquid equilibrium, defined as the absence of net vapour mass flux across the interface, is not established immediately at the macroscopic level. Such a delay underscores the importance of investigating the transient nature of vapour-driven Marangoni flows. Furthermore, due to the continuous vapour absorption and the finite solubility of vapour in the droplet (Ryu, Ko & Kim Reference Ryu, Ko and Kim2022), the vapour mass flux along the droplet interface evolves over time: from an absorption-dominated regime to a state where absorption and evaporation coexist. This leads to the local establishment of vapour–liquid equilibrium, particularly at the contact line. Before this equilibrium is reached, the vapour concentration at the droplet interface is not directly proportional to the vapour distribution in the air, but they nearly align once equilibrium is achieved. Consequently, the vapour concentration gradient along the interface and corresponding Marangoni flow undergo a significant change, particularly around the onset of the local equilibrium. This phenomenon has been indirectly captured as the critical increase in flow magnitude after reaching a quasi-steady state in previous studies (Hegde et al. Reference Hegde, Chakraborty, Kabi and Basu2018; Park et al. Reference Park, Ryu, Sung and Kim2020; Kabi et al. Reference Kabi, Pal and Basu2020). However, while various phenomena using vapour have been reported, there is a lack of systematic research on internal changes occurring within the evaporating droplet in various vapour environments. To further understand the vapour-driven Marangoni flow throughout the entire droplet lifetime, a detailed investigation into the time-varying flow structure and vapour mass flux across the droplet interface is necessary.

Figure 1. Schematic of time-dependent vapour-driven Marangoni flows inside evaporating droplets before and after the vapour–liquid equilibrium at the droplet contact line, depending on the external vapour concentration gradient along the droplet radius. The greenish shadows represent the vapour concentration distribution in air. The reddish contours indicate axisymmetric vapour mass flux at the interface. The purple arrows represent the typical resulting flow patterns inside the droplet. Here,
$\gamma _{\textit{a}}$
and
$\gamma _{\textit{c}}$
are the surface tension at the apex and contact line, respectively,
$\gamma _{\textit{s}}$
is the surface tension at the stagnation point where two opposing interfacial flows meet, and
$\Delta \gamma$
indicates the surface tension difference (e.g.
$\Delta \gamma _{{a-c}}=\gamma _{{a}} - \gamma _{{c}}$
), determining the direction of Marangoni stress described by the blue arrows.
In this study, we observe the time-dependent flow structures using both experiment and numerical simulation throughout the entire lifetime of a droplet. To examine the effect of vapour concentration distribution around the droplet, we test two distinct conditions: low and high external vapour concentration gradients, representing uniform and localised vapour distributions, respectively. As depicted in figure 1, a representative flow transition leading to significant changes in flow structure is captured under each vapour concentration: a transition from inward–outward competing flows to outward interfacial flows in the high vapour concentration gradient, and a transition from inward interfacial flows to multiple vortices in the low vapour concentration gradient. By observing ethanol mass flux along the droplet interface using our numerical simulation, it can be noticed that the flow transition generally occurs when the ethanol mass flux becomes zero at the droplet contact line, indicating the local establishment of vapour–liquid equilibrium at the contact line. Focusing on the vapour mass flux and the resulting surface tension gradient, particularly before and after vapour–liquid equilibrium, we elucidate the reason and timing of the flow transition within a droplet in vapour environments.
In § 2, we describe the experimental set-up and the flow visualisation technique used to measure the internal flow patterns. We also present our numerical simulations to investigate the vapour mass flux across the interface, and provide a quantitative comparison with the experimental results. The flow transition and corresponding vapour mass flux at the interface for the low and high external vapour concentration gradients are discussed in §§ 3.1 and 3.2, respectively. In § 4, we focus on the key feature, namely a dramatic and somehow surprising increase in flow magnitude during the flow transition period, when the droplet is surrounded by the high external vapour concentration gradient. The paper ends with conclusions in § 5.
2. Methods
2.1. Experimental set-up
For the experiments, we deposited a
$2\pm 0.05$
$\unicode {x03BC}$
l water droplet (deionised water, 18.2 M
$\Omega$
cm Merck Millipore, Direct-Q3) on the cover glass (SciLab Korea). To reproduce a sessile droplet with the same diameter (
$R = 1.5 \pm 0.1$
mm) and height (
$H = 0.54 \pm 0.05$
mm), we used a PDMS mask and a high-frequency generator BD-10AS (Electro-Technic Products, USA) to create a circular hydrophilic area 3 mm in diameter on the hydrophobic glass substrate. From this substrate treatment, the droplet radius was maintained during the whole evaporation process. The ethanol vapour was exposed on top of the sessile droplet from the capillary tube, and the location of the tube was controlled by a motorised stage as illustrated in figure 2(a). The distance between the droplet and the tube was set at 2 and 10 mm to make the external vapour concentration gradient high and low, respectively. To isolate the external effects, they were covered with an acrylic hood to block the external flows. The temperature and relative humidity were maintained at 20
$^\circ$
C and 55 %, respectively.

Figure 2. Experimental set-up. (a) Ethanol in the capillary tube was located above the droplet. Here, we set the distance
$d$
between the droplet and the tip of the tube as 2 mm (high vapour concentration gradient) and 10 mm (low vapour concentration gradient) by the motorised stage. The inner diameter
$a$
of the capillary tube was 1 mm. The droplet radius (
$R$
) and height (
$H$
) were
$1.5 \pm 0.1$
mm and
$0.54 \pm 0.05$
mm, respectively. (b) To implement PIV, fluorescent particles in the droplet were excited by the laser (
$\lambda=532$
nm), and particles’ images were captured by the high-speed camera through the optical filter (
$\lambda \gt540$
nm). The images were acquired at the different focal planes, e.g.
$z=50$
and 200
$\unicode {x03BC}$
m.
2.2. Particle image velocimetry
To investigate the Marangoni flow transition, flow patterns inside the droplet were measured by two-dimensional (2-D) particle image velocimetry (PIV) with the set-up described in figure 2(b). To estimate internal flow patterns in the droplet, we observed the flow pattern at two different heights,
$z=50$
and 200
$\unicode {x03BC}$
m, and the experiments were repeated at least three times for each focal plane. As a tracer particle, we used 1.9
$\unicode {x03BC}$
m diameter fluorescent particles (PS-FluoRed-Fi320, microparticles GmbH, Germany) with concentration
$1.25\times10^{-5}$
wt %, which were excited by a 530 nm wavelength laser and emitted a 607 nm fluorescence signal. To record the particle movements, we used a high-speed camera (FASTCAM Mini AX100, Photron, USA) at 1000 frames per second. Here, the Stokes number
$St=\Delta \rho\, d^2_{{p}}U/(18\mu _{{w}} R)$
, which indicates the ratio between the particle relaxation time and the fluid relaxation time, is in the range
$10^{-10}$
to
$10^{-8}$
, where
$\Delta \rho=\rho _{{p}} - \rho _{{w}}$
indicates the density difference between the particle (
$\rho _{{p}}=1050$
kg m
$^{-3}$
) and the water (
$\rho _{{w}}=1000$
kg m−3),
$d_{{p}}$
is the particle diameter (1.9
$\unicode {x03BC}$
m),
$\mu _{{w}}$
is the dynamic viscosity of water (
$10^{-3}$
Pa s), and
$U$
is the typical flow velocity measured from the PIV experiment, which is
$0.01{-}1$
mm s−1. Therefore, it shows that the particles well follow the fluid flow because
$St$
is much smaller than unity (Kim et al. Reference Kim, Lee, Kim and Kim2015).
To acquire the velocity vectors, we conducted iterative 2-D cross-correlations with various widths of interrogation windows: 50
$\%$
overlapped 64
$\times$
64 pixels for a coarse grid, and 50
$\%$
overlapped 32
$\times$
32 pixels for a fine grid (Stamhuis & Thielicke Reference Stamhuis and Thielicke2014). Furthermore, a high-pass filter was used to remove the blurred particles out of a focal plane, and contrast-limited adaptive histogram equalisation (CLAHE) was applied to enhance the image quality. From the mass conservation for the incompressible flow (
$\boldsymbol{\nabla }\boldsymbol{\cdot }u = 0$
), where
$u = (\tilde {u},\tilde {v})$
is the 2-D velocity vector in the
$x$
–
$y$
Cartesian coordinates, the uncertainty of PIV measurements can be derived from the estimate

where
$\sigma _{\Delta x}$
is the overall error amplitude of the particle displacement during the time difference
$\Delta t$
, and
$D_{{I}}$
is the size of the interrogation window. Here, the error amplitude ranged from 0.001 to 0.056 pixels, which was within the previously reported uncertainty of
$O(0.1\ \text{pixels})$
(Adrian & Westerweel Reference Adrian and Westerweel2011; Sciacchitano Reference Sciacchitano2019). Thus the reliability of the PIV results was validated. Further validation of our PIV set-up through coffee-ring flow observations is provided in supplementary material at https://doi.org/10.1017/jfm.2025.10567.
2.3. Numerical simulation
Experimentally investigating the time-dependent change in the mass flux of ethanol vapour across the interface is almost impossible, so we carried out numerical simulations. The time-dependent internal Marangoni flow patterns were also numerically and experimentally compared and examined. The numerical simulations were based on a sharp-interface finite element method in an axisymmetric arbitrary Eulerian–Lagrangian description (Diddens Reference Diddens2017).
For the gas phase, it was sufficient to solve the quasi-stationary vapour diffusion equations, i.e. the Laplace equations
${\nabla} ^2c_\alpha =0$
, for both vapour partial mass concentrations
$c_{{e}}$
and
$c_{{w}}$
for ethanol and water, respectively. Due to the low concentrations of ethanol and water (below 1 wt %), Fick’s diffusion was a good approximation of the Maxwell–Stefan diffusion in the ternary gas mixture of water, ethanol and air (Krishna & Wesselingh Reference Krishna and Wesselingh1997). Additional test simulations with the consideration of the full advection–diffusion equations, along with Stefan flow and natural convection, confirmed the validity of the quasi-stationary diffusion-limited approximation in the gas phase. As far field conditions, vanishing ethanol vapour and ambient water vapour based on the experimentally determined relative humidity (
$\approx55$
%) were imposed. At the ethanol vapour source, saturated ethanol vapour concentration was imposed, and vanishing water vapour flux was assumed. Full simulations, which included the evaporation/condensation dynamics of ethanol and water, and considered the flow inside the vapour source capillary, did not give considerably better agreement with the experiments. The present simulations assume diffusion-limited evaporation, which is justified for ethanol–water droplets at room temperature, where Stefan flow, natural convection and vapour-phase contributions to Marangoni stress are negligible (Diddens et al. Reference Diddens2017b
). This is further supported by our own test simulations including full vapour-phase dynamics, which showed no significant deviation from the diffusion-limited results (see supplementary material and movie 3). Nevertheless, we note that in evaporation regimes where vapour transport is not diffusion-limited – for example, when advection or vapour asymmetry becomes dominant – the resulting flow and evaporation dynamics may differ. Understanding such regimes remains an important open question. At the liquid–gas interface of the droplet, the vapour concentrations were imposed based on the vapour–liquid equilibrium according to Raoult’s law, including activity coefficients to account for the non-ideality of ethanol–water mixtures. More details and the used activity coefficients can be found in Diddens et al. (Reference Diddens, Kuerten, Van der Geld and Wijshoff2017a
).
From the vapour fields, the evaporation rates
$j_\alpha =-D_\alpha \boldsymbol{\nabla }c_\alpha \boldsymbol{\cdot }\boldsymbol{n}$
were obtained, where
$\alpha$
is
$e$
for ethanol or
$w$
for water,
$D_\alpha$
is the diffusion coefficient of the component,
$\boldsymbol{n}$
is the outward normal unit vector at the liquid–gas interface, and negative values indicate vapour absorption. At the liquid–gas interface, the mass transfer led to an augmented kinematic boundary condition
$\rho (\boldsymbol{u}-\boldsymbol{u}_{{I}})\boldsymbol{\cdot }\boldsymbol{n}=j_{{w}}+j_{{e}}$
, connecting the liquid velocity
$\boldsymbol{u}$
with the interface motion
$\boldsymbol{u}_{{I}}$
. Furthermore, the different evaporation/absorption rates included a change of the local ethanol mass fraction
$w_{\textit{e}}$
in the droplet by the rate
$J=-\rho D\,\boldsymbol{\nabla }w_{\textit{e}}\boldsymbol{\cdot }\boldsymbol{n}=w_{\textit{e}}j_{{w}}+(1-w_{\textit{e}})j_{\textit{e}}$
. The ethanol mass fraction
$w_{\textit{e}}$
was solved fully coupled with the droplet flow by the following system (Diddens et al. Reference Diddens2017b
):



with the stress tensor
${\boldsymbol{T}}=-p{\boldsymbol{I}}+\mu (\boldsymbol{\nabla }\boldsymbol{u}+\boldsymbol{\nabla }\boldsymbol{u}^{\textrm{T}})$
. Besides the kinematic boundary condition and the compositional change rate, the Laplace pressure and the Marangoni stress were considered at the liquid–gas interface,

with the surface tension
$\gamma$
, curvature
$\kappa$
, and surface gradient operator
$\boldsymbol{\nabla} _S$
. Note that density
$\rho$
, dynamic viscosity
$\mu$
,
$\gamma$
, and the mass diffusivity
$D$
of the ethanol–water mixture in the droplet are functions of the local ethanol mass fraction
$w_{\textit{e}}$
, which were obtained by experimental data (Vazquez, Alvarez & Navaza Reference Vazquez, Alvarez and Navaza1995; González et al. Reference González, Calvar, Gómez and Domínguez2007; Pařez et al. Reference Pařez, Guevara-Carrion, Hasse and Vrabec2013). Thermal influences due to the latent heat of evaporation/absorption and change in solution have been disregarded since they were inferior compared to the solutal dependency of these properties. In this study, thermal Marangoni effects were neglected based on the dominance of solutal effects in binary mixtures of water and ethanol under ambient conditions. Previous studies have demonstrated that the solutal Marangoni number typically exceeds the thermal Marangoni number in such systems (Bennacer & Sefiane Reference Bennacer and Sefiane2014; Diddens et al. Reference Diddens2017b
). Moreover, we confirmed that the simulation results – including thermal effects such as latent heat of evaporation, thermal Marangoni flow, and thermally-driven natural convection – exhibit trends similar to those observed under isothermal conditions (see supplementary material and movie 3).
At the substrate, a Navier-slip condition with slip length 1
$\unicode {x03BC}$
m was used to resolve the incompatibility of the kinematic boundary condition with mass transfer at the pinned contact line. The explicit choice of the slip length did not influence the macroscopic results, as long as it was sufficiently small – specifically, smaller than 1
$\%$
of the initial droplet height. In this range, the slip length has no noticeable effect on the droplet dynamics, either qualitatively or quantitatively. The effect of the slip length plays a role in the finite element adjacent to the contact line, where it is relevant to regularise the incompatibility of the kinematic boundary condition with mass transfer and a vanishing velocity at the substrate. We therefore take a slip length corresponding to the typical size of the element at the contact line. More detailed explanations of the numerical method, including the used weak formulation can be found in the supplementary information of Jalaal et al. (Reference Jalaal, ten Hagen, Diddens, Lohse and Marin2022).

Figure 3. Comparison of the spatially averaged time-dependent radial velocity
$\bar {U}_{r,A}(t)$
depending on the vapour exposure conditions (a)
$d=10$
mm and (b)
$d=2$
mm, between the experimental (triangle symbols) and numerical (grey solid lines) results at different focal planes
$z=50$
and 200
$\unicode {x03BC}$
m. Here,
$\bar {U}_{r,A}(t)$
is the averaged radial velocity in the 2-D measurement area
$A$
, i.e.
$(\int _AU_r\,\text{d}A)/A$
, where
$U_r$
is the
$r$
-direction velocity at each focal plane. Error bars are obtained from three independent measurements.
To check the numerical results, the averaged radial velocity inside the droplet was compared with the experimental results, which were calculated as
$\bar {U}_{r,A} (t)= (\int _AU_r\,\text{d}A)/A$
, from 2-D PIV measurements for two different source positions (
$d=2$
and 10 mm) and two different focal planes (
$z=50$
and 200
$\unicode {x03BC}$
m). As shown in figure 3, there is relatively good agreement in flow magnitude and trend of flow transition, indicating that the numerical results can resemble the experimental results. The analysis of the detailed flow structures and the comparison between numerics and experiments will be discussed in the following sections.
3. Droplet interaction with external vapour
We will explore a water droplet evaporating in two distinct ethanol vapour environments characterised by low and high external vapour concentration gradients, which are distinguished by the vapour concentration difference between the droplet apex and the contact line. Figures 4 and 6 show both experimental and numerical observations of the transitions from inward interfacial flow to complicated asymmetric vortices and from inward–outward competing flow to outward interfacial flow for low and high concentration gradients, respectively. By comparing the flow patterns, we confirm that our numerical simulations agree reasonably well with experiments. To elucidate the reason for the flow transitions, we numerically investigate the distribution of the ethanol mass flux along the droplet interface. In figures 5 and 7, we primarily discuss the ethanol mass fraction (
$w_{\textit{e}}$
) determining the local surface tension and the tangential velocity (
$V_{\textit{t}}$
) profiles at the liquid–gas interface because the surface tension gradient along the interface and the induced interfacial flow are predominant.
3.1. Droplet surrounded by the low external vapour concentration gradient

Figure 4. Experimental and numerical results of the transient flow patterns within the droplet under the far-field source (
$d=10$
mm). The 2-D flow fields are observed using (a) PIV and (b) numerical simulation at the focal planes (
$a1, b1$
)
$z=200$
$\unicode {x03BC}$
m and (
$a2, b2$
)
$z=50$
$\unicode {x03BC}$
m. The black arrows are velocity vectors, and the white arrows represent typical flow patterns. The contour
$U_r$
is the
$r$
-direction velocity at each focal plane. We set
$t = 0$
for the moment when the ethanol vapour source is placed above the droplet apex,
$t_{{eq}}$
indicates the time when the vapour–liquid equilibrium is achieved, and
$t_{{binary}}$
indicates the moment at which the internal flow shows a complicated mixing flow, which is similar to the case of evaporating binary droplets. (c) The ethanol mass flux (green arrows at the droplet interface) and typical flow pattern at early and late stages reconstructed from numerical simulation (see supplementary movie 1). Blue and red arrows indicate the clockwise and anticlockwise directions, respectively.
For the low external vapour concentration gradient, we experimentally observe three stages based on the flow structure: (i) inward interfacial flow during the early stage, (ii) no flow during the transition stage, and (iii) asymmetric vortices in the late stage. In the first stage, radial outward and inward flows are observed experimentally at
$z=50$
and 200
$\unicode {x03BC}$
m, respectively, as shown in figure 4(a). These observations align with the numerical results presented in figure 4(b). The inward interfacial flow at the droplet side view, as depicted in the early stage schematic in figure 4(c), is consistent with previous studies (Majumder et al. Reference Majumder2012; Wang et al. Reference Wang, Karapetsas, Valluri, Sefiane, Williams and Takata2021). When the surrounding vapour concentration is relatively uniform, the vapour concentration difference across the droplet interface is also nearly uniform along the interface. Consequently, based on the diffusion equation, the absorption flux distribution of surrounding vapour is mainly determined by the geometry of the droplet interface rather than by the concentration difference. For a sessile droplet with contact angle less than 90
$^\circ$
, the evaporation flux at the contact line is the highest due to the geometry (Deegan et al. Reference Deegan, Bakajin, Dupont, Huber, Nagel and Witten1997). For a uniform distribution of surrounding vapour, the vapour concentration gradient across the interface that induces absorption would be opposite to the gradient caused by droplet evaporation. Thus, as the green arrows illustrate in figure 4(c), the absorption flux of surrounding vapour would be the highest at the contact line, as demonstrated by Wang et al. (Reference Wang, Karapetsas, Valluri, Sefiane, Williams and Takata2021), resulting in the inward interfacial flow.
As time goes on, the magnitude of the inward interfacial flow gradually decreases and reaches zero at the transition stage, as shown in figure 4(a). Although the inward circulating flow is still observed in the numerical results, as presented in figure 4(b), it is comparable to the no-flow state observed in the experimental results because the flow magnitude is close to zero. The disappearance of flow indicates a negligible ethanol mass flux across the entire droplet interface, signifying that the vapour–liquid equilibrium has been achieved. We define this moment as
$t_{{eq}}$
.
Despite reaching the vapour–liquid equilibrium, the no-flow state is not sustained until the droplet fully evaporates. Asymmetric vortices suddenly appear in the late stage of the experiment, as depicted in figure 4(a). While the flow structure differs, an outward interfacial flow is also created in the late stage of the simulation, as shown in figures 4(b) and 4(c), indicating that the vapour–liquid equilibrium is not maintained. This can be attributed to water evaporation, as the relative humidity in our setup is approximately 55 %, not 100 %. After the vapour–liquid equilibrium, the ethanol vapour no longer absorbs into the droplet, but the water continues to evaporate into the air. Consequently, the ethanol concentration in the droplet increases, causing the evaporation of ethanol from the droplet. This triggers solutal Marangoni flow, resembling the behaviour observed in the evaporation of binary mixture droplets (Christy et al. Reference Christy, Hamamoto and Sefiane2011; Bennacer & Sefiane Reference Bennacer and Sefiane2014). We set
$t=t_{{binary}}$
when the additional Marangoni flow is initiated. The discrepancy between the experimental and numerical results in the flow pattern would be due to the breaking of axial symmetry caused by Marangoni instability (Diddens et al. Reference Diddens2017b
). We acknowledge that the current numerical simulation was conducted in axisymmetry, limiting its ability to fully capture the three-dimensional features observed in the experimental results. Nevertheless, the current numerical model shows that the outward interfacial flow can correspond to multiple vortices observed in experiments on binary mixture evaporation (Diddens et al. Reference Diddens2017b
). Therefore, we consider the simulated and experimental results to represent the same underlying flow behaviour. Moreover, the 2-D model shows good agreement with the experiments in terms of the timing and trend of flow transitions. Based on these points, we believe that the numerical results demonstrate qualitatively good agreement with the experiments.

Figure 5. Numerical analysis for the low external vapour concentration gradient, namely the far-field source case (
$d=10$
mm). (a) The schematic of a droplet under the source, where
$\varPhi$
is the angle from the bottom substrate,
$V_{\textit{t}}$
is the tangential velocity, and
$w_{\textit{e}}$
indicates ethanol mass fraction just below the interface, illustrated as a dark grey area. (b) Normalised ethanol mass fraction profiles
$\bar {w_{\textit{e}}}$
(
$=w_{\textit{e}}(t, \varPhi )/w_{\textit{e}}(t, \unicode{x03C0} /2)$
) along the droplet interface over time. The inset shows the profiles near the droplet contact line. (c) Time-dependent tangential velocity profiles
$V_{\textit{t}}(t,\varPhi )$
along the droplet interface. The inset shows the profiles near the droplet contact line during the transition stage (
$t \approx412\pm6$
s). (d) The time evolution of the distribution of
$\varPhi ^*$
, an angle indicating locations of stagnation in which two opposing flows meet near the droplet interface. The lower inset shows the stagnation points at the transition stage (
$t \approx412\pm6$
s). In the upper inset,
$\gamma _{{a}}$
and
$\gamma _{{c}}$
are the surface tensions at the apex and the contact line, respectively, and
$\Delta \gamma _{{a-c}}$
indicates the surface tension difference between the apex and the contact line (
$\gamma _{{a}} - \gamma _{{c}}$
).
Then we examine the flow transitions in detail using the numerical results, focusing on the ethanol mass fraction (
$w_{\textit{e}}$
) and tangential velocity (
$V_{\textit{t}}$
), as illustrated in figure 5(a). To analyse the relative distribution of ethanol mass fraction along the droplet interface that determines surface tension gradient, we observe the normalised ethanol mass fraction
$\bar {w_{\textit{e}}} = w_{\textit{e}}(t, \varPhi )/w_{\textit{e}}(t, \unicode{x03C0} /2)$
, as depicted in figure 5(b), where
$\varPhi$
is the angle from the bottom substrate indicating the location at the droplet interface. As shown in the profile of
$\bar {w_{\textit{e}}}$
at
$t = 300$
s (see solid line), the ethanol mass fraction is relatively higher near the contact line in the early stage, inducing inward interfacial flow, which can be explained as negative tangential velocity in figure 5(c). Then the dotted lines in figures 5(b) and 5(c) show that the ethanol mass fraction along the droplet interface becomes nearly homogeneous at approximately
$t = 412$
s, resulting in tangential velocity approaching zero. This can be explained by many random stagnation points near
$t = 412$
s in figure 5(d), indicating that the surface tension gradients along the droplet interface fluctuate near zero. Since the ethanol mass fraction distribution becomes proportional to the uniform vapour concentration distribution, it is evident that the vapour–liquid equilibrium is achieved throughout the entire droplet interface at
$t \approx 412$
s. After the vapour–liquid equilibrium, the relative ethanol mass fraction at the contact line decreases over time, causing an outward interfacial flow in the late stage, as depicted by the dash-double-dotted lines in figures 5(b) and 5(c). This phenomenon can be attributed to the higher evaporation rate of ethanol at the contact line, similar to the evaporation in an ethanol–water binary mixture droplet, as observed in the experiment. Therefore, in a uniform vapour environment, the flow direction completely changes from inward to outward before and after vapour–liquid equilibrium, respectively.
3.2. Droplet surrounded by the high external vapour concentration gradient

Figure 6. Experimental and numerical results of the transient flow patterns within the droplet under the near-field source (
$d=2$
mm). The 2-D flow fields are observed using (a) PIV and (b) numerical simulation at the focal planes (
$a1, b1$
)
$z=200$
$\unicode {x03BC}$
m and (
$a2, b2$
)
$z=50$
$\unicode {x03BC}$
m. The black arrows are velocity vectors, and the white arrows represent typical flow patterns. The contour
$U_r$
is the
$r$
-direction velocity at each focal plane. We set
$t = 0$
for the moment when the ethanol vapour source is placed above the droplet apex, and
$t_{{eq}}$
indicates the time when the vapour–liquid equilibrium is achieved at the droplet contact line. (c) The ethanol mass flux (green arrows at the droplet interface) and typical flow structures at early and late stages reconstructed from numerical simulation (see supplementary movie 2). Blue and red arrows indicate the clockwise and anticlockwise directions, respectively.
Under the high external vapour concentration gradient, we experimentally observe two distinct flow patterns: (i) inward–outward competing flow in the early stage, and (ii) outward interfacial flow in the late stage. As shown in figure 6(a), radial inward flows form near the droplet centre at
$z=50$
and 200
$\unicode {x03BC}$
m, but the different direction flows occur near the droplet interface. The numerical results also showed a similar trend, as shown in figure 6(b), demonstrating that the inward and outward interfacial flows compete in the early stage (see figure 6
c). Especially in the experiment, the flow at
$z=200$
$\unicode {x03BC}$
m is non-axisymmetric and complicated because of the competition between flows in different directions. The generation of two opposing flows arises from the higher absorption flux at the contact line and the apex compared to the middle (see green arrows in figure 6
c). The absorptive flux near the apex is relatively high due to exposure to the high external vapour concentration, while near the contact line, the flux is elevated because of the low contact angle geometry (Wang et al. Reference Wang, Karapetsas, Valluri, Sefiane, Williams and Takata2021). This results in a reversed surface tension gradient near the contact line, as compared to near the apex.
In the late stage, we observe the radial inward flow at
$z=50$
$\unicode {x03BC}$
m, and the radial outward flow at
$z=200$
$\unicode {x03BC}$
m near the droplet interface in both experiment and simulation, as depicted in figures 6(a) and 6(b). This indicates that the inward interfacial flow near the contact line diminishes over time, and only the outward interfacial flow is left, which has been typically observed in previous literature (Malinowski et al. Reference Malinowski, Volpe, Parkin and Volpe2018; Park et al. Reference Park, Ryu, Sung and Kim2020). This flow transition occurs when the absorption flux at the contact line becomes almost zero, as illustrated in the late stage schematic in figure 6(c). Under the near-field source, the vapour concentration in air at the contact line is the lowest, so the required ethanol concentration in water to achieve the vapour–liquid equilibrium is the lowest according to Henry’s law. However, due to the droplet geometry, the absorption flux at the contact line is relatively high, inducing the most rapid establishment of vapour–liquid equilibrium at the contact line. After the equilibrium at the contact line, the reversed surface tension gradient near the contact line disappears, leading to the transitions from competing flow to outward interfacial flow. Then we define this moment when the vapour–liquid equilibrium reaches at the contact line as
$t_{{eq}}$
.

Figure 7. Numerical analysis for the high external vapour concentration gradient, namely the near-field source case (
$d=2$
mm). (a) A schematic of a droplet under the source, where
$\varPhi$
is the angle from the bottom substrate,
$V_{{t}}$
is the tangential velocity, and
$w_{{e}}$
indicates the ethanol mass fraction just below the interface, illustrated as a dark grey area. (b) Normalised ethanol mass fraction profiles
$\bar {w_{\textit{e}}}$
(
$=w_{\textit{e}}(t, \varPhi )/w_{\textit{e}}(t, \unicode{x03C0} /2)$
) along the droplet interface over time. The inset shows the profiles near the droplet contact line. (c) Time-dependent tangential velocity profiles
$V_{{t}}(t,\varPhi )$
along the droplet interface. The inset shows the profiles near the droplet contact line. (d) The time evolution of the distribution of
$\varPhi ^*$
, an angle indicating the locations of stagnation in which two opposing flows meet near the droplet interface. In the inset,
$\gamma _{{a}}, \gamma _{{c}},\, {\rm and}\,\gamma_s$
are the surface tensions at the apex, contact line and stagnation point, respectively, and
$\Delta \gamma$
indicates the surface tension difference (e.g.
$\Delta \gamma _{{a-s}}=\gamma _{{a}} - \gamma _{{s}}$
).
Numerical analysis reveals that in the early stage, the ethanol mass fractions at the droplet contact line and the apex are higher than at the middle between the apex and the contact line, as depicted in figure 7(b) (see solid and dashed lines). The ethanol concentration gradient is opposite near the contact line and the apex, causing a reversed surface tension gradient near the contact line compared to near the apex. This results in opposite-direction flows, supported by the negative and positive tangential velocities near the contact line and apex, respectively, as shown in figure 7(c). Over time, as illustrated in figure 7(d), the stagnation point – representing the location where inward and outward interfacial flows compete – gradually moves closer to the contact line. This indicates that the inward interfacial flow (characterised by negative tangential velocity) diminishes over time. At
$t \approx 145$
s, the inward flow completely disappears, as reflected in the ethanol mass fraction profile, which shows a gradual decrease from the apex to the contact line, and the tangential velocity profile, which displays only positive values during the later stages (see dotted and dash-dotted lines in figures 7
b,c). This flow transition indicates that the reversed surface tension gradient near the contact line in the early stage resolves into a relatively homogeneous surface tension gradient in the late stage. Notably, in the concentrated vapour field, vapour–liquid equilibrium is never fully established across the entire droplet interface. The transition from inward–outward competing flows to an outward interfacial flow is driven by the local establishment of vapour–liquid equilibrium near the contact line.
4. An unexpected significant surge in flow magnitude during the transition phase
Unlike the uniform vapour field, the concentrated vapour field exhibits a distinct flow pattern: a dramatic increase in flow magnitude immediately after the transition. While a similar flow trend – the development of flow magnitude over time – was previously observed by Park et al. (Reference Park, Ryu, Sung and Kim2020) and Kabi et al. (Reference Kabi, Pal and Basu2020), the underlying mechanism remains unexplored, and the generation of competing flows has not been clearly addressed yet. In this section, we aim to elucidate this phenomenon based on the numerical results and a first-order scaling model.

Figure 8. (a) Snapshots of the ethanol mass flux (
$J_{{e}}$
) along the interface and the resulting flow field in the droplet under the near-field source (
$d=2$
mm) from the numerical simulation at
$t=100$
, 180 and 250 s. The red arrows indicate the ethanol mass flux at the droplet interface. Inside the droplet, the ethanol mass fraction (left) and the flow velocity field (right) are described with different contours. The black arrows are the velocity vectors. See supplementary movie 2 for the simulation over the entire droplet lifetime. (b) Ethanol mass flux profiles along the droplet interface at
$t=100$
, 180 and 250 s. The dashed line means
$J_{{e}}=0$
. (c) Time evolution of the ethanol mass flux at the droplet contact line (
$J_{{e,c}}$
, black line) and the averaged tangential velocity (
$\bar {V}_{{t}}$
, red line). The blue dot and dashed line indicate the moment when the ethanol mass flux at the contact line reaches zero (
$t=180$
s). The red dot and dashed line represent the point at which the flow magnitude dramatically increases. (d) Time evolution of the theoretical velocities from (4.1) (
$V_{{th}}$
, green line) and (4.2) (
$V_{\textit{th}}^{\prime}$
, black dashed line), and the averaged tangential velocity (
$\bar {V_{{t}}}$
, red line).
As depicted in figure 8(a), during the flow transition, when the inward interfacial flow near the contact line (see
$t=100$
s) disappears at
$t=180$
s, a strong Marangoni flow suddenly emerges near the contact line, with the velocity magnitude increasing to approximately 10–
$10^2$
times that of the flow at
$t=100$
s. After that, the strong outward interfacial flow becomes dominant in the late stage (see
$t=250$
s). The generation of inward interfacial flow can be explained by the slightly increasing flux near the contact line at
$t=100$
s, as depicted in figure 8(b), indicating a reversed surface tension gradient near the contact line compared to that near the apex. Subsequently, the strong Marangoni flow forms when this reversed surface tension gradient diminishes after
$t=180$
s. Similarly, as depicted in figure 8(c), a notable increase in the slope of
$\bar {V}_{\textit{t}}$
(marked by a red dot) occurs at approximately
$t=180$
s when the ethanol mass flux at the contact line reaches zero (marked by a blue dot), leading to continuous decrease of the ethanol mass flux from the apex to the contact line. These results indicate that the sudden velocity increase follows the zero ethanol mass flux at the contact line right after the reversed surface tension gradient is resolved. This suggests that the competition between different surface tension gradients retards the flow magnitude in the early stage.
To clarify this, we derive a theoretical model for predicting the flow speed, with and without considering the competing surface tension gradients. For the sake of simplicity, we assume a quasi-steady state, and a droplet has a relatively flat shape (
$\epsilon = H/R \ll 1$
), so that
$Re\,\epsilon = \rho_w U R \epsilon / \mu_w$
is smaller than unity. Within a lubrication approximation, the Marangoni force along the droplet surface (
$\sim \Delta \gamma\, R$
), which drives the flow, can be balanced with the viscous force along the droplet height (
$\sim \mu _{{w}} V_{{th}}R^2/H$
) (Malinowski et al. Reference Malinowski, Volpe, Parkin and Volpe2018; Park et al. Reference Park, Ryu, Sung and Kim2020). As in previous studies, without accounting for the presence of surface tension gradients in different directions,
$\Delta \gamma$
is simply defined as the magnitude of the surface tension difference between the contact line (
$\gamma _{{c}}$
) and the apex (
$\gamma _{{a}}$
) of the droplet. The theoretical interfacial velocity is

where
$V_{{M}}=(\gamma _{\textit{c}}-\gamma _{\textit{a}})/\mu _{{w}}$
is the Marangoni velocity. In this case, the theoretical velocity
$V_{{th}}$
, estimated using the surface tension value obtained from the numerical simulation (the green line in figure 8
d), overestimates the averaged tangential velocity
$\bar {V}_{\textit{t}}$
(the red line in figure 8
d) in the early stage. This indicates that two surface tension gradients in different directions, arising from the competing flows, should be considered separately. To address this, we account for the location of the stagnation point
$\varPhi ^*$
to distinguish the surface areas where each surface tension gradient is applied. The surface areas are defined based on the ratio of angles captured from the inset of figure 7(d):
${\varPhi ^*}/({\unicode{x03C0} /2})$
and
$ 1-{\varPhi ^*}/({\unicode{x03C0} /2} )$
. Then the modified theoretical velocity could be expressed as

where
$\gamma _{{s}}$
is a surface tension at the stagnation point. The black dashed line in figure 8(d), representing
$V_{{th}}^{\prime}$
, could estimate the trend of the flow retardation in the early stage of
$\bar {V}_{{t}}$
. This equation suggests that the coexistence of two surface tension gradients in the early stage causes a relatively weak flow due to the reduction of Marangoni force. The flow retardation resolves when
$\varPhi ^*$
disappears, indicating the achievement of vapour–liquid equilibrium at the droplet contact line. As a result, monitoring the time-varying vapour mass flux and the resulting surface tension gradients is indispensable for estimating not only the flow structure but also the flow magnitude, which varies significantly before and after vapour–liquid equilibrium.
5. Conclusions
In this study, we conducted both experimental and numerical investigations into the time-dependent ethanol vapour mass transport across the liquid–gas interface and the resulting internal flow patterns of the water droplets under different external ethanol vapour distributions. We observed two distinct flow transitions in the experiments under low and high external vapour concentration gradients, which have not been previously reported. For a low vapour concentration gradient, we observed not only a typical inward interfacial flow as reported by Majumder et al. (Reference Majumder2012) and Kim (Reference Kim2022), but also the formation of multiple vortices. In the case of a high external vapour concentration gradient, we noticed a competition between inward and outward interfacial flows along the droplet surface before generating the typical outward interfacial flow previously observed by Malinowski et al. (Reference Malinowski, Volpe, Parkin and Volpe2018) and Ryu et al. (Reference Ryu, Ko and Kim2022). These flow transitions mainly occur due to the achievement of vapour–liquid equilibrium at the droplet contact line. In the early stage, absorption flux near the droplet contact line is relatively high due to the droplet geometry (contact angle is less than 90
$^\circ$
). This high absorption flux, compared to the low vapour concentration in air (indicating that the required vapour concentration in water to achieve equilibrium is low), leads to the fast achievement of vapour–liquid equilibrium at the droplet contact line. It causes a shift in the direction of the surface tension gradient near the droplet contact line in the late stage, inducing the transition in the flow structure.
Especially in the high external vapour concentration gradient, a dramatic increase in flow speed occurs during the flow transition. Even though this velocity increase was observed in previous studies (Park et al. Reference Park, Ryu, Sung and Kim2020; Kabi et al. Reference Kabi, Pal and Basu2020), a reason behind this has not been discussed before due to the lack of observation of vapour mass flux. According to our findings, the generation of the inward and outward competing flows retards the flow magnitude in the early stage due to the presence of opposing surface tension gradients at the interface. In the late stage, as the inward interfacial flow diminishes and the flow direction changes near the droplet contact line, the flow magnitude then significantly increases. By estimating the Marangoni stresses applied on the droplet interface considering the competing surface tension gradients, we reveal that the flow magnitude can be significantly retarded by the opposing surface tension gradients in the small area near the contact line. Consequently, we can predict when the flow retardation ends by capturing the vapour–liquid equilibrium at the droplet contact line.
Our results demonstrate that in any vapour environment, the flow structure and the magnitude vary significantly depending on the achievement of vapour–liquid equilibrium at the droplet contact line. For applications that utilise flow control via the vapour-driven solutal Marangoni effect, such as mixing and particle deposition, it is crucial to estimate the flow structure during the whole lifetime of a droplet. Therefore, understanding the time-varying vapour mass flux along the droplet interface and the associated transition in flow patterns is critical. Our experimental observation and numerical simulation would be helpful for facilitating better utilisation of vapour-driven Marangoni flow.
Supplementary material and movies
Supplementary material and movies are available at https://doi.org/10.1017/jfm.2025.10567.
Funding
This work was supported by the Basic Science Research Programme through the National Research Foundation of Korea funded by the Ministry of Science (NRF-2021R1A2C2007835). This paper is partially based on research that has been conducted as part of Samsung Electronics (IO201216–08212-01). We also acknowledge support from: an Industrial Partnership Programme, High Tech Systems and Materials (HTSM) of the Netherlands Organisation for Scientific Research (NWO); a funding for public-private partnerships (PPS) of the Netherlands Enterprise Agency (RVO) and the Ministry of Economic Affairs (EZ); Canon Production Printing Netherlands BV; IamFluidics BV, TNO Holst Centre, and the University of Twente, Eindhoven University of Technology and Utrecht University.
Declaration of interests
The authors report no conflict of interest.
Data availability statement
The data that support the findings in this study are available from the authors upon reasonable request.
Author contributions
J.R. and H.K designed the research, J.R. performed the experiments, and C.D. performed the simulations. All authors contributed equally to analysing data and reaching conclusions, and in writing the paper.