# Intricate Role of Thermal Properties and Volatility in Droplet Spreading: A
Generalization to Tanner’s Law

Zhenying Wang12<EMAIL_ADDRESS>George Karapetsas3 Prashant
Valluri4  Chihiro Inoue1  1Department of Aeronautics and Astronautics, Kyushu
University, Nishi-Ku, Motooka 744, Fukuoka 819-0395, Japan 2International
Institute for Carbon-Neutral Energy Research (WPI-I2CNER), Kyushu University,
Nishi-Ku, Motooka 744, Fukuoka 819-0395, Japan 3Department of Chemical
Engineering, Aristotle University of Thessaloniki, Thessaloniki 54124, Greece
4Institute of Multiscale Thermofluids, School of Engineering, University of
Edinburgh, Edinburgh EH9 3JL, UK

###### Abstract

Droplet spreading is ubiquitous and plays a significant role in liquid-based
energy systems, thermal management devices, and microfluidics. While the
spreading of non-volatile droplets is quantitatively understood, the spreading
and flow transition in volatile droplets remains elusive due to the complexity
added by interfacial phase change and non-equilibrium thermal transport. Here
we show, using both mathematical modeling and experiments, that the wetting
dynamics of volatile droplets can be scaled by the spatial-temporal interplay
between capillary, evaporation, and thermal Marangoni effects. We elucidate
and quantify these complex interactions using phase diagrams based on
systematic theoretical and experimental investigations. A spreading law of
evaporative droplets is derived by generalizing Tanner’s law (valid for non-
volatile liquids) to a full range of liquids with saturation vapor pressure
spanning from $10^{1}$ to $10^{4}$ Pa and on substrates with thermal
conductivity from $10^{-1}$ to $10^{3}$ W/m/K. Besides its importance in
fluid-based industries, the conclusions also enable a unifying explanation to
a series of individual works including the criterion of flow reversal and the
state of dynamic wetting, making it possible to control liquid transport in
diverse application scenarios.

###### keywords:

Authors should not enter keywords on the manuscript, as these must be chosen
by the author during the online submission process and will then be added
during the typesetting process (see
http://journals.cambridge.org/data/relatedlink/jfm-keywords.pdf for the full
list)

##  1 Introduction

Wetting presents the basic state when a finite volume of liquid contacts a
solid substrate (De Gennes (1985), Bonn et al. (2009)). The requirement for
wetting differs for different applications, and hence the design principles.
For example, good wetting enables efficient liquid supply and avoids dry-out
in phase change thermal management, _e.g._ , electronics cooling (Mudawar
(2001)). While in condensation, a quick wetting transition will increase the
substrate’s thermal resistance and worsen the condensation performance (Rose
(2002)). A full understanding of the interplay between wetting and phase
change is thus necessary for in-demand design of various energy systems (Tu et
al. (2018), Breitenbach et al. (2018)) as well as in other microfluids based
industries (Zheng et al. (2010), Teh et al. (2008)).

The wetting dynamics of nonvolatile liquids are known to be governed by
capillary force and viscous dissipation. Specifically, the spreading rate of a
nonvolatile Newtonian droplet in perfect wetting cases can be quantitatively
described by Tanner’s law, $R(t)\sim t^{1/10}$ (Tanner (1979), Lelah & Marmur
(1981), Blake (2006)), which can be derived mathematically by employing the
lubrication approximation and balancing the effect of capillary forces with
the resisting viscous forces generated by the droplet motion (Carlson (2018),
Pang & Ó Náraigh (2022), Voinov (1976)). The spreading law was then validated
by successive experiments (Marmur (1983)) and further extended to non-
Newtonian liquids (Rafaï et al. (2004), Jalaal et al. (2021)), gravity-
influenced cases (Cazabat & Stuart (1986)), near-critical cases (Saiseau et
al. (2022)), liquids on non-rigid solid and on thin liquid films (Carré et al.
(1996), Cormier et al. (2012)), as well as under-liquid systems (Mitra & Mitra
(2016), Goossens et al. (2011)).

Compared to nonvolatile liquids, the spreading of volatile liquids is more
complex and remains elusive owing to the complexity added by interfacial phase
change and non-equilibrium thermal transport (Jambon-Puillet et al. (2018),
Wang et al. (2022), Berteloot et al. (2008), Wang et al. (2019)).
Specifically, the non-uniform mass flux tends to recede the liquid-air
interface especially near the three phase contact line. The effect of
evaporation cooling and the non-uniform pathway for heat dissipation induces
temperature gradient and therefore thermal Marangoni flow, which complexes the
flow pattern inside the droplet and further affects the spreading dynamics.

A full elucidation of these interacting mechanisms requires a quantitative
description of each physical process, as well as a thorough understanding of
all influencing parameters. In this research, we combine experimental and
numerical approaches to address these issues and revisit the milestone
descriptions of flow interaction in liquid spreading. The explorations enable
us to quantify the strength of dominating physics that varies spatially and
temporally, which on the other hand provides a unifying criterion for flow
reversal (Hu & Larson (2005), Xu & Luo (2007), Ristenpart et al. (2007), Xu et
al. (2010)), pattern formation (Hu & Larson (2006), Deegan et al. (1997), Li
et al. (2015)), and Marangoni-capillary interactions (Tsoumpas et al. (2015),
Shiri et al. (2021), Yang et al. (2022)) in drying droplets.

Specifically, we investigate the flow state and the resulting wetting dynamics
of volatile liquids on thermal conductive substrates. In the numerical aspect,
previous lubrication-type models usually consider boundary conditions with
uniform heating (Anderson & Davis (1995), Ajaev (2005)) or preset substrate
temperature/temperature gradient (Karapetsas et al. (2013), Dai et al. (2016),
Charitatos & Kumar (2020), Wang et al. (2021)), which failed to illustrate the
non-negligible role of non-isothermal heat conduction in the overall process.
Additionally, the numerical and experimental investigations have at large been
developing in parallel without strict comparisons that consider the specific
application scenarios. In this work we take account of the mass, momentum, and
energy transport in the liquid phase, as well as the heat conduction through
the solid substrate. We apply the assumption of a precursor film (ultrathin
adsorbed liquid layer, Kavehpour et al. (2003), Hoang & Kavehpour (2011),
Ajaev (2005)) existing in front of the three-phase contact line, which
eliminates the stress singularity and avoids experimental inputs on the
contact line motion, _e.g._ , preset advancing/receding contact angles or a
finite slip. A strict one-to-one comparison is then conducted between
experimental and numerical results for a broad range of liquid-solid
combinations. The reliability of the theoretical models further give rise to a
general principle of flow transition and non-isothermal spreading of volatile
droplets.

The rest of this paper is organized as follows. In §2, we introduce the
formulated mathematical model, the scaling and simplification, and the
numerical method. §3 introduces the material preparation, the detailed
experimental setup and procedures. §4 presents the numerical and experimental
results including the dominating mechanisms, the transition of flow pattern,
and their relation to the spreading dynamics. §5 summarizes the key findings
of this work, and share our perspectives on the spreading of Newtonian/non-
Newtonian fluids with/without interfacial phase change.

##  2 Formulation of Numerical Model

We consider a thin liquid droplet (aspect ratio
$\varepsilon=\hat{H}_{0}/\hat{R}_{0}\ll 1$, incompressible Newtonian fluid)
sitting on a thermal conductive solid substrate with finite thickness and
surrounded by a mixture of its own vapor and noncondensing gas. The low aspect
ratio allows for the application of lubrication theory which simplifies the
Navier-Stokes equations while effectively describes the dominating physical
processes. The lubrication approximation is applicable for droplets with small
and moderate contact angles typically less than 40°, which is suitable for
current research on complete wetting and partial wetting cases. A precursor
film is assumed to exist at the solid surface in front of the three phase
contact line. This avoids the stress singularity that may arise due to the
inherent contradiction in simultaneously assuming the no-slip boundary
condition and expecting a displacement between liquid and gas at the moving
contact line.

Figure 1: (1) A single component droplet sitting on a solid substrate in
contact with its own vapor and noncondensing gas: The aspect ratio,
$\hat{H}_{0}/\hat{R}_{0}$, is assumed as very small, $c_{\mathrm{sat}}$
denotes the saturation vapor concentration at the liquid-air interface,
$c_{\mathrm{g}}$ refers to the vapor concentration in the gas phase and keeps
constant, $\hat{H}_{\mathrm{solid}}$ denotes the substrate thickness and
$\hat{k}_{\mathrm{solid}}$ denotes the thermal conductivity. $\vec{n}$ and
$\vec{t}$ denote the outward units vectors acting in normal and tangential
directions to the interface. The center of droplet base, _O_ , is defined as
the origin of the coordinate. (2) Governing equations in the gas, liquid and
solid phases along with the boundary equations at the central axis, the
liquid-gas interface and the liquid-solid interface.

Governing Equations A cylindrical coordinate system,
$(\hat{r},\theta,\hat{z})$, is applied to solve the velocity field,
$\hat{u}=(\hat{u},\hat{v},\hat{w})$, where $\hat{u}$, $\hat{v}$, and $\hat{w}$
correspond to the horizontal, azimuthal and vertical components of the
velocity field respectively. The liquid-vapor interface locates at
$\hat{z}=\hat{h}(\hat{r},\hat{t})$ and the liquid-solid and solid-vapor
interfaces locate at $\hat{z}=0$. The liquid phase is governed by the
incompressible mass, momentum and energy equations, Eqs. (1, 2, 3),

$\hat{\nabla}\cdot\hat{\boldsymbol{u}}=0,$ (1)
$\hat{\rho}(\frac{\partial\hat{\boldsymbol{u}}}{\partial\hat{t}}+\hat{\boldsymbol{u}}\cdot\hat{\nabla}\hat{\boldsymbol{u}})=\hat{\nabla}\cdot\hat{\boldsymbol{E}},$
(2)
$\hat{\rho}\hat{c}_{p}(\frac{\partial\hat{T}}{\partial\hat{t}}+\hat{\boldsymbol{u}}\cdot\hat{\nabla}\hat{T})=\hat{\nabla}\cdot\hat{k}\hat{\nabla}\hat{T},$
(3)

where liquid density $\hat{\rho}$, specific heat capacity $\hat{c}_{p}$, and
thermal conducitivity $\hat{k}$ are set as constant, $\hat{T}$ denotes
temperature, $\hat{t}$ denotes time, and $\hat{\boldsymbol{E}}$ denotes the
total stress tensor at the liquid side,

$\hat{\boldsymbol{E}}=-\hat{p}\boldsymbol{I}+\hat{\mu}(\hat{\nabla}\hat{\boldsymbol{u}}+\hat{\nabla}\hat{\boldsymbol{u}}^{T}),$
(4)

where $\hat{\mu}$ is the dynamic viscosity (linear function of liquid
temperature), and $\boldsymbol{I}$ is the identity tensor.

At the liquid-gas interface, the liquid velocity and geometry satisfy a force
balance both in the normal and tangential directions. The normal stress
balance relates the normal stress jump, surface tension, mean curvature and
van der Waals interactions,

$-\hat{p}+\boldsymbol{n}\cdot\hat{\boldsymbol{\tau}}\cdot\boldsymbol{n}=2\hat{\kappa}\hat{\sigma}+\hat{\Pi}-\hat{p}_{g},$
(5)

where $\hat{p}$ is the pressure of the liquid phase, $\hat{p}_{g}$ is the
total pressure of the gas phase, $\hat{\boldsymbol{\tau}}$ is the shear stress
tensor of the liquid phase.
$2\hat{\kappa}=-\hat{\nabla}_{S}\cdot\boldsymbol{n}$ is twice the mean
curvature of the free surface and
$\hat{\nabla}_{S}=(\boldsymbol{I}-\boldsymbol{n}\boldsymbol{n})\cdot\hat{\nabla}$
is the surface gradient operator. $\hat{\Pi}$ denotes the disjoining pressure
accounting for intermolecular interactions near the contact line (De Gennes et
al. (1990)),

$\hat{\Pi}=\frac{\hat{A}}{6\pi\hat{h}^{3}},$ (6)

with $\hat{A}$ being the dimensional Hamaker constant.

The tangential stress balance relates the shear stress jump and the surface
tension gradient. By ignoring the shear stress from the gas phase due to small
viscosity, we arrive at,

$\boldsymbol{n}\cdot\hat{\boldsymbol{E}}\cdot\boldsymbol{t}=\hat{\nabla}_{S}\hat{\sigma}\cdot\boldsymbol{t},$
(7)

where $\boldsymbol{t}$ refers to the outward vectors which is tangential to
the interface.

The motion of free surface can be described with the kinematic boundary
condition, expressed as,

$\frac{\partial\hat{h}}{\partial\hat{t}}+\hat{u}_{S}\frac{\partial\hat{h}}{\partial\hat{r}}+\frac{\hat{v}_{S}}{\hat{r}}\frac{\partial\hat{h}}{\partial\theta}=\hat{w}_{S}.$
(8)

Along the droplet interface, the boundary equation of mass balance can be
expressed by the relationship between the velocity of the liquid solution,
$\boldsymbol{\hat{u}}$, and the velocity of the interface,
$\boldsymbol{\hat{u}}_{s}$,

$(\boldsymbol{\hat{u}}-\boldsymbol{\hat{u}}_{S})\cdot
n=\frac{\hat{J}}{\hat{\rho}},$ (9)

where $\hat{J}$ denotes the interfacial mass flux of water vapor. The jump
energy balance can then be derived taking account of the heat release at the
liquid-vapor interface and ignoring the heat conduction into the gas phase,

$\hat{J}\hat{L}+\hat{k}\hat{\nabla}\hat{T}\cdot\boldsymbol{n}=0,$ (10)

where $\hat{L}$ denotes the latent heat of vaporization. The Hertz–Knudsen
equation is commonly used for predicting the mass flux induced by evaporation
or condensation towards a liquid–vapor interface. At thermodynamic
equilibrium, the chemical potentials of the gas phase and liquid phase across
the liquid–air interface reach balance, and the following relation can be
derived combining the Hertz–Knudsen equation (Plesset & Prosperetti (1976),
Moosman & Homsy (1980)) and the chemical balance relations (Atkins et al.
(2014)).

$\hat{J}=\hat{p}_{\mathrm{v,sat}}\sqrt{\frac{\hat{M}}{2\pi\hat{R}_{g}\hat{T}_{g}}}(\frac{\hat{M}}{\hat{\rho}\hat{R}_{g}\hat{T}_{g}}(\hat{p}-\hat{p}_{g})+\frac{\hat{M}\hat{L}}{\hat{R}_{g}\hat{T}_{g}^{2}}(\hat{T}_{i}-\hat{T}_{g})+ln(\frac{1}{\chi_{\mathrm{vapor}}})),$
(11)

where $\hat{M}$ is the molar mass of water, $\hat{R}_{g}$ is the gas constant,
$\hat{T}_{g}$ is the temperature of gas phase, $\hat{T}_{i}$ is the
temperature at the liquid-gas interface, $\hat{p}_{\mathrm{v,sat}}$ is the
saturation vapor pressure in the gas phase, and is the relative vapor
concentration, defined as the ratio of vapor pressure to the saturation vapor
pressure in the gas phase.

Scaling The key parameters are scaled as follows,

$\begin{split}\hat{r}=\hat{R}_{0}r,\ \hat{z}=\hat{H}_{0}z,\
(\hat{u},\hat{w})=(\hat{u}^{*}u,\frac{\hat{H}_{0}}{\hat{R}_{0}}\hat{u}^{*}w),\\\
\hat{\mu}=\hat{\mu}_{0}\hat{\mu},\
\hat{T}=\hat{T}_{\mathrm{ref}}+T\Delta\hat{T},\
\hat{\sigma}=\hat{\sigma}_{0}\sigma=\hat{\sigma}_{0}(1+\frac{\hat{\eta}_{\sigma}\Delta\hat{T}}{\hat{\sigma}_{0}}T),\\\
\hat{c}_{p}=\hat{c}_{p,0}c_{p},\ \hat{k}=\hat{k}_{0}k,\
\hat{t}=\frac{\hat{R}_{0}}{\hat{u}^{*}}t,\
\hat{p}=\hat{p}_{g}+\frac{\hat{\mu}_{0}\hat{u}^{*}\hat{R}_{0}}{\hat{H}_{0}^{2}}p,\
\hat{J}=\frac{\hat{k}_{0}\Delta\hat{T}}{\hat{H}_{0}\hat{L}}J.\end{split}$ (12)

We define characteristic velocity as
$\hat{u}^{*}=\frac{\varepsilon\hat{\eta}_{\sigma}\Delta\hat{T}}{\hat{\mu}_{0}}$
and derive the dimensionless group, including the Reynolds number,
$\mathrm{Re}=\frac{\hat{\rho}_{0}\hat{u}^{*}\hat{H}_{0}}{\hat{\mu}_{0}}$,
Stokes number,
$\mathrm{St}=\frac{\hat{\rho}_{0}\hat{g}_{G}\hat{H}_{0}^{3}}{\hat{\mu}_{0}\hat{u}^{*}\hat{R}_{0}}$,
Prandtl number, $\mathrm{Pr}=\frac{\hat{\mu}_{0}\hat{c}_{p,0}}{\hat{k}_{0}}$,
evaporation number,
$E=\frac{\hat{k}_{0}\hat{R}_{0}\Delta\hat{T}}{\hat{\rho}_{0}\hat{L}\hat{H}_{0}^{2}\hat{u}^{*}}$,
Marangoni number,
$\mathrm{Ma}=\frac{\hat{\mu}_{0}\hat{u}^{*}}{\varepsilon\hat{\sigma}_{0}}=\frac{\hat{\eta}_{\sigma}\Delta\hat{T}}{\hat{\sigma}_{0}}$,
Hamaker constant,
$\mathcal{A}=\frac{\hat{\mathcal{A}}}{6\pi\hat{\mu}_{0}\hat{u}^{*}\hat{R}_{0}\hat{H}_{0}}$,
the coefficients in the expression of mass flux,
$\delta=\frac{\hat{M}\hat{\eta}_{\sigma}\Delta\hat{T}}{\hat{\rho}_{0}\hat{R}_{g}\hat{T}_{g}\hat{H}_{0}}$
(measure of Kelvin effect),
$\psi=\frac{\hat{M}\hat{L}\Delta\hat{T}}{\hat{R}_{g}\hat{T}_{g}^{2}}$ (effect
of local temperature difference on the mass flux), and Jakob number,
$\mathrm{Ja}=\frac{\hat{k}_{0}\Delta\hat{T}}{\hat{H}_{0}\hat{L}\hat{p}_{v,sat}}\sqrt{\frac{2\pi\hat{R}_{g}\hat{T}_{g}}{\hat{M}}}$,
measuring the importance of kinetic effects at the interface.

For sessile droplets with slow interfacial phase change, the inertia effect
can be neglected with very small Reynolds number (usually Re $<10^{-1}$ for
droplet wetting and evaporation), therefore we set the Reynolds number to zero
in the formulation. By further ignoring the parameter variation in the
azimuthal direction, we arrive at the dimensionless governing equations,
including the conservation equations of mass, Eq. (13), momentum, Eq. (14) and
energy, Eq. (15), as well as the dimensionless boundary equations, including
the normal and tangential stress boundary balance, Eq. (16), the kinematic
boundary condition, Eq. (17), the jump energy balance, Eq. (18), and the
expression of interfacial mass flux, Eq. (19).

$\frac{1}{r}\frac{\partial(ru)}{\partial r}+\frac{\partial w}{\partial z}=0,$
(13) $\frac{\partial p}{\partial r}=\frac{\partial}{\partial
z}(\mu\frac{\partial u}{\partial z}),$ (14) $\frac{\partial}{\partial
z}(k\frac{\partial T}{\partial z})=0,$ (15)
$p=-\frac{\varepsilon^{2}\sigma}{Ma}(\frac{1}{r}\frac{\partial}{\partial
r}(r\frac{\partial h}{\partial r}))-\frac{\mathcal{A}}{h^{3}},\
\tau_{zr}=\frac{1}{Ma}\frac{\partial\sigma}{\partial r},$ (16) $\frac{\partial
h}{\partial t}+u\frac{\partial h}{\partial r}-w+EJ=0,$ (17) $J+k\frac{\partial
T}{\partial z}=0,$ (18) $\mathrm{Ja}J=\delta
p+\psi(T_{i}-T_{g})+\mathrm{ln}(\frac{1}{\chi_{\mathrm{vapor}}}),$ (19)

The above equations are integrated over the vertical direction, with the
integration of velocity and temperature expressed as $f=\int_{0}^{h}udz,\
\Theta=\int_{0}^{h}Tdz$, and the governing equations are subsequently
discretized using a finite element/Galerkin method; the variables are
approximated using quadratic Lagrangian basis functions $\Phi_{i}$. After
applying the divergence theorem, we arrive at the weak form of five governing
equations, Eqs. (20-24), with five independent variables, $h$, $p$, $f$,
$\Theta$, and $J$.

$R=\int(\frac{\partial h}{\partial
t}\Phi_{i}+EJ\Phi_{i}-f\frac{\partial\Phi_{i}}{\partial
r})rdr+(rf\Phi_{i})_{r=0}^{r=r_{\infty}},$ (20) $R=\int(h\frac{\partial
p}{\partial r}-(\mu\frac{\partial u}{\partial z})_{0}^{h})\Phi_{i}rdr,$ (21)
$R=\int(k\frac{\partial T}{\partial z})_{0}^{h}\Phi_{i}rdr,$ (22)
$R=\int((p+\frac{\mathcal{A}}{h^{3}})\frac{Ma}{\varepsilon^{2}\sigma}\Phi_{i}-\frac{\partial
h}{\partial r}\frac{\partial\Phi_{i}}{\partial r})rdr+(r\frac{\partial
h}{\partial r}\Phi_{i})_{r=0}^{r=r_{\infty}},$ (23)
$R=\int(\mathrm{Ja}J-\delta
p-\psi(T|_{z=h}-T_{g})+\mathrm{ln}\chi_{\mathrm{vapor}})\Phi_{i}rdr.$ (24)

Finally, we consider an Ansatz distribution function of liquid temperature _T_
and velocity _u_ in _z_ axis to complete the mathematical formulation. We fit
the flow velocity and the temperature distribution with a parabolic
expression, _i.e._ $y=c_{3}+c_{2}z+c_{1}z^{2}$, then derive their expressions
by applying the boundary conditions at $z=0$ and $z=h$.

At the liquid-gas interface, $z=h$, the flow velocity meets the tangential
stress balance Eq. (16), which relates the surface tension gradient to the
shear stress (multiplication of dynamic viscosity and velocity gradient). At
the liquid-solid interface, $z=0$, a no slip boundary condition applies. The
flow velocity is thereafter derived as,

$u=(\frac{3f}{h^{2}}-\frac{1}{2\mu\mathrm{Ma}}\frac{\partial\sigma}{\partial
r})z+(\frac{3}{4\mu h\mathrm{Ma}}\frac{\partial\sigma}{\partial
r}-\frac{3f}{2h^{3}})z^{2}.$ (25)

For the temperature distribution, the temperature meets the jump energy
balance at the droplet surface, $z=h$, Eq. 18, which relates the temperature
gradient in the liquid side to the interfacial heat flux due to phase change.
At the liquid-solid interface, we consider substrates both with constant
temperature and with downward heat dissipation. For cases where the substrate
temperature is kept constant with temperature controller, the temperature
distribution in the liquid side is derived as,

$T=T_{W}+(\frac{3\Theta}{h^{2}}-\frac{3T_{W}}{h}+\frac{J}{2k})z+(\frac{3T_{W}}{2h^{2}}-\frac{3\Theta}{2h^{3}}-\frac{3J}{4kh})z^{2}.$
(26)

For a thermal conductive substrate with finite thickness $H_{\rm solid}$ and
thermal conductivity $k_{\rm solid}$, the process can be simplified into a
quasi-steady one-dimensional problem. Similar approximations were also adopted
in the studies of Ristenpart et al. (2007), Xu et al. (2010) and Dunn et al.
(2009). At small Reynold numbers, the heat transfer in the vertical direction
is dominated by heat conduction rather than convection. In such cases, the
boundary condition at $z=0$ in this situation can be described as,

$k\frac{\partial T}{\partial z}=\frac{T-T_{W}}{\rm Bi},$ (27)

where Biot number, $\mathrm{Bi}=H_{\mathrm{solid}}\textfractionsolidus
k_{\mathrm{solid}}=\hat{H}_{\mathrm{solid}}\hat{k}_{0}\textfractionsolidus\hat{H}_{0}\hat{k}_{\mathrm{solid}}$,
denotes the ratio of thermal resistances inside the solid phase and the liquid
phase. The distribution of liquid temperature in such cases is thereafter
derived as,

$\begin{split}T=c_{3}+c_{2}z+c_{1}z^{2},\mathrm{where}\ \\\
c_{1}=\frac{1}{3\mathrm{Bi}k+h}(-\frac{3\mathrm{Bi}J}{2h}-\frac{3J}{4k}-\frac{3\Theta}{2h^{2}}+\frac{3T_{W}}{2h}),\\\
c_{2}=\frac{1}{3\mathrm{Bi}k+h}(\frac{hJ}{2k}+\frac{3\Theta}{h}-3T_{W}),\\\
c_{3}=\frac{1}{3\mathrm{Bi}k+h}(\frac{\mathrm{Bi}hJ}{2}+\frac{3\mathrm{Bi}k\Theta}{h}+hT_{W}).\end{split}$
(28)

Initial and Boundary Conditions In the region of $0\leq r\leq 1$, the initial
conditions, $t=0$, are set as,

$h(r,0)=h_{\infty}+1-r^{2},\ f(r,0)=0,\ \Theta(r,0)=h(r,0)T_{W},$ (29)

where $h_{\infty}$ is the thickness of the precursor film, $T_{W}$ is the
initial temperature of the solid substrate. For zero mass flux at the
precursor region, the initial value of $h_{\infty}$ is set as
$\sqrt[3]{\frac{\mathcal{A}\delta}{\psi(T_{W}-T_{g})-\mathrm{ln}(\chi_{\mathrm{vapor}})}}$.
The parameters at the droplet center, $r=0$, satisfy the symmetric boundary
conditions, expressed as,

$\frac{\partial h}{\partial r}(0,t)=0,\ f(0,t)=0,\
\frac{\partial\Theta}{\partial r}(0,t)=0.$ (30)

In the region of precursor film, $r>1$, the parameters are set as,

$h(r,0)=h_{\infty},\ f(r,0)=0,\ \Theta(r,0)=h_{\infty}T_{W}.$ (31)

We apply finite element Galerkin method to solve the system of governing
equations. One dimensional mesh is constructed along the _r_ -direction, and
the length of the computational domain is set properly according to the actual
situation of each condition to ensure the free development of the droplet
morphology. The domain is discretized along _r_ from 0 to $r_{\infty}$ into
equally spaced nodes, the total number being $N_{\mathrm{tot}}$, with grid
independence test of the solution. At each time step, we apply the Newton-
Raphson method to obtain the solution across the computational domain
progressively. The solution evolves forward with an adaptive time interval
which adjusts according to the maximum residual errors of the governing
equations from the previous time step.

##  3 Experimental Method

We use high speed camera (Photron FASTCAM Mini AX200 type with frame rate of
3,000 fps) to trace the fast spreading of droplet upon its contact with a
solid surface. To avoid the influence of initial impinging force of the
droplet, we set the distance between the syringe and the substrate to a
similar value with the initial diameter of the spherical droplet ($V_{0}=3\pm
0.2\mu L$ for 1-Butanol and 2-Propanol, and $0.6\pm 0.05\mu L$ for FC-72), so
that the droplet can touch the solid surface right after it is generated from
the syringe tip. We selected 1-Butanol, 2-Propanol (IPA) and FC-72 as the
working fluids, each with one order of magnitude larger saturation vapor
pressure (Table 1). Three types of solid materials are utilized, _i.e._ ,
copper, SUS430 and glass, each with one order of magnitude smaller thermal
conductivity. We also chose copper substrates of thickness 0.3mm, 1mm and 3mm
to demonstrate the effect of substrate thickness. The three types of
substrates exhibit a highly wetting state for 1-Butanol, 2-Propanol and FC-72
which matches the requirement of low contact angle (less than 40°) for
lubrication approximation. For each condition, we repeated the experiments for
at least 5 times. The average values of the spreading rates along with
standard deviations are derived for comparison with the numerical predictions
(Table 3). The key thermophysical properties of different liquids and solids
are given in Table 1 and 2.

Liquid |  $\hat{\rho}$ ($kg/m^{3}$) |  $\hat{\mu}$ ($mPa\cdot s$) |  $\hat{\sigma}$ ($mN/m$) |  $\hat{\eta_{\sigma}}$ ($mN/m/K$) |  $\hat{k}$ ($W/m/K$) |  $\hat{c}_{p}$ ($J/kg/K$) |  $\hat{p}_{\mathrm{sat}}$ (Pa) |  $\hat{L}$ ($kJ/kg$)
---|---|---|---|---|---|---|---|---
Water | 998 | 1.005 | 72.75 | -0.15 | 0.60 | 4180 | 2341 | 2257
1-Butanol | 810 | 2.573 | 24.57 | -0.07 | 0.1553 | 2386 | 580 | 584
2-Propanol | 786 | 2.012 | 23.0 | -0.09 | 0.135 | 1527 | 4420 | 663
FC-72 | 1680 | 0.64 | 10 | - | 0.057 | 1100 | 30900 | 88
Methanol | 792 | 0.59 | 22.7 | -0.0773 | 0.204 | 2460 | 10300 | 1191
Ethanol | 789 | 1.144 | 22 | -0.0832 | 0.166 | 2570 | 5786 | 1030
Chloroform | 1490 | 0.563 | 27.5 | -0.1295 | 0.13 | 957 | 25900 | 248
Heptane | 684 | 0.389 | 20.2 | -0.098 | 0.128 | 2220 | 4621 | 3553
Hexane | 659 | 0.301 | 18.4 | -0.1022 | 0.125 | 2243 | 19140 | 382
Table 1: Thermophysical properties of representative liquids at 293 K except FC-72 (298K). 1-Butanol, 2-Propanol and FC-72 are utilized in present experiments. Water, 2-Propanol, Methanol, Ethanol and Chloroform are utilized in previous studies on circulation reversal (Ristenpart et al. (2007), Xu et al. (2010), Li et al. (2015)) as marked in Fig. 4. Alkanes (e.g., Heptane and Hexane) are utilized in the study of Shiri et al. (2021) on the influence of thermal Marangoni flow on droplet geometry as marked in Fig. 6. Material Type | Copper | SUS430 | Glass | Polydimethylsiloxane (PDMS)
---|---|---|---|---
$k_{\mathrm{solid}}$ ($W/m/K$) | 385 | 25 | 1.2 | 0.23
Table 2: Thermal conductivity of representative substrate materials. Copper,
SUS430 and glass are utilized in present study. Glass is utilized in the work
of Xu et al. (2010) and Hu & Larson (2006). PDMS is utilized in the work of
Ristenpart et al. (2007).

To remove any oxide residuals and ensure the uniformity and smoothness of the
solid surface, we fabricate substrate surfaces with mirror finishing. The
surface roughness of glass, SUS430 and copper is characterized to be
$\mathrm{Sa}_{\mathrm{glass}}=0.258\mu m$,
$\mathrm{Sa}_{\mathrm{SUS430}}=0.134\mu m$,
$\mathrm{Sa}_{\mathrm{copper}}=0.291\mu m$ with 3D laser scanning microscope
(OLYMPUS/LEXT OLS4000). Before each set of experiments, we prepare the solid
substrate by 15 min ultrasonic bath with 99.8$\%$ ethanol, then flush with
large quantity of deionized water, and dry the surface with high speed clean
air flow. We start the high speed camera right before the deposition of the
droplet and stop recording as the droplet fully spreads. The videos are
further processed with imageJ (Schneider et al. (2012)) and self-developed
Matlab codes, where the start of the measurement is set as the moment when a
droplet fully detaches from the syringe tip and forms a spherical cap at the
solid surface.

Microparticle image velocimetry ($\mu$PIV) is implemented to measure the flow
field near the bottom of an evaporating droplet. The measurements are based on
an inverted fluorescent microscope (Nikon Ts2-FL). Specifically, a green laser
(468nm) is utilized to excite the fluorescent micro tracers (Fluoro-Max; green
fluorescent polymer microspheres: excitation/emission 468/508 nm, diameter 1.0
$\mu$m). A CCD camera (STC-MCS510U3V) connected to the microlens captures the
fluorescent signal emitted by the tracer particles. We use butanol and IPA as
the test fluids, which have similar surface tension but different volatility
(Table 1). Experiments are conducted on slide glass for microscope with
thickness of $1.1\pm 0.2$ mm. The environmental condition is controlled as
25°C and 50$\%$ RH during the experiments. The droplet size for microPIV
experiments is set as $0.5\pm 0.002\mu L$, and is focused with 20x lens from
the bottom. Each set of experiments is repeated for more than five times. The
representative videos are presented as in Supplementary Movie 3, where the
video of Butanol is speeded up by 20 times and the video of IPA is slowed down
by 3 times with ImageJ.

##  4 Results and Discussion

###  4.1 Transition of flow states – physical mechanisms

Fig. 2 shows the evolution of flow field in a representative case, which can
be classified into three stages. At the initial stage ($\scriptsize{1}⃝$) when
a liquid droplet contacts a solid surface, the capillary effect drives a
continuous outward flow, leading to rapid droplet spreading. In the second
stage ($\scriptsize{2}⃝$), the evaporation cooling effect and the nonuniform
heat dissipation results in temperature gradient across the droplet surface,
and induces thermal Marangoni flow from the droplet edge toward the apex. The
opposing thermal Marangoni and capillary forces lead to the formation of two
stagnation points (SP) at the droplet surface. The upper stagnation point
moves towards the droplet apex as the Marangoni flow gets stronger and
disappears with a recirculation vortex that engulfs the whole droplet. The
lower stagnation point remains near the contact line ($\scriptsize{3}⃝$) as a
result of capillary-Marangoni interaction (Xu & Luo (2007), Ristenpart et al.
(2007), Wang & Harris (2018)).

Figure 2: Evolution of flow state for an evaporative droplet on a thermal
conductive substrate. The three stages are classified according to the state
of internal flow, where continuous outward capillary flow dominates droplet
spreading at the initial stage, the Marangoni flow develops at the second
stage, and balances with the capillary flow at the third stage. The key
dimensionless numbers in this case are set as $\rm Ja=5\times 10^{-3}$, $\rm
Bi=10^{-4}$, $\rm Ma=10^{-2}$.

Yet two questions arise as our analysis goes on. One is ‘does the stagnation
point universally exist?’ Another is ‘does thermal Marangoni flow always play
a role in wetting dynamics?’ These two questions are intricately related as
the relative strength of thermal Marangoni flow and capillary flow decides the
flow state and the formation of the stagnation point, and therefore the
wetting dynamics.

To elucidate this, we conducted detailed parametric analyses on all possible
combinations of liquids and solids. Two dominant dimensionless numbers are
extracted which govern the interacting physics, namely, the Jakob number, Ja,
which evaluates the balance between interfacial phase change and heat transfer
into the liquid phase, and the liquid-solid Biot number, Bi, which indicates
the ratio of liquid phase to solid phase thermal resistance.

$\mathrm{Ja}=\frac{\hat{k}_{0}\Delta\hat{T}}{\hat{H}_{0}\hat{L}\hat{p}_{\text{v,sat}}}\sqrt{\frac{2\pi\hat{R}_{\text{g}}\hat{T}_{\text{g}}}{\hat{M}}},\
\mathrm{Bi}=\frac{\hat{H}_{\text{solid}}\hat{k}_{0}}{\hat{H}_{0}\hat{k}_{\text{solid}}},$
(32)

where $\hat{k}_{0}$ is the thermal conductivity of the liquid, $\Delta T$ is
the maximum temperature difference across the droplet, $\hat{H}_{0}$ is the
initial droplet height, $\hat{L}$ is the latent heat of the liquid, $\hat{M}$
is the molar mass of the liquid, $\hat{p}_{\mathrm{v,sat}}$ is the saturation
vapor pressure of the liquid, $\hat{R}_{\mathrm{g}}$ is the gas constant,
$\hat{T}_{\mathrm{g}}$ is the gas phase temperature,
$\hat{H}_{\mathrm{solid}}$ is the substrate thickness, and
$\hat{k}_{\mathrm{solid}}$ is the substrate thermal conductivity.

With decreasing Bi, heat dissipation into the substrate becomes efficient. In
such cases, thermal resistance in the liquid phase dominates, thus the
temperature gradient across the droplet surface is large from the droplet apex
to the edge (large difference in thermal resistance across the liquid). The
large temperature gradient strengthens the Marangoni flow and retards droplet
spreading.

The Jakob number, mostly related to liquid volatility, decides the strength of
evaporation cooling. For moderately volatile liquids, the evaporation mass
flux takes heat away from the droplet surface, while the heat supply from the
substrate is also considerable. Near the contact line, the liquid film is
thin, enabling good heat supply and therefore a higher temperature in
comparison to the droplet apex; this is the case for many volatile liquids
such as water and alcohol. At very high liquid volatility (small Ja, _e.g._
flash evaporation), heat supply from the substrate can no longer meet the heat
loss from the droplet surface, as a result, the droplet surface gets
‘uniformly cold’, leading to weak thermal Marangoni flow. For liquids that
hardly evaporate (low volatility and large Ja), the effect of evaporation
cooling is negligible, and the process is near isothermal, which also leads to
weak Marangoni effect.

Taking an overview of the dominating mechanisms, we can see that Ja
(interfacial thermal effect) and Bi (heat conduction into the substrate)
determine the temperature field in the droplet-solid system. Once the
temperature gradient establishes, the Marangoni number determines the strength
of the Marangoni flow, defined as $\rm
Ma=\frac{\hat{\mu}_{0}\hat{\emph{u}}^{*}}{\varepsilon\hat{\sigma}_{0}}=\frac{\hat{\eta}_{\sigma}\Delta\hat{\emph{T}}}{\hat{\sigma}_{0}}$,
where the characteristic velocity $\hat{u}^{*}$ is defined as
$\hat{u}^{*}=\frac{\varepsilon\hat{\eta}_{\sigma}\Delta\hat{T}}{\hat{\mu}_{0}}$,
$\hat{\mu}_{0}$ denotes the liquid viscosity, $\hat{\sigma}_{0}$ is the
liquid-gas surface tension at standard conditions, $\varepsilon$ is the aspect
ratio of the droplet, $\hat{\eta}_{\sigma}$ is the coefficient of surface
tension to temperature, and $\Delta\hat{T}$ is the maximum temperature
difference across the droplet. The Marangoni flow interacts with the capillary
flow, which together affect the flow field inside the droplet, and further
determine the droplet dynamics along with the nonuniform distribution of
evaporation mass flux, as indicated by Fig. 3.

Figure 3: Decomposition of the dominating mechanisms in wetting dynamics of
evaporative droplets.

###  4.2 Conditions for flow reversal

To quantify the role of different effects on droplet dynamics, we decompose
the dimensionless interfacial flow velocity $u_{\mathrm{s}}$ into thermal
Marangoni velocity, $u_{\mathrm{tg}}$, resulting from interfacial temperature
gradient, and the capillary velocity, $u_{\mathrm{ca}}$, related to the local
surface curvature.

$\begin{split}u_{\text{s}}=-\frac{h^{2}}{2\mu}\frac{\partial p}{\partial
r}-\frac{h}{\mu}\frac{\partial T_{\text{s}}}{\partial r},\
u_{\text{tg}}=\frac{h}{\mu}\frac{\partial T_{s}}{\partial r},\
u_{\text{ca}}=-\frac{h^{2}}{2\mu}\frac{\partial p}{\partial r},\end{split}$
(33)

where _h_ is the local height of the liquid film, $\mu$ is the liquid
viscosity, _p_ is the pressure at the liquid phase, _r_ is the coordinate in
the radial direction, $T_{\mathrm{s}}$ is the interfacial temperature.

Fig. 4 maps the relative strength of thermal Marangoni flow and capillary
flow, $|u_{\mathrm{tg,max}}|/|u_{\mathrm{ca,max}}|$, for a wide range of
combinations of Ja and Bi, along with decomposed interfacial velocity
(a)$\sim$(e). With decreasing Bi from (d) to (b) to (e), the thermal Marangoni
flow enhances and outweighs the capillary effect. As a result, the direction
of interfacial flow reverses at a position close to the contact line, forming
a stagnation point as marked in Fig. 4(b) and (e).

Figure 4: Phase diagram for the relative strength of thermal Marangoni flow to
capillary flow. Figures (a) to (e) show the decomposition of surface velocity
$u_{\mathrm{s}}$ to thermal Marangoni velocity $u_{\mathrm{tg}}$ and capillary
velocity $u_{\mathrm{ca}}$, corresponding to letters (a) to (e) in the
colormap. Sketches I, II, III indicate the dominating mechanism at different
regimes of the phase diagram. The experimental data by Xu et al. (2010) (1a:
water on glass of thickness 1 mm, 1b: water on glass of thickness 0.15 mm, 1c:
IPA on glass of thickness 1 mm), Ristenpart et al. (2007) (2 $\mu$L droplets
of 2a: Methanol, 2b: Ethanol, 2c: IPA, 2d: Chloroform, on PDMS of thickness
4mm) and Li et al. (2015) (3a to 3b with leftward dotted arrows and
representative deposition patterns: 2.5 $\mu$L droplets on glass substrate of
30, 50, 65, and 80 ${}^{\circ}\rm C$) are marked out, which indicate the flow
reversal inside the droplet with changing environmental settings. Note that
the dotted line presents the criterion for the formation of stagnation point,
and does not necessarily correspond to the contour line of
$|u_{\mathrm{tg,max}}|/|u_{\mathrm{ca,max}}|=1$.

From (a) to (b) to (c), the value of Ja increases, corresponding to decreasing
liquid volatility. At a small value of Ja, severe liquid-vapor phase change
takes place, preferentially at the contact line region, leading to the
concentrated interplay of Marangoni and capillary effects near the contact
line (Fig. 4(a)). The severe interfacial phase change also leads to uniformly
cold droplet surface, corresponding to weak Marangoni effect. At large value
of Ja (Fig. 4(c)), the whole process is close to isothermal due to weak
evaporation. In such cases, the thermal Marangoni flow is also weak, and the
droplet spreading is dominated by capillary flow, close to the case described
by Tanner’s law.

Overall, the thermal Marangoni effect becomes apparent and even outweighs the
capillary effect at limited combinations of Ja and Bi. In regime II circled by
dot lines in Fig. 4, the thermal Marangoni effect is strong enough to reverse
the flow direction and form a stagnation point at the droplet surface. The
droplet dynamics within this regime is therefore a joint result of the thermal
Marangoni effect, evaporation effect, and capillary effect, in contrast to the
high-Ja regime where capillary effect dominates (Regime I) and the low-Ja
regime where liquid depletion due to evaporation near the contact line plays a
significant role (Regime III).

Previous studies revealed the circulation reversal in evaporating droplets
with varying substrate conductivity. We recalculated the experimental data in
these pioneer work (Ristenpart et al. (2007), Xu et al. (2010)), which are
demonstrated in the phase diagram and show good correspondence with our
theoretical predictions. Due to technical limitations, liquids in those
studies are limited to water and alcohols where the evaporation is moderate
with sufficient time for observation, and the substrates are mostly glass for
microPIV visualization. With the criterion proposed in Fig. 4, we are able to
quickly find out the flow state for liquids with varying thermal conductivity
and volatility, and on substrates with varying thickness and thermal
properties, which could not have been possible by experimental observations.

Corresponding to the flow state inside the droplet, the deposition patterns
from drying droplets differ, which provides a possible way for controllable
pattern formation in applications including printing, coating and thin film
fabrication. As indicated in Fig. 4, the deposition patterns transit from
coffee ring (Deegan et al. (1997)) to coffee eye with increasing liquid
volatility (Hu & Larson (2006)) and increasing substrate temperature (Li et
al. (2015)) where the thermal Marangoni flow strengthens, reverses the flow
direction, and leads to the formation of coffee eye at the center as indicated
by diamond marks 3a to 3d in Fig. 4. Further visualizations with microPIV also
demonstrate the flow reversal as the liquid changes from low volatility
Butanol to high volatility IPA on slide glass (Fig. 5, Supplementary Movie 3).
In correspondence to the flow direction, the deposition transits from coffee-
ring to more uniform patterns as indicated by the deposited fluorescent
tracing particles suspended in (1) Butanol and (2) IPA droplets in Fig. 4.

Figure 5: Flow pattern near three-phase contact line visualized by MicroPIV
where the tracing tracks of fluorescent particles are overlap of 100
continuous frames: (a) a 0.5$\mu L$ Butanol droplet on slide glass, and (b) a
0.5$\mu L$ IPA droplet on slide glass.

###  4.3 Spreading law of volatile droplets

The interfacial phase change and interacting flows further decide the
spreading dynamics of the droplet. Specifically, the speed of contact line can
be quantified with a power law in the form of $R(t)=t^{n}$, where n is the
spreading exponent. We then check the flow state and evaluate the
corresponding spreading rate for a wide range of Ja and Bi. A phase diagram is
thereafter derived with representative cases shown in Fig. 6. Overall, the
spreading exponent decreases from 1/10 (Tanner’s law) to smaller values (1/11,
1/12, … 1/25) with decreasing Ja and decreasing Bi. As Bi decreases, heat
transport across the substrate gets efficient and the thermal resistance in
the liquid phase becomes dominant. This strengthens both the evaporation mass
flux and the thermal Marangoni flow (enlarged temperature gradient).

With decreasing Ja, the evaporation effect enhances, while the thermal
Marangoni flow strengthens first (from isothermal state to large temperature
gradient) and then weakens (as the droplet surface gets uniformly cold due to
strong evaporation cooling effect). Overall, the influence of substrate
thermal properties on the spreading rate is weak at small liquid volatility
and gets significant at high liquid volatility as the evaporation cooling
effect and the thermal transport process becomes eminent. These conclusions
allow for optimal design of the thermal properties of working fluids and
substrates, where a balance between liquid wetting and thermal efficiency is
needed for specific application scenarios, _e.g._ , in avoiding transition
from dropwise to filmwise condensation, or in avoiding dry-out in cooling of
high-power-density electronics.

Figure 6: Phase diagram for the spreading rate of droplets along with
demonstrations of flow and temperature fields in representative conditions.
Contour lines of the spreading exponent are marked out with corresponding
values. Experimental results by Shiri et al. (2021) ((1) Alkanes on glass
substrates with decreasing thickness) and Kumar et al. (2022) ((2) HFE 7500 to
7100 with increasing volatility) are presented with yellow dotted arrows.
Experimental results in our research are marked with stars numbered as
$\scriptsize{1}⃝$ Butanol on 1.1 mm glass, $\scriptsize{2}⃝$ Butanol on 1.0 mm
copper, $\scriptsize{3}⃝$ IPA on 1.1mm glass, $\scriptsize{4}⃝$ IPA on 1.0mm
SUS430, $\scriptsize{5}⃝$ IPA on 1.0mm copper, $\scriptsize{6}⃝$ FC-72 on
1.1mm glass, $\scriptsize{7}⃝$ FC-72 on 1.0mm SUS430, $\scriptsize{8}⃝$ FC-72
on 1.0mm copper ($\scriptsize{1}⃝\sim\scriptsize{8}⃝$ are cases with substrate
temperature $T_{\rm sub}$= 20${}^{\circ}\rm C$), $\scriptsize{9}⃝$ FC-72 on
1.0mm copper ($T_{\rm sub}$ = 25${}^{\circ}\rm C$), $\scriptsize{10}⃝$ FC-72
on 1.0mm copper ($T_{\rm sub}$= 35${}^{\circ}\rm C$) and $\scriptsize{11}⃝$
FC-72 on 1.0mm copper ($T_{\rm sub}$ = 45${}^{\circ}\rm C$).

To further validate our theoretical findings, a series of experiments are
conducted with selected liquids and solid substrates. To ensure uniform
spreading and minimize hysteresis, we prepare substrates with mirror finishing
(roughness characterized as $\rm Sa_{glass}=0.258\ \mu m$, $\rm
Sa_{SUS430}=0.134\ \mu m$, $\rm Sa_{copper}=0.291\ \mu m$), and the substrate
temperature is controlled as $20\pm 0.2\ ^{\circ}\rm C$ with a Peltier module.
Experiments are conducted with liquids of different magnitudes of volatility,
_i.e._ 1-Butanol ($p_{\rm sat,20^{\circ}C}=580\ \rm Pa$), 2-Propanol (IPA,
$p_{\rm sat,20^{\circ}C}=4420\ \rm Pa$) and FluorinertR (FC-72, $p_{\rm
sat,25^{\circ}C}=30900\ \rm Pa$) and on substrates with stepwise varying
thermal conductivity and thickness, _i.e._ glass ($k_{\rm glass}=0.8\ \rm
W/m/K$, $h_{\rm glass}=1.1\ \rm mm$), SUS430 ($k_{\rm SUS430}=25\ \rm W/m/K$,
$h_{\rm SUS430}=1.0\ \rm mm$), copper ($k_{\rm copper}=398\ \rm W/m/K$,
$h_{\rm copper}=0.1,\ 1.0,\ 3.0\ mm$). The results are marked with numbered
stars in Fig. 6, which indicate robust decreasing trend of spreading rate with
increasing liquid volatility
($\scriptsize{1}⃝$→$\scriptsize{3}⃝$→$\scriptsize{6}⃝$,
$\scriptsize{2}⃝$→$\scriptsize{5}⃝$→$\scriptsize{8}⃝$), increasing substrate
conductivity ($\scriptsize{3}⃝$→$\scriptsize{4}⃝$→$\scriptsize{5}⃝$,
$\scriptsize{6}⃝$→$\scriptsize{7}⃝$→$\scriptsize{8}⃝$), and decreasing
substrate thickness (from 3 mm copper to 0.3 mm copper).

We further apply the calculated dimensionless numbers to the proposed model,
and get the spreading component at the steady spreading stage. Very nice
quantitative correspondence is found between experimental and theoretical
values as summarized in Table 3. Here, the experimental liquids and substrates
cover most of the cases that span from non-volatile to most volatile liquids
in practical applications and on substrates with a full range of thermal
conductivity, thus allowing for reliable prediction of the liquid behaviors by
quickly locating its Ja-Bi coordinate in the phase diagram.

No. | Exp. | Num. | No. | Exp. | Num. | No. | Exp. | Num.
---|---|---|---|---|---|---|---|---
$\scriptsize{1}⃝$ | 0.0995±0.0028 | 0.103 | $\scriptsize{2}⃝$ | 0.092±0.0028 | 0.0936 | $\scriptsize{3}⃝$ | 0.0894±0.0026 | 0.0910
$\scriptsize{4}⃝$ | 0.0869±0.0026 | 0.0813 | $\scriptsize{5}⃝$ | 0.0803±0.0053 | 0.0788 | $\scriptsize{6}⃝$ | 0.0634±0.0024 | 0.0752
$\scriptsize{7}⃝$ | 0.0601±0.0011 | 0.0578 | $\scriptsize{8}⃝$ | 0.0568±0.0052 | 0.0504 | $\scriptsize{9}⃝$ | 0.0481±0.0036 | 0.0490
$\scriptsize{10}⃝$ | 0.0473±0.0030 | 0.0446 | $\scriptsize{11}⃝$ | 0.0429±0.0010 | 0.0435 |  |  |
Table 3:  Experimental (Exp.) spreading exponent _n_ in comparison to the
numerical (Num.) value by matching its position (Ja, Bi) in the phase diagram.
Numbers $\scriptsize{1}⃝\sim\scriptsize{11}⃝$ correspond with the cases marked
in 6. Deviations between experimental and numerical results exist mainly in
the high-volatility and low-substrate-conductivity region (small Ja and large
Bi, _e.g._ , case $\scriptsize{6}⃝$), where the heat dissipation into the gas
phase is no longer negligible.

Besides the experiments in present research, our conclusions on the effect of
substrate thermal properties correspond with the work of Shiri et al. (2021)
where alkane droplets indicate increasing contact angle and retarded spreading
as the thickness of the glass substrate decreases stepwise from 6 mm to 0.063
mm. The influence of liquid volatility, we recently found, also meet with a
very recent work by Kumar et al. (2022) (APS DFD abstract) where the spreading
rate of HFE-7500 to 7100 droplets decreases with decreasing chain length, that
is, increasing liquid volatility.

Looking at the overall life span of a deposited droplet, the experiments
additionally indicate shortening spreading time of an evaporating droplet with
increasing substrate conductivity, increasing liquid volatility, and
increasing substrate temperature. For highly volatile liquids (lower left
corner in the phase diagram), the spreading stage even disappears
(indistinguishable) as the strong evaporation effect counteracts the capillary
spreading, in which cases the contact line recedes right after the droplet
deposition. This is exactly the case according to our observation of FC-72
droplets on heated substrates, where spreading slows down and the spreading
period shortens as the substrate temperature increases stepwise from
$25^{\circ}\rm C$ to $45^{\circ}\rm C$ as exampled by
$\scriptsize{9}⃝$→$\scriptsize{10}⃝$→$\scriptsize{11}⃝$ in Fig. 6. When the
substrate temperature exceeds $45^{\circ}\rm C$ ($T_{\rm
FC-72,boiling}=56^{\circ}\rm C$ ), the spreading stage becomes
indistinguishable, _i.e._ the contact line starts to recede right after
deposition, and the droplet ($\rm 0.6\ \mu L$) fully evaporates within 3
seconds, which represents a general scenario in flash evaporation.

##  5 Conclusions and Perspectives

In this research, we revisit the classical problem of droplet wetting with
interfacial phase change, which is of fundamental importance to a wide range
of industrial processes that involve three phases. The objective originates
from the urgent need for a unified explanation to a couple of basic scientific
problems including the flow transition for different liquid-solid
combinations, the existence of stagnation point near the three phase contact
line, the relative strength of capillary and Marangoni flows, the deposition
patterns from colloidal suspensions, as well as the corresponding spreading
dynamics in diverse scenarios. Due to the requirement for substrate
transparency and the available frame rate of existing visualization
techniques, explorations on these scientific issues are limited to substrates
like glass, and common liquids like water and alcohol, while important liquids
and substrates in the vast industrial fields are less explored, such as
refrigerants and metallic substrates. This requires a close collaboration
between experimental and numerical efforts, where information that cannot be
revealed from one side can be disclosed by the other.

Here, we develop a lubrication type model that considers the non-isothermal
heat conduction into the solid substrate. We further conducted experiments on
the flow field and the spreading rate for a wide range of combinations of
liquids and solids. With one-to-one comparison between the numerical and
experimental results, we elucidate the key dimensionless numbers that come
into effect, _i.e._ , the Jacok number, Ja, the liquid-solid Biot number, Bi,
and the Marangoni number, Ma. A unified criteria is proposed for the
transition of flow state, and provide unbiased explanations to a series of
individual works on deposition patterns from drying droplets. A phase diagram
for spreading dynamics is further established along with experimental
demonstrations, which generalizes Tanner’s law (for non-volatile liquids in an
isotherm state) to a full range of liquids with saturation vapor pressure
spanning from $10^{1}$ to $10^{4}$ Pa and on substrates with all range of
thermal conductivity spanning from $10^{-1}$ to $10^{3}$ W/m/K.

In the theoretical aspect, this work illustrates the relative strength of
interacting physical mechanisms as Regime I, II, III in Fig. 4, where the
capillary effect dominates at low liquid volatility in Regime I, the Marangoni
effect, capillary effect and evaporation effect interact in Regime II, and the
strong evaporation effect (preferential mass loss near contact line) dominates
in Regime III at very high liquid volatility. Such decomposition of the
mechanisms allows for quick location of the influencing factors, and enables
efficient control of liquid behaviors in various phase change energy devices,
printing-coating processes, and microfluidic systems.

Overall, Tanner’s law describes the spreading of non-volatile Newtonian fluids
on smooth solid substrates. Here ’non-volatile’, ’Newtonian’ and ’smooth’ are
necessary conditions for the universality of the spreading exponent. With
interfacial phase change (volatile liquids), the process becomes non-
isothermal in which case the non-uniform distribution of mass flux and the
thermal Marangoni effect come into play along with the capillary and viscous
effect, as quantified in the present work. For non-Newtonian fluids, the
viscosity largely depends on the shear stress. Theories for the spreading
dynamics of dilatant fluids (shear thickening) and pseudo-plastic fluids
(shear thinning) were established in previous studies based on the lubrication
approximation (Rafaï et al. (2004), Wang et al. (2007a)) and experiments (Wang
et al. (2007b)). For topological surfaces with micro structures, the spreading
dynamics can be more complex. A systematic understanding is yet reached due to
the diverse surface topology and the complex liquid-solid interaction, where
the small-scale surface structure affects the dynamics of the contact line and
thereafter the overall droplet behaviors. By taking advantage of the
asymmetric effect of surface topology, it is also possible to realize
directional fluid motion, which indeed has been widely exploited in recent
years with promising applications in directional liquid transport (Li et al.
(2017), Li et al. (2020)), water harvesting (Zheng et al. (2010), Dai et al.
(2018)), and dry-out prevention (Li et al. (2021)).

The authors gratefully acknowledge the support received from Japanese Society
for the Promotion of Science (JSPS), and ThermaSMART project of European
Commission (Grant No. EC-H2020-RISE-ThermaSMART-778104). ZW and CI acknowledge
the support from JSPS KAKENHI (Grant Nos. JP21K14097, JP21H01251, JP22K18771).

Declaration of Interests The authors report no conflict of interest.

##  Appendix A Detailed Experimental Results

Fig. 7 and Fig. 8 present the variation of droplet contact radius along with
time for representative combinations of liquids and solids (logarithmic
coordinates). A summary and comparison of the spreading rate through
experimental and numerical approaches are given in Table 3. Taking an overview
of all the variations curves, we can see that the effect of liquid volatility
is much stronger than the substrate conductivity or thickness on the spreading
dynamics of a droplet.

For IPA and butanol droplets, spreading starts right after the droplet gets
contact with the solid surface and lasts for $10^{0}\sim 10^{1}$ seconds
before the contact line slows down and recedes. By changing the substrate
material from glass to SUS430 and to copper as in Fig. 7, the slope of the
variation curve decreases, indicating a decreasing trend of spreading rate
with increasing substrate conductivity ($\scriptsize{1}⃝$→$\scriptsize{2}⃝$,
$\scriptsize{3}⃝$→$\scriptsize{4}⃝$→$\scriptsize{5}⃝$).

Figure 7: Spreading behaviors of droplets at substrates with different thermal
conductivity ($T_{\mathrm{sub}}$ = 20°C). $\scriptsize{1}⃝$: IPA on 1.1 mm
glass, $\scriptsize{2}⃝$: IPA on 1.0 mm SUS430, $\scriptsize{3}⃝$: IPA on 1.0
mm copper, $\scriptsize{4}⃝$: Butanol on 1.1 mm glass, $\scriptsize{5}⃝$:
Butanol on 1 mm copper. These cases correspond with the mark numbers in Fig.
6.

Fig. 8 indicates the effect of liquid volatility on droplet spreading. For
FC-72 droplets, a slight fluctuation (spread and contract) takes place at the
very initial moment ($<10\ ms$) when a droplet contacts the solid surface as a
result of the inertia effect (even through the distance between the syringe
tip and the solid surface is controlled as small as possible, the inertia
effect can cause a slight fluctuation when the liquid surface tension is
small, which can be evaluated by the Weber number). The droplet
($\scriptsize{8}⃝$ in Fig. 8) then spreads, slows down and recedes after
several hundred milliseconds. By further increasing the substrate temperature,
the period of spreading shortens along with decreasing spreading rate as
indicated by
$\scriptsize{8}⃝$→$\scriptsize{9}⃝$→$\scriptsize{10}⃝$→$\scriptsize{11}⃝$ in
Fig. 8. For substrate temperature higher than 45 °C, the spreading stage is no
longer distinguishable, in which case the evaporation effect is strong enough
to counteract the capillary spreading. As a result, no apparent spreading can
be observed, _i.e._ , the contact line recedes right after the droplet
deposition due to significant mass loss at the contact line region.

Figure 8: Spreading behaviors of droplets with different volatility on 1.0 mm
copper substrate. $\scriptsize{2}⃝$: Butanol at $T_{\mathrm{sub}}$ = 20°C,
$\scriptsize{5}⃝$: IPA at $T_{\mathrm{sub}}$ = 20°C,
$\scriptsize{8}⃝\sim\scriptsize{10}⃝$: FC-72 at $T_{\mathrm{sub}}$ = 20°C,
25°C, 35°C, 45°C. These cases correspond with the mark numbers in Fig. 6.

In the mathematical model, we formulate the expression of interfacial mass
flux by combining the Hertz-Knudsen equation, the chemical balance relation
across the liquid-air interface, and the ideal gas assumption. The Hertz-
Knudsen equation is often utilized to describe the evaporation rates based on
the statistics of vapor molecules adsorbed onto the liquid surface. The
evaporation and condensation coefficients contained in the equation are thus
important while they have inexplicably been found to span 3 orders of
magnitude when comparing with the actual experimental data as reported in a
number of studies (Persad & Ward (2016), Young (1991), Fang & Ward (1999)). In
the formulation, we set the evaporation and condensation coefficients as 1
with the assumption that the phase change process is near equilibrium. A
modified Ja number, Ja*, is utilized for numerical calculations to compensate
for the overprediction of the Hertz-Knudsen equation, where liquids with
higher volatility (higher molecular activity) tend to be more overpredicted.
Specifically, Ja* is calculated by magnifying Ja for
$500\frac{\mathrm{log}_{\mathrm{max}}-\mathrm{log}Ja}{\mathrm{log}_{\mathrm{max}}-\mathrm{log}_{\mathrm{min}}}$
times, where $\mathrm{log}_{\mathrm{max}}$ is the maximal order of magnitude
of Ja for all possible types of liquids and $\mathrm{log}_{\mathrm{min}}$ is
the minimal order of magnitude ($\mathrm{log}_{\mathrm{max}}=0$,
$\mathrm{log}_{\mathrm{min}}=-4$ in our present work). The calculated values
of Ja* and Bi for cases $\scriptsize{1}⃝$ to $\scriptsize{11}⃝$ are summarized
in Table 4, which correspond with the data points marked in phase diagram Fig.
6.

Case | $\scriptsize{1}⃝$ | $\scriptsize{2}⃝$ | $\scriptsize{3}⃝$ | $\scriptsize{4}⃝$ | $\scriptsize{5}⃝$ | $\scriptsize{6}⃝$ | $\scriptsize{7}⃝$ | $\scriptsize{8}⃝$ | $\scriptsize{9}⃝$ | $\scriptsize{10}⃝$ | $\scriptsize{11}⃝$
---|---|---|---|---|---|---|---|---|---|---|---
Bi | 2.8e-1 | 8.1e-4 | 2.5e-1 | 1.1e-2 | 7.0e-4 | 1.7e-1 | 7.6e-3 | 4.9e-4 | 4.9e-4 | 4.9e-4 | 4.9e-4
Ja* | 7.0e-1 | 7.0e-1 | 1.1e-1 | 1.1e-1 | 1.1e-1 | 1.9e-2 | 1.9e-2 | 1.9e-2 | 1.6e-2 | 1.1e-2 | 7.9e-3
Table 4: Values of Ja* and Bi for cases $\scriptsize{1}⃝$ to
$\scriptsize{11}⃝$ (expressed in scientific notation).

Additional PIV experiments indicate the transition of flow field as the liquid
volatility increases. For low-volatility Butanol, continuous outward flow is
observed driven by the capillary force, which forms a ‘coffee ring’ as the
droplet dries out (Fig. 5(a), Supplementary Movie 3, left). As the liquid
volatility increases (IPA), the PIV particles tend to change its direction at
a position near the three-phase contact line (stagnation point) as shown by
the tracing track of fluorescent particles in Fig. 5(b). The direction change
avoids particle accumulation at the periphery, which leads to more uniform
deposition patterns as the droplet fully dries out. The transition of flow
direction and the formation of deposition patterns correspond with our
numerical predictions by locating the Ja-Bi positions for (1) Butanol-slide
glass and (2) IPA-slide glass in Fig. 4.

## References

  * Ajaev (2005) Ajaev, Vladimir S 2005 Spreading of thin volatile liquid droplets on uniformly heated surfaces. Journal of Fluid Mechanics 528, 279–296.
  * Anderson & Davis (1995) Anderson, DM & Davis, SH 1995 The spreading of volatile liquid droplets on heated surfaces. Physics of Fluids 7 (2), 248–265.
  * Atkins et al. (2014) Atkins, Peter, Atkins, Peter William & de Paula, Julio 2014 Atkins’ physical chemistry. Oxford university press.
  * Berteloot et al. (2008) Berteloot, Guillaume, Pham, C-T, Daerr, A, Lequeux, François & Limat, Laurent 2008 Evaporation-induced flow near a contact line: Consequences on coating and contact angle. Europhysics Letters 83 (1), 14003.
  * Blake (2006) Blake, Terence D 2006 The physics of moving wetting lines. Journal of Colloid and Interface Science 299 (1), 1–13.
  * Bonn et al. (2009) Bonn, Daniel, Eggers, Jens, Indekeu, Joseph, Meunier, Jacques & Rolley, Etienne 2009 Wetting and spreading. Reviews of Modern Physics 81 (2), 739.
  * Breitenbach et al. (2018) Breitenbach, Jan, Roisman, Ilia V & Tropea, Cameron 2018 From drop impact physics to spray cooling models: a critical review. Experiments in Fluids 59, 1–21.
  * Carlson (2018) Carlson, Andreas 2018 Fluctuation assisted spreading of a fluid filled elastic blister. Journal of Fluid Mechanics 846, 1076–1087.
  * Carré et al. (1996) Carré, Alain, Gastel, Jean-Claude & Shanahan, Martin ER 1996 Viscoelastic effects in the spreading of liquids. Nature 379 (6564), 432–434.
  * Cazabat & Stuart (1986) Cazabat, AM & Stuart, MA Cohen 1986 Dynamics of wetting: effects of surface roughness. The Journal of Physical Chemistry 90 (22), 5845–5849.
  * Charitatos & Kumar (2020) Charitatos, Vasileios & Kumar, Satish 2020 A thin-film model for droplet spreading on soft solid substrates. Soft Matter 16 (35), 8284–8298.
  * Cormier et al. (2012) Cormier, Sara L, McGraw, Joshua D, Salez, Thomas, Raphaël, Elie & Dalnoki-Veress, Kari 2012 Beyond tanner’s law: Crossover between spreading regimes of a viscous droplet on an identical film. Physical Review Letters 109 (15), 154501\.
  * Dai et al. (2016) Dai, Qingwen, Khonsari, MM, Shen, Cong, Huang, Wei & Wang, Xiaolei 2016 Thermocapillary migration of liquid droplets induced by a unidirectional thermal gradient. Langmuir 32 (30), 7485–7492.
  * Dai et al. (2018) Dai, Xianming, Sun, Nan, Nielsen, Steven O, Stogin, Birgitt Boschitsch, Wang, Jing, Yang, Shikuan & Wong, Tak-Sing 2018 Hydrophilic directional slippery rough surfaces for water harvesting. Science advances 4 (3), eaaq0919.
  * De Gennes et al. (1990) De Gennes, PG, Hua, X & Levinson, P 1990 Dynamics of wetting: local contact angles. Journal of fluid mechanics 212, 55–63.
  * De Gennes (1985) De Gennes, Pierre-Gilles 1985 Wetting: statics and dynamics. Reviews of Modern Physics 57 (3), 827.
  * Deegan et al. (1997) Deegan, Robert D, Bakajin, Olgica, Dupont, Todd F, Huber, Greb, Nagel, Sidney R & Witten, Thomas A 1997 Capillary flow as the cause of ring stains from dried liquid drops. Nature 389 (6653), 827–829.
  * Dunn et al. (2009) Dunn, GJ, Wilson, SK, Duffy, BR, David, S & Sefiane, Khellil 2009 The strong influence of substrate conductivity on droplet evaporation. Journal of Fluid Mechanics 623, 329–351.
  * Fang & Ward (1999) Fang, G & Ward, CA 1999 Temperature measured close to the interface of an evaporating liquid. Physical Review E 59 (1), 417.
  * Goossens et al. (2011) Goossens, Sarah, Seveno, David, Rioboo, Romain, Vaillant, Alexandre, Conti, J & De Coninck, Joël 2011 Can we predict the spreading of a two-liquid system from the spreading of the corresponding liquid–air systems? Langmuir 27 (16), 9866–9872.
  * Hoang & Kavehpour (2011) Hoang, A & Kavehpour, HP 2011 Dynamics of nanoscale precursor film near a moving contact line of spreading drops. Physical review letters 106 (25), 254501.
  * Hu & Larson (2005) Hu, Hua & Larson, Ronald G 2005 Analysis of the effects of marangoni stresses on the microflow in an evaporating sessile droplet. Langmuir 21 (9), 3972–3980.
  * Hu & Larson (2006) Hu, Hua & Larson, Ronald G 2006 Marangoni effect reverses coffee-ring depositions. The Journal of Physical Chemistry B 110 (14), 7090–7094.
  * Jalaal et al. (2021) Jalaal, Maziyar, Stoeber, Boris & Balmforth, Neil J 2021 Spreading of viscoplastic droplets. Journal of Fluid Mechanics 914, A21.
  * Jambon-Puillet et al. (2018) Jambon-Puillet, Etienne, Carrier, Odile, Shahidzadeh, Noushine, Brutin, David, Eggers, Jens & Bonn, Daniel 2018 Spreading dynamics and contact angle of completely wetting volatile drops. Journal of Fluid Mechanics 844, 817–830.
  * Karapetsas et al. (2013) Karapetsas, George, Sahu, Kirti Chandra & Matar, Omar K 2013 Effect of contact line dynamics on the thermocapillary motion of a droplet on an inclined plate. Langmuir 29 (28), 8892–8906.
  * Kavehpour et al. (2003) Kavehpour, H Pirouz, Ovryn, Ben & McKinley, Gareth H 2003 Microscopic and macroscopic structure of the precursor layer in spreading viscous drops. Physical review letters 91 (19), 196104\.
  * Kumar et al. (2022) Kumar, Manish, Parimalanathan, Senthil, Rednikov, Alexey & Colinet, Pierre 2022 Wetting dynamics of volatile hydrofluoroether (hfe) liquid droplets. Bulletin of the American Physical Society .
  * Lelah & Marmur (1981) Lelah, Michael D & Marmur, Abraham 1981 Spreading kinetics of drops on glass. Journal of Colloid and Interface Science 82 (2), 518–525.
  * Li et al. (2017) Li, Jiaqian, Zhou, Xiaofeng, Li, Jing, Che, Lufeng, Yao, Jun, McHale, Glen, Chaudhury, Manoj K & Wang, Zuankai 2017 Topological liquid diode. Science advances 3 (10), eaao3530.
  * Li et al. (2021) Li, Jiaqian, Zhou, Xiaofeng, Tao, Ran, Zheng, Huanxi & Wang, Zuankai 2021 Directional liquid transport from the cold region to the hot region on a topological surface. Langmuir 37 (16), 5059–5065.
  * Li et al. (2020) Li, Xing, Li, Jiaqian & Dong, Guangneng 2020 Bioinspired topological surface for directional oil lubrication. ACS applied materials & interfaces 12 (4), 5113–5119.
  * Li et al. (2015) Li, Yanshen, Lv, Cunjing, Li, Zhaohan, Quéré, David & Zheng, Quanshui 2015 From coffee rings to coffee eyes. Soft Matter 11 (23), 4669–4673.
  * Marmur (1983) Marmur, Abraham 1983 Equilibrium and spreading of liquids on solid surfaces. Advances in Colloid and Interface Science 19 (1-2), 75–102.
  * Mitra & Mitra (2016) Mitra, Surjyasish & Mitra, Sushanta K 2016 Understanding the early regime of drop spreading. Langmuir 32 (35), 8843–8848.
  * Moosman & Homsy (1980) Moosman, Steven & Homsy, GM 1980 Evaporating menisci of wetting fluids. Journal of Colloid and Interface Science 73 (1), 212–223.
  * Mudawar (2001) Mudawar, Issam 2001 Assessment of high-heat-flux thermal management schemes. IEEE Transactions on Components and Packaging Technologies 24 (2), 122–141.
  * Pang & Ó Náraigh (2022) Pang, Khang Ee & Ó Náraigh, Lennon 2022 A mathematical model and mesh-free numerical method for contact-line motion in lubrication theory. Environmental Fluid Mechanics 22 (2-3), 301–336.
  * Persad & Ward (2016) Persad, Aaron H & Ward, Charles A 2016 Expressions for the evaporation and condensation coefficients in the hertz-knudsen relation. Chemical reviews 116 (14), 7727–7767.
  * Plesset & Prosperetti (1976) Plesset, Milton S & Prosperetti, Andrea 1976 Flow of vapour in a liquid enclosure. Journal of Fluid Mechanics 78 (3), 433–444.
  * Rafaï et al. (2004) Rafaï, Salima, Bonn, Daniel & Boudaoud, Arezki 2004 Spreading of non-newtonian fluids on hydrophilic surfaces. Journal of Fluid Mechanics 513, 77–85.
  * Ristenpart et al. (2007) Ristenpart, WD, Kim, PG, Domingues, C, Wan, J & Stone, Howard A 2007 Influence of substrate conductivity on circulation reversal in evaporating drops. Physical Review Letters 99 (23), 234502.
  * Rose (2002) Rose, JW 2002 Dropwise condensation theory and experiment: a review. Proceedings of the Institution of Mechanical Engineers, Part A: Journal of Power and Energy 216 (2), 115–128.
  * Saiseau et al. (2022) Saiseau, Raphael, Pedersen, Christian, Benjana, Anwar, Carlson, Andreas, Delabre, Ulysse, Salez, Thomas & Delville, Jean-Pierre 2022 Near-critical spreading of droplets. Nature Communications 13 (1), 7442.
  * Schneider et al. (2012) Schneider, Caroline A, Rasband, Wayne S & Eliceiri, Kevin W 2012 Nih image to imagej: 25 years of image analysis. Nature methods 9 (7), 671–675.
  * Shiri et al. (2021) Shiri, Samira, Sinha, Shayandev, Baumgartner, Dieter A & Cira, Nate J 2021 Thermal marangoni flow impacts the shape of single component volatile droplets on thin, completely wetting substrates. Physical Review Letters 127 (2), 024502.
  * Tanner (1979) Tanner, LH 1979 The spreading of silicone oil drops on horizontal surfaces. Journal of Physics D: Applied Physics 12 (9), 1473.
  * Teh et al. (2008) Teh, Shia-Yen, Lin, Robert, Hung, Lung-Hsin & Lee, Abraham P 2008 Droplet microfluidics. Lab on a Chip 8 (2), 198–220.
  * Tsoumpas et al. (2015) Tsoumpas, Yannis, Dehaeck, Sam, Rednikov, Alexey & Colinet, Pierre 2015 Effect of marangoni flows on the shape of thin sessile droplets evaporating into air. Langmuir 31 (49), 13334–13340.
  * Tu et al. (2018) Tu, Yaodong, Wang, Ruzhu, Zhang, Yannan & Wang, Jiayun 2018 Progress and expectation of atmospheric water harvesting. Joule 2 (8), 1452–1475.
  * Voinov (1976) Voinov, OV 1976 Hydrodynamics of wetting. Fluid Dynamics 11 (5), 714–721.
  * Wang & Harris (2018) Wang, Lihui & Harris, Michael T 2018 Stagnation point of surface flow during drop evaporation. Langmuir 34 (20), 5918–5925.
  * Wang et al. (2007a) Wang, XD, Lee, DJ, Peng, XF & Lai, JY 2007a Spreading dynamics and dynamic contact angle of non-newtonian fluids. Langmuir 23 (15), 8042–8047.
  * Wang et al. (2007b) Wang, XD, Zhang, Y, Lee, DJ & Peng, XF 2007b Spreading of completely wetting or partially wetting power-law fluid on solid surface. Langmuir 23 (18), 9258–9262.
  * Wang et al. (2021) Wang, Zhenying, Karapetsas, George, Valluri, Prashant, Sefiane, Khellil, Williams, Adam & Takata, Yasuyuki 2021 Dynamics of hygroscopic aqueous solution droplets undergoing evaporation or vapour absorption. Journal of Fluid Mechanics 912, A2.
  * Wang et al. (2019) Wang, Zhenying, Orejon, Daniel, Sefiane, Khellil & Takata, Yasuyuki 2019 Coupled thermal transport and mass diffusion during vapor absorption into hygroscopic liquid desiccant droplets. International Journal of Heat and Mass Transfer 134, 1014–1023.
  * Wang et al. (2022) Wang, Zhenying, Orejon, Daniel, Takata, Yasuyuki & Sefiane, Khellil 2022 Wetting and evaporation of multicomponent droplets. Physics Reports 960, 1–37.
  * Xu & Luo (2007) Xu, Xuefeng & Luo, Jianbin 2007 Marangoni flow in an evaporating water droplet. Applied Physics Letters 91 (12), 124102\.
  * Xu et al. (2010) Xu, Xuefeng, Luo, Jianbin & Guo, Dan 2010 Criterion for reversal of thermal marangoni flow in drying drops. Langmuir 26 (3), 1918–1922.
  * Yang et al. (2022) Yang, Xiuyuan, Wu, MengMeng, Doi, Masao & Man, Xingkun 2022 Evaporation dynamics of sessile droplets: The intricate coupling of capillary, evaporation, and marangoni flow. Langmuir 38 (16), 4887–4893.
  * Young (1991) Young, JB 1991 The condensation and evaporation of liquid droplets in a pure vapour at arbitrary knudsen number. International journal of heat and mass transfer 34 (7), 1649–1661.
  * Zheng et al. (2010) Zheng, Yongmei, Bai, Hao, Huang, Zhongbing, Tian, Xuelin, Nie, Fu-Qiang, Zhao, Yong, Zhai, Jin & Jiang, Lei 2010 Directional water collection on wetted spider silk. Nature 463 (7281), 640–643.