1. Introduction
The drying of water droplets on window glass after it has rained is commonly observed in daily life. Droplets drying looks like a simple problem, but shows surprisingly complicated fundamental science, including evaporation, particle convection, contact line (CL) motion, etc. When the droplets are laden with nonvolatile solutes and evaporate on a solid surface, this phenomenon becomes even more complex because versatile deposition patterns of solutes are left on the substrate. To understand the dynamics of the drying of liquid droplets and to tailor or control the final deposition pattern has not only held a special interest in scientific research, but important applications in inkjet printing [1–3], microelectronic device manufacturing [4] and even medicinal diagnostics [5, 6]. As Larson pointed out in his review published in Nature that ‘In modern times, a world of physical chemistry can similarly be observed by watching a droplet dry' [7], the drying of liquid droplets has spawned thousands of publications, which lie at the crossroad of physics, chemistry, materials science and even biology. In this topical review, we will give brief descriptions of the scientific advances in the understanding of droplets drying. We will first review the versatility of the deposition patterns, mainly about ring-like patterns, of drying sessile droplets that were observed in experiments. Then, we show the recent development of the theoretical studies of the drying of liquid droplets based on the Onsager variational principle (OVP) theory.
1.1. Ring-like patterns in experiments
In 1997, Deegan et al reported that a spilt droplet of coffee leaves a circular ring, rather than a uniform spot after drying on a coffee mug [8]. This is known as the famous coffee-ring phenomenon. Inspired by this seminal work, other ring-like patterns have been subsequently reported, including mountain-like [9–11], volcano-like [12, 13] and even multi-ring patterns [14, 15]. Besides ring-like patterns, a large number of complicated deposition patterns have recently been found in experiments, such as eye-like [16], finger-like [17], spoke-like [18] and spiral patterns [19], as shown in figure 1. As a special topical review, we hereafter review the formation mechanism of the ring-like deposition patterns in terms of both experimental and theoretical aspects.
Figure 1. Schematic diagram of experimentally observed deposition patterns of drying droplets: (a) coffee-ring; (b) mountain-like; (c) volcano-like; (d) multi-ring; (e) eye-like; (f) spoke-like, which are the density distribution of solutes left on the substrate after droplet drying. |
The coffee-ring pattern, figure 1(a), is usually observed in our daily life, where most solutes are deposited along the initial perimeter of the droplet. Deegan [8, 20] first pointed out that the formation mechanism of coffee-ring is the combined effects of the pinning of the CL and the evaporation-induced fluid flow inside the droplets. When the CL is pinned, the direction of the height-averaged fluid flow inside the droplet is from the droplet center to the edge. Such convection flux causes the nonvolatile solutes to accumulate and then deposit in the area around the CL, forming the coffee-ring pattern. In recent years, the utilization and suppression of the formation of coffee-ring has been applied to many aspects of industrial production [1, 21] and for more details one can refer to another recent review [22].
On the other hand, a freely moving CL leads to a mountain-like deposition pattern, where most of the solutes deposit at the center of the droplet [figure 1(b)]. Pauchard and Allain [23] studied the drying of a water/polysaccharid Dextran droplet on a clean glass slide. They found that when the polymer concentration is around 0.20 g/ml and the relative humidity about ±5%, a mountain-like deposition pattern is left on the slide. Later, Willmer et al [9] studied aqueous poly(ethylene oxide) solution droplets drying on glass substrates. They also obtained the mountain-like deposition pattern and found that the radius of the deposition pattern continuously increased to the size of the droplet's initial radius by increasing the concentration of polymer. Besides the solute concentration, Li et al [11] found that the wetting characteristics of the droplet solution with respect to the substrate, i.e. the contact angle hysteresis (CAH), is another important factor in obtaining the mountain-like pattern. With a weak CAH substrate (such as a silica glass or polycarbonate substrate for the water droplet), the CL freely recedes and forms the mountain-like deposition pattern.
When the moving ability of the CL is in between the pinned and freely moving cases, volcano-like patterns are observed, of which the peak of the density distribution of the deposited solutes is located between the droplet center and its initial perimeter, as shown in figure 1(c). Kajiya et al [12] studied the drying of a water-poly(N, N-dimethylacrylamide) PDMA droplet on a glass substrate and observed a volcano-like deposition pattern. They showed that the location of the peak approaches the initial perimeter as the equilibrium contact angle of the droplet increases by using different glass slides as the substrates. They explained the volcano-like pattern based on the mobility of the CL. They showed that both the CL and solutes move from the droplet edge to the center. However, the CL moves faster than the solutes, forming the volcano-like deposition pattern.
When the CL has a stick-slip motion, the deposits form concentric rings called a multi-ring deposition pattern, as shown in figure 1(d). In these cases, the droplet CL is first pinned, generating a flow from the droplet center to its edge. This flow convects solutes to the boundary and deposits them near the CL. This process is similar to the formation of the coffee-ring. The difference is that as the droplet volume decreases, the contact angle decreases and creates an inward unbalanced force acting on the CL. When the contact angle becomes less than the receding contact angle, the CL starts to slip and move quickly towards the center accompanied by the increase of the contact angle. The CL keeps receding until the contact angle increases to its equilibrium value. The repetition of this stick-slip motion of the CL forms the multi-ring deposition pattern.
In 1995, even earlier than the observation of the coffee-ring, Adachi et al [24] observed the multi-ring pattern in the drying of water-polystyrene particle droplets on borosilicate glass plates. Later, Deegan et al [25] also found such interesting deposition patterns in the drying of the same droplets on mica, but with a different solute size of the polystyrene particle (0.1 μm). Figure 2(a) is the experimental results of the multi-ring pattern found in the evaporation of deionized water-DNA droplets on glass surfaces. An interesting phenomenon is that the multi-ring pattern is made of a solid circle surrounded by many concentric rings. Besides sessile droplets, confined evaporation systems can also generate the multi-ring deposition pattern. In figure 2(b), Xu et al [14] obtained these patterns when a solution evaporated on the sphere-on-flat geometry, where the solution is made by poly[2-methoxy-5-(2-ethylhexyloxy)-1,4-phenylenevinylene] (MEH-PPV) (MW = 50–300 kg/mol) and toluene.
Figure 2. Examples of multi-ring patterns through the evaporation of droplets. (a) Typical DNA stain patterns with multi-ring formation. Reprinted figure with permission from [15]. Copyright (2020) the American Physical Society. (b) MEH-PPV ring patterns in a sphere-on-flat configuration. Reprinted figure with permission from [14]. Copyright (2020) the American Physical Society. |
It is clear that the final deposition pattern is mainly determined by the motion of the CL. Many experimental conditions, including the evaporation rate, viscosity of the solution, wetting property of the substrate, etc, can affect this motion. Various combinations of these factors result in versatile deposition patterns. In fact, we are far from discovering all the possible droplet deposition patterns [26].
1.2. Conventional hydrodynamic theory of drying droplets
In this subsection, we briefly review the development of the conventional hydrodynamic theory of drying droplets. We will review the basic theoretical framework. Some simple reduced scale theories are also reviewed here, to give fundamental understanding of this complex problem.
The deposition pattern is the density distribution of nonvolatile solutes left on the substrate after droplets have dried. Figure 3(a) describes the process of the deposition. The solutes are convected to the CL by the evaporation-induced fluid flow and then left behind as the CL recedes. Therefore, the fluid flow inside the droplet and moving ability of the CL are two crucial factors in determining the final deposition pattern. Calculating the flow velocity and motion of the CL during evaporation is the main task of theoretical models of drying.
Figure 3. (a) Sketch of the essential core part of the geometry of every deposition process where material is left behind by a moving CL of suspension or solution with a volatile solvent. Reprinted from publication [27]. Copyright (2020), with permission from Elsevier. (b) Schematic of a droplet with axis-symmetry in a cylindrical coordinate system (side view). Relevant parameters are the radius of the CL R(t), height of the droplet at the center H(t), contact angle θ(t), evaporation rate J(r, t) and profile of liquid-vapor interface h(r, t). |
For simplicity, common assumptions used in the theory are that a sessile droplet has axis-symmetry in a cylindrical coordinate, as shown in figure 3(b), and the solutes have the same height-averaged radial velocity as the solvent [25, 28]. The latter assumption is valid for droplets where the size of the laden particles is large enough to neglect its diffusion, and on the other hand, it should be small enough to neglect the sedimentation (vertical motion) [25] . Once we have the velocity field of the evaporation-induced fluid flow, we have the motion of all solutes.
It is known that the fluid flow inside the droplet is determined by the shape evolution of the droplet, which is often described by the liquid/vapor interface profile h(r, t). According to the conservation of mass [equation (1 )], the time evolution of h(r, t) can be obtained if the radial velocity field ur(r, t) and the evaporation rate J(r, t) are given.
$\begin{eqnarray}\displaystyle \frac{\partial h}{\partial t}=-\displaystyle \frac{1}{r}\displaystyle \frac{\partial ({{rhu}}_{r})}{\partial r}-J(r,t)\sqrt{1+{\left(\displaystyle \frac{\partial h}{\partial r}\right)}^{2}}.\end{eqnarray}$
Let us consider the J(r, t) first. The evaporation process of droplets has been studied both experimentally and theoretically for many years. Various evaporation models have been proposed in theory [20, 27, 29], including J = J0, $J={J}_{0}$ and $J={J}_{0}\exp (-{{Ar}}^{2})$, where J0 and A are the phenomenological constant. Usually, the evaporation rate can be obtained by the vapor-liquid diffusion model [30, 31].
First, the distribution of vapor concentration c(r, t) can be solved by the diffusion equation,2 ) reduces to ∇2c = 0, with the following boundary conditions,
$\begin{eqnarray}\displaystyle \frac{\partial c}{\partial t}=D{{\rm{\nabla }}}^{2}c.\end{eqnarray}$
The time required for the water concentration in the vapor phase to respond to the droplet shape evolution is of the order of R2/D, where R is the droplet contact radius and D is the diffusion constant of the vapor in the air. Hu et al [31] pointed out that when the ratio of R2/Dtf = cv(1 − H)/ρ is small (where, tf is the lifetime of the droplet, and cv is the saturated vapor concentration), the evaporation process can be regarded as a quasi-steady process. Then, equation ( $\begin{eqnarray}\begin{array}{rcl}r\lt R\,z & = & h(r)\quad c={c}_{v},\\ r\gt R\,z & = & 0\quad J=0,\\ r=\infty \,z & = & \infty \quad c={H}_{v}{c}_{v},\end{array}\end{eqnarray}$
where Hv is the relative humidity. The finite element method (FEM) [31] has been used to solve this problem to obtain the vapor concentration at time t, c(r, t). Figure 4 is a typical solution of the density distribution of c(r, t). Once we have the c(r, t), the evaporation rate J(r, t) can be obtained by the equation, J(r, t) = D∇c.
Figure 4. Contour plot of the vapor concentration distribution above a droplet of radius R = 1 mm and height h0 = 0.364 mm. Parameters used in the FEM method are vapor diffusivity D = 26.1 mm2/s, relative humidity H = 0.40 and saturated vapor concentration on the droplet surface cv = 2.32 × 10−8 g/mm3. Gray bars represent the vapor concentration in g/mm3. Reprinted with permission from [31]. Copyright (2020) the American Chemical Society. |
Besides the FEM calculations of the evaporation rate, Deegan et al [20, 31] analyzed the evaporation of spherical-cap droplets by mimicking the problem to the calculation of the electrostatic potential of a charged conductor with a similar droplet shape. They obtained an analytical expression of the evaporation rate,4 ) are,
$\begin{eqnarray}J(r,t)={J}_{0}(\theta ){\left(1-{\tilde{r}}^{2}\right)}^{-\lambda (\theta )},\end{eqnarray}$
where the rescaled $\tilde{r}$ is defined as r/R and θ(t) is the droplet contact angle. Both Deegan et al and Hu et al have compared this analytical evaporation rate with the full numerical results and showed that this can be a good approximation of J(r, t) when the droplet contact angle is between 0° and 90°. The detailed forms of the parameters used in equation ( $\begin{eqnarray}\lambda (\theta )=\displaystyle \frac{1}{2}-\displaystyle \frac{\theta }{\pi },\end{eqnarray}$
$\begin{eqnarray}\begin{array}{l}{J}_{0}(\theta )=\displaystyle \frac{{{Dc}}_{v}(1-{H}_{v})}{R}(0.27{\theta }^{2}+1.30)\\ \quad \ \times \ \left(0.6381-0.2239{\left(\theta -\displaystyle \frac{\pi }{4}\right)}^{2}\right).\end{array}\end{eqnarray}$
Given the evaporation rate, the velocity field of the fluid flow inside the droplet can be obtained by solving the continuity and Navier–Stokes equations. For the quasi-steady process, the inertia term can be neglected. Then, the equations in a cylindrical coordinate are,
$\begin{eqnarray}\displaystyle \frac{1}{r}\displaystyle \frac{\partial ({{ru}}_{r})}{\partial r}+\displaystyle \frac{\partial {u}_{z}}{\partial z}=0,\end{eqnarray}$
$\begin{eqnarray}\eta \left(\displaystyle \frac{\partial }{\partial r}\left(\displaystyle \frac{1}{r}\displaystyle \frac{\partial }{\partial r}({{ru}}_{r})\right)+\displaystyle \frac{{\partial }^{2}{u}_{r}}{\partial {z}^{2}}\right)=\displaystyle \frac{\partial P}{\partial r},\end{eqnarray}$
$\begin{eqnarray}\eta \left(\displaystyle \frac{1}{r}\displaystyle \frac{\partial }{\partial r}\left(r\displaystyle \frac{\partial {u}_{z}}{\partial r}\right)+\displaystyle \frac{{\partial }^{2}{u}_{z}}{\partial {z}^{2}}\right)=\displaystyle \frac{\partial P}{\partial z},\end{eqnarray}$
where uz is the vertical velocity, η is the liquid viscosity and P is the pressure inside the droplet. Hu and Larson [32] numerically solved the set of equations and obtained the full velocity field including vertical and radial velocities, as shown in figure 5, providing an intuitive image of the flow field inside the droplet.
Figure 5. Streamline plots of the evaporation-induced fluid flow inside the droplets from the FEM for droplets with a contact angle of (a) 40° and (b) 10°. Reprinted with permission from [32]. Copyright (2020) the American Chemical Society. |
It is difficult to obtain an analytical solution of the full model mentioned above. Assumptions have been applied to the full model for special cases to obtain analytical solutions. When the droplet contact angle is small, the lubrication approximation can been applied to reduce the model by ignoring the vertical velocity uz. The reduced model was given as [33, 34]:1 ) and (10 ). As the solutes are assumed to have the same velocity as the solvent, the density distribution of the solute, φ(r, t), can be obtained by the mass conservation equation,
$\begin{eqnarray}\displaystyle \frac{\partial h(r,t)}{\partial t}=\displaystyle \frac{1}{3\eta }\displaystyle \frac{1}{r}\displaystyle \frac{\partial }{\partial r}\left({{rh}}^{3}(r,t)\displaystyle \frac{\partial P(r,t)}{\partial r}\right)-J(r,t),\end{eqnarray}$
where the ∂P/∂r is the gradient of the pressure, which corresponds to the form of h(r, t). Then, the time evolution of the profile h(r, t) and the radial velocity field u(r, t) can be obtained by the combination of equations ( $\begin{eqnarray}\displaystyle \frac{\partial (\phi h)}{\partial t}=-\displaystyle \frac{1}{r}\displaystyle \frac{1}{\partial r}({rh}\phi {u}_{r}).\end{eqnarray}$
Experiments observed two special evaporation modes of drying droplets: the constant contact angle (CCA) and the constant contact radius (CCR) modes, where either the contact angle or contact radius is nearly unchanged during evaporation. For these two cases, analytical solutions of the radial velocity can be obtained. Deegan et al [20] combined the global mass of conservation [equation (12 )] and the spherical-cap assumption [equation (13 )],14 ) can be inserted in equation (11 ) directly and further analysis can be conducted.
$\begin{eqnarray}\begin{array}{l}\displaystyle \frac{{\rm{d}}V}{{\rm{d}}t}={\displaystyle \int }_{0}^{R}{\rm{d}}r^{\prime} 2\pi r^{\prime} h(r^{\prime} ,t)={\displaystyle \int }_{0}^{R}{\rm{d}}r^{\prime} 2\pi r^{\prime} J(r^{\prime} ,t)\\ \quad \,\ \times \ \sqrt{1+{\left(\displaystyle \frac{\partial }{\partial r^{\prime} }h(r^{\prime} ,t)\right)}^{2}},\end{array}\end{eqnarray}$
$\begin{eqnarray}h(r,t)=\sqrt{{\left(\displaystyle \frac{h{\left(0,t\right)}^{2}+{R}^{2}}{2h(0,t)}\right)}^{2}-{r}^{2}}-\displaystyle \frac{{R}^{2}-h{\left(0,t\right)}^{2}}{2h(0,t)}.\end{eqnarray}$
Then, the dimensionless height-averaged radial velocity ${\tilde{u}}_{r}$ is derived as: $\begin{eqnarray}{\tilde{u}}_{r}=\frac{1}{4}\frac{1}{1-\tilde{t}}\frac{1}{\tilde{r}}\left[{\left(1-{\tilde{r}}^{2}\right)}^{-\lambda (\theta )}-\left(1-{\tilde{r}}^{2}\right)\right].\end{eqnarray}$
Equation (In conclusion, the conventional hydrodynamic theory has been mostly based on the conservation law and N-S equations. Under the basic framework, the main task is to solve the high-order partial differential equations (PDEs) [equation (10 )] either analytically or numerically. It is worth mentioning that even more complicated conditions such as drying on a porous substrate [33], drying with a receding CL [27, 28] and drying with a non-uniform temperature system [35, 36] have also been studied in recent years.
Besides theory, simulations such as the Monte Carlo (MC) [37–43], molecular dynamics (MD) [44–46], density functional theory (DFT) [43, 47–52] and the lattice Boltzmann method (LBM) [53–59] are also used to study the drying of liquid droplets. Simulations such as experiments can provide some insight into these problems, which is very useful to construct the theoretical model. A detailed review of the simulation of drying droplets was given by Thiele et al [43].
2. OVP theory of drying droplets
In this section, we will review the OVP theory of the evaporation process and the formation of the deposition pattern of drying droplets. We show how such a theoretical framework can be formalized and the essential physics underlying this drying problem. We further show that such a simple theory can be used as the basis of MC simulation to predict deposition patterns in high dimension and even for cases without radial symmetry.
For a soft matter system, we can usually choose a set of parameters x = (x1, x2, ⋯ , xn) to describe the state of the system. If it is a force free system, i.e. the dissipation force is balanced by the potential force, then the time evolution of the system can be obtained as:15 ) represents the potential force driving the system to its equilibrium state, while the first term is the frictional force resisting this change.
$\begin{eqnarray}-\displaystyle \sum _{j}{\zeta }_{{ij}}(x){\dot{x}}_{j}-\displaystyle \frac{\partial A(x)}{\partial {x}_{i}}=0,\end{eqnarray}$
where A(x) is the free energy of the system and ζij(x) the associated frictional coefficient. The second term of equation (In 1931, Onsager noticed that many phenomenological equations for the time evolution of non-equilibrium systems have the same structure as those of particle-fluid systems [60, 61]. Then, with the Onsager's reciprocal relation ζij = ζji, we can construct time-evolution equation (16 ), by minimizing the Rayleighian function,
$\begin{eqnarray}{\mathfrak{R}}(\dot{x};x)=\displaystyle \frac{1}{2}\displaystyle \sum _{i,j}{\zeta }_{{ij}}{\dot{x}}_{i}{\dot{x}}_{j}+\displaystyle \sum _{i}\displaystyle \frac{\partial A}{\partial {x}_{i}}{\dot{x}}_{i},\end{eqnarray}$
where $\displaystyle \frac{1}{2}\sum {\zeta }_{{ij}}{\dot{x}}_{i}{\dot{x}}_{j}$ is the dissipation function Φ, which is equal to half of the energy dissipation taking place in the system per unit time. Therefore, for a given soft matter dynamic problem, we can first construct a Rayleighian function, $\begin{eqnarray}{\mathfrak{R}}(\dot{x};x)={\rm{\Phi }}(\dot{x};x)+\dot{A}(\dot{x};x),\end{eqnarray}$
which contains the energy dissipation function Φ and the time change rate of the free energy $\dot{A}$ of the system. Then, the evolution equations are obtained by, $\begin{eqnarray}\partial {\mathfrak{R}}/\partial {\dot{x}}_{i}=0,\end{eqnarray}$
where i = 1, 2, …, f stands for the f system characterizing variables.This theoretical framework is useful to construct not only conventional hydrodynamic equations, but some unknown evolution equations for complex soft matter problems [62, 63]. These dynamic equations are usually high-order partial differential equations, which require complicated numerical solutions. In recent years, Doi and co-workers found that the OVP theory can be used as a tool of approximation to derive first-order evolution equations of non-equilibrium systems. Such simple models simplify the numerical calculations or can even be solved analytically, bringing insightful fundamental understanding of non-equilibrium systems.
We use the sessile droplet as an example of how we use the OVP theory to study non-equilibrium problems. For a liquid droplet evaporating on the substrate, we consider R(t) as the droplet contact radius and H(t) as the height in the center. The surface profile is described by h(r, t) at time t, which is the same as figure 3(b).
Applying the lubrication approximation, the energy dissipation due to the fluid flow inside the droplet can be written as:
$\begin{eqnarray}{\rm{\Phi }}=\displaystyle \frac{1}{2}{\int }_{0}^{R}2\pi r\displaystyle \frac{3\eta }{h}{v}^{2}{\rm{d}}r,\end{eqnarray}$
where η is the liquid viscosity and v(r, t) is the height-averaged velocity of the fluid at position r inside the droplet.The free energy F is the sum of the interfacial energy:
$\begin{eqnarray}\begin{array}{rcl}F & = & \left({\gamma }_{\mathrm{LS}}-{\gamma }_{\mathrm{SV}}\right)\pi {R}^{2}+{\gamma }_{\mathrm{LV}}{\displaystyle \int }_{0}^{R}2\pi r\sqrt{1\,+\,h{{\prime} }^{2}}{\rm{d}}r\\ & \approx & \ \displaystyle \frac{{\gamma }_{\mathrm{LV}}\pi {R}^{2}{\theta }_{e}^{2}}{2}+{\gamma }_{\mathrm{LV}}{\displaystyle \int }_{0}^{R}\pi {rh}{{\prime} }^{2}{\rm{d}}r,\end{array}\end{eqnarray}$
where $h^{\prime} =\partial h/\partial r$ and $| h^{\prime} | \ll 1$ when droplets have a small contact angle. The surface tension of liquid/vapor, substrate/liquid and substrate/vapor are γLV, γSL and γSV, respectively.We first show that minimizing the Rayleighian function without any assumption is a systematic way of constructing the corresponding conventional hydrodynamic equation.
The height-averaged velocity satisfies the conservation law,
$\begin{eqnarray}\dot{h}=-\displaystyle \frac{1}{r}\displaystyle \frac{\partial }{\partial r}({rvh})-J.\end{eqnarray}$
Then, the time change of the free energy can be obtained from equations (20 ) and (21 ) as:19 ) and (22 ). We minimized it with respect to v, resulting in:
$\begin{eqnarray}\dot{F}=-2{\gamma }_{\mathrm{LV}}\pi {\int }_{0}^{R}\displaystyle \frac{\partial }{\partial r}\left[\displaystyle \frac{1}{r}\displaystyle \frac{\partial }{\partial r}\left(r\displaystyle \frac{\partial h}{\partial r}\right)\right]\left({rhv}\right){\rm{d}}r,\end{eqnarray}$
where we have dropped all terms that do not depend on v, since we only need to minimize it later with respect to v. Then, the Rayleighian function is given by the addition of equations ( $\begin{eqnarray}v=\displaystyle \frac{{\gamma }_{\mathrm{LV}}{h}^{2}}{3\eta }\displaystyle \frac{\partial }{\partial r}\left[\displaystyle \frac{1}{r}\displaystyle \frac{\partial }{\partial r}\left(r\displaystyle \frac{\partial h}{\partial r}\right)\right].\end{eqnarray}$
Inserting the expression v(r) into equation (21 ), the mass conservation equation becomes,10 ).
$\begin{eqnarray}\dot{h}=-\displaystyle \frac{{\gamma }_{\mathrm{LV}}}{3\eta r}\displaystyle \frac{\partial }{\partial r}\left[{h}^{3}r\displaystyle \frac{\partial }{\partial r}\left(\displaystyle \frac{1}{r}\displaystyle \frac{\partial }{\partial r}(r\displaystyle \frac{\partial h}{\partial r})\right)\right]-J.\end{eqnarray}$
We therefore obtained the droplet shape evolution equation, which is the same as the conventional hydrodynamic equation (Actually, we can follow the same procedure to construct most conventional hydrodynamic equations and even for problems when we do not know their dynamic equations. These hydrodynamic equations are usually high-order partial differential equations that require complicated numerical solutions. The OVP theory can be further used as an approximation to reduce high-order PDEs to first-order dynamic equations. Then, the numerical solutions become much simpler and even analytical solutions can be obtained.
The droplet drying model equation (24 ) can be simplified by introducing an assumption that the droplet has a parabolic profile when the contact angle θ is small. The surface profile now becomes,
$\begin{eqnarray}h(r,t)=H(t)\left(1-\displaystyle \frac{{r}^{2}}{{R}^{2}}\right).\end{eqnarray}$
With this assumption, the volume of the droplet has the form $V=\displaystyle \frac{1}{2}\pi {{HR}}^{2}$ and the contact angle θ of the droplet can be written as θ = 4V/πR3.For the entire evaporation process, the volume change rate $\dot{V}(t)$ is known to be proportional to the contact radius R(t),
$\begin{eqnarray}\dot{V}(t)={\dot{V}}_{0}\displaystyle \frac{R(t)}{{R}_{0}},\end{eqnarray}$
where V0 and R0 are the initial values of V(t) and R(t). We assume a homogeneous evaporation rate [64], $\begin{eqnarray}J=-\displaystyle \frac{\dot{V}}{\pi {R}^{2}}.\end{eqnarray}$
It is worth noting that a realistic evaporation rate should be a spatial function of r and can even be divergent at the CL. Although this homogeneous evaporation rate is used to simplify the model, the series of works reviewed here showed that ignoring the position dependence of the evaporation rate is reasonable in the study of the ring-like deposition patterns of drying droplets. However, if one wants to study more complicated deposition patterns or focus on the shape evolution of drying droplets, a spatially dependent evaporation rate has to be used.
Combining equations (21 ), (25 ), (26 ) and (27 ), we can solve v(r) as:
$\begin{eqnarray}v(r)=r\left(\displaystyle \frac{\dot{R}}{R}-\displaystyle \frac{\dot{V}}{4V}\right)=r\left(\displaystyle \frac{\dot{R}}{2R}-\displaystyle \frac{\dot{H}}{4H}\right).\end{eqnarray}$
Inserting the expression of h(r, t) and v(r) into the dissipation function, we have: $\begin{eqnarray}{\rm{\Phi }}=\pi R{\xi }_{\mathrm{hydro}}{\left(\dot{R}-\displaystyle \frac{R\dot{V}}{4V}\right)}^{2},\end{eqnarray}$
where the hydrodynamic friction parameter ξhydro = 3C/ηθ, $C=\mathrm{ln}\,(R/2\epsilon )-1$ and ε is the molecular cutoff length.On real surfaces, there is an extra energy dissipation due to the CL motion, which should be considered in the dissipation function,29 ) and (30 ), the total dissipation function is written as:
$\begin{eqnarray}{{\rm{\Phi }}}_{\mathrm{cl}}=\displaystyle \frac{1}{2}{\xi }_{\mathrm{cl}}2\pi R{\dot{R}}^{2},\end{eqnarray}$
where ξcl is a phenomenological parameter representing extra friction due to the surface roughness and other sources of dissipation. Combining equations ( $\begin{eqnarray}{\rm{\Phi }}=2\pi R\left[\displaystyle \frac{1}{2}{\xi }_{\mathrm{hydro}}{\left(\dot{R}-\displaystyle \frac{R\dot{V}}{4V}\right)}^{2}+\displaystyle \frac{1}{2}{\xi }_{\mathrm{cl}}{\dot{R}}^{2}\right].\end{eqnarray}$
We now turn to the free energy of the system. As we choose R and V as the slow variables, the free energy will be written as a function of the two slow variables. By assuming a parabolic form of the droplet interface profile, the free energy now has an explicit expression:
$\begin{eqnarray}F={\gamma }_{\mathrm{LV}}\left(\displaystyle \frac{4{V}^{2}}{\pi {R}^{4}}+\displaystyle \frac{1}{2}\pi {R}^{2}{\theta }_{e}^{2}\right).\end{eqnarray}$
Accordingly, the time change rate of the free energy becomes, $\begin{eqnarray}\dot{F}={\gamma }_{\mathrm{LV}}\left(-\displaystyle \frac{16{V}^{2}}{\pi {R}^{5}}+\pi {\theta }_{e}^{2}R\right)\dot{R}+{\gamma }_{\mathrm{LV}}\displaystyle \frac{8V\dot{V}}{\pi {R}^{4}}.\end{eqnarray}$
The Rayleighian function now becomes, $\begin{eqnarray}\begin{array}{l}{\mathfrak{R}}=2\pi R\left(\displaystyle \frac{1}{2}{\xi }_{\mathrm{hydro}}{\left(\dot{R}-\displaystyle \frac{R\dot{V}}{4V}\right)}^{2}+\displaystyle \frac{1}{2}{\xi }_{\mathrm{cl}}{\dot{R}}^{2}\right)\\ \quad \ +\ {\gamma }_{\mathrm{LV}}\left(-\displaystyle \frac{16{V}^{2}}{\pi {R}^{5}}+\pi {\theta }_{e}^{2}R\right)\dot{R}+{\gamma }_{\mathrm{LV}}\displaystyle \frac{8V\dot{V}}{\pi {R}^{4}}.\end{array}\end{eqnarray}$
With this Rayleighian function, we can obtain the droplet shape evolution and deposition patterns.2.1. Shape evolution equations of evaporating droplets
The droplet shape evolution is described by $\dot{R}$ and $\dot{\theta }$, which can be obtained from the above Rayleighian function. The time evolution of contact radius $\dot{R}(t)$ is obtained by $\partial {\mathfrak{R}}/\partial \dot{R}=0$. Combined with the evaporation rate, equation (26 ), we can write down the set of evolution equations of the system in a scaled form [64]:
$\begin{eqnarray}\left(1+{k}_{\mathrm{cl}}\right){\tau }_{\mathrm{ev}}\dot{R}=-\displaystyle \frac{{R}^{2}{V}_{0}}{4{R}_{0}V}+\displaystyle \frac{{V}_{0}^{\tfrac{1}{3}}\theta \left({\theta }^{2}-{\theta }_{e}^{2}\right)}{6{{Ck}}_{\mathrm{ev}}{\theta }_{e}^{3}},\end{eqnarray}$
$\begin{eqnarray}{\tau }_{\mathrm{ev}}\dot{V}=-{V}_{0}\displaystyle \frac{R(t)}{{R}_{0}},\end{eqnarray}$
$\begin{eqnarray}\theta (t)=\displaystyle \frac{4V(t)}{\pi {R}^{3}(t)},\end{eqnarray}$
where ${\tau }_{\mathrm{ev}}=-{V}_{0}/\dot{{V}_{0}}$ is the evaporation characteristic time and ${\tau }_{\mathrm{re}}=\eta {V}_{0}^{\tfrac{1}{3}}/{\gamma }_{\mathrm{LV}}{\theta }_{e}^{3}$ is the relaxation time.Here, we define ${k}_{\mathrm{ev}}={\tau }_{\mathrm{re}}/{\tau }_{\mathrm{ev}}$ to represent the evaporation rate and kcl = ξclθ/3Cη to describe the moving ability of the CL. It is clear that the set of equations (36 ), (37 ) and (38 ) are just first-order ordinary differential equations. It will be much easier for the numerical solutions than the conventional hydrodynamic equations, i.e. high-order PDE, equation (10 ).
We notice that the parabolic assumption [equation (25 )] is only valid when θe is small, and cannot be used to study droplets with a large contact angle. In a recent work [65], Jiang et al used a spherical shape of the droplet liquid-vapor interface and set the radius of sphere $\tilde{R}$ and contact angle θ as the slow variables. By using the OVP theory, they obtained the equation of $\dot{\tilde{R}}$ and $\dot{\theta }$. Therefore, this model can be used to study the droplet shape evolution during evaporation for droplets with equilibrium contact angle θe from 0 to π/2. (figure 6 shows the droplet shape evolution of a droplet with θe = 0.8). Based on the full model, various evaporation modes can be obtained. It is interesting to note that for special cases such as CCA and CCR, analytical results can be derived for $\dot{R}$ and $\dot{\theta }$.
Figure 6. Droplet shape evolution for evaporation on frictionless substrate for three situations: (a) θ0 < θe, (b) θ0 = θe and (c) θ0 > θe. Corresponding evolution of the CL r/r0 is shown in panel (d) and the contact angle θ is shown in panel (e). For all calculations, θe = 0.8 and kev = 0.001. Reproduced with permission from [65]. |
2.2. Deposition patterns
After obtaining the droplet shape evolution, we now discuss how to calculate the deposition patterns. Moreover, we will show that the dynamic model from the OVP theory can also be used for the MC simulations of deposition patterns for high-dimensional and two-droplet systems.
We assume that the velocity of the solute is the same as the solvent given by equation (28 ). As we now know the velocity of solutes, we can trace the position of solutes first located at r0 of any time t, $\tilde{r}({r}_{0},t)$. If at time td the solute is out of the droplet, $\tilde{r}=R({t}_{d})$, the solute deposits at this time and at the location of $\tilde{r}({r}_{0},{t}_{d})$. According to the mass conservation of particles, the concentration μ can be calculated as:
$\begin{eqnarray}\mu \left[\tilde{r}({r}_{0},t)\right]=h({r}_{0},0){\phi }_{0}\frac{2\pi {r}_{0}{\rm{d}}{r}_{0}}{2\pi \tilde{r}{\rm{d}}\tilde{r}},\end{eqnarray}$
where φ0 is the homogeneous concentration at the initial time.Therefore, $\mu \left[\tilde{r}({r}_{0},{t}_{d})\right]$ gives the final density distribution of the deposited solutes on the substrate. In a short summary that when $\tilde{r}({r}_{0},{t}_{d})=R({t}_{d})$, solute with an initial position r0 deposits at time td and at the position of R(td). The final deposition pattern is described by $\mu \left[\tilde{r}({r}_{0},{t}_{d})\right]$.
A series of works about the study of deposition patterns of drying droplets have been conducted by using the above-mentioned theoretical framework. In [64], Man and Doi clarified how the CL friction (described by kcl) and the evaporation rate (described by kev) affect the final deposition pattern. Figure 7 shows that the deposition pattern changes continuously from a coffee-ring to volcano-like and to mountain-like as the kcl increases. When kcl = 100, the CL hardly moves from the initial position and the coffee-ring pattern appears. On the other hand, when kcl = 0, the CL recedes freely and a mountain-like pattern appears. Especially when the evaporation rate is fast (kev ≫ 1), the model reduces to a simple form and an analytical expression of the peak position of the pattern rp can be obtained as:
$\begin{eqnarray}{r}_{p}={R}_{0}{\left(\displaystyle \frac{{k}_{\mathrm{cl}}}{2{k}_{\mathrm{cl}}+1}\right)}^{\tfrac{1}{2(1+{k}_{\mathrm{cl}})}}.\end{eqnarray}$
This theoretical prediction has been confirmed by a later experiment [66].
Figure 7. (a) and (b) are the profile of the deposits left on the substrate when the drying is completed. Transition from coffee-ring to volcano-like and then mountain-like pattern induced by changing the value of kcl from 100 to 0. (a) The case of a fast evaporation rate characterized by kev = 1 and (b) the case where the evaporation rate is small with kev = 10−3. For both cases, θe = θ0 = 0.2 and Δt/τev = 10−5. (c) Sketch map of the transition from coffee-ring to mountain-like pattern. Reprinted figure with permission from [64]. Copyright (2020) the American Physical Society. |
Instead of taking kcl as constant, Wu et al [67] assumed that it has two values, a low value for freely moving CL and a high value for pinned CL, written as:
$\begin{eqnarray}{k}_{\mathrm{cl}}=\left\{\begin{array}{cc}0 & \mathrm{for}\ \theta \leqslant {\theta }_{R}\ \mathrm{or}\ \dot{\theta }\gt 0,\\ \alpha & \mathrm{for}\ \theta \gt {\theta }_{R}\ \mathrm{or}\ \dot{\theta }\leqslant 0,\end{array}\right.\end{eqnarray}$
where θR is the receding contact angle (the angle below which the CL starts to recede) and α is a constant representing the CL moving ability. Then a stick-slip motion of the CL can be obtained by carefully choosing the value of α. They obtained various multi-ring deposition patterns of a drying droplet and pointed out that the coffee-ring and mountain-like deposition pattern can be regarded as a multi-ring pattern. The coffee-ring pattern is a multi-ring pattern with an infinite large inter-ring spacing, while a mountain-like pattern is a multi-ring pattern with an 0 inter-ring spacing, as shown in figure 8. In addition, [67] shows that as the ring radius decreases and below a critical value (the innermost ring radius) the multi-ring is replaced by a solid-circle pattern. An analytical expression of the radius of this solid circle in terms of experimental parameters has been given as: $\begin{eqnarray}{R}_{S}=\left\{\begin{array}{cc}\displaystyle \frac{2C{\theta }_{e}^{3}{V}_{0}^{2/3}}{\pi {R}_{0}}\displaystyle \frac{{k}_{\mathrm{ev}}}{{\theta }_{R}^{2}({\theta }_{e}^{2}-{\theta }_{R}^{2})} & {\theta }_{R}\gt {\theta }_{c},\\ \displaystyle \frac{8{{CV}}_{0}^{2/3}{k}_{\mathrm{ev}}}{\pi {R}_{0}{\theta }_{e}} & {\theta }_{R}\lt {\theta }_{c},\end{array}\right.\end{eqnarray}$
where ${\theta }_{c}={\theta }_{e}/\sqrt{2}$ in that model.
Figure 8. Various multi-ring deposition patterns of drying droplets with different values of θR: (a) 0.00, (b) 0.08, (c) 0.16, (d)0.24, (e) 0.32 and (f) 0.4. The larger θR is, the higher the moving ability of the CL. For all calculations, θ0 = θe = 0.4, R0 = 4 and kev = 10−3. Reprinted with permission from [67]. Copyright (2020) the American Chemical Society. |
Later, Wu et al [68] proposed an extended model based on the OVP theory for the drying of liquid droplets of surfactant solutions. Surfactants can change the droplet liquid/vapor surface tension,
$\begin{eqnarray}{\gamma }_{\mathrm{LV}}={\gamma }_{\mathrm{LV}}^{0}-{k}_{B}T{{\rm{\Gamma }}}_{\infty }\mathrm{ln}\left(1+\displaystyle \frac{c}{b}\right),\end{eqnarray}$
where ${\gamma }_{\mathrm{LV}}^{0}$ is the surface tension of the solvent without surfactants, kB is the Boltzmann constant, T is the temperature, Γ∞ is the maximum surfactant surface concentration and b−1 is the Langmuir constant. As the surface tension now is a function of the surfactant concentration, the added surfactant can change the motion of the droplet CL, either hindering droplet receding or facilitating droplet spreading. Correspondingly, the deposition pattern can change from mountain-like to volcano-like to coffee-ring by increasing the surfactant initial concentration, as shown in figure 9(a), providing another way to control the deposition pattern. They also showed that when the CL motion undergoes a clear receding-advancing transition, a double-ring pattern is formed, as shown in figure 9(c). It is worth noting that this two-ring pattern is induced by adding surfactants into the droplet, which is different to the stick-slip mechanism of the formation of the multi-ring pattern mentioned above.
Figure 9. (a) shows the effect of a surfactant in the deposition patterns. When the initial surfactant concentration β is increased, the deposition pattern changes from mountain-like to volcano-like and then to coffee-ring pattern. All parameters are kcl = 0, kev = 0.002 and α = 1.0. (b) and (c) are the formation of the two-ring deposition pattern; (b) is the evolution of the profile of the surface and (c) is the evolution of the density distribution of solutes left on the substrate, μ(r, t)/(H0φ0). For (b) and (c), kev = 10−4, kcl = 0, α = 0.3 and β = 0.1. Reprinted with permission from [68]. Copyright (2020) the American Chemical Society. |
The model is not limited to the single-droplet system and can be extended for the study of the drying of two droplets. Hu et al [69] studied the drying of two neighboring droplets, focusing on the corresponding deposition patterns. The two droplets evaporate simultaneously and each affects the evaporation of the other. Therefore, the evaporation rate does not have axis-symmetry anymore and an asymmetric evaporation rate J(r, t) has to be used,
$\begin{eqnarray}J(\vec{r},t)=-\displaystyle \frac{\dot{V}(t)}{\pi {R}^{2}(t)}+{J}_{e}\left(\displaystyle \frac{x-{x}_{c}}{{x}_{c}}\right),\end{eqnarray}$
where Je is an input parameter characterizing the asymmetry of the evaporation rate and xc is the location of the droplet center. Due to the asymmetric evaporation rate, the evaporation-induced fluid flow does not have axis-symmetry as the single-droplet case anymore. Instead, an x-direction flow is induced, as sketched in figure 10(a). Correspondingly, new deposition patterns including fan-like and eclipse-like patterns are obtained by the authors, as shown in figures 10(b)–(d). In this work, as the axis-symmetry of the fluid flow is broken, one has to use an MC simulation to calculate the deposition patterns, of which the sample evolution is governed by the dynamic equation derived from the OVP theory.
Figure 10. (a) is the sketch map of the influence of the existence of asymmetric evaporation rate Je in the evaporation-induced fluid flow inside the droplet. (b), (c) and (d) are the asymmetric deposition patterns calculated by the MC simulation. (b) Fan-like deposition is obtained when kcl = 0 and an asymmetric volcano-like deposition pattern is obtained when kcl = 10. (c) Eclipse-like deposition pattern is obtained when kcl = 50. For all the calculations θe = θ0 = 0.2, Δt/τev = 10−5, kev = 0.001 and Je = 0.3. Adapted with permission from [69]. Copyright (2020) the American Chemical Society. |
3. Perspectives and conclusion
We gave detailed descriptions of how to use the OVP theory to construct a conventional hydrodynamic equation of drying droplets. Moreover, we also showed how to use the OVP theory as an approximation tool to reduce a high-order PDE to the first-order ordinary equation of drying droplets. The OVP theory can be treated as a systematic way to study non-equilibrium problems in soft matter. How to choose the proper slow variables of a given problem is one of the key challenging issues in applying this theory. Some extension of the model given in this review can be used to tackle these various problems, like the drying of droplets that contain different shapes of solute particles and the drying of thin film. Wang et al [70] have taken into account the contact angle hysteresis (CAH) based on the model in [64] and analyzed the effects of the CAH on the coffee-ring deposition. Xu et al [71] employed the OVP theory to derive the evolution equations of a mixture of binary fluid thin film. We believe that besides droplet systems, the drying of thin films and the corresponding deposition patterns can be well studied by the OVP theory reviewed here.
In 2011, Yodh et al [72] pointed out that the coffee-ring effect can be suppressed by shape-dependent capillary interactions. In the experiments, they dried water drops containing a suspension of micrometer-sized polystyrene spheres stretched asymmetrically to different aspect ratios. Results showed that droplets laden with sphere particles form coffee-ring patterns, but with ellipsoidal particles can form uniform thin film, as shown in figure 11(a). It is clear that the shape of the solute particles provides a new way to control the deposition patterns. Besides particle shape, Patil et al [73] proved that particle size can also change the patterns. When ambient temperature is Ts = 27°, smaller monodispersed colloidal polystyrene particles (dP = 0.1, 0.46 and 1.1 μm) form mountain-like deposits on the wet-oxidized silicon surface after the water dried, whereas larger particles (dP = 3 μm) left ring formation. Indeed, it is similar to mountain-like deposition where particles concentrate in the center. When the radius of droplets is beyond the capillary length, we have to consider the effect of gravity on the deposition patterns [74, 75]. This effect is especially important when droplets are on inclined surfaces, which has been discussed both theoretically and experimentally [76, 77].
Figure 11. Images of the final distributions of ellipsoids (a) and spheres (b) after evaporation. Coffee-ring pattern changes to a uniform thin film by replacing the spherical solutes with ellipsoidal solutes, while all other experimental conditions remain the same. Reprinted by permission from [72], Copyright (2020). |
From recent experimental studies, we can generally identify two new factors in tailoring the final deposition pattern of drying droplets. First, the properties of the solute can largely affect the final deposition pattern. We may include these effects by adjusting the diffusion coefficient in the model. The second one is the effect of gravity. This effect can be inserted into the model by introducing an extra term into the free energy. Then, one can analyze the effect of gravity on the fluid flow inside the droplets and see this effect in the final deposition patterns. The development of the OVP theory in soft matter is still at the beginning stage. Various successful examples of using this theory to deal with dynamic problems are still needed to enrich the understanding of the OVP theory.
In this review, we have shown the theoretical framework of the formation of the deposition pattern of drying droplets. By using the OVP theory proposed by Lars Onsager, the high-order conventional hydrodynamics equations can be reduced to first-order evolution equations, which provides us with a new way to understand these complex problems and makes the physical picture clearer. Based on the evolution framework, the formation mechanism of different ring-like patterns has been analyzed and explained. The purpose of this review is to introduce the background of the research on deposition patterns and also review the recent development of the OVP theory in drying droplets. We hope that this short review can help readers who are interested in this problem and especially help readers to build their own models of soft matter dynamic problems of interest.