Sanjeev Kumar

Marc Madele

David Brutin

Numerical model for sessile drop evaporation on heated substrate under microgravity

Keywords:

Sessile drop have simple geometry, but the physics involved in the evaporation process is complex due to the numerous intricate interactions and the fluid nature of the sessile drop. An accurate quantitative model of the evaporation process can lead to a greater understanding and control over it. The numerical model is developed for sessile drop evaporation on a heated substrate under microgravity, based on the rocket and parabolic experiments understanding to gauge the drop's internal dynamics. The quantitative validation of our model is performed via comparison with experiments. Then, we suggested a correlation for the evaporation rate of Ethanol sessile drop. Simultaneously, we are analysing the flow motion with 3D DNS computations of evaporating sessile drop. This provides an insight of the Marangoni effect in the dynamics of evaporation process and the occurrence of secondary instabilities. We are able to capture first time the fine influences of the secondary instabilities appearance on the evaporation rate. Our numerical model is valid in the absence of convection in vapour phase thus lead to interface of evaporating drop in under fully saturated vapour that typically exist in microgavity conditions.

Introduction

Despite the significant improvement in numerical solving capabilities in recent years, a number of fluid mechanics problems coupled with thermal aspects remain partially solved. Even for very simple geometries, such as a spherical cap drop anchored on a flat surface and maintained at a constant substrate temperature, the evaporation kinetics are not fully understood. In the process of the sessile drop evaporation, hydrodynamics within in the liquid domain and ambient of sessile drop play pivotal role in the process. Several experimental and numerical investigated have been conducted in the literature to characterize the internal flow for pure drop [START_REF] Girard | Evaporation and marangoni driven convection in small heated water droplets[END_REF][START_REF] Ye | Experimental investigation on evaporation dynamics of sessile ethanol droplets on a heated substrate[END_REF][START_REF] Zhang | Numerical investigation on flow instability of sessile ethanol droplets evaporating in its pure vapor at low pressure[END_REF][START_REF] Josyula | Insights into the evolution of the thermal field in evaporating sessile pure water drops[END_REF] and natural convection under gravity conditions [START_REF] Shahidzadeh-Bonn | Evaporating droplets[END_REF][START_REF] Navaz | Convective evaporation model of sessile droplets in a turbulent flow-comparison with wind tunnel data[END_REF][START_REF] Ljung | Evaporation of a sessile water droplet subjected to forced convection in humid environment[END_REF].

The steady-state diffusion model under predicting evaporation rate of pure ethanol drop of initial volume 5 µL at 𝑇 𝑠 = 60 • C with lifetime 12 ± 1 s [START_REF] Gurrala | Evaporation of ethanol-water sessile droplet of different compositions at an elevated substrate temperature[END_REF]), whereas Gurrala et al, 2019 [START_REF] Gurrala | Evaporation of ethanol-water sessile droplet of different compositions at an elevated substrate temperature[END_REF] reported that even though steady state diffusion model remains satisfactory for evaporation of water droplet of lifetime 190 ± 5 s with same initial volume and substrate temperature. Interestingly, Josyula et al, 2021 [START_REF] Josyula | Insights into the evolution of the thermal field in evaporating sessile pure water drops[END_REF] reported that water drop evaporation time with an initial volume of 5.6 ± 0.2 µL on a hydrophilic substrate with initial contact angle of 67 • maintained at 50 • C is over predicted by diffusion model.

The under prediction by diffusion model at substrate T s = 60 • C for pure ethanol for the experiment performed by [START_REF] Gurrala | Evaporation of ethanol-water sessile droplet of different compositions at an elevated substrate temperature[END_REF] can be understood by dimensionless heat equilibrium time however, for the water with an initial volume of 5.6 ± 0.2 µL with substrate temperature 50 • C enough though it is still under limit of the assumption of quasi-steady sanjeev0611@yahoo.com (S. KUMAR); marc.medale@univ-amu.fr (M. Madele); david.brutin@univ-amu.fr (D. BRUTIN) ORCID(s): 0000-0002-5884-0753 (S. KUMAR); 0000-0003-2661-883X (M. Madele); 0000-0002-0006-6822 (D. BRUTIN) evaporation i.e 𝑡 ℎ ∕𝑡 𝐹 < 0.1, the higher evaporation rate is the effect of natural convection under gravity conditions that are mostly not accounted in diffusion model [START_REF] Shahidzadeh-Bonn | Evaporating droplets[END_REF][START_REF] Navaz | Convective evaporation model of sessile droplets in a turbulent flow-comparison with wind tunnel data[END_REF].

The non-isothermal models which are recently developed such as Nguyen et al, 2018, [START_REF] Nguyen | Analytical model for diffusive evaporation of sessile droplets coupled with interfacial cooling effect[END_REF] analytical Model for diffusive evaporation of sessile drops coupled with interfacial cooling effect with the two-way coupling model of quasi-steady thermal diffusion within the droplet and quasi-steady diffusion-controlled droplet evaporation. This coupled problem of diffusive vapor and heat conduction also has been solved numerically by Xu and Ma [START_REF] Xu | Analysis of the effects of evaporative cooling on the evaporation of liquid droplets using a combined field approach[END_REF], and Dunn et al. [START_REF] Dunn | A mathematical model for the evaporation of a thin sessile liquid droplet: Comparison between experiment and theory[END_REF]. Pan et al. [START_REF] Pan | Assessment of water droplet evaporation mechanisms on hydrophobic and superhydrophobic substrates[END_REF] also carried out numerical modelling to study the competing effects of external natural convection and evaporating cooling. It needs to be considered that the above models are still not applicable for heated substrate.

Similarly, as we increase the substrate temperature, the coexistence of tangential temperature gradient at the gasliquid interface and the normal temperature gradient inside the drop will induce various thermal convections, such as, thermocapillary flow and buoyancy convection. The intense and flow direction of thermal convections are influenced by liquid properties [START_REF] Sefiane | Self-excited hydrothermal waves in evaporating sessile drops[END_REF][START_REF] Mamalis | Nonisothermal spreading dynamics of self-rewetting droplets[END_REF], drop geometry [START_REF] Josyula | Insights into the evolution of the thermal field in evaporating sessile pure water drops[END_REF][START_REF] Chen | Marangoni flow induced evaporation enhancement on binary sessile drops[END_REF], and other factors such as natural convection in the ambient. The heated substrate also changes the energy transfer during evaporation thus the evaporation rate [START_REF] Mahmud | Experimental investigation of interfacial energy transport in an evaporating sessile droplet for evaporative cooling applications[END_REF][START_REF] Ghasemi | Energy transport by thermocapillary convection during sessile-water-droplet evaporation[END_REF].

Recently, Ye et al [START_REF] Ye | Experimental investigation on evaporation dynamics of sessile ethanol droplets on a heated substrate[END_REF] experimentally investigated evaporation dynamics of sessile ethanol drops on a heated substrate with contact radius of the drop on substrate is fixed at 2.5 mm and the droplet height varies from 0.4 to 1.2 mm. The experiment indicates that with the decrease of the drop height, the temperature distribution near the drop centre becomes non-uniform, evaporation rate first decreases, and then increases with the decrease of drop height.

In parallel, Sefiane et al. [START_REF] Sefiane | Self-excited hydrothermal waves in evaporating sessile drops[END_REF] reported secondary instability in sessile drop evaporating under ambient conditions and later Brutin et al. [START_REF] Brutin | Infrared visualization of thermal motion inside a sessile drop deposited onto a heated surface[END_REF] presented a series of experiments on heated droplets exhibiting similar patterns. Multiple researcher reported similar instability under various conditions experimentally [START_REF] Ghasemi | Energy transport by thermocapillary convection during sessile-water-droplet evaporation[END_REF][START_REF] Wang | Transition of marangoni convection instability patterns during evaporation of sessile droplet at constant contact line mode[END_REF][START_REF] Shi | Marangoni convection instability in a sessile droplet with low volatility on heated substrate[END_REF][START_REF] Sobac | Thermocapillary instabilities in an evaporating drop deposited onto a heated substrate[END_REF] and numerical simulation are performed to resolved the secondary instability [START_REF] Zang | Evaporation of a droplet: From physics to applications[END_REF][START_REF] Semenov | 3d unsteady computations of evaporative instabilities in a sessile drop of ethanol on a heated substrate[END_REF][START_REF] Zhu | Bénard-marangoni instability in sessile droplet evaporating at constant contact angle mode on heated substrate[END_REF][START_REF] Gavrilina | Modeling unsteady bénard-marangoni instabilities in drying volatile droplets on a heated substrate[END_REF].

However, Karapetsas et al. [START_REF] Karapetsas | Convective rolls and hydrothermal waves in evaporating sessile drops[END_REF] first to attempt numerically to demonstrate the physical mechanisms that drive instability and pattern formation with "one-sided" twodimensional (2D) numerical model and using linear stability analysis. Another closely related "one-sided" 3D numerical modeling study conducted by Semenov et al. [START_REF] Semenov | 3d unsteady computations of evaporative instabilities in a sessile drop of ethanol on a heated substrate[END_REF] to observed the flow in drops. Semenov et al. with "one-sided" 3D numerical model provide no direct comparison of the sessile drop evaporation rate, but only qualitative comparison with the top view of experiment drop. Semenov et al. concluded with the model don't use any fitting parameters, however, equation 9 certainly use fitting parameter namely 𝐵 and Δℎ that were previously introduced [START_REF] Semenov | Boundary conditions for a one-sided numerical model of evaporative instabilities in sessile drops of ethanol on heated substrates[END_REF] and required 2D simulation to acquired the values. Semenov et al. [START_REF] Semenov | 3d unsteady computations of evaporative instabilities in a sessile drop of ethanol on a heated substrate[END_REF] "onesided" 3D numerical model due to fitting and some other inconsistent is restricted to achieve its objective. However, Karapetsas et al. [START_REF] Karapetsas | Convective rolls and hydrothermal waves in evaporating sessile drops[END_REF] and Semenov et al. [START_REF] Semenov | 3d unsteady computations of evaporative instabilities in a sessile drop of ethanol on a heated substrate[END_REF] model remains fundamental and inspiring numerical model.

The vapour distribution around a sessile droplet was recently studied using experimental [START_REF] Kelly-Zion | Vapor distribution above an evaporating sessile drop[END_REF] and numerical [START_REF] Pan | Assessment of water droplet evaporation mechanisms on hydrophobic and superhydrophobic substrates[END_REF][START_REF] Saada | Evaporation of a sessile drop with pinned or receding contact line on a substrate with different thermophysical properties[END_REF] approaches for purely diffusive evaporation. Kelly-Zion et al. [START_REF] Kelly-Zion | Vapor distribution above an evaporating sessile drop[END_REF] have recently shown using infrared spectroscopy and computed tomography that vapour emitted by sessile drops at room temperature behaves differently compared to the commonly accepted diffusion-limited model. Models taking into account convection have been developed to predict evaporation time and obtain the best reagent concentration for combustion [START_REF] Bertoli | A finite conductivity model for diesel spray evaporation computations[END_REF][START_REF] Sazhin | Advanced models of fuel droplet heating and evaporation[END_REF]. The most commonly used model, the Spalding evaporation model, takes into account convection using the calculation of mass and heat balance. This evaporation could be more complicated in presence of external forces such as electric field and force convection [START_REF] Navaz | Convective evaporation model of sessile droplets in a turbulent flow-comparison with wind tunnel data[END_REF][START_REF] Bin | Evaporation of methanol droplet on the teflon surface under different air velocities[END_REF] especially with heat substrate. Clearly, all cases mentioned above do not favour producing a neat and reliable benchmark case, and therefore, in ambient, however, complexity can be reduced by removing the gravity component under microgravity experiment.

As expected evaporation under microgravity is lower than terrestrial gravity [START_REF] Carle | Hydrothermal waves on ethanol droplets evaporating under terrestrial and reduced gravity levels[END_REF][START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF] and the deviation continue to increase with the drop liquid volatility and substrate temperature as observed in the case of Ethanol and HFE7100. Also, the comparison of the average of the evaporation rate in microgravity (numerical) and terrestrial gravity (experiments) for water at substrate temperature 𝑇 = 313 K, evaporation rate under microgravity is 8.85% lower than the under terrestrial gravity condition but difference increase further up to 18.14% at 𝑇 = 343 K [START_REF] Sanjeev | Numerical and experimental studies of thermo-convective instabilities in sessile drop in the course of evaporation[END_REF][START_REF] Sáenz | Evaporation of sessile drops: a three-dimensional approach[END_REF]. Thus, to obtain correct evaporation rate under gravity condition, a nonsaturated vapour pressure boundary condition is required at the sessile drop interface for one-way coupled diffusive model. Alternative, two-way coupled model to account for the nonlinear coupling between the vapour phase and the liquid one. Indeed in the gravity field, one has to account for the convection in the gas phase induce by thermalsolutal gradients. Our validation of saturated vapor pressure boundary condition under microgravity with current model is an opportunity for reconsideration of vapor pressure influence on interface under gravity conditions for the one-sided diffusive models.

A microgravity experiment not only helps to classify the effect of natural convection but also the absence of gravity typically decreases the computation times by an order of magnitude, which can be made possible to compare with expensive 3D computations and parametric studies at a reasonable cost. So, we proposed new numerical model for microgravity conditions and validated it with experiment data from the sounding rocket and quantitatively compare with parabolic flight, for more detail on experiments refer to results published by Kumar et al.,2020 [34] on sounding rocket experiments and Carle et al. 2013 [START_REF] Carle | Experimental evidence of the atmospheric convective transport contribution to sessile droplet evaporation[END_REF] on parabolic flight experiments. In this communication we aims to numerically investigate the thermo-hydrodynamics of sessile drop evaporation by focusing on evaporation dynamics of the drops on heated substrates, the flow motion inside volatile drops, and the effect of secondary Marangoni instability by leveraging microgravity conditions.

Theoretical Basis

The evaporation under microgravity conditions without heated substrate is primarily quasi-steady diffusive in nature [START_REF] Carle | Experimental evidence of the atmospheric convective transport contribution to sessile droplet evaporation[END_REF]. Evaporation of sessile drop in microgravity is assumed under saturated vapour pressure as evaporation only takes place on the surface of the liquid, and the surface of the liquid is under saturated vapour density due to the collection of vapour around the sessile drop in the absence of gravity. So the liquid interface is saturated with vapour particles, equilibrium occurs on the surface of the liquid. However, due to heated substrate the forward change (liquid to vapour) that is endothermic can be sustained and heat needed to convert the liquid into the vapour is available at the interface. The average energy of the particles leaving from liquid interface is governed by the interface temperature that is sustained by the Marangoni convection.

The validity of proposition that sessile drop is under saturated vapour pressure in microgravity during evaporation put question on the consideration of various diffusion models for evaporation of sessile drop in gravity. As some diffusion models (mostly one-way coupled [START_REF] Gurrala | Evaporation of ethanol-water sessile droplet of different compositions at an elevated substrate temperature[END_REF][START_REF] Xu | Analysis of the effects of evaporative cooling on the evaporation of liquid droplets using a combined field approach[END_REF][START_REF] Dunn | A mathematical model for the evaporation of a thin sessile liquid droplet: Comparison between experiment and theory[END_REF][START_REF] Chen | Marangoni flow induced evaporation enhancement on binary sessile drops[END_REF][START_REF] Semenov | Evaporation of sessile water droplets: Universal behaviour in presence of contact angle hysteresis[END_REF][START_REF] Baghel | Modeling of heat transfer through a liquid droplet[END_REF][START_REF] Hu | Effect of liquid-vapor interface area on the evaporation rate of small sessile droplets[END_REF]) also considered evaporation of sessile drop under saturated vapour pressure with gravity. However under same conditions the evaporation rate under gravity is higher compared to microgravity [START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF][START_REF] Carle | Experimental evidence of the atmospheric convective transport contribution to sessile droplet evaporation[END_REF]. In gravity, natural convection influence (reduces) the vapour pressure, so vapour pressure should be not be considered equal to saturated vapour pressure. That the reason, diffusive models with saturated vapour pressure assumption in gravity predicting lower evaporation rate compare to experiments even for the non-heated substrate at the ambient conditions.

Governing Equations

We are considering sessile drop evaporation under microgravity conditions. The sessile drop is a spherical cap with a constant contact area with substrate (pinned mode). During the evaporation process, sessile drop evaporation rate is mainly controlled by the vapour diffusion in the surrounding atmosphere. The energy transfer from substrate to the interface of sessile drop under microgravity conditions is considered mainly due to Marangoni convection [START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF]. Marangoni convection occurs inside the sessile drop as an effect of surface tension gradient due to its temperature dependency. We neglected the energy transfer into gas phase via radiative and conductive mode. Due to the absence of gravity in space, there is no convection in vapour phase that plays an important role in vapor for mass and energy transport that leads to higher evaporation rate on the Earth conditions [START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF][START_REF] Carle | Experimental evidence of the atmospheric convective transport contribution to sessile droplet evaporation[END_REF] and we can also neglect the Stefan flow in microgravity conditions [START_REF] Carle | Experimental evidence of the atmospheric convective transport contribution to sessile droplet evaporation[END_REF][START_REF] Carle | Droplet Wetting and Evaporation[END_REF]. So, most importantly, we assume in microgravity conditions, sessile drop is evaporating under fully saturated vapour pressure condition based on local interface temperature, i.e interface of the sessile drop is fully covered with vapour. In summary, we consider the fully coupled fluid flow-heat transfer problem inside the domains and at the interface we successfully managed to accommodate the effect of vapor (ambient conditions) via our surface boundary conditions. Surface boundary conditions are only valid in the absence of external forces such as an external electric field that could affect the vapour transportation around the sessile drop interface thus the evaporation rate [START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF]. For surface boundary conditions we are going to take help from Pickentt & Bexon's analytical equation [START_REF] Picknett | The evaporation of sessile or pendant drops in still air[END_REF] which is mainly used for the diffusion base evaporation at the ambient temperature for capillary length 𝐿 𝑐 < √ 𝜎∕(𝜌𝑔), to translate into the local evaporation rate at the sessile drop interface. The local evaporation rate would later be used to calculate and define the global evaporation rate (𝐸 𝑣 ) for the sessile drop on the heated substrate under microgravity conditions.

Geometry and environment conditions

The model considers sessile drop as a three -dimensional spherical cap with a constant base radius, (𝑅 = 𝑐𝑜𝑛𝑠𝑡𝑎𝑛𝑡) (i.e CR mode) which is also known as the pinned sessile drop, and moving interface (Γ) with substrate over time as volume changes during the evaporation process. Figure 1 shows the details for sessile drop geometries. This is valid for the 90% of the drop lifetime of the evaporation process. Substrate is considered isothermal (heated substrate) with constant substrate temperature (𝑇 sub ).

Numerical model formulation

The model solves the in-compressible Navier-Stokes equations for conservation of momentum and continuity equations for conservation of mass in the drop liquid domain along with heat transfer equations for single-phase fluid flow which in their most general form are:

𝜌 ( 𝜕 u 𝜕 𝑡 | | | |𝑓 + [( u -u 𝑚 ) ⋅𝛁 ] u ) = 𝛁 ⋅ σ, (1)
𝛁⋅u = 0, (2)
𝜌𝑐 𝑝 ( 𝜕 𝑇 𝜕 𝑡 | | | |𝑓 + (u -u 𝑚 ) ⋅𝛁 𝑇 ) + 𝛁 ⋅ 𝐪 = 0, (3)
𝐪 = -𝑘∇𝑇 (4)
Here, u is the velocity vector (m ⋅ s -1 ), u 𝑚 mesh velocity (m ⋅ s -1 ), 𝑇 is the temperature field (K), and 𝜌(𝑇 ), 𝑐 𝑝 (𝑇 ) and 𝑘(𝑇 ) are temperature dependent density ( kg ⋅ m -3 ) , specific heat capacity at constant pressure (J⋅(kg ⋅ K) -1 ) and thermal conductivity (W⋅(m ⋅ K) -1 ) of the drop respectively, subscript 𝑓 means that time derivative is taken at a fixed mesh node (fixed mesh coordinates), total stress tensor is σ = -𝑝I + π where 𝑝 is hydrodynamic pressure (Pa), I is the identity matrix and π = 𝜇 ( 𝛁u + (𝛁u) T ) is viscous stress tensor (Pa), 𝜇 dynamic viscosity, (Pa ⋅ s), q is conductive heat flux, that satisfies the Fourier's law of heat conduction. k is the thermal conductivity (W ⋅ (m ⋅ K) -1 ).

Initial conditions

Since, the physical problem of sessile drop that we are solving is an unsteady, non-isothermal evaporation of pinned sessile drop, so, the initial conditions are as follows: at time 𝑡 = 0, the drop is a spherical cap of contact angle 𝜃 𝑡 = 𝜃 0 with substrate and base radius 𝑅; sessile drop initial temperature 𝑇 | 𝑡=0 = 𝑇 𝑠𝑢𝑏 , the fluid flow inside sessile drop is at rest ( u| 𝑡=0 = 0). The initial pressure filed within the drop is set to the Laplace pressure (Δ𝑝) plus the atmospheric pressure (𝑃 𝑎𝑡𝑚 = 𝑃 ), i.e 𝑝| 𝑡=0 = 𝑃 𝑎𝑡𝑚 + 2𝛾 sin 𝜃 0 ∕𝑅, where 𝛾 is the liquid-gas interfacial tension and 𝜃 0 is initial contact angle. The initial temperature in the whole domain is set to the substrate temperature (𝑇 𝑠𝑢𝑏 ).

Surface boundary conditions

The above governing equations are subject to boundary conditions along liquid-gas interface and the substrate. So, now we layout details about the constraints applied to the sessile drop along liquid-gas interface and the substrate. In our work, the substrate is considered isothermal at 𝑇 𝑠𝑢𝑏 . The fluid flow boundary condition at the substrate is no slip. Now, let us consider the interface (Free-surface boundary) conditions Γ at the sessile drop moving interface. To begin with, the temperature at the interface is non-uniform based on the latent heat of evaporation, provided that local evaporation rate is known (see eq. 8), one can define the boundary conditions at Γ as follows:

-𝑘𝛁𝑇 ⋅ n = 𝐽 𝑒𝑣𝑎𝑝 Γ 𝐻 𝑒𝑣𝑎𝑝 𝜋(𝑅 2 + 𝐻 2 𝑡 ) , ( 5
)
𝜌 ( u ⋅ n -𝑢 Γ ) = 𝐽 𝑒𝑣𝑎𝑝 Γ 𝜋(𝑅 2 + 𝐻 2 𝑡 ) , (6)
In equation 5, n is the unit normal vector, 𝐽 𝑒𝑣𝑎𝑝 Γ is local evaporation rate, 𝐻 𝑒𝑣𝑎𝑝 is latent heat of evaporation, and 𝐻 𝑡 is the height of sessile drop at time. In equation 6 conservation of mass across the interface Γ imposes the velocity for deforming mesh. So, the moving interface Γ is considered as a boundary, moving along normal vector n with velocity 𝑢 Γ . Local interface velocity 𝑢 Γ and local evaporation rate 𝐽 𝑒𝑣𝑎𝑝 Γ are defined later. Similarly, the stress balance at the moving liquid interface Γ is as follows:

σ⋅n ⏟ ⏟ ⏟ Total stress = -𝛾 (𝛁⋅n) Γ n ⏟⏞⏞⏞⏞ ⏟⏞⏞⏞⏞ ⏟ Normal stress + 𝑑𝛾 𝑑𝑇 𝛁 Γ 𝑇 ⏟⏞⏟⏞⏟ Tangential stress , (7)
The normal and tangential stresses at the liquid boundary should be balanced across and along the interface, where (𝛁⋅n) Γ is the divergence of vector n at Γ, that is curvature of Γ, which is equal to 2 sin 𝜃∕𝑅, 𝛁 Γ 𝑇 is the surface gradient of temperature at the interface Γ and 𝑑𝛾 𝑑𝑇 is temperature coefficient of surface tension.

Evaporation rate

In order to compare evaporation rates of sessile drops measured in microgravity experiments, one can refer to the analytical model for evaporation limited by diffusion adapted from Picknett and Bexon [START_REF] Picknett | The evaporation of sessile or pendant drops in still air[END_REF] for constant contact area (up to de-pinning) and spherical cap shape. However, the analytical model derived by Picknett and Bexon is for the sessile drop in the ambient temperature on the surface (no heating). In order to adapt and to calculate the correct evaporation rate based on the interface temperature, we modified the equation 13 derived by Picknett and Bexon [START_REF] Picknett | The evaporation of sessile or pendant drops in still air[END_REF]. The modified equation can be read as follows:

𝐽 𝑒𝑣𝑎𝑝 Γ ( 𝑅, 𝑃 , 𝑇 𝑙,Γ , 𝜃 𝑡 ) = 2𝜋𝐷 𝑐𝑜𝑒𝑓 𝑓 (𝑃 , 𝑇 𝑙,Γ )𝐶 𝑠𝑎𝑡 (𝑇 𝑙,Γ )𝑅𝐹 (𝜃 𝑡 ) (8)
where, 𝐽 𝑒𝑣𝑎𝑝 Γ (

𝑅, 𝑃 , 𝑇 𝑙,Γ , 𝜃 𝑡 )

Γ is local evaporation rate on the Γ instead of the total evaporation rate derived by Picknett and Bexon, and 𝑥, 𝑦 and 𝑧 are the coordinates of computational domain in the direction of the 𝑥-axis, 𝑦axis, 𝑧-axis respectively, 𝑇 𝑙,Γ is the interface temperature at the mesh node in Kelvin, 𝑅 is the base radius, 𝐶 𝑠𝑎𝑡 is the saturated vapour concentration, and 𝐷 𝑐𝑜𝑒𝑓 𝑓 is the diffusion coefficient of liquid drop in ambient environment. 𝐷 𝑐𝑜𝑒𝑓 𝑓 is calculated according to the Fuller -Schettler -Giddings equation [START_REF] Fuller | Diffusion of halogenated hydrocarbons in helium. the effect of structure on collision cross sections[END_REF] as follows:

𝐷 𝑐𝑜𝑒𝑓 𝑓 (𝑃 , 𝑇 𝑙,Γ ) = 0.01013𝑇 1.75 𝑙,Γ ( 1 𝑀 𝑎 + 1 𝑀 𝑙 ) 0.5 𝑃 [ ( ∑ 𝑣 𝑎 ) 1∕3 + ( ∑ 𝑣 𝑙 ) 1∕3 ] 2 , (9)
where 𝑀 𝑎 and 𝑀 𝑙 are the molecular weights of the ambient fluid and liquid of sessile drop in g ⋅ mol -1 respectively; 𝑃 is the ambient pressure in Pa, and 𝑣 𝑎 and 𝑣 𝑙 are the atomic diffusion volumes.

The concentration of vapor at drop interface Γ, 𝐶 𝑠𝑎𝑡 can be calculated by assuming it as ideal gas and can be written as:

𝐶 𝑠𝑎𝑡 (𝑇 𝑙,Γ ) = 𝑃 𝑠𝑎𝑡 𝑀 𝑙 𝑅 𝑔𝑎𝑠 𝑇 𝑙,Γ (10)
In equation 10, 𝑅 𝑔𝑎𝑠 is the universal gas constant, 𝑃 𝑠𝑎𝑡 is the saturation vapour pressure based on local temperature 𝑇 𝑙,Γ , as, we assume that due to absence of the microgravity there is no convection transport in the vapour, 𝑃 𝑠𝑎𝑡 can be determined according to the Clausuis-Clapeyron equation, defined by

𝑃 𝑠𝑎𝑡 (𝑇 ) = 𝑃 𝑟𝑒𝑓 𝑙,𝑠𝑎𝑡 exp [ - (𝐻 𝑣𝑎𝑝 𝑀 𝑙 ) 𝑅 𝑔𝑎𝑠 ( 1 𝑇 𝑙,Γ - 1 𝑇 𝑟𝑒𝑓 𝑠𝑎𝑡 )] (11)
In the Clausuis-Clapeyron equation 𝑃 𝑟𝑒𝑓 𝑙,𝑠𝑎𝑡 is the reference saturation vapour pressure based on the reference saturation temperature 𝑇 𝑟𝑒𝑓 𝑠𝑎𝑡 .

𝐹 (𝜃

𝑡 ) = ⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩ ( 0.6366 𝜃 + 0.09591 𝜃 2 -0.06144 𝜃 3 ) ∕ sin 𝜃 for 0 ≤ 𝜃 < π∕18, ( 0.00008957 + 0.6333 𝜃 + 0.116 𝜃 2 -0.08878 𝜃 3 + 0.01033 𝜃 4 ) ∕ sin 𝜃 for π∕18 ≤ 𝜃 ≤ π. ( 12
)
𝐹 (𝜃 𝑡 ) is the function of contact angle with substrate derived by Picknett and Bexon [START_REF] Picknett | The evaporation of sessile or pendant drops in still air[END_REF] using the Snow's series solution [START_REF] Snow | Potential problems and capacitance for a conductor bounded by two intersecting spheres[END_REF]. Finally, the local evaporation rate 𝐽 evap Γ at Γ can be derived by using the equation 8 along with the equation 9 10, 11 and 12 :

𝐽 𝑒𝑣𝑎𝑝 Γ ( 𝑅, 𝑃 , 𝑇 𝑙,Γ , 𝜃 𝑡 ) = 2𝜋𝑅 0.01013𝑇 1.75 𝑙,Γ ( 1 𝑀 𝑎 + 1 𝑀 𝑙 ) 0.5 𝑃 [ ( ∑ 𝑣 𝑎 ) 1∕3 + ( ∑ 𝑣 𝑙 ) 1∕3 ] 2 × 𝑀 𝑙 𝑅 𝑔𝑎𝑠 𝑇 𝑙,Γ × 𝑃 𝑟𝑒𝑓 𝑙,sat exp [ - ( 𝐻 vap 𝑀 𝑙 ) 𝑅 ( 1 𝑇 𝑙,Γ - 1 𝑇 𝑟𝑒𝑓 sat )] 𝐹 (𝜃 𝑡 ) (13)
As 𝐽 𝑒𝑣𝑎𝑝 Γ is the local evaporation rate at the interface Γ, in other words, is the spatial dependent evaporation rate, however to validate and compare with experiment we define the average evaporation rate of sessile drop at any given 𝜃 (or time), that could be defined as:

- 𝑑𝑚 𝑑𝑡 = 𝐽 𝑒𝑣𝑎𝑝 𝑎𝑣𝑔,𝜃 = 1 ∫ Γ 𝑑Γ ∫ Γ 𝐽 𝑒𝑣𝑎𝑝 Γ ( 𝑅, 𝑃 , 𝑇 𝑙,Γ , 𝜃 𝑡 ) 𝑑Γ (14)
And then finally, the evaporation rate independent of base radius 𝑅 that we are going to use in our further discussion can be written as:

𝐸 𝑣 = - 1 𝑅 ⋅ 𝑑𝑚 𝑑𝑡 = 𝐽 𝑒𝑣𝑎𝑝 𝑎𝑣𝑔,𝜃 𝑅 ( 15
)
Finally the interface velocity 𝑢 Γ at Γ is derived based on average evaporation rate and that can be calculated by using the equation 14, 𝐽 𝑡𝑜𝑡 = 𝐽 𝑒𝑣𝑎𝑝 𝑎𝑣𝑔,𝜃 ; with two assumptions: pinned contact line (𝑅 = 𝑐𝑜𝑛𝑠𝑡) and spherical-cap shape of the drop in the course of evaporation [START_REF] Semenov | 3d unsteady computations of evaporative instabilities in a sessile drop of ethanol on a heated substrate[END_REF].

𝑢 Γ = - 𝐽 𝑡𝑜𝑡 π𝜌𝑅 3 ( 1 + cos 𝜃 𝑡 ) 2 sin 𝜃 𝑡 𝑧. ( 16
)
In an effort to obtain observable comparisons between numerical and experimental results, we post-processed the computed temperature field into simplified IR images. Computation is done by summing the IR radiation coming from substrate surface (assumed to be a grey body) with the integral of IR radiation sources distributed across the thickness of a semi-transparent drop [START_REF] Semenov | 3d unsteady computations of evaporative instabilities in a sessile drop of ethanol on a heated substrate[END_REF]. The equation 31 of the [START_REF] Semenov | 3d unsteady computations of evaporative instabilities in a sessile drop of ethanol on a heated substrate[END_REF] is expressed as:

𝑇 4 𝐼𝑅,𝑛𝑢𝑚 (𝑥, 𝑦) =𝑇 4 𝑠𝑢𝑏 𝑒 -𝑘 𝑆𝑊 𝑓 (𝑥,𝑦) + ∫ 𝑓 (𝑥,𝑦) 0 𝑇 4 (𝑥, 𝑦, 𝑧)𝑘 𝑆𝑊 𝑒 -𝑘 𝑆𝑊 (𝑓 (𝑥,𝑦)-𝑧) 𝑑𝑧 (17)
where, 𝑇 𝐼𝑅,𝑛𝑢𝑚 is IR numerical, 𝑘 𝑆𝑊 (mm -1 ) absorption coefficient, f(𝑥,𝑦) represents the local drop thickness and 𝑇 𝑠𝑢𝑏 is substrate temperature. The equation 31 of [START_REF] Semenov | 3d unsteady computations of evaporative instabilities in a sessile drop of ethanol on a heated substrate[END_REF] is corrected as the 𝜀 * is not the emissivity but it is IR camera correction.

Top View of the mesh.

Cut side view of the mesh in vertical median plane.

Numerical Modelling

Spatial and time discretizations

A simulation campaign was dedicated to determining the appropriate mesh that can provide numerical outputs independent of the mesh used. We tested various types of spatial discretizations, from fully unstructured mesh to mixed with an aim to achieve evaporation and Marangoni instabilities occurrence independent of grid sensitivity. Figure 3 & 2 presents the final configuration used for numerical simulations. Generated meshes are axisymmetric with the combination of hexahedra, prisms, quads and triangle elements. Mesh elements are controlled by defining the distribution of the elements in vertical, radial and circumferential directions, which are adjusted according to the size of the sessile drop, substrate temperature and fluid properties.

To perform the solution of our time-dependent problem we selected the implicit time-stepping method BDF at order two for solving the set of partial differential equations. BDF methods have been used for a long time and are known for their stability. The maximum time-step (Δ𝑡) is set to 1 s in case of primary Marangoni instability. However, during the secondary instability phase of sessile drop evaporation Δ𝑡 is between 10 -2 to 10 -3 s for Ethanol and in the range of 10 -3 to 10 -4 s for HFE7100.

Time-dependent moving interface problem where the sessile drop interface Γ is changing its shape due to evaporation, the Moving Mesh interface based on an Arbitrary Lagrangian-Eulerian (ALE) formulation is used to control the spatial frame at Γ. The Moving Mesh (ALE) interface is that which defines the deformation of the material frame relative to the geometry frame To model the sessile drop evaporation process we are using the mixed-order discretization for Laminar Flow, i.e P2 + P1 elements -that is, tri-quadratic interpolations for the velocity fields and piecewise linear interpolation for pressure. Discretization for the temperature field for heat transfer in fluids is tri-quadratic Lagrange. The geometry shape order is using the quadratic shape functions. The governing equations have been implemented and solved with the commercial COMSOL Multiphysics® software, which is based on the finite element method (FEM).

Solution strategies

At time 𝑡 = 0, the initial conditions are imposed on the computational domain (ref 2.4), the terms 𝐽 𝑒𝑣𝑎𝑝 Γ and 𝑢 Γ from equations 13 and 16 are respectively calculated. Computation then proceed to next iteration 𝑡 𝑖+1 = 𝑡 𝑖 + Δ𝑡 𝑖 , with the solution of equations 2 3 and 4 in the domain and equations 5, 6 and 7 at the interface Γ. The unknown variables 𝑢, 𝑣, 𝑤, 𝑝, and 𝑇 are calculated, then accordingly the local evaporation rate ( 13) is calculated in turn to define the interface and mesh velocities. During the solution procedure at every iteration, the local evaporation rate (13) depends only on the local temperature 𝑇 𝑙,Γ , so to calculate the global evaporation rate we need to perform some post-processing (refer equation 14).

In the cases of high substrate temperatures we introduce some ramping for the substrate temperature up to 𝑡 = 1 s that allows to increase substrate temperature in steps and prevent numerical oscillations inside the domain. In some cases, this could even effect the primary instability and secondary instability phase. Similarly in the case of highly volatile fluid (or low ambient pressure) such as HFE7100 we use another strategy that follows in three steps. Firstly, we perform only thermal diffusion for 𝑡 = 1 s to prevent numerical instabilities without moving the mesh. Secondly, the last results from the first step are used to initialize the coupled computations, i.e fluid flow and thermal equation, without moving mesh for 0.1 s that stabilizes the thermal distribution inside the drop. Finally, results from the latter are used as the initial solution for the standard simulation. With the above strategies, we can reach a stable primary instability phase that otherwise could effect the evaporation rate and subsequent evaporation process.

To achieve our objective we used our main model Coupled Fluid flow and Heat transfer Diffusion Model for microgravity conditions. In the following discussion we will refer to it as FHDM, however in some cases, for comparison, we also use the HDM that designates Heat transfer (by conduction only in the drop) Diffusion Model is just to understand how does HDM models perform without accounting the fluid flow for the sessile drop evaporation on the heated substrate with lower computational cost. The main difference between the FHDM and HDM comes from number of PDEs equation solved, FHDM solved the coupled fluid flow and heat transfer inside the sessile domain (equations 2, 3 and 4 inside the domain, equations 5, 6 and 7 at interface) and when in the same model, we switch off the fluid flow module and only solve heat transfer (equation 4 inside the domain, equation 5 at the interface) from heated substrate to sessile drop interface via conduction (is called HDM). Thus these modelling choices strongly affect the computational cost of the solution and accuracy of flow dynamic inside the sessile domain which in turn reflects on global results. However, in some cases, HDM method could be advantageous over the FHDM with some tolerance in the accuracy. Later, we will also discuss the HDM method's limitations and domains of its applicability.

Model validation 4.1. Physical properties

The thermo-physical properties considered for HFE7100 and Ethanol are presented in table 1. Universal gas constant 𝑅 𝑔𝑎𝑠 is 8.314 J ⋅ (mol ⋅ K) -1 , molar mass of ambient air is 28.97 g ⋅ mol -1 , atomic diffusion volumes of air is 19.7, molar mass of ambient nitrogen 28.014 g ⋅ mol -1 and atomic diffusion volumes nitrogen 18.5 [START_REF] Fuller | Diffusion of halogenated hydrocarbons in helium. the effect of structure on collision cross sections[END_REF]. The reference saturated vapour pressure 𝑃 𝑟𝑒𝑓 𝑙,sat is 26900 Pa and 7833.2 Pa, for HFE7100 and ethanol respectively at reference temperature

𝑇 𝑟𝑒𝑓 sat = 298.15 K [48].

Numerical validation with sounding rocket experiments HFE7100

For simulation, sessile drop's initial geometrical parameter, substrate and ambient conditions corresponding to drop 7DPµ𝑔 are from sounding rocket experiment [START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF], which are as follows: the drop base radius 𝑅 = 2 mm and initial contact angle 𝜃 0 = 60 • , whereas ambient gas is nitrogen at pressure 105400 Pa and substrate temperature 𝑇 sub = 301.15 K (refer [START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF] for more on sounding rocket experiment).

The domain is discretized by the hybrid mesh elements (90% elements Hexahedra elements, see fig. 3). The discretization of sessile drop domain consist 18 elements in radial direction (RD), 17 elements in vertical direction (VD) and 64 circumferential elements (CD)) that lead to generate total number of 21074 Mesh elements. The FHDM method uses 1240906 number of degrees of freedom (DOF) with quadratic geometry shape order whereas for the same shape order and mesh, HDM method solved 756530 DOF.

In figure 4, volume evolution with time compared between numerical simulation and experiment [START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF]. The experiment data extrapolated further to mark the expected critical volume for the appearance of secondary instability. More on the HFE7100 instabilities and extrapolation provided in the PhD of S. Kumar [START_REF] Sanjeev | Numerical and experimental studies of thermo-convective instabilities in sessile drop in the course of evaporation[END_REF]. The comparison between experiment and numerical estimated within 4.6% deviation at 𝑡 = 30 s from experiment, that also within the limit of error bar in the experiment (see fig 5). In figure 5, we compare evaporation rate from simulations, experiments and analytical model. The evaporation rate under earth conditions is nearly twice as compared to microgravity conditions due to convective transport in ambient and vapour phase [START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF][START_REF] Carle | Experimental evidence of the atmospheric convective transport contribution to sessile droplet evaporation[END_REF]. Evaporation HFE7100 [START_REF] Rausch | Density, surface tension, and kinematic viscosity of hydrofluoroethers hfe-7000, hfe-7100, hfe-7200, hfe-7300, and hfe-7500[END_REF]46] Ethanol [47, rate under microgravity experiments and simulations (both FHDM and HDM) are comparable which validates our assumption (and numerical model) for sessile drop evaporation under microgravity that it is purely diffusive in nature and sessile drop interface is fully under saturated vapour pressure (of HFE7100 vapour). This can be explained due to the absence of gravity and convective transport in vapour, HFE7100 vapour generated during evaporation, saturates ambient nitrogen near the sessile drop interface that results in the reduction of the evaporation rate of sessile drop in microgravity which also validates that gravity dominant effects are actually located in vapor phase (vapour density gradient induces thermal-solutal convection). So, it necessitates using a 2-sided model for gravity conditions. As expected analytical model predicts higher evaporation rate in comparison to microgravity experiments by considering sessile drop domain as isothermal, at substrate temperature for analytical model. It is also noteworthy that isothermal boundary conditions for substrate temperature 𝑇 𝑠𝑢𝑏 could be seen as a strong boundary condition as in experiments, the difference between substrate centre and edge temperatures is in the range of 0.15 -0.21 • C.

It is worth to mention computational cost of FHDM and HDM simulations for HFE7100 sessile drop on a heated substrate. The FHDM simulation with 1240906 DOF, total computational time on 48 cores cpu for 𝑡 = 75.3 s is 70 days, 17 hours! in comparison with 756530 DOF, HDM method on 32 cores cpu for 𝑡 = 86.6 s is 22 hours, 40 min only. The cost of simulation with FHDM is very high for resolving secondary instability i.e 57.5 to 75.3 s in 62 days, 19 hours. The computational cost for the FHDM is 73 times of HDM method for HFE7100 sessile drop due to its volatility, in same duration HDM under predicted 15% evaporation to the FHDM.

Results and Discussion

Quantitative comparison with Ethanol

The physical properties of Ethanol are presented in table 1. For simulation, sessile drop's initial geometrical parameter, substrate and ambient conditions are corresponding to drop parabolic flight experiment [START_REF] Carle | Experimental evidence of the atmospheric convective transport contribution to sessile droplet evaporation[END_REF] which are as follows: the sessile drop base radius 𝑅 is of different range from 1.5 to 3.6 mm, However for most of the comparisons the effect of diameter is neutralized by dividing the global evaporation rate by radius of sessile drop. The ethanol sessile drop evaporating in ambient air which is at cabin pressure

Microgravity Condition

Substrate Temperature,T sub (K) of elements in radial, vertical and circumferential direction changes with substrate temperature.

In figure 6, we presented quantitative comparison of Ethanol sessile drop evaporation rate (𝐸 𝑣 ) from numerical simulation with parabolic flight experiments and Pickentt & Bexon analytical equation [START_REF] Picknett | The evaporation of sessile or pendant drops in still air[END_REF]. For low substrate temperatures (comparable to ambient temperature), numerical results from FHDM, HDM and analytical model are comparable to experiments. As the temperature of substrate increased, evaporation rate increases following a power law and a divergence between analytical, FHDM and HDM evaporation rates emerges clearly. For comparison, values of evaporation rate is 4.61 × 10 -5 , 4.02 × 10 -5 and 5.29 × 10 -5 kg ⋅ (s ⋅ m) -1 for FHDM, HDM and analytical model respectively at the 𝑇 sub = 333.15 K, contact angle 𝜃 0 = 23 • and cabin pressure 𝑃 = 0.835 bar. As described by Carle et al, 2013 [START_REF] Carle | Experimental evidence of the atmospheric convective transport contribution to sessile droplet evaporation[END_REF], experimental results have a dispersion of data due to perturbations caused by the aircraft flight, for example, by vibrations. However, in our opinion lower values of experimental evaporation rates to numerical ones also result from the combination of two different factors. The first one is an indirect method (Heat flux method) used by Carle et al, 2013 [START_REF] Carle | Experimental evidence of the atmospheric convective transport contribution to sessile droplet evaporation[END_REF] to evaluate experimental evaporation rate and the second one is numerical modelling one assuming isothermal substrate boundary conditions.

However, as shown in figure 8 (top) divergence in evaporation rate is much broader for higher contact angles. Evaporation rate for FHDM at initial contact angle 𝜃 0 = 87 • is the highest contrary to HDM where it is the lowest during the lifetime of sessile drop with initial contact angle 𝜃 0 = 87 • . The divergence between the two not only reduces with decrease of contact angle with time, but also divergence reduced for low substrate temperature even for higher contact angles. This divergence represents the important role played by the Marangoni convection in transport of energy from substrate to sessile drop interface Γ thus in elevating interface temperature 𝑇 𝑙,Γ and evaporation rate. For more clarity and to follow full lifetime of sessile drop evaporation with FHDM and HDM at 𝑇 sub = 318 K and 323.15 K.

Figure 7, shows the plot of average interface temperatures 𝑇 𝑎𝑣𝑔,Γ with respect to substrate temperatures. The 𝑇 𝑎𝑣𝑔,Γ from HDM evolves with contact angles meanwhile 𝑇 𝑎𝑣𝑔,Γ from FHDM does not (see figure 8 (bottom)) up to secondary instability appear. So, according to the FHDM, 𝑇 𝑎𝑣𝑔,Γ is approximately independent of the contact angle and evaporation rate reduced with contact angle due to the reduction of the evaporating surface area Γ or in another words f(𝜃)). The insert of fig. 7) shows the obtained fitting of 𝑇 𝑎𝑣𝑔,Γ with respect to 𝑇 𝑠𝑢𝑏 , where,

𝑇 𝑎𝑣𝑔,Γ = 1 Γ ∫ Γ 𝑇 𝑙,Γ 𝑑Γ ( 18
)
for the Ethanol sessile drop, from fitting

𝑇 𝑎𝑣𝑔,Γ = 𝑎 𝐸 × 𝑇 𝑏 𝐸 𝑠𝑢𝑏 ( 19
)
with a E = 1.41148±0.047 and b E = 0.93887±0.0058. Now we can replace 𝑇 𝑙,Γ in the equation 13 with 𝑇 𝑎𝑣𝑔,Γ :

𝐸 𝑣,𝐸𝑡ℎ = 2𝜋 0.01013𝑇 1.75 𝑎𝑣𝑔,Γ ( 1 𝑀 𝑎 + 1 𝑀 𝑙 ) 0.5 𝑃 [ ( ∑ 𝑣 𝑎 ) 1∕3 + ( ∑ 𝑣 𝑙 ) 1∕3 ] 2 𝑀 𝑙 𝑅 𝑔𝑎𝑠 𝑇 𝑎𝑣𝑔,Γ × 𝑃 𝑟𝑒𝑓 𝑙,sat exp [ - ( 𝐻 vap 𝑀 𝑙 ) 𝑅 ( 1 𝑇 𝑎𝑣𝑔,Γ - 1 𝑇 𝑟𝑒𝑓 sat )] 𝐹 (𝜃) (20)
The correlation equation 20 is valid for Ethanol pinned sessile drop of spherical shape under microgravity conditions (under saturation vapour density on the interface Γ) on a heated substrate for any given pressure, with neglecting the effect of secondary Marangoni instability for the lower contact angle.

Evaporation of Ethanol sessile drops for two substrate temperature (whole cycle)

In figure 8, plotted the comparison of FHDM and HDM for the substrate temperatures 𝑇 sub = 318.15 K and 323.15 K. Both computations start with an initial contact angle of 𝜃 0 = 90 • and with a uniform initial temperature field (at substrate heating temperature) in the drop and an imposed isothermal heating temperature at the substrate surface. The liquid domain initial discretization consists hybrid mesh of 21312 elements. The computational domain consists of 18 elements in the radial direction (RD), 15 elements in the vertical direction (VD) and 64 circumferential elements (CD). The time discretization remains the same as discussed earlier.

For FHDM, the current simulation is performed in two stages. In the first stage simulation starting from 𝜃 = 90 • to 27 • , then the computational domain is re-meshed. In the second stage simulation started from 𝜃 = 33 • to 9.45 • . Then simulations were later merged to get a complete simulation for 𝜃 = 90 • to 9.45 • . For computation cost and resource used for the simulation refer to [START_REF] Sanjeev | Numerical and experimental studies of thermo-convective instabilities in sessile drop in the course of evaporation[END_REF].

By comparing the first and the second stage for FHDM, we can clearly see the high cost for simulation after the appearance of secondary instability which are further costlier for high substrate temperature. However, the computational Contact angle, 𝜃 (deg.) cost of HDM is far lower than the FHDM method but with its own disadvantages as inaccuracy is high to medium range contact angle. As we can see for 𝑇 𝑠𝑢𝑏 = 318.15 K, evaporation starting from 𝜃 = 90 • to 9.45 • , took 𝑡 = 132.8 s and 143.6 s for FHDM and HDM respectively. In figure 8, at 𝑡 = 0, computations start with a uniform initial temperature (at substrate heating temperature) in sessile drop liquid domain and imposed heating temperature at substrate surface. Thus, measured evaporation rate is the highest (equivalent to the analytical model) at 𝑡 = 0, due to latent heat of vaporization, drop cools down from interface; meanwhile, temperature near the contact line remains higher due to heat conduction from the substrate through a thin layer of Ethanol. This creates a vertical temperature gradient in drop bulk and a tangential one along the drop surface near the contact line. These temperature gradients promote the development of tangential stress that lead to development of Marangoni convection, after initial intense unsteady phenomena finished (at about 𝑡 = 5 s), then due to Marangoni convection, the drop reaches thermal equilibrium. The important role of Marangoni convection can be seen clearly in the figure 8 by comparing the evaporation rate of FHDM and HDM at contact angle 𝜃 = 80 • which is for FHDM, the evaporation rate is 𝐸 𝑣 = 3.1 × 10 -5 kg ⋅ (s ⋅ m) -1 corresponding average interface temperature is 𝑇 𝑎𝑣𝑔,Γ = 315.40 K, while in comparison, for HDM the evaporation rate and average interface temperature are only 2.1 × 10 -5 kg ⋅ (s ⋅ m) -1 and 307.25 K respectively. Similar trends followed by evaporation at 𝑇 𝑠𝑢𝑏 = 323.15 K (see fig [START_REF] Gurrala | Evaporation of ethanol-water sessile droplet of different compositions at an elevated substrate temperature[END_REF], comparison shows that evaporation rate for the HDM is nearly 32% lower than that of FHDM at contact angle 𝜃 = 80 • .

In figure 8 (top), FHDM has the highest evaporation rate initially, which starts reducing with contact angle, meanwhile, HDM has the lowest initial evaporation rate and then it starts increasing with contact angle. The difference could be understood from figure 8 (bottom), indeed with FHDM, Marangoni convection increases interface temperature which remains nearly constant until the secondary Marangoni instability appear, while with HDM the energy transfer is only by thermal diffusion in liquid drop, interface temperature increases as time goes on and contact angle reduces, narrowing the substrate to interface distance. The difference between FHDM and HDM evaporation rate reduces as time goes on and contact angle decreases. So, the higher is the interface temperature due to Marangoni convection, the higher is the supplied energy to account for latent heat of vaporisation resulting in higher evaporation rate.

The effect of secondary instability is shown in insert and circle of figure 8, where we can see transition from the base primary Marangoni instability (axisymmetric torus) to the secondary Marangoni instability characterized by fully 3D unsteady fluid flows. After this transition one can observe a dynamic multi-cellular thermo-convective pattern similar to shown in figure 9. The secondary Marangoni instability cells clearly appear for sessile drop radius 𝑅 = 1.5 mm and pressure 𝑃 = 0.835 bar at contact angle 𝜃 = 17.7 • and 19.2 • for 𝑇 𝑠𝑢𝑏 = 318.15 K and 323.15 K, respectively. The evaporation rate suddenly rises and becomes constant instead of decreasing while interface temperature starts to increase with the decrease in contact angle (see figure 8b).

So, the change in evaporation rate due to secondary instability are 2.49% and 2.53% for 𝑇 𝑠𝑢𝑏 = 318.15 K and 323.15 K, respectively. It shows that despite primary Marangoni instability has a strong effect on sessile drop evaporation under microgravity, the secondary Marangoni instability have comparatively a low effect on the evaporation rate of sessile drop. Similar results are seen for HFE7100 experiments, both on the ground as well as in microgravity experiments [START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF]. Figure 3 of [START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF] shows no change in the slope of the volume vs time plot after the secondary instability appear. They seem to play a minimal role in evaporation of pure liquid sessile drop, meanwhile, it plays important role in particle deposition [START_REF] Chatterjee | Evaporation of initially heated sessile droplets and the resultant dried colloidal deposits on substrates held at ambient temperature[END_REF].

The effect of the secondary Marangoni instability on evaporation rate is approximately 2 -3%, but the computational cost is approximately 6 to 8 times higher due to finer mesh and time step requirements. The cost can increase drastically with volatile fluids as in case of HFE7100, for which the computational cost for the same physical time is approximately 30 times higher. As in most of the cases secondary instability appear at low contact angle so it is advisable to use HDM in those latter cases. HDM is a good option for low temperature and low contact angle as in case of ethanol it could be used below 𝜃 = 20 • and 𝑇 𝑠𝑢𝑏 = 323 K [START_REF] Sanjeev | Numerical and experimental studies of thermo-convective instabilities in sessile drop in the course of evaporation[END_REF]. For comparison, the evaporation rate at 𝜃 = 19.2 • and 𝑇 𝑠𝑢𝑏 = 323 K is 𝐸 𝑣 = 3 × 10 -5 kg ⋅ (s ⋅ m) -1 and 2.8 × 10 -5 kg ⋅ (s ⋅ m) -1 for FHDM and HDM, respectively. So, HDM has approximately 7% lower evaporation rate compared to FHDM and this difference would reduce with decreasing contact angle.

Sessile drop dynamics and cells patterns evolution with time

For our analysis, we take a sessile drop of ethanol at substrate temperature 𝑇 𝑠𝑢𝑏 = 318.15 K, pressure 𝑃 = 0.835 bar and radius 𝑅 = 2 mm. The total time taken for evaporation from initial contact angle 𝜃 = 33 • up to 15 • is 𝑡 = 35.7 s. The average evaporation rate during the primary instability phase for sessile drop is approximately 𝐸 𝑣 = 2.4 ± 0.05 × 10 -5 kg ⋅ (s ⋅ m) -1 and average interface temperature is around 𝑇 𝑎𝑣𝑔,Γ = 316 K.

In figure 9, the illustration presents the evolution of the sessile drop (top view) with contact angle (left to right in the sequence of the appearance) and interface surface temperature, top view (IR) using equation 17, local evaporation rate using equation 13, surface velocity and local normal total flux (on the substrate, opposite to isothermal boundary surface layer). Similarly, figure 10 displays streamlines (colored with temperature) with velocity vector (back) in quarter geometry. Whereas figure 11 is top view of the vorticity vector field on the surface of sessile drop and figure 12 is temperature (top), velocity vector (middle) and velocity magnitude (bottom) at the cut plane (half), inside the domain of Ethanol sessile drops. In Figure 13 At 𝑡 = 0 s, the evaporation rate is maximum due to the effect of the initial boundary conditions, i.e initial drop temperature is equal to the substrate temperature. With time at 𝑡 = 2 s thermal energy is adjusted with fluid flow and thermal distribution in domain according to boundary conditions, and it arrives in equilibrium state. After that, from 𝜃 = 32 to 19 • ) a primary instability is developed having toroidal shape with axis at the center of sessile drop (see fig. 10a), with the lowest temperature around the apex, accordingly the lowest evaporation rate at the apex that increases along the radial position towards the contact line (see fig. 13, follow the legend 𝜃 = 20.6 • ). Concurrently, while fluid flowing from contact line toward apex at interface, velocity magnitude is the highest near surrounding of the apex and reduces towards the contact line similar pattern followed by normal total flux at substrate (see fig. 14a. In Figure 14a, velocity field 𝑥-component (𝑢), and velocity field 𝑧-component (𝑤) shows that flow moving toward the apex without any velocity field 𝑦-component (𝑣) representing no azimuthal velocity in flow which is also affirm by the vorticity fields in figure 11a. The cut plane in domain shows temperature and velocity profile corresponding to the primary instability phase (see fig. 10a, and figure 12a). The occurrence of secondary instability cells begin with the break of axis-symmetric torus roll at 𝑡 = 27.4 s (𝜃 = 19.4 • ), instability cells appear near the middle of sessile drop radius. Instability cells can be seen in figure 9 (a1) in the direction moving from substrate to sessile drop interface while moving towards the contact lines, having approximately 1 K higher temperature than the surrounding liquid (see fig. 13). Velocity field of sessile drop interface in fig 14b shows sudden change and fluctuation in the velocity field that corresponds to secondary instability cell. Figure 10b,c & d, also show that cells are moving in the opposite direction of main flow and it is confirmed by the velocity field 𝑢 and 𝑣 (see fig. 14b). Also, figure 11b shows increase in the In general, the velocity magnitude is lower in comparison to initial primary instability phase except around the instability cells (see fig. 14c). As seen in figures 9d2 and 10c, velocity is comparatively higher around the periphery of the cells. The velocity field 𝑥-component (𝑢), and velocity field 𝑧-component (𝑤) along with velocity field 𝑦-component (𝑣) in figure 14c shows that flow is moving towards the apex while the secondary instability cells are moving towards the contact line contrary to primary instability phase now velocity field 𝑦-component (𝑣) is also inducing the azimuthal velocity in the flow which is also affirmed by the vorticity fields in figure 11c. Interestingly, vorticity fields in figure 11b,c & d shows that the rotation of the flow within the secondary instability cells is higher and opposite to the sessile drop primary flow. As the secondary instability have minimal contribution in the evaporation so is the vorticity however, vorticity could be interesting for the particle transportation and deposition studies.

With time (or reducing contact angle) number of cells increased as we can see in figures 10d and the cells are also developed closer to the center of sessile drop. The secondary instability cells have an intrinsic wavelength that roughly depends on the local liquid thickness along the radial direction of the sessile drop. The secondary instability cells wavelength decreases along with the aspect ratio. The wave number strongly depends on the aspect ratio, substrate temperature and volatility of the fluid. However, evaporation rate nearly remains constant and average interface temperature continues to increase slowly (similar to fig. 8) whereas, Δ𝑇 continues to reduce with the increasing average interface temperature. The local evaporation rate along the radial position of sessile drop interface also elevates from apex to the position approximately at |𝑅| = 1.7 mm, after that local evaporation rate ( or temperature)remains the same as it was in primary instability phase (see fig. 13). Globally, we can see in the figure 9 that the local distribution pattern of evaporation rate is following the temperature profile where the heat flux from the substrate is following the velocity profile. We can also mention that the density of the number of cells neither heavily influences the evaporation rate nor the velocity of the sessile drop.

Conclusion

In microgravity conditions, the driving instability takes place in the liquid phase, therefore we have developed a one-sided model to account for the gas phase in the simple manner, thanks to ad-hoc boundary conditions. Despite the fully coupled fluid flow-heat transfer problem inside the domain and at the interface we successfully managed to accommodate the effect of vapor in microgravity, thus, the high spatial resolution in drop.

Quantitative validation of our model is performed by comparing it with the experimental results. Both, the numerical and the parabolic flight results over the range of substrate temperature for Ethanol sessile drop show that evaporation rate increases following a power law with substrate temperature. So, a correlation is proposed for evaporation rate which is valid for Ethanol pinned sessile drop of spherical shape under microgravity conditions (or saturation vapour density) on a heated substrate for any given pressure, with a neglecting effect of secondary instability for the lower contact angle.

The three relevant questions related to the dynamics of sessile drop evaporation that have been addressed are as following: i) Global and local evolution of evaporation rate and temperature. ii) Instability structure and its occurrence. iii) Secondary instability and its influence on evaporation.

The evaporation rate starts to first decrease with the contact angle, then it suddenly increases approximately by 2 -3% with the occurrence of the secondary Marangoni instability at the critical height, however, then it becomes constant for Ethanol drop but for HFE7100 it starts decreasing again. The important contribution of Marangoni convection is transport of energy from substrate to interface thus elevating interface temperature and evaporation rate. This elevated interface temperature remains nearly constant and is approximately independent of the contact angle until the secondary instability occurs. Thus, the evaporation rate reduces with contact angle due to the reduction of the evaporating surface area. The local distribution pattern of evaporation rate is following the temperature pattern whereas the heat flux from the substrate is following the velocity field over the interface. It can also be added that in reporting the evaporation rate one should express it in versus with initial and final contact angles instead of time. As, the evaporation rate at given contact angle is independent of the previous history of sessile drop.

Evaporation starts with a Marangoni convection (primary instability) that results in an axisymmetric torus roll, with the lowest temperature and accordingly the lowest evaporation rate around the apex, which increases along the radial position towards the contact line. The velocity magnitude is the highest along with the free surface surrounding the apex and then flow dives toward the substrate close to the symmetry axis. This torus roll appears as soon as it exists temperature gradient over the liquid drop surface.

Then, a transition to a secondary instability occurs once a critical aspect ratio is overcome and a fully 3D unsteady fluid flow is developed. The axisymmetric torus roll of the primary instability breaks into cells, core cells get too bigger, and split in turn. Then dynamic multicellular thermoconvective flower-like patterns emerge. The secondary instability cells have an intrinsic wavelength that roughly depends on the local liquid thickness along the radial direction of the sessile drop. The secondary instability cells wavelength decreases along with the aspect ratio. The wave number strongly depends on the aspect ratio, substrate temperature and volatility of the fluid. The radially outward moving instability cells create an illusion of the main flow moving from apex to contact line, however, velocity vectors confirm that instability cells are moving in the opposite direction of the main thermo-convective flow. Despite primary instability has a strong effect on sessile drop evaporation rate, secondary instability does not influence the evaporation rate heavily (merely 2 -3%), even though the evaporation rate and average interface temperature trend change with its occurrence.

Figure 1 :

 1 Figure 1: Schematic of the sessile drop. 𝑅 is the base radius of the pinned sessile drop, 𝑇 sub is the substrate temperature, 𝜃 𝑡 is the contact of sessile drop, 𝐻 𝑡 maximum height of the interface from the substrate, Γ is the interface of the sessile drop, 𝑇 𝑙,Γ (𝑥, 𝑦, 𝑧) is the temperature at the mesh node, 𝑢 Γ is the interface Γ velocity due to evaporation, and 𝐽 𝑒𝑣𝑎𝑝 Γ is the local evaporation rate.

Figure 2 :

 2 Figure 2: Spatial discretization: Meshing of the sessile drop.

Figure 3 :

 3 Figure 3: Meshed sessile drop of the HFE7100 with initial contact angle 𝜃 0 = 60 • and base radius 𝑅 = 2 mm.

Figure 6 :

 6 Figure 6: Experimental and numerical simulation comparison at different substrate temperatures. Numerical simulations of pure Ethanol evaporation rate with microgravity experimental results [37] at contact angle 𝜃 = 23 • and cabin pressure (𝑃 = 0.835 bar) for different substrate temperature.

Figure 8 :

 8 Figure 8: Comparison of FHDM and HDM for ethanol sessile drops at 𝑇 𝑠𝑢𝑏 = 318.15 and 323.15 K with base radius 𝑅 = 1.5 mm and pressure 𝑃 = 0.835 bar. a) Evaporation rate (top); b) Average interface temperature (bottom).

  , shown the interface temperatures along the radial position of sessile drop interface in X-Y plane at different contact angles. The divided illustrations in figure 10, 11 and 12 are representing the four phases of Ethanol sessile drop's evaporation: (a) at contact angle 𝜃 = 20.6 • : primary instability phase; (b) at contact angle 𝜃 = 19.4 • : transition begins from primary to secondary instability; (c) at contact angle 𝜃 = 19 • : complete appearance of secondary instability cells; (d) at contact angle 𝜃 = 15.4 • : fully developed flow with secondary Marangoni instability cells at the low contact angle.

Figure 9 :

 9 Figure 9: Ethanol sessile drops with surface temperature, IR temperature, local evaporation rate distribution, normal total flux (substrate) and velocity and at the substrate temperature 𝑇 𝑠𝑢𝑏 = 318.15 K, pressure 𝑃 = 0.835 bar and radius 𝑅 = 2 mm. (a-e) at contact angle 𝜃 = 20.6 • showing in primary instability phase; (a1 -e1) at contact angle 𝜃 = 19.4 • , starting of the secondary instability cells; (a2 & e2) at contact angle 𝜃 = 19 • fully developed stage of the secondary instability with flower like pattern; (a3 & e3) at contact angle 𝜃 = 15.4 • secondary instability core cells are getting bigger. Refer movies

  vorticity magnitude. Concurrently, figure 12 (b) clearly show the dynamics inside the cell in liquid domain, it also shows high temperature fluid moving from substrate to interface and going towards the contact line. While figures 9 a1, b1, c1, d1, & e1 and figure 10b show the top view of sessile drop with appearance of the instability cells. The non-axis-symmetric flower like pattern is developed corresponding to secondary instability phase at contact angle 𝜃 𝑐𝑟𝑖𝑡 = 19 • ( height 𝐻 𝑡 = 3.35 × 10 -4 m). The change in evaporation rate due to secondary instability is around 𝐸 Δ𝐽 = 2.25%, at critical aspect ratio (AR crit , 𝐻 𝑡 ∕𝑅) = 0.17. With secondary instability development the reducing evaporation rate with contact angle became constant and

igure 10 :Figure 11 :

 1011 Figure 11: Vorticity field vector on the surface of Ethanol sessile drop on heated substrate. Numerical simulation and subplot conditions are same as fig. 10. Readers are advised to refer to web version of this figure for better display.

Figure 13 :

 13 Figure 13: Comparison of interface temperature 𝑇 Γ along the radial position of sessile drop interface. The smooth profile of interface temperature at 𝜃 = 20.6 • corresponding to primary instability phase and rest remaining interface temperature along the radial position with fluctuation represent the secondary instability. Numerical condition are same as fig. 9.

Figure 14 :

 14 Figure 14: Velocity magnitude and velocity field component along the radial position of the sessile drop. Numerical parameters are same as fig. 9. Refer movie.

Figure 15 :

 15 Figure 15: The probe point at initial location (𝑥 = 0.75 mm, 𝑦 = 0 mm, 𝑧 = 0.51 mm) on the sessile drop interface. (Top) Evaporation rate and temperature, (Bottom) Velocity magnitude and velocity field component. Numerical parameters are same as Fig. 9.

Table 1

 1 Physical properties of pure liquid HFE7100, and Ethanol at 298.15 K and 1 atm.

	48]

  Numerical simulation for HFE7100 on heated substrate under microgravity condition. Comparison of volume evolution with respect to time between numerical simulation and experiment (with experiment[START_REF] Kumar | Sessile volatile drop evaporation under microgravity[END_REF]). Numerical simulation time shifted according to match with experiment maximum initial volume. The vertical bar on plot at 𝑡 = 65.2 s and Vol.= 2.2 µL marks the appearance of secondary instability. For time evolution refer to the Supplementary Materials for complete movies .

	9
	8						HFE-7100 Microgravity
	7											Experiment (Kumar et al, 2020) Experiment extrapolated
	6											Numerical simulation
	5
	4
	3												(65.2, 2.2)
	2
	1
	0 0	6	12	18	24	30 35 Time(s)	38	44	50	56	62	68	75	78
	0.10 0.15 0.20 Figure 4: 0.198								Exp. data (Carle et al., 2013) FHDM HDM Picknett & Bexon at Tsub.
	E
	0.05		0.087	0.091	0.083	0.11
									298 303 308 313 318 323 328 333 338
	0.00
		Earth's	Microgravity	FHDM	HDM		Analytical
		gravity Exp.	Exp
	Figure 5: Comparison of HFE7100 sessile drop evaporation
	rate from experiment and numerical simulations. For Earth
	gravity exp.: average evaporation rates of sessile drops 1DP1g
	and 3DP1g, Microgravity Exp.: evaporation rates of sessile
	drops 7DPµg [34] (refer table 1 of [34]). The parameters of
	sessile drop 7DPµg were used for the numerical simulation
	FHDM and HDM and also for the calculation of analytical
	diffusion-limited evaporation rate without an electric field [42].
	The evaporation rate is calculated from initial contact angle
	𝜃 0 = 57 • to 𝜃 0 = 46 • .

𝑃 = 83500 Pa, substrate temperature 𝑇 sub from 298.15 K to 333.15 K, and at contact angle 𝜃 = 23 • . The spatial discretization used in numerical simulation for the Ethanol sessile drop is shown in figure

2

. The mesh structure for all the substrate temperatures remains the same only the number

  Figure 7: Numerical simulation's average interface temperature (𝑇 𝑎𝑣𝑔,Γ ) vs substrate temperatures. 𝑇 𝑎𝑣𝑔,Γ of ethanol sessile drop between 𝜃 0 = 87 • to 𝜃 0 = 15 • with FHDM and HDM. Insert shows fitting for the 𝑇 𝑎𝑣𝑔,Γ vs 𝑇 𝑠𝑢𝑏 only for FHDM. Fitting details (for FHDM): 𝑇 𝑎𝑣𝑔,Γ = 𝑎 𝐸 ×𝑇 𝑏 𝐸 𝑠𝑢𝑏 with a E = 1.41148±0.047 and b E = 0.93887 ± 0.0058. Readers are advised to refer to web version of this figure for better display.

	335
			335
	295 300 305 310 315 320 325 Avg. g nterface remperature, T avg,𝛤 (K) 330	Avg.Interface Temperature ( ) (K) sub,Γ	295 300 305 310 315 320 325 330 335 290 Substrate Temperature ( ) (K) sub T 295 300 305 310 315 320 330 325	a b Reduced Chi-Sqr R-Square (COD) Adj. R-Square odel Equation Plot	1.41148 ± 0.0473 0.93887 ± 0.00582 0.19969 0.99828 0.99824 Allometric1 y = a*x^b Concatenated Data	FHDM θ : FHDM θ: FHDM θ: FHDM θ: FHDM θ: FHDM θ: HDM θ: HDM θ: HDM θ : HDM θ: HDM θ : HDM θ:	15°2 5°2 0°3 0°5 0°7 0°8 7°M 0°3 0°5 0°7 0°8 7°1
	290
	295	300	305	310	315	320	325	330	335
				Substrate Temperature, T sub (K)

Acknowledgement

We would like to acknowledge the financial support of the Centre National d'Études Spatiales (CNES) and the European Space Agency (ESA). Also, Centre de Calcul Intensif d'Aix-Marseille is acknowledged for granting access to its high performance computing resources.

Conflict of interest

The authors declared that there is no conflict of interest.

CRediT authorship contribution statement