International Journal of Chemical Engineering

Volume 2019, Article ID 7954965, 12 pages

https://doi.org/10.1155/2019/7954965

## Deposition of Colloidal Particles during the Evaporation of Sessile Drops: Dilute Colloidal Dispersions

Davidson School of Chemical Engineering, Purdue University, West Lafayette, Indiana, USA

Correspondence should be addressed to Michael T. Harris; ude.eudrup@sirrahtm

Received 17 May 2019; Revised 31 July 2019; Accepted 5 September 2019; Published 3 October 2019

Academic Editor: Antonio Brasiello

Copyright © 2019 Pei-Fang Sung et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

The deposition of colloidal silica particles during the evaporation of sessile drops on a smooth substrate has been modeled by the simultaneous solution of the Navier–Stokes equations, the convective-diffusive equation for particles, and the diffusion equation for evaporated vapor in the gas phase. Isothermal conditions were assumed. A mapping was created to show the conditions for various deposition patterns for very dilute suspensions. Based on values of the Peclet (Pe) number and Damkholer numbers (Da and Da_{−1}), the effects of adsorption and desorption were discussed according to the map. Simulations were also done for suspensions with a high particle concentration to form a solid phase during the evaporation by using a packing criterion. The simulations predicted the height and width of the ring deposit near the contact line, and the results compared favorably to experimental particle deposition patterns.

#### 1. Introduction

Using the evaporating sessile drop could be a simple way to make uniform patterns for printing sensors and pharmaceutical drugs and applying novel film coatings [1–4]. The deposition pattern of colloidal silica particles during the evaporation of sessile drops on smooth substrates has been investigated experimentally and computationally. Deegan et al. proposed that the capillary flow inside an evaporating droplet carries the solute radially to the contact line, which forms a ring-patterned deposit [5]. This ring pattern was also observed by other groups for the deposition of polystyrene latex colloids on glass substrates [6, 7]. While some studies focus on the influence of microfluidic flow on the deposition [8, 9], the deposition patterns can also vary from a ring to a uniform layer according to the substrate roughness [10], particle morphology [11], the evaporation speed [12] and mode [13], and the composition of the solvent [14]. Sommer compared the deposition pattern in an evaporating sessile drop between polystyrene latex particles and the hydroxyapatite particles [15]. While both are depositing onto a titanium disk, a ring formed in the polystyrene-laden solution while a uniform layer was obtained from the hydroxyapatite solution.

Bhardwaj et al. studied the role of the Derjaguin–Landau–Verwey–Overbeek (DLVO) interactions on the final deposition patterns [16]. The particle (titania) surface potential and the substrate (glass) surface potential were measured under different pH values ranging from 1.4 to 11.7 by adding HCL (0.1 M) or NaOH (0.01 M). The competition between the hydrodynamic radial flow and the particle/surface interaction has a significant effect on the final deposition. If the attractive force between the particles and the substrate is large enough and the Peclet number is low enough, the particles tend to stick to the surface before they are swept to the contact line by the radial flow and finally form a uniform deposition pattern.

To calculate the particle deposition growth, Deegan found that the mass of the ring grows versus time based on the power law by assuming depth-average velocity and particle concentration, which is valid only in very thin droplets [5, 17]. Fischer solved the particle deposition profiles under three types of evaporation flux distribution along the droplet surface using lubrication theory, also assuming a very thin droplet and neglecting particle diffusion [18]. Hu and Larson used Brownian dynamics simulation to calculate the particle deposition profile based on the analytical form of the velocity in both radial and vertical components [19].

The adsorption and desorption of particles to the substrate can be modeled as a first-order reaction [20]. Bhardwaj et al. solved the vapor concentration, fluid profile, and particle concentration distribution inside the droplet using the finite element method [21]. The depinning phenomena of the droplet were also considered. In their simulation, once the particle concentration reaches the maximum packing *X* = 0.7, the particle deposit starts growing and the depinning happens. Widjaja and Harris solved the vapor concentration, fluid profile, and particle concentration distribution inside the droplet simultaneously by solving the Stokes equations of fluid motion and the convective-diffusive equations using the finite element method [22]. The particle deposition profile was calculated assuming that the flux of the particle deposition follows a first-order reaction. A phase diagram of the final deposition profile under different Peclet (Pe) and Damkholer (Da) numbers was obtained.

Desorption of particles can happen since there are flows inside the drop that can move the particles that have deposited. The fact that flows facilitate desorption of particles from the fluid-solid interface has been revealed previously [23, 24]. Yet limited research has incorporated the desorption behavior of particles to the studies of drop evaporation. The lack of desorption in the model may lead to incorrect prediction of the deposition profiles.

In this paper, the vapor, fluid, and particle phases were also solved simultaneously using the finite element technique. In addition to solving the particle concentration in the drop phase, the equation for particle deposition at the fluid/substrate interface was included, and the deposition profile was calculated from the interface concentration. The addition of the interface concentration makes it possible to include not only the adsorption effect but the desorption effect. The results are interpreted in two parts. First, a parametric study is performed with different convection, adsorption, and desorption effects, i.e., final deposition profiles were obtained under different Pe, Da, and Da_{−1} (desorption Damkohler number). Second, the maximum packing criterion was applied to account for the formation of particle accumulation near the contact line. Different from Bhardwaj’s setup [21], this paper works on a drop whose contact line does not move, like a water drop with suspended silica particles on a glass substrate.

#### 2. Numerical Simulation

##### 2.1. System Description

The finite element technique is used to solve the equation governing the vapor concentration, fluid velocity, and the particle concentration inside the droplet. The droplet was assumed to be axisymmetric, and the fluid properties such as viscosity, density, and the surface tension are constant throughout the entire evaporation process. While the volume of the droplet is decreasing, the contact line stays pinned, and the simulations were stopped when the contact angle decreased to about 5°. The governing equations and boundary conditions for each of the three phases will be shown in the following section.

##### 2.2. Governing Equations and Boundary Conditions

The system contains three phases: vapor, liquid, and solid (particle) phases. The evaporation of the droplet is assumed to undergo a pseudosteady state evaporation process, and the vapor is stagnant; therefore, the governing equation simplifies to the Laplace equation as follows (dimensionless form):

The boundary conditions for the vapor phase are as follows:where *c* is the dimensionless vapor concentration, *r* and *z* are the dimensionless radial and height distance from the center of the drop base, *H* is the relative humidity, and **n** is the unit normal vector. Explanation of the boundary condition (3) is in the SI. The characteristic length for the entire system was chosen to be the drop base radius *R*_{c}; the characteristic vapor concentration is the saturated vapor concentration of the liquid in the evaporating droplet. Once the vapor concentration is solved, the dimensionless evaporation flux can be calculated according to Fick’s Law, , where **n**_{s} is the unit normal vector at the drop-air interface.

To solve the fluid flow inside the evaporating droplet, we assume the drop to be incompressible and Newtonian and solve the continuity equation (5) and the Stokes equation (6). Here, we consider the case where the droplet is small enough so that the shape of the drop is part of the spherical cap where gravity has a negligible effect on the drop shape. The fluid properties (e.g., density, viscosity, and surface tension) are assumed to be constants during the evaporation process. The contact line is assumed to stay pinned until the contact angle is about 5°. The governing equations and the boundary conditions are listed as follows (dimensionless form):where **V** is the dimensionless velocity vector, .

The boundary conditions are (dimensionless form)where is the nodal velocity along the drop surface, 2*H*_{S} is twice the local mean curvature of the surface, and Ca is the capillary number, which is defined by Ca = *μv*_{c}/*σ*, where *μ* is the fluid viscosity and *σ* is the drop-air interface surface tension. is the characteristic velocity which is defined by *J*_{c}*/ρ* where *ρ* is the liquid density and *J*_{c} is characteristic evaporation flux defined by where is the diffusivity of water/air and *c*_{0} is the saturated vapor concentration. For water diffusing into air, and the saturation concentration for water, .

The convective-diffusive equation was used to solve the particle concentration inside the droplet as follows (dimensionless form):where the Peclet number, is the diffusivity of the silica particles in water with the value of 4.37 × 10^{−12}. The characteristic particle concentration is *C*_{P,0}, which is the initial particle concentration and is assumed to be initially uniform in the drop phase. With this model, the particles are assumed small enough so that Brownian forces dominate over gravitational forces, and thus, gravitational effects on the particles are negligible.

We assumed that there are no particles flowing across the axisymmetric line. The boundary condition along the droplet surface is the particle mass conservation at the interface during the evaporation, which is

The boundary condition on the substrate is the deposition of the particles onto the substrate, in which the deposition flux *J*_{p} is the sum of the adsorption and desorption [25] as follows:

Here, is the Damkholer number which is a ratio of particle adsorption over particle diffusion; is the desorption Damkholer number which is a ratio of particle desorption over particle convection, where is the characteristic time. *k* and *k*_{−1} are the adsorption and desorption rate constants. It is assumed that the rates of adsorption/desorption of particles onto the substrate are the same as the particle adsorption/desorption onto particles that have already deposited, and they both take the form of a first-order reaction. Details on the derivation and explanation of this equation are discussed in the SI. Γ_{S} is the dimensionless surface particle concentration and is also one of the variables to be solved. The governing equation to solve along the nodes on the substrate is as follows (dimensionless form):

Equations (8)∼(11) were used to solve for *C*_{P} and Γ_{S}. In this study, we investigated the effect of the dimensionless groups Pe, Da, and Da_{−1} on particle final deposition patterns. In order to compare the final deposition patterns for different Pe numbers, the final deposition will be plotted as Γ_{S}/Pe versus position *r*.

Furthermore, the case with high particle concentrations was simulated. When particle concentration is high, packing of particles happen near the contact line, as what will be discussed in Section 3.3. Therefore, the contact line region was treated as follows.

The particle concentration for each element was checked at each time step to see if the critical concentration *C*_{P,max} is achieved. If the concentration in any element achieves the critical concentration, the convective-diffusive equation in this “solid” element will not be solved anymore, and the nodal position was “fixed.”

Once the “solid” elements are determined, the *C*_{P} values on the nodes inside the solid element will be equal to the critical value, noted as *C*_{P,max} = 130.86. The common nodes between the solid element and the neighboring fluid elements are fixed and are applied under no-flux boundary condition of particles as in equation (10) in the Method section.

We assumed that the particles in the “solid” elements are still wetted by the liquid. The fluid dynamics were still being solved, i.e., evaporation is still considered on the surface of the solid phase and drives the liquid in the drop phase to flow toward the contact line. A sketch is shown in the SI. In the solid phase, the equations of fluid dynamics keep the same as the liquid phase. A “superficial velocity” resulting from the fluid dynamics equations describes the averaged velocity in the solid phase.

##### 2.3. Code Validation

After the vapor concentration profile was calculated, the evaporation flux along the droplet surface and the total evaporation rate were computed. The total evaporation rate as a function of the contact angle computed from our FEM algorithm was compared to the analytical results from Picknett and Bexon [26] as shown in Figure 1.