Gravity induced shape effects on the time-dependent evaporation of pendant drops

The paper presents a method to model the time-dependent evaporation of pendant drops taking into account the effect of drop deformation induced by gravity. The model is based on the solution to the time-dependent drop mass and energy conservation equations, where the mass and energy fluxes through the gas mixture are numerically evaluated for a range of Bond numbers and contact angles. The evaporation characteristics of pendant and sessile drops on hydrophobic and hydrophilic substrates are compared in terms of evaporation times and evaporative cooling, for both constant contact angle and constant contact radius modes.


Introduction
The interest of the scientific community on the evaporation of a drop in contact with a solid surface is witnessed by the huge literature available on this topic.This is mainly due to its importance in many applications (like, printing technology [1], cooling systems [2], surface patterning [3], biomedical processes, etc.; see [4] for a recent review), but also to the complexity of the phenomena occurring in this process.Pendant drops, which are drops hanging from a surface, are often used for experimental studies on drop evaporation [5,6], and they are important for surface tension measurements of air-liquid surfaces [7,8], in heat transfer processes [9] or in large volume protein single-crystals growth [10].Differently from sessile drops, pendant drops become unstable when their size increases, falling from the substrate, a fact that limits the size of the drops.The shape of a drop in contact with a solid surface has been deeply studied since the pioneering work of [11] and a copious literature on this subject is available (see for example [12] and the cited literature therein).
The effect of gravity on the processes involving pendant drops has been often related to the natural convection caused by buoyancy in the gas surrounding the drop [13].Gravity also influences the flow inside the drop, leading to fluid structures with interesting practical implications (see for example [14,15]).The first study on the specific effect of gravity induced deformation on the evaporation of a pendant drop is possibly that of [16], who found experimentally an increment in the evaporation rate up to about 30%.In an experimental study presented in [17], it was shown that a pendant drop on super-hydrophobic substrates evaporates faster than the sessile one of equal volume, suggesting that this may be partially ascribed to the different drop shapes, caused by the opposite effect of gravity, considered the nonnegligible value of the Bond number (Bo = 0.2).
Pendant and sessile drop evaporation can be considered a typical example of a diffusion-controlled phenomenon.Quasi-steadiness is often assumed, considering that the drop life-time is orders of magnitude larger than the transport time-scale.In the absence of external convection, and assuming constant gas mixture thermophysical properties, the analytical problem reduces to the solution to a Laplace equation [18,19].For drops with a spherical-cap shape (i.e. with negligible effect of gravity on their shape), the analytical solution given by [20] in terms of Mehler-Fock transform is almost universally used (see [21] for some different types of solution).Even when thermophysical properties are temperature dependent, the phenomenon can be modelled through the solution of a Laplace problem (see [22]).The analytical-numerical model reported in the next sections considers the evolution of the shape and temperature of evaporating pendant drops, and quantifies the effect of the actual drop shape and of the evaporative cooling for different evaporation modes and substrate wettability.

Mathematical modelling
In this paper, a mixed analytical-numerical approach will be used to model the diffusion-limited processes taking place during the evaporation of a pendant drop deformed by the effect of gravity.

Modelling the energy and mass fluxes through the gas phase
In the following, the gas phase is assumed to be a binary mixture of the vapourising species and the ambient gas, behaving as a single component ideal gas.Under quasi-steady conditions, the species conservation equations can be written as.
where p=0,1 for the ambient gas and the vapour, respectively, and the molar fluxes are is the total molar flux, c the molar density, D10 the binary diffusion coefficient, Mm (p) the molar mass of the p-component and y (p) is the molar fraction of the p-component.The steady-state energy conservation equation is where kmix is the gaseous mixture conductivity, ( ) c is the vapour specific heat capacity, ρ is the mass density of the mixture and Uj is the mass average velocity.In equation (2) minor terms (like viscous dissipation), assumed to be irrelevant for the present analysis, are neglected (refer to [24] p. 465, or [25] p. 589 for a more complete form).The effect of external convection is neglected, but the Stefan flow is fully considered (Uj is the Stefan flow velocity).Assuming negligible absorption of the gas in the liquid drop and uniform Dirichlet boundary conditions on the drop surface, an analytical solution to the PDE system (1, 2) can be found, also when the gas mixture thermophysical properties (ρ, ( ) c and kmix) depend on temperature (see [22] and Appendix), and the solution is given in terms of an auxiliary function, Φ, solution of the following Laplace problem ) In analytical modelling, pendant and sessile drops are often assumed to have the shape of a spherical cap, since the drop size is usually smaller than the capillary length, ( ) . However, it was shown [26] that the evaporation rate of a sessile drop may be perceptibly influenced by the gravity induced drop deformation even when the drop size is smaller than Lc.The analytical solution of the problem (3) for spherical caps was given in [20] (see also [21] for other types of solutions), but when the drop is deformed by gravity, analytical solutions are not available and numerical tools are needed (see next subsection).The shape of sessile and pendant drops in the gravitational field is fully defined by its volume and the substrate wettability [26], that is by the Bond number Bo R L = (where ( ) = and V is the drop volume,) and the contact angle, θc.Since the evaporation and the heat rates can be found from the integral n A dA ∇ Φ  (see Appendix), both rates depend on Bo and θc and, to find such dependence, a parametric analysis is performed, using the numerical methods described in the following subsections.

Energy transfer through the solid substrate
The heat transfer from the solid substrate can be modelled assuming quasi-steadiness, and assuming that the substrate is much thicker than the drop, the problem is equivalent to the heat conduction from a semiinfinite solid to a disc of radius equal to the drop basis radius, Rc, laying on the solid surface.This problem can be solved using the conductive shape factor method obtaining ( ) where S is the conductive shape factor.For the given configuration (a disk on a semi-infinite body) the shape factor is equal to S=4Rc [27].This simple method can be considered acceptable when the drop temperature is almost constant, and this happens just after the first period of temperature variation, which is usually very short compared to the whole evaporation period.

Numerical methods for the solution of the Laplace problem
The shape of a pendant drop in a gravitational field is found solving numerically the classical Bashford-Adams equation [11] by a Runge-Kutta method with variable time steps.A double iteration cycle was necessary to yield the equilibrium drop shape for given Bo and θc, and the relative accuracy on the drop volume was better than 2×10⁻⁴, while the average accuracy on the contact angle was better than 0.4deg.The obtained drop shape was compared to the numerical solution of the integral method proposed in [28] (see [26] for details), and relative differences of the drop profile coordinates were less than 10⁻⁴.The drop profiles where then used to set the numerical domain to solve the problem (3) using the commercial code COMSOL Multiphysics®. Figure 1 shows some samples of the pendant drop shape on two substrates.The upper limit of the Bond number for the two contact angles is due to the stability constrains for pendant drops (see section 3).A finite element discretisation approach is used on 2D axis-symmetric grids.The computational domain is a half-sphere with a radius R∞=30Req (see figure 2 for a schematic representation); this choice was suggested by previous works [29,30].The Dirichlet condition at infinity is substituted by an equivalent condition on the outer boundary, R= R∞, through an iterative procedure that accounts for the finiteness of the domain size (more details can be found in [26,29,30]).The mesh was unstructured, made by triangular tessellation, with grid refinements around the triple line.An analysis of grid independence was performed (more details can be found in [29,31]) and in the present study grids made of about 10 5 cells were used.The evaporation and the sensible heat rates are evaluated by the integral ɺ is the corresponding value for an equivalent volume spherical cap, which is given by the Picknett-Bexon correlation [23].Figure 2 shows the parameter G as a function of Bo for pendant drops on twelve substrates, with contact angle ranging from 20° up to 110°.For a fixed Bond number, the parameter G increases with the contact angle.For a fixed contact angle θc the parameter G reaches a minimum value at a specific Bond number, depending on the substrate wettability.On hydrophobic substrates, pendant drops deformed by gravity show higher evaporation rate than the corresponding iso-volumic spherical caps.For hydrophilic substrates, the minimum value of G is obtained with deformed pendant drops.The smaller the contact angle, the larger the value of Bo corresponding to the minimum value of G. Within the range of the selected operating conditions, neglecting the effect of deformation can lead to an underestimation of the evaporation rate from pendant drops on hydrophobic substrates up to about 7%, while for hydrophilic substrates the overestimation can reach about 16%.Validation of the numerical method was done on drop geometries where analytical solutions to the problem (3) exist [20], and the agreement on the values of the evaporation rates was of the order of 10⁻³ for hydrophilic substrates and 10⁻⁴ for hydrophobic substrates.This numerical tool was then used to build up a database for the quantity n A dA ∇ Φ  that will be used to solve the mass and energy balance on the pendant drop, as described in the next subsection.

The mass and energy balance on the drop
The mass and energy conservation equations for the pendant drop, neglecting the non-uniformity of the inner temperature distribution, are ev dM m dt = − ɺ ; ( ) where M is the drop mass, cL is the liquid specific heat, h'LV is the modified latent heat of vaporisation [32,33], ev m ɺ is the evaporation rate, s Q ɺ is the sensible heat rate exchanged through the drop/gas interface and c Q ɺ is that exchanged with the substrate.s Q ɺ and ev m ɺ are calculated from the quantity n A dA ∇ Φ  (see Appendix), which can be found from the above-mentioned database by a bivariate interpolation, given the instantaneous value of Bo and θc.An implicit Euler method was used to numerically solve equations (6).

Results and discussion
Time-dependent evaporation of water pendant and sessile drops is analysed, to point out the effect of drop deformation on the evaporation characteristics.The two main evaporation modes of sessile and pendant drops are the so-called constant contact radius (CCR) and constant contact angle (CCA).The CCR mode, often observed for hydrophilic substrates, is characterised by the pinning of the triple line and the constancy of the drop base radius.The CCA is often observed for hydrophobic substrates, the contact angle remains constant while the drop shrinks and the triple line slides on the substrate surface.However, when a drop evaporates on real substrates, both modes may appear and also a combination of the two is observed.With the aim of evidencing the sole effect of the drop shape on the evaporation, simplified scenarios will be used, assuming that the drop evaporates in one of the two modes during the whole observation period.One hydrophobic (θc=100°) and one hydrophilic (θc=50°) substrates are considered, and pendant, sessile and spherical-cap shaped drops are assumed to evaporate under both modes.Figure 3 shows the shape and the stability map (see [12] for a clear explanation) of sessile and pendant drops, respectively, and the paths followed by drops on hydrophobic and hydrophilic substrates are reported for both evaporation modes (thick black lines).
The initial Bond number is 0.6 for the hydrophobic substrate and 1.6 for the hydrophilic substrate.When the CCA evaporation mode is set, the evaporation process is recorded until the Bond number reduces to 0.2, while with the CCR evaporation mode the process is analysed until the contact angle reduces to 20° for the hydrophilic substrate (θc,i=50°) and until the Bond number reduces to 0.25 for the hydrophobic substrate (θc,i=100°).The evaporative cooling effect is considered, setting the wall thermal conductivity equal to 1.05 W/mK (glass substrate) and it is further analysed later in this section.The initial drop temperature is fixed equal to 20°C in equilibrium with the wall and the surrounding gas.
Figure 4 shows, for the CCA evaporation mode and θc,i=50°, the transient profiles of the nondimensional equivalent drop radius squared versus the non-dimensional time, , where D10 is the mass diffusion coefficient of the vaporizing species in air, calculated at the initial conditions, for sessile, pendant and a spherical cap drops, with the same initial volume.The evaporation process follows the well-known D 2 -law for all the drop shapes (with the exception of the short initial evaporative cooling period), with the sessile drop vaporizing faster.Figure 4 reports the corresponding profiles of the drop temperature, which evidences the initial cooling (about 0.5 degrees) followed by an almost constant temperature period.The time to reach Bo=0.2 is 3.3% and about 7% greater for the spherical and the pendant drops, respectively.For the CCA evaporation mode and θc,i=100° the differences in the time needed to reach the Bond number equal to 0.2 differ by less than 0.9%.Figure 5a shows the transient profiles of the heat rates as function of τ for the three drop shapes vaporizing in CCR mode on a hydrophilic substrate.When the drop contact angle reduces to 20°, the final Bond number is equal to 0.874, 0.673 and 0.8 for the sessile, the pendant and the spherical drops, respectively.For the CCR mode too, the pendant drop vaporizes slower, always experiencing smaller evaporation and heat rates.Differently from the CCA cases, the evaporation and heat rates of the spherical-cap drop are initially lower than those of the sessile drop; at about τ =5000, the spherical-cap drop starts to show higher evaporation and heat rates.The time to reach Bo=0.9 is about 12% longer for the pendant drop compared to the sessile drop, and about 0.7% shorter for the spherical drop.The drop temperature varies in a range of less than 0.5 degrees for all the three drop shapes (not shown here).On a hydrophobic substrate (θc,i=100°), the effect of drop shape is rather limited, as for the corresponding CCA case.When the Bond number reduces to 0.25, θc is equal to 40.7°, 67.9° and 47.8° for the sessile, the pendant and the spherical drops, respectively.The lifetime to evaporate up to this point is 1.5% and 0.2% higher for the pendant and the spherical drops, respectively, compared to the sessile one.The effect of wall thermal conductivity on the cooling and evaporation characteristics for the CCR mode and θc=50° is briefly analysed, for three values of the wall thermal conductivity, namely 0.033 W/mK (polyurethane), 1.05 W/mK (glass) and 12.1 W/mK (stainless steel).Figure 5b shows the profiles of drop temperature as function of the non-dimensional time within the initial evaporative cooling period.The cooling effect is important for the polyurethane substrate, with drop temperature reducing of about 4.4°C.On more conductive walls the cooling effect is less evident and the evaporation becomes faster.The effect of drop shape is similar for all the scenarios with the three selected thermal conductivities.The larger the thermal conductivity, the smaller the time when the spherical-cap shaped  ( )( ) , where y (1) is the vapour molar fraction and ( ) ( ) ( ) where K₀, H₀, Φ₀ and Φ₁ are constants to be obtained from the boundary conditions, ( ) ( ) The heat and evaporation rates are then found by integrating the normal component of the fluxes over the drop surface ( ) 40th UIT International Heat Transfer Conference (UIT 2023) Journal of Physics: Conference Series 2685 (2024) 012065 IOP Publishing doi:10.1088/1742-6596/2685/1/0120653
surface (see Appendix), and also on a spherical surface containing the drop, as an additional check.To quantify the effect of drop deformation on evaporation, a parameter G rate from a deformed drop with a given Bond number and

Figure 2 .
Figure 2. Schematic of the numerical domain and sample of the numerical mesh (size not to scale) (left).Non-dimensional evaporation rate G (equation 5) versus the Bond number for pendant drops on twelve hydrophilic and hydrophobic substrates (right).

Figure 3 .
Figure 3.Stability map for pendant drops (left) and shape map for sessile drops (right).The top grey curve in the graph on the left is the stability limit for pendant drops.The dashed lines are iso-Bo, the continuous lines are iso-θc, the thick black curves show the paths of CCA and CCR evaporating drops.The meaning of the geometrical parameters h and Rc is explained in figure2.

Figure 4 .Figure 5 .
Figure 4. Non-dimensional drop radius squared (left) and drop temperature (right) versus the nondimensional time, for a sessile, a pendant and a spherical cap drop, on a hydrophilic substrate.
is the free stream temperature and the pedex ref refers to an arbitrary reference condition.The mixture thermal conductivity is obtained by the Wassilieva[34] relation, and the exponents m, qp and b depend on the species involved.Under uniform Dirichlet boundary conditions on the drop surface and at infinity, the solution to the conservation equations for the vapour concentration and the temperature fields can be written in a closed form when the harmonic auxiliary function Φ (see equations 3) is known.Defining the function ( ) heat flux and the vapour mass flux on the drop surface can then be calculated in terms of the gradient of Φ as Heat Transfer Conference (UIT 2023) Journal of Physics: Conference Series 2685 (2024) 012065 IOP Publishing doi:10.1088/1742-6596/2685/1/0120658