Vinícius Braga Leite

Denize Kalempa

Irina Martin Graur

Kinetic modelling of evaporation and condensation phenomena around a spherical droplet

Introduction

The evaporation and condensation phenomena are of great interest in many fields such as environmental sciences, micro and nanotechnology, security of nuclear plants, etc. For instance, evaporation and condensation processes in mini and microchannels have been widely used in air-cooled condensers, in heat pipes and other devices for thermal control of different systems. In particular, microchannel condensers [START_REF] Awad | A critical review on condensation heat transfer in microchannels and minichannels[END_REF] are used to increase performance in heat transfer, to reduce components size and improve energy efficiency. The authors of Ref. [START_REF] Kanklikar | Heat transfer in microchannels -2012 status and research needs[END_REF] state that innovative cooling technologies must necessarily use the processes of boiling, evaporation and condensation in novel reduced geometrical configurations intended for the heat transfer enhancement. Moreover, according to Ref. [START_REF] Bar-Cohen | Thermal management of on-chip hot spot[END_REF], a big gap still exists between the chip-level heat generation and the system-level heat removal in high performance electronic systems, such as computers, electronics and solid-state lasers. Thus, one can say that evaporation and condensation phenomena play an important role in the development of saving energy and environmentally friendly new technologies.

The theory of evaporation and condensation has its foundations in the classical kinetic theory of gases. Hertz [START_REF] Hertz | Uber die Verdungstung der Flussigkeiten, inbesondere des Quecksilbers, im luftleeren Raume[END_REF] and Knudsen [5] provided the molecular basis for understanding these phenomena in the free molecular flow regime, where the intermolecular collisions are negligible. Later, Pao [START_REF] Pao | Temperature and density jumps in the kinetic theory of gases and vapors[END_REF][START_REF] Pao | Application of kinetic theory to the problem of evaporation and condensation[END_REF] and Thomas et al. [START_REF] Thomas | Reverse temperature gradient in the kinetic theory of evaporation[END_REF] took into account the effect of the intermolecular collisions. These authors predicted theoretically the discontinuity of the pressure and temperature profiles through a liquid-vapor interface and explained this phenomenon by the non-equilibrium state of a gas inside the thin layer around the interface, the Knudsen layer. Labuntsov and coworkers also contributed for the development of the kinetic theory of evaporation and condensation inside the Knudsen layer, see e.g. Refs. [START_REF] Labuntsov | Analysis of the evaporation and condensation processes[END_REF][START_REF] Labuntsov | Analysis of intensive evaporation and condensation[END_REF], and proposed useful equations to describe both weak and strong evaporation and condensation processes. Experimentally, the authors of Refs. [START_REF] Shankar | On the temperature distribution in liquid-vapor phase change between plane liquid surfaces[END_REF][START_REF] Ward | Expression for predicting liquid evaporation flux: Statistical rate theory approach[END_REF][START_REF] Ward | Interfacial conditions during evaporation or condensation of water[END_REF] were the first who have measured the temperature jump at the liquid-vapor interface. Later, various campaigns were undertaken to measure the temperature discontinuity at the interface under different physical conditions, e.g. Refs. [START_REF] Badam | Experimental and theoretical investigations on interfacial temerature jumps during evaporation[END_REF][START_REF] Gatapova | The temperature and pressure jumps at the vapor-liquid interface: Application to a two-phase cooling system[END_REF][START_REF] Gatapova | The temperature jump at water-air interface during evaporation[END_REF][START_REF] Kazemi | Experimental and numerical study of the evaporation of water at low pressures[END_REF][START_REF] Jafari | Evaporation mass flux: a predictive model and experiments[END_REF]. However, all experimental data on the temperature jump are very scattered, ranging from 15K to 0.1K.

In practice, the Navier-Stokes (NS) equations, usually valid only under equilibrium conditions, are widely used to describe phase transition phenomena. For instance, according to Ref. [START_REF] Kryukov | Evaporation of diesel fuel droplets: kinetic versus hydrodynamic models[END_REF], the hydrodynamic approach, which is usually used in engineering models of diesel fuel droplet evaporation and incorporated in relevant computational fluid dynamic codes, provides results in disagreement with realistic diesel engines, so more general kinetic models need to be employed.

Theoretical approaches for the modelling of the evaporation and condensation phenomena, based on the gas kinetic theory, e.g. Refs. [START_REF] Cipolla | Kinetic theory of condensation and evaporation[END_REF][START_REF] Sone | Kinetic theory of evaporation and condensation -Hydrodynamic equation and slip boundary condition[END_REF][START_REF] Sone | Kinetic theory of evaporation and condensation[END_REF][START_REF] Sone | Evaporation and condensation on a plane condensed phase: numerical analysis of the linearized Boltzmann equation for hard-sphere molecules[END_REF][START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF][START_REF] Young | The condensation and evaporation of liquid droplets in a pure vapour at arbitrary Knudsen number[END_REF], confirm that the phase transition occurs only under non-equilibrium conditions and that the continuity of the thermodynamic variables, such as temperature and pressure, cannot be ensured. Most of these studies focus on the evaporation-condensation in plane geometry of one or two condensed phases and only few authors considered the evaporation from a spherical droplet [START_REF] Sone | Kinetic theory of evaporation and condensation -Hydrodynamic equation and slip boundary condition[END_REF][START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF][START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF]. In Ref. [START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF] the problem of evaporation-condensation on a spherical droplet was studied by implementing a linearized approach for a kinetic model to the Boltzmann equation, while in Ref. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF] the same problem was solved via the linearized Boltzmann equation. The asymptotic analysis for small Knudsen number was proposed in Ref. [START_REF] Sone | Kinetic theory of evaporation and condensation -Hydrodynamic equation and slip boundary condition[END_REF]. The author of Ref. [START_REF] Young | The condensation and evaporation of liquid droplets in a pure vapour at arbitrary Knudsen number[END_REF] derived several useful correlations for the evaporation rate and in-terfacial heat transfer and compared them with the numerical results of [START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF]. The question on the limits of the applicability of the linearized kinetic type equation for the simulation of the evaporation and condensation phenomena are discussed only in the case of the plane geometry [START_REF] Aoki | Gas flows caused by evaporation and condensation on two parallel condensed phases and the negative temperature gradient: numerical analysis by using a nonlinear kinetic equation[END_REF][START_REF] Polikarpov | Heat and mass transfer in a rarefied gas confined between its two parallel condensed phases[END_REF]. In addition, the influence of the non-complete evaporation-condensation was very often ignored.

The objective of the present work is the modelling of the steady evaporation and condensation phenomena which take place on the surface of a spherical droplet of the condensed phase of a vapor on the basis of the kinetic theory. The paper is organized as follows. First, the model proposed by Shakhov [START_REF] Shakhov | Approximation kinetic equations in rarefied gas theory[END_REF] for the linearized Boltzmann equation, referenced here as S-model, is solved numerically by the Discrete Velocity Method (DVM) implemented to take into account the discontinuity of the distribution function of molecular velocities around a convex body, see e.g. Refs. [START_REF] Sone | Thermal creep in rarefied gas[END_REF][START_REF] Sone | Discontinuity of the velocity distribution function in a rarefied gas around a convex body and the S layer at the bottom of the Knudsen layer[END_REF]. Second, the analytical solution is derived in the case of free molecular and hydrodynamic regimes. In the hydrodynamic regime, the thermodynamic analysis is used to obtain the temperature and pressure jump boundary conditions at the liquid-vapor interface. Third, the non-linear S-model kinetic equation is solved numerically by using the DVM in the velocity space and the Total Variation Diminishing (TVD) type numerical scheme [START_REF] Ho | Heat transfer through rarefied gas confined between two concentric spheres[END_REF] for the approximation of the spatial derivatives, which does not require a special treatment for the discontinuity of the distribution function on the convex droplet surface. The mass and energy flow rates as well as the macroscopic characteristics of the gas flow around the droplet due to the phase transition were calculated in a range of the gas rarefaction which covers the free molecular, transition and hydrodynamic regimes. Moreover, complete and non-complete evaporationcondensation on the liquid-vapor interface were considered by varying the evaporation-condensation coefficient. The linearized and non-linear results were compared to find the limits of applicability of the linearized approach.

Statement of the problem

Let us consider a gas flow around its spherical condensed-phase droplet of radius R 0 caused by steady evaporation and condensation at the interface. The temperature of the liquid droplet is supposed to be constant and equal to T d . The corresponding saturation pressure, p d = p sat (T d ), could be calculated by the Clausius-Clapeyron equation [START_REF] Landau | Fluid Mechanics[END_REF]. Far from the droplet the gas is at thermodynamic equilibrium with pressure p 0 and temperature T 0 , which are different from the pressure and temperature of the droplet.

It is worth noting that the evaporation and condensation processes are time dependent so that the droplet radius R 0 and temperature T d are func-tions of time. However, according to the literature, see e.g. [START_REF] Fuchs | Evaporation and Droplet Growth in Gaseous Media[END_REF][START_REF] Hidy | The Dynamics of Aerocolloidal Systems[END_REF], the steady state solution is valid under the following assumptions: (i) the difference between the vapor density far from the droplet and at the interface are very small so that the droplet condensation growth and evaporation depletion can be considered as quasi-stationary processes and steady state solutions can be applied; (ii) the droplet surface temperature is maintained constant by an inner heat source within the droplet volume.

The gas rarefaction is characterized by the rarefaction parameter, δ, which is inversely proportional to the well known Knudsen number, and defined here as

δ = R 0 0 , 0 = µ 0 v 0 p 0 , v 0 = 2k B T 0 m , ( 1
)
where 0 is the equivalent molecular free path, v 0 is the most probable molecular velocity, k B is the Boltzmann constant and m is the molecular mass of the gas particles. In case of rigid-spheres model, the viscosity of a gas, µ, is related to the viscosity µ 0 at the equilibrium temperature T 0 as follows

µ(T ) = µ 0 T T 0 . (2)
Note that the limit δ → 0 corresponds to the free molecular regime, while the limit δ → ∞ corresponds to the hydrodynamic or continuum limit.

In following we apply the S-model kinetic equation [START_REF] Shakhov | Approximation kinetic equations in rarefied gas theory[END_REF] to simulate the behavior of the spherical condensed-phase droplet due to the evaporation and condensation at the spherical surface. The kinetic equation is solved numerically in its linearized and non-linear forms and the obtained results are compared to establish the limits of applicability of the linearized approach. The analytical solution obtained from the linearized equations of continuum mechanics with appropriate temperature and pressure jump boundary conditions is also presented and compared to the numerical solution of the linearized model kinetic equation.

Kinetic equation

Considering the symmetry of the problem relatively to the center point, the kinetic model proposed by Shakhov [START_REF] Shakhov | Approximation kinetic equations in rarefied gas theory[END_REF] for the Boltzmann equation in the absence of external forces reads

∂f ∂t + v r ∂f ∂r - v t r ∂f ∂θ = ν(f S -f ), (3)
where f (t , r , v) is the distribution function of molecular velocities, t is the time, r is the radial coordinate, and v=(v r , v Θ , v Φ ) is the molecular velocity vector, where v r , v Θ , and v Φ are the radial, polar and azimuthal components in spherical coordinates (v p , θ, φ). Thus, in Eq. ( 3), v r =v p cos θ and v t =v p sin θ denote, respectively, the radial and tangential components of the molecular velocity vector v. Note that, v p = v 2 r + v 2 t is the magnitude of the molecular velocity vector, while θ and φ are its orientational angles. The quantity ν is the molecular collision frequency which is assumed to be independent of the molecular velocity and evaluated as

ν = p µ , (4)
where µ, given in (2), is the viscosity of the gas at an arbitrary temperature T . The Shakhov equilibrium distribution function f S reads

f S (n, T, V, q ) = f M 1 + 2mVq 15n(k B T ) 2 mV 2 2k B T - 5 2 , f M (n, T, V) = n m 2πk B T 3/2 exp - mV 2 2k B T , (5)
where f M is the local Maxwellian distribution function; V = vu is the peculiar velocity vector, u = (u r , 0, 0) is the bulk velocity vector; q = (q r , 0, 0) is the heat flux vector. The macroscopic quantities corresponding to the gas number density, the temperature, the bulk velocity and heat flux vectors, are calculated in terms of the distribution function of molecular velocities, see e.g. Refs. [START_REF] Kogan | Rarefied Gas Dynamics[END_REF][START_REF] Ferziger | Mathematical Theory of Transport Processes in Gases[END_REF], as follows

n(t , r ) = ∞ 0 π 0 2π 0 f (t , r , v) dv, (6)
T (t , r ) = m 3nk B ∞ 0 π 0 2π 0 V 2 f (t , r , v) dv, (7)
u (t , r ) = 1 n ∞ 0 π 0 2π 0 vf (t , r , v) dv, (8)
q (t , r ) = m 2 ∞ 0 π 0 2π 0 V 2 vf (t , r , v) dv, (9)
where dv = v 2 p sin θdv p dθdφ. The mass and energy flow rates, which are the quantities of practical interest, are also calculated in terms of the distribution function of molecular velocities as

Ṁ = ∞ R 0 dA ∞ 0 π 0 2π 0 mvf (t , r , v) dv, (10)
Ė = ∞ R 0 dA ∞ 0 π 0 2π 0 v mv 2 2 f (t , r , v) dv, (11)
where dA =8πr dr is an area element of a spherical surface with radius r .

On the droplet surface we assume that the fraction σ of incident molecules condensates at the surface and then evaporates from it with the equilibrium Maxwellian distribution function, while the fraction 1σ is reflected under the assumption of diffuse scattering. Therefore, the boundary condition at the vapor-liquid interface is written as

f (t, R 0 , v) = σf s + (1 -σ)f r at v r > 0, (12)
where f s and f r are the equilibrium Maxwellian functions given by

f i = n i m 2πk B T s 3/2 exp - mv 2 p 2k B T s , i = s, r, (13)
where n d is the number density calculated from equation of state using the droplet temperature and pressure as n d = p d /(k B T d ) and n r is the number density calculated from the impermeability condition for the molecules reflected diffusively from the surface. Far from the sphere the gas is assumed to be in thermodynamic equilibrium so that

f ∞ = lim r →∞ f (t , r , v) = f M 0 , (14)
where the Maxwellian distribution function is given by

f M 0 = n 0 m 2πk B T 0 3/2 exp - mv 2 p 2k B T 0 , (15)
where the number density n 0 is calculated from the equation of state n 0 = p 0 /(k B T 0 ). In the numerical calculations, a maximum radial distance from the centre of the spherical droplet is set to ensure the thermodynamic equilibrium far from the sphere.

It is worth mentioning that, although our problem concerns the steady evaporation and condensation at the vapor-liquid interface, the time dependence of the distribution function was maintained in this section, essentially because the numerical method used for the solution of the non-linear kinetic equation is based on the iteration in time, see Section 5.2. The implicit numerical scheme is adopted here to decrease the computational time. However, the explicit scheme could be also applied to quantify the evaporation rate in time.

Linearized approach

Linearized kinetic equation

To simulate the evaporation condensation phenomena, two thermodynamic driven forces are introduced as follows

X n = |n d -n 0 | n 0 , X T = |T d -T 0 | T 0 . ( 16
)
Since the saturation pressure depends on the saturation temperature, these two thermodynamic forces are not independent quantities and could be related by the Clausius-Clapeyron equation. However, in the frame of this study, these forces are assumed as being independent from each other, to see clearly the impact of each force on the evaporation process.

When the system is weakly disturbed from thermodynamic equilibrium, the thermodynamic forces given in [START_REF] Gatapova | The temperature jump at water-air interface during evaporation[END_REF] are very small, i.e. X n 1 and X T 1. Therefore, the kinetic equation ( 3) can be linearized by representing the distribution function of molecular velocities as

f (r , v) = f M 0 1 + h (n) (r , v)X n + h (T ) (r , v)X T , (17)
where the Maxwellian distribution function corresponding to the gas equilibrium state far from the spherical droplet is given by [START_REF] Gatapova | The temperature and pressure jumps at the vapor-liquid interface: Application to a two-phase cooling system[END_REF]. The functions h (n) and h (T ) are the perturbation functions due to the density and the temperature deviations, X n and X T , respectively. For further derivations, the dimensionless quantities corresponding to the radial coordinate r, molecular velocity vector c=(c r , c Θ , c Φ ), bulk velocity vector u=(u r , 0, 0) and heat flux vector q=(q r , 0, 0) are introduced as

r = r R 0 , c = v v 0 , u = u v 0 , q = q p 0 v 0 . ( 18
)
Thus, the representation (17) allows us to write the macroscopic characteristics of the gas flow given in ( 6)-( 9) as follows

n = n 0 (1 + ν (n) X n + ν (T ) X T ), T = T 0 (1 + τ (n) X n + τ (T ) X T ), u = u (n) X n + u (T ) X T , q = q (n) X n + q (T ) X T . (19)
Moreover, the substitution of the representation (17) into the kinetic equation (3) allows us to obtain two independent linearized kinetic equations for each thermodynamic force which read

c r ∂h (i) ∂r - c t r ∂h (i) ∂θ = δ ν (i) + c 2 - 3 2 τ (i) +2c r u (i) r + 4 15 c r c 2 - 5 2 q (i) r -h (i) , i = n, T, (20)
where c r =c p cos θ and c t =c p sin θ denote the dimensionless radial and tangential components of the molecular velocity vector c. The dimensionless moments appearing on the right hand side of Eq. ( 20) corresponding to each thermodynamic force are obtained from ( 19) and ( 6)-( 9) as

ν (i) (r) = 2 √ π ∞ 0 π 0 h (i) (r, θ, c p )e -c 2 p dc, (21)
τ (i) (r) = 4 3 √ π ∞ 0 π 0 c 2 - 3 2 h (i) (r, θ, c p )e -c 2 p dc, (22)
u (i) r (r) = 2 √ π ∞ 0 π 0 c p cos θh (i) (r, θ, c p )e -c 2 p dc, (23)
q (i) r (r) = 2 √ π ∞ 0 π 0 c p cos θ c 2 - 5 2 h (i) (r, θ, c p )e -c 2 p dc, (24)
where dc=c 2 p sin θ dc p dθ. For convenience, hereafter, the subscript "r" denoting the radial components of the gas bulk velocity and heat flux is omitted.

From the boundary condition ( 12), the representation (17) leads to the following linearized boundary conditions for the perturbation functions due to X n and X T h

(n) (1, θ, c p ) = σ + (1 -σ)n (n) r , at c r > 0, (25)
h (T ) (1, θ, c p ) = c 2 - 3 2 + (1 -σ)n (T ) r at c r > 0, (26)
where the parameters n

(n) r
and n

(T ) r are obtained from the impermeability condition for the gas at the interface as

n (n) r = -4 ∞ 0 π π/2 c 3 p h (n) (1, θ, c p ) sin θ cos θ dc p dθ (27)
and

n (T ) r = - 1 2 -4 ∞ 0 π π/2 c 3 p h (T ) (1, θ, c p ) sin θ cos θ dc p dθ. ( 28
)
Far from the sphere, the substitution of the representation ( 17) into ( 14) leads to

h (i) ∞ = lim r→∞ h (i) = 0, i = n, T. (29)
The infinite interval [0, ∞) for the radial distance r is substituted by a finite interval [0, r max ], where the maximum distance r max must be chosen to guarantee the thermodynamic equilibrium far from the spherical droplet.

Thermodynamic analysis

From the mass and energy conservation laws we can show that the dimensionless mass and energy flow rates read

Ṁ = Ṁ 4πR 2 0 mn 0 v 0 = r 2 u, (30)
Ė = Ė 4πR 2 0 v 0 p 0 = r 2 q + 5 2 u , (31)
where Ṁ and Ė are defined in [START_REF] Labuntsov | Analysis of intensive evaporation and condensation[END_REF] and [START_REF] Shankar | On the temperature distribution in liquid-vapor phase change between plane liquid surfaces[END_REF]. Note that the quantities on the right hand side of ( 30) and ( 31) are constant, i.e. they do not depend on the radial coordinate.

According to the thermodynamics of irreversible processes [START_REF] De Groot | Non-Equilibrium Thermodynamics[END_REF], the entropy production S can be written as

S = k B [J n X n + (J n + J T )X T ], (32)
where the thermodynamic fluxes are defined as

J n = n 0 v 0 u dA , J T = n 0 v 0 q dA , (33)
and the thermodynamic forces X n and X T are defined in [START_REF] Gatapova | The temperature jump at water-air interface during evaporation[END_REF]. These thermodynamic fluxes satisfy the linear relations

J n = Λ nn X n + Λ nT X T , (34)
J n + J T = Λ T n X n + Λ T T X T , (35)
where Λ ij are the so called kinetic coefficients which satisfy the Onsager reciprocity relation Λ nT =Λ T n . Thus, the mass and energy flow rates defined in ( 10) and ( 11) can be written in terms of the thermodynamic fluxex J n and J T as

Ṁ = mJ n , Ė = k B T J T + 5 2 J n . ( 36
)
The advantage of using J n and J T instead of Ṁ and Ė is that the number of independent kinetic coefficients is reduced from four to three due to the reciprocity relation. Moreover, the reciprocity relation represents an important criterion to verify the accuracy of the numerical calculations. For convenience, the dimensionless kinetic coefficients are introduced as

Λ ij = Λ ij 4πR 2 0 n 0 v 0 , i, j = n, T (37)
and are obtained from [START_REF] Labuntsov | Analysis of intensive evaporation and condensation[END_REF], [START_REF] Shankar | On the temperature distribution in liquid-vapor phase change between plane liquid surfaces[END_REF], [START_REF] Fuchs | Evaporation and Droplet Growth in Gaseous Media[END_REF], [START_REF] Hidy | The Dynamics of Aerocolloidal Systems[END_REF] and [START_REF] Kogan | Rarefied Gas Dynamics[END_REF] as

Λ nn = r 2 u (n) , Λ nT = r 2 u (T ) , (38)
Λ T n = r 2 (u (n) + q (n) ), Λ T T = r 2 (u (T ) + q (T ) ). ( 39
)
Therefore, the reciprocity relation for the problem in question reads

u (T ) (r) = u (n) (r) + q (n) (r). ( 40
)

Analytic solution in the free molecular regime

When δ 1, the intermolecular collisions can be neglected and, consequently, the linearized kinetic equations given in [START_REF] Cipolla | Kinetic theory of condensation and evaporation[END_REF] are reduced to

c r ∂h (i) ∂r - c t r ∂h (i) ∂θ = 0, i = n, T. (41)
Moreover, in this flow regime, the perturbation function of the incident particles at the interface corresponds to that far from the sphere, given by [START_REF] Shakhov | Approximation kinetic equations in rarefied gas theory[END_REF]. As a consequence, the boundary conditions given by ( 25) and ( 26) are reduced to

h (n) (1, θ, c p ) = σ at c r > 0, (42)
h (T ) (1, θ, c p ) = c 2 -2 + σ 2 at c r > 0. ( 43
)
The kinetic equations given in ( 41) can be solved analytically by employing the known method of the characteristics to solve a differential equation. Thus, the solution of Eq. ( 41) subject to the boundary condition corresponding to each thermodynamic force, given by ( 42) and ( 43), reads

h (n) (r, θ, c p ) = σ, 0 ≤ θ ≤ θ 0 , 0, θ 0 < θ ≤ π, (44)
h (T ) (r, θ, c p ) = c 2 -2 + σ 2 , 0 ≤ θ ≤ θ 0 , 0, θ 0 < θ ≤ π, (45)
where the angle θ 0 is given by

θ 0 = arcsin 1 r . ( 46
)
After substituting the solutions ( 44) and ( 45) into ( 21)-( 24), the following expressions for the macroscopic characteristics of the gas flow are obtained

ν (n) (r) = σ 2   1 -1 - 1 r 2   , τ (n) (r) = 0, (47)
ν (T ) (r) = - (1 -σ) 4   1 -1 - 1 r 2   , τ (T ) (r) = 1 2   1 -1 - 1 r 2   , (48)
u (n) (r) = σ 2 √ π 1 r 2 , u (T ) (r) = σ 4 √ π 1 r 2 , (49)
q (n) (r) = - σ 4 √ π 1 r 2 , q (T ) (r) = 1 √ π 1 - σ 8 1 r 2 . ( 50
)
The dimensionless mass flow rate at the interface is obtained from (30) as

Ṁ = σ 2 √ π X n + σ 4 √ π X T , (51)
and corresponds to widely used Hertz-Knudsen formula for the mass flow rate of vapor from the droplet. The energy flow rate is obtained from [START_REF] Sone | Discontinuity of the velocity distribution function in a rarefied gas around a convex body and the S layer at the bottom of the Knudsen layer[END_REF] as

Ė = σ √ π X n + 1 √ π 1 + σ 2 X T . (52)

Numerical solution for arbitrary gas rarefaction

The linearized kinetic equations given in [START_REF] Cipolla | Kinetic theory of condensation and evaporation[END_REF] subject to the corresponding boundary condition for each thermodynamic force are solved numerically by applying the DVM, whose details can be found in the literature, see e.g. Ref. [START_REF] Sharipov | Rarefied Gas Dynamics. Fundamentals for Research and Practice[END_REF]. Moreover, the split method proposed in Ref. [START_REF] Naris | The driven cavity flow over the whole range of the Knudsen number[END_REF] to deal with the problem of the discontinuity of the distribution function of molecular velocities on the spherical interface is employed. In rarefied gas dynamics, the problem of the discontinuity of the distribution function is a peculiarity inherent to gas flows around convex bodies, see e.g. Ref. [START_REF] Sone | Molecular Gas Dynamics. Theory, Techniques and Applications[END_REF], and must be treated carefully when a finite difference scheme is used. According to the analytic solutions in the free molecular regime, given in ( 44) and ( 45), the discontinuity occurs when θ=θ 0 . For small values of rarefaction parameter this behavior still remains for θ=θ 0 and, as a consequence, the numerical solution based on the finite difference scheme to calculate the derivative ∂h/∂θ leads to a significant numerical error. The split method consists on the decomposition of the perturbation function into two parts as

h (i) (r, θ, c p ) = h (i) 0 (r, θ, c p ) + h(i) (r, θ, c p ), i = n, T, (53)
where h

(i) 0 satisfy the following partial differential equation c r ∂h (i) 0 ∂r - c t r ∂h (i) 0 ∂θ + δh (i) 0 = 0, (54)
subject to the boundary conditions

h (n) 0 (1, θ, c p ) = σ at c r > 0, (55)
h (T ) 0 (1, θ, c p ) = c 2 -2 + σ 2 at c r > 0. ( 56
)
The functions h(n) and h(T) satisfy the kinetic equations given by ( 20), just changing h (i) by h(i) (i = n, T ), with boundary conditions

h(i) (1, θ, c p ) = (1 -σ)ñ (i) r at c r > 0, (57)
where

ñ(i) r = ∞ 0 π π/2 c 3 h(i) (1, θ, c p ) sin θ cos θ dc p dθ. ( 58
)
The advantage of this method is that the function h

(i)
0 is found analytically by employing the method of the characteristics, while the function h(i) is sufficiently smooth so that a finite difference scheme can be used and it leads to a smaller numerical error. For convenience, the analytic solution h (i) 0 (i = n, T ), for each thermodynamic force is presented in Appendix A. It is worth noting that the moments given in ( 21)-( 24) are also decomposed into two parts due to the representation (53). Thus, they are written as

ν (i) (r) = ν (i) 0 (r) + ν(i) (r), τ (i) (r) = τ (i) 0 (r) + τ (i) (r), u (i) (r) = u (i) 0 (r) + ũ(i) (r), q (i) (r) = q (i) 0 (r) + q(i) (r), (59)
where the quantities with tilde are calculated via the expressions given in ( 21)-( 24), just changing h (i) by h(i) , while the quantities with subscript "0" are calculated in terms of the function h (i) 0 and their expressions are provided in Appendix A. To show the typical behavior of the perturbation functions h (n) and h (T ) depending on the angle θ for small values of δ when the DVM with and without the split technique is employed, some figures are provided as online supplementary material.

The kinetic equation for each thermodynamic force subject to the corresponding boundary condition and asymptotic behavior was solved via the DVM with an accuracy of 0.1% for the gas flow macroscopic characteristics at the vapor-liquid interface. The calculations were carried out in a range of the gas rarefaction which covers the free molecular, transition and hydrodynamic flow regimes. The Gaussian-Hermite quadrature was used to discretize the molecular velocity space and calculate the moments of the perturbation function. The numerical technique to calculate the nodes and weights is described in Ref. [START_REF] Krylov | Approximate Calculation of Integrals[END_REF]. A central finite difference scheme was used to approximate the derivatives appearing in the kinetic equation. The accuracy was estimated by varying the grid parameters N r , N θ and N c , corresponding to the number of nodes in the radial coordinate r, angle θ and magnitude of the molecular velocity, c p . The values of these parameters, N c and N θ , were fixed at 12 and 200, respectively, while N r varied according to the maximum value of radial coordinate, r max , so that the increment ∆r ∼ 10 -3 . The distance r max was chosen so that the dimensional distance r max varied from 10 0 to 100 0 when the rarefaction parameter varied from 0.01 to 10. Note that, from Eqs. ( 1) and ( 18), the distance r max can be written in terms of the equivalent mean free path 0 as r max = r max δ 0 . In the limits δ → 0 and δ → ∞, the numerical solution was compared to those obtained analytically in the free molecular and hydrodynamic regimes. The reciprocal relation [START_REF] Naris | The driven cavity flow over the whole range of the Knudsen number[END_REF] was verified within the numerical error.

Hydrodynamic solution

When δ 1, the problem can be solved via the linearized continuum equations with appropriate pressure and temperature jump conditions at the vapor-liquid interface. Although the modelling via the kinetic equation provides more accurate information concerning the evaporation-condensation process, in this flow regime the convergence of the numerical scheme based on the discrete velocity method is very slow due to the large number of the grid points necessary to guarantee the accuracy of the calculations. Then, in this flow regime, it is more convenient to use the continuum equations. In our notation, the linearized balance equations of mass, momentum and energy read

d(r 2 u) dr = 0, (60)
dp dr = 0, (61)
d dr r 2 dT dr = 0, ( 62
)
where the pressure is expressed via the equation of state of an ideal gas. The solution of these equations requires the use of appropriate boundary conditions which take into account the pressure and temperature jumps at the vapor-liquid interface. Since the processes of heat and mass transfer at the interface can be described by the non-equilibrium thermodynamics, the boundary conditions for the pressure and temperature at the interface can be obtained from the inverted form of the flux-force relations. Thus, according to Ref. [START_REF] Fujikawa | Vapor-liquid interface, bubbles and droplets[END_REF], one can write

p = p s (1 -ζ P J n -ζ T J T ) , at r = 1, (63)
T = T s (1 -ζ T J n -ζ TP J T ) , at r = 1, (64)
where the dimensionless thermodynamic fluxes J i (i = n, T ) are introduced as

J i = J i 4πR 2 0 n 0 v 0 , (65)
and J i are defined in [START_REF] Landau | Fluid Mechanics[END_REF]. The jump coefficients are obtained from the kinetic theory, see e.g. Ref. [START_REF] Cipolla | Kinetic description of cylindrical heat conduction in a polyatomic gas[END_REF], and in case of monoatomic gas under the assumption of diffuse reflection of the fraction of particles which do not condense at the interface read

ζ P = 2 √ π 1 σ -0.40044 , ζ T = 0.446658, ζ TP = 1.042203. ( 66
)
From the balance equations ( 60)-( 62) and the vapor-liquid interface conditions (63) and (64), the bulk velocity and temperature deviation of the gas flow from equilibrium are obtained as follows

u(r) = 1 r 2 [X n + (1 -ζ T ξ)X T ] (ζ P -ζ 2 T ξ) , (67)
τ (r) = - 1 r(1 + κ) [ζ T X n -(ζ P -ζ T )X T ] (ζ P -ζ 2 T ξ) , (68)
where

κ = ζ TP 2δ γ (γ -1) 1 Pr , ξ = κ ζ TP 1 (1 + κ) . (69)
In the previous expression, γ = c p /c v is the ratio of the specific heat capacities and Pr is the Prandtl number. These parameters for a monatomic gas are equal to 5/3 and 2/3, respectively. According to the balance equation (61), the pressure is a constant quantity. Therefore, since the temperature profile is known from (68), the density deviation from equilibrium reads

ν(r) = -τ (r). ( 70
)
The heat flux is obtained from the Fourier law as

q(r) = - 1 r 2 κ ζ TP (1 + κ) [ζ T X n -(ζ P -ζ T )X T ] (ζ P -ζ 2 T ξ) . ( 71
)
Note that, in the limit δ → ∞, the expressions (67), ( 68) and (71) tend to

u(r) = 1 ζ P r 2 (X n + X T ), (72)
τ (r) = - 1 r ζ T ζ P X n -1 - ζ T ζ P X T , (73)
q(r) = - 1 r 2 κ ζ TP ζ T ζ P X n -1 - ζ T ζ P X T . (74)
The dimensionless mass and energy flow rates at the vapor-liquid interface are calculated from [START_REF] Sone | Thermal creep in rarefied gas[END_REF] and [START_REF] Sone | Discontinuity of the velocity distribution function in a rarefied gas around a convex body and the S layer at the bottom of the Knudsen layer[END_REF]. It is worth to underline that, in our notation, the widely used Schrage relation [START_REF] Schrage | A Theoretical Study of Interphase Mass Transfer[END_REF] for the mass flow rate under the assumption of very small driving forces reads

ṀSch = σ (2 -σ) √ π X n + 1 2 X T . (75)
Thus, in order to compare the present results with those obtained from the Schrage relation, the dimensionless mass flow rate at the interface obtained from ( 30) and ( 72) is used and it reads

Ṁ = 1 ζ P (X n + X T ), (76)
The comparison between the results obtained from the continuum approach (76), the kinetic modelling and Schrage expression (75) is presented in Section 6.

Non-linear approach

Non-linear kinetic equation and boundary conditions

To establish the limits of applicability of the linearized approach, the non-linear S-model kinetic equation ( 3) is solved numerically. Additionally to the dimensionless quantities introduced in (18), the defined dimensionless variables are used in the following

f = f v 3 0 n 0 , t = t v 0 R 0 , n = n n 0 , T = T T 0 , μ = µ µ 0 . (77)
Thus, the dimensionless form of the kinetic equation ( 3) reads

∂ f ∂t + c r ∂ f ∂r - c t r ∂ f ∂θ = δn T ( f S -f ). ( 78
)
The dimensionless boundary condition for solving the previous equation is obtained from [START_REF] Ward | Expression for predicting liquid evaporation flux: Statistical rate theory approach[END_REF] as

f (t, 1, θ, c p ) = σ fs + (1 -σ) fr , c r > 0. ( 79
)
Moreover, far from the spherical droplet, the normalized distribution function satisfies the initial condition obtained from [START_REF] Badam | Experimental and theoretical investigations on interfacial temerature jumps during evaporation[END_REF] as

f (0, r, θ, c p ) = f∞ = lim r→∞ f (t, r, θ, c p ) = 1 π 3/2 e -c 2 p . (80)
The flow fields given in ( 6)-( 9) are written in terms of the dimensionless distribution function introduced here as

n(t, r) = 2π ∞ 0 π 0 f (t, r, θ, c p )dc, (81)
T (t, r) = 2π 3n ∞ 0 π 0 [(c r -u r ) 2 + c 2 t ] f (t, r, θ, c p )dc, (82) ûr
(t, r) = 2π n ∞ 0 π 0 c r f (t, r, θ, c p )dc, ( 83
) qr (t, r) = 2π ∞ 0 π 0 c r [(c r -u r ) 2 + c 2 t ] f (t, r, θ, c p )dc. ( 84
)
Once again, for convenience, the subscript "r" denoting the radial components of the gas bulk velocity and heat flux is omitted as well as the super-scriptˆfor dimensionless variables. The dimensionless mass and energy flow rates are given by the previous expressions [START_REF] Sone | Thermal creep in rarefied gas[END_REF] and (31).

Method of solution

To obtain the flow parameters in the transitional regime, Eq. ( 78) is solved numerically by applying the technique developed previously in Ref. [START_REF] Ho | Heat transfer through rarefied gas confined between two concentric spheres[END_REF]. First, the DVM is applied to split the continuum molecular velocity space c p in the governing equation (78) into discrete velocity set c pm , in which m = 1, 2, .., N c P . The Gaus-Hermite quadrature is used to calculate the integrals in the velocity space c p appearing in the calculation of the gas flow fields, given in (81)-(84). As a consequence, the velocity magnitudes c pm are taken to be the roots of the Hermite polynomial of order N c P and associated with their corresponding weights. Then the set of N c P kinetic equations, corresponding to N c P values of discrete velocity c pm , is discretized in time and space by finite difference method (FDM). The range of orientation angle θ (0 ≤ θ ≤ π) is divided into N θ uniform segments defined by θ l , l = 0, 1, .., N θ , and the integrals in the angle θ appearing in (81)-( 84) are calculated by the Simpson 3/8 rule.

A maximum radial coordinate, r max , is set and the gap between r max and r 0 is divided into N r equal intervals characterized by r i , i = 0, 1, .., N r . Let us denote the numerical solution of Eq. ( 78) at the time level k as

f k i,l,m = f t k , r i , θ l , c pm and define ∆t k = t k+1 -t k , ∆r i = r i+1 -r i , ∆θ l = θ l+1 -θ l , ∆f k = f k+1 i,l,m -f k i,l,m .
In this paper we are interested only in the steady-state solution, therefore the fully time-implicit Godunov-type scheme is applied [START_REF] Yang | Rarefied flow computations using nonlinear model Boltzmann equations[END_REF][START_REF] Titarev | Numerical study of unsteady evaporation and heat transfer from a spherical surface[END_REF] so that the discretized form of the kinetic equation reads

1 ∆t k + c p cos θ ∂ ∂r - c p sin θ r ∂ ∂θ ∆f k = RHS k , RHS k = nδT 1-ω f S k -f k -c p cos θ ∂f k ∂r + c p sin θ r ∂f k ∂θ . (85)
The left hand side of Eq. ( 85) is the implicit part, where the spacial derivatives are approximated by the first order upwind scheme. The advantage of this discretization is that the solution ∆f k can be obtained directly, i.e. without calculation of the inverse matrix, just by marching in appropriate direction. The right hand side of Eq. ( 85) RHS k is the explicit part, where the spacial derivatives are approximated by the second order TVDlike scheme. Minute details concerning the numerical scheme are provided in Ref. [START_REF] Graur | Simulation of the transient heat transfer between two coaxial cylinders[END_REF]. The calculations were carried out for three pairs of temperature T s /T 0 and density n s /n 0 ratios and for the same values of rarefaction parameter and evaporation-condensation coefficient previously considered when the linear approach was implemented. The grid parameters were chosen to ensure an accuracy of 0.1% for the moments of the distribution function. While the values of the parameters N c and N θ were fixed at 12 and 400, respectively, the parameter N r varied according to the value of the maximum radial coordinate, r max , so that the increment ∆r i ∼ 10 -2 (i=0, ..., N r ). For the three pairs of T s /T 0 and n s /n 0 considered in the calculations, r max was chosen so that the dimensional quantity r max =r max δ 0 varied from 10 0 to 100 0 as δ varied from 0.01 to 10.

It is worth mentioning that the numerical method used to solve the nonlinear kinetic equation did not require a special treatment to deal with the problem of the discontinuity of the perturbation function which appeared when the DVM was employed to solve the linearized kinetic equation. To show the behavior of the normalized distribution function f on the angle θ and fixed r and c p some figures are provided as online supplementary material for small values of δ.

Results and discussion

We present the results in the following manner. First, the results obtained numerically from the linearized S-model kinetic equation and the analytical solution of the linearized Navier-Stokes equations subject to the pressure and temperature jumps boundary conditions at the interface are presented in order to show the limits of the applicability of this analytical solution in terms of the rarefaction parameter. The cases of complete and non-complete evaporation-condensation are considered. Then, the linearized and non-linear numerical solutions of the S-model kinetic equation are compared to show the limits of the linearized approach application.

Linearized approach

A comparison between the results of the present work, calculated using the discrete velocity method, and those provided in Ref. [START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF], obtained via the integral-moment method, is presented in Tables 1 and2 for complete evaporation-condensation at the liquid-vapor interface (σ = 1) and in a wide range of the gas rarefaction. From these tables one can see a good agreement between the results for the bulk velocities, u (n) and u (T ) , and heat fluxes, q (n) and q (T ) , at the interface. The difference between the results is less than 2% for u (n) , u (T ) and q (T ) , and less than 10% for q (n) .

The results for the mass flow rate and heat flux from the sphere were also compared to the results presented in Ref. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF], which were obtained from the solution of the linearized Boltzmann equation for hard-spheres molecules and complete evaporation-condensation at the interface. The comparison is shown in Tables 3 and4. For convenience, the results are presented in terms of the temperature X T and pressure X P =X n +X T driving forces as it was done in Ref. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF]. Moreover, the results are given as functions of the parameter k=( √ π/2)Kn introduced in Ref. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF], which is proportional to the Knudsen number Kn and can be defined either in terms of the gas viscosity or heat conductivity. As one can see from these tables, there is a good agreement between the present results obtained via the numerical solution of the linearized S-model kinetic equation and the data reported in Ref. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF]. The difference between the results is smaller than 2% in the whole range of the parameter k. The figures showing the comparison between the present results and those from both Refs. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF][START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF] for the mass flow rate and heat flux is available as online suplementary material. The present results are in a better agreement with those provided in Ref. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF] than those obtained by the authors of Ref. [START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF]. It is worth noting that the integral-moment method employed in Ref. [START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF] consists on obtaining a set of integral equations for the moments of the distribution function and its advantage is that only the physical space must be discretized. Nevertheless, this method requires much more computational memory and CPU time than that required by the discrete velocity method. Moreover, even with the computational infrastructure available nowadays, the numerical solution of the Boltzmann equation is still a hard task. Therefore, the implementation of the model kinetic equation allows us to solve the problem with good accuracy and moderate computational effort. Since the results of the present work and those provided in Ref. [START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF] were obtained by the solution of the same kinetic equation, one can conclude that the better agreement with the results from the Boltzmann equation is due to the numerical scheme used in this work.

Figures 123show the profiles of the bulk velocities, u (n) and u (T ) , and heat fluxes, q (n) and q (T ) , as functions of the dimensionless radial coordinate r, for the rarefaction parameter, δ, equal to 0.1, 1 and 10, respectively. For each value of rarefaction parameter, four values of the evaporation-condensation coefficient are considered, σ=0.1, 0.4, 0.6 and 1. In these figures, the dimensionless radial distance r=10 corresponds to a dimensional radial distance r =10R 0 =10δ 0 . According to these figures, the larger the rarefaction parameter the larger the bulk velocities in the vicinity of the sphere due to the driving forces X n and X T . The opposite behavior occurs for the heat fluxes, i.e. the larger the rarefaction parameter the smaller the magnitude of the heat fluxes in the vicinity of the spherical droplet. Note that q (n) is negative, while q (T ) is positive. These figures also show us that the influence of the evaporation-condensation coefficient, σ, on the bulk velocities u (n) and u (T ) , and heat flux q (n) is significant, while the influence on the heat flux q (T ) can be neglected whatever the value of the gas rarefaction. In fact, the heat flux q (T ) depends only on the temperature difference between the droplet surface and the gas in equilibrium far from the droplet. The same qualitative behavior was found in Ref. [START_REF] Polikarpov | Heat and mass transfer in a rarefied gas confined between its two parallel condensed phases[END_REF], in which the problem of noncomplete evaporation-condensation between two parallel condensed phases was simulated. For instance, in case of δ=10, the numerical results obtained from the kinetic equation at the interface, given in Tables 5 and6, show a difference of about 93% for the bulk velocities u (n) and u (T ) , and heat flux q (n) , when the evaporation-condensation coefficient varies from 1 to 0.1. On the other hand, the difference for the heat flux q (T ) is less than 20%. These results show us that the evaporation-condensation on the spherical droplet depends not only on the gas rarefaction degree, but depends strongly on the evaporation-condensation coefficient at the interface.

The comparison between the solution of the linearized kinetic equation and the solution obtained from the Navier-Stokes equations given in (72) and ( 71) is shown in Figure 4 for δ=10 and complete evaporation-condensation at the interface. According to this figure, there is a good agreement between both results. Thus, in order to do just an estimation of the evaporation parameters from the spherical droplet, as well as of the flow fields around it, the analytical solution obtained from the linearized Navier-Stokes equations subjected to the temperature and pressure jump boundary conditions can be used instead of the solution of the linearized kinetic equation. Tables 5 and6 show the comparison in case of complete and non-complete evaporationcondensation at the interface when δ=10. Note that, just for comparison, the results obtained in the limit of δ → ∞, given in Eqs. ( 72) and (74), were also included in the tables. Thus, as one can see from Tables 5 and6, there is a good agreement between the resuls, and such an agreement will be improved with the increase of the rarefaction parameter. Moreover, a better agreement depends on the use of more precise values of the pressure and temperature jump coefficients. In the present work, the temperature and pressure jump coefficients at the interface provided in Ref. [START_REF] Cipolla | Kinetic theory of condensation and evaporation[END_REF] were used, but these coefficients were calculated via the solution of the BGK-model kinetic equation [START_REF] Bhatnagar | A model for collision processes in gases[END_REF] and the variational method.

The comparison of the widely used Schrage formula (76) with the analytical solution obtained from the linearized Navier-Stokes equations with the temperature and pressure jump boundary conditions (75) at the interface and the numerical solution of the linearized kinetic equation for δ = 10 are presented in Table 7. From this table it is clear that the analytical solution derived from the Navier-Stokes equations shows a better agreement with the numerical solution of the kinetic equation compared to the Schrage expression. As it was mentioned previously, the discrepancy between the solution of the kinetic equation for δ=10 and the analytical solution obtained from the Navier-Stokes equations could be explained by the accuracy of the coefficients involved in the temperature and pressure jump conditions.

Comparison between linearized and non-linear S-model solution

The comparison of the results obtained from the linearized approach and that from the solution of the non-linear kinetic equation are presented in Tables 8910for different values of the rarefaction parameter and evaporation coefficients. This comparison can only be done when the driving forces are specified to calculate the macroscopic parameters from the numerical data given by the linearized approach. Three pairs of the temperature ratios, T s /T 0 , and corresponding number density ratios, n s /n 0 , are chosen for this comparison. According to Table 8, for very small deviation from the equilibrium values, T s /T 0 =1.01 and n s /n 0 =1.01, both approaches provide very similar results, with a maximal deviation of less than 1% for the velocity and around 2.5% for the heat flux. However, in practice, following the Clausius-Clayperon equation [START_REF] Landau | Fluid Mechanics[END_REF], the small deviation in temperature leads to a deviation in pressure and, therefore, in the number density, ten times larger. Thus, the results for the ratios T s /T 0 =1.01 and n s /n 0 =1.1 are presented in Table 9 for analysis. As it is clear from Table 9, in this case the difference between the bulk velocity and heat flux calculated by linearized and non-linear approaches is larger, of the order of 4% for the velocity and 6% for the heat flux. Finally, the larger temperature deviation is considered, corresponding to 5% from its equilibrium value. In this case, according to the Clapeyron data [START_REF] Landau | Fluid Mechanics[END_REF], the density deviation is ten times larger, i.e. n s /n 0 =1.5. As it is clear from Table 10, a large discrepancy between linearized and non-linear results are observed, of the order of 25% for both quantities, the macroscopic velocity and the heat flux. It is to note that for all considered values of the evaporation coefficient the difference between linearized and non-linear results is maximum for the complete evaporation, σ = 1. Thus, according to the results, the linearized approach could be applied only when the temperature deviation is less than 5%.

Conclusion

The steady weak evaporation from or condensation onto a spherical condensed phase in its vapor was studied numerically by using the linearized and non-linear S-model kinetic equations. When using the linearized approach, the obtained results have been found in good agreement with those obtained from the linearized Boltzmann equation (maximum difference of 2%) and with those obtained from the solution of the S-model via the integral-moment method (maximum difference of 10%) in the case of complete evaporationcondensation at the interface. In the free molecular flow regime, the analytical solution was obtained in case of complete and non-complete evaporationcondensation at the droplet surface. In the hydrodynamic regime, the lin-earized Navier-Stokes equations subjected to the temperature and pressure jump boundary conditions were solved analytically and the obtained results are in good agreement to those calculated from the numerical solution of the the kinetic equation in the case of δ=10. The improvement of this agreement could be probably reached when using the more precise values for the pressure and temperature jumps at the interface. The widely used Schrage expression for the mass flow rate was compared to that derived from the analytical solution of the Navier-Stokes equations and also with the numerical results obtained from the linearized kinetic equation when σ=0.1, 0.5 and 1. It was found that the analytical solution of the Navier-Stokes equations with the jump conditions provides the values of the mass flow rate closer to the values obtained from the numerical solution of the S-model kinetic equation than tha values predicted by the Schrage expression.

The non-linear S-model kinetic equation was solved by the implicit TVD type numerical method, which does not need any special treatment for the discontinuity of the molecular velocity distribution function at the interface. The results obtained from the linearized and non-linear approaches were compared for several typical ratios of temperatures and the number densities. It was found that the linearized approach could be applied when the temperature difference between the droplet surface and the gas in equilibrium far from the doplet does not exceed 5%. If this difference is larger than 5%, the non-linear approach must to be applied.

The results obtained in the present work show us a strong dependence of the flow fields on the evaporation-condensation coefficient so that one cannot always rely on the results for complete evaporation-condensation at the liquid-vapor interface. In fact, a better understanding of the evaporationcondensation at the interface of a spherical droplet relies on a better description of the phase transition at the interface of the droplet and further investigation on this research topic must be encouraged. Moreover, since the results available in the literature for non-complete evaporation-condensation at the interface are still scarce, the results obtained in the present work represent a significant contribution towards a better understanding of the evaporation-condensation phenomena at the interface of spherical droplets. (A.12)

The functions I m are defined as

I m (δS) = ∞ 0 c m p exp -c 2 p - δS c p dc p , (A.13)
and are calculated by using the power series representation given in Ref. [START_REF] Abramowitz | Handbook of Mathematical Functions with formulas, graphs and mathematical tables[END_REF]. 49), free molecular regime. b Eq. ( 72), hydrodynamic regime.

Table 2: Heat fluxes, q (n) and q (T) , at the liquid-vapor interface (r=1) in case of complete evaporation-condensation (σ=1). Comparison with the numerical data of Ref. [START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF].

-q (n) q (T) δ present work Ref. [START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF] 50), free molecular regime. b Eq. (71), hydrodynamic regime.

Table 3: Mass flow rates from the sphere in case of complete evaporation-condensation (σ=1). Comparison with the numerical data of Ref. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF]. Here, the parameter k is related to the rarefaction parameter δ as k=1/(γ 1 δ) with γ 1 =1.270042.

r 2 u (P ) -r 2 u (T ) k present work Ref. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF] Comparison with the numerical data of Ref. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF]. Here, the parameter k is related to the rarefaction parameter δ as k=3/(2γ 2 δ) with γ 2 =1.92284.

-r 2 q (P ) r 2 q (T ) k present work Ref. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF] Table 5: Comparison between the results at the interface liquid-vapor due to the driven force X n : kinetic equation versus Navier-Stokes when δ=10.

u (n) -q (n) σ Eq. (67) Eq. ( 72) kin. eq. Eq. (71) Eq. ( 74 u (T ) q (T ) σ Eq. (67) Eq. ( 72) kin. eq. Eq. (71) Eq. ( 74 5.573×10 -4 1.01

1.1 3.232×10 -3 3.118×10 -3 3.185×10 -3 1.05 1.5 1.616×10 -2 1.559×10 -2 1.593×10 -2 σ=0.5 1.01 1.01 3.527×10 -3 2.821×10 -3

3.242×10 -3 1.01

1.1 1.940×10 -2 1.975×10 -2 1.858×10 -2 1.05 1.5 9.700×10 -2 9.873×10 -2 9.270×10 -2 σ=1

1.01 1.01 9.410×10 -3 8.463×10 -3 8.157×10 -3 1.01 1.1 5.175×10 -2 5.924×10 -2 4.664×10 -2 1.05 1.5 2.590×10 -1 2.960×10 -1 2.332×10 -1

  ν

  4)I 3 (δS) dθ.

Figure 4 :

 4 Figure 1: Rarefaction parameter δ=0.1

Table 1 :

 1 Bulk velocities, u (n) and u (T ) , at the liquid-vapor interface (r=1) in case of complete evaporation-condensation (σ=1). Comparison with the numerical data of Ref.[START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF].

		u (n)		u (T)
	δ	present work Ref. [26] present work Ref. [26]
	0 a	0.282	-	0.141	-
	0.01	0.283	0.284	0.141	0.142
	0.1	0.289	0.293	0.147	0.148
	0.5	0.312	0.317	0.172	0.176
	1	0.330	0.336	0.203	0.206
	5	0.401	0.402	0.333	0.338
	7	0.415	0.415	0.361	0.366
	10	0.424	0.426	0.385	0.390
	∞ b	0.470	-	0.470	-

a Eq. (

Table 4 :

 4 Heat fluxes from the sphere in case of complete evaporation-condensation (σ=1).

	present work Ref. [24]

Table 6 :

 6 Comparison between the results at the interface liquid-vapor due to the driven force X T : kinetic equation versus Navier-Stokes when δ=10.

	) kin. eq.

Table 7 :

 7 Comparison between the results obtained for the mass flow rate in the hydrodydnamic regime with those provided by the Schrage relation

	) kin. eq.

Acknowledgments

Vinicius Braga Leite and Denize Kalempa acknowledge the Fundação de Amparo à Pesquisa do Estado de São Paulo (FAPESP), grants 2016/20716-9, 2016/09118-2 and 2015/20650-5, for the support of the research.

Supplementary material

Supplementary material is available. Figures 1 and2 show the comparison between our results and those provided in Refs. [START_REF] Takata | Evaporation from or condensation onto a sphere: numerical analysis of the Boltzmann equation for hard-sphere molecules[END_REF][START_REF] Chernyak | The kinetic theory of heat and mass transfer from a spherical particle in a rarefied gas[END_REF] for the mass and heat flux from the sphere as functions of the rarefaction parameter and in case of complete evaporationcondensation at the interface.

Figure 3 shows the profiles of the perturbation functions h (n) and h (T ) on the angle θ at r=1.975 and molecular speed c p ≈ 2.3 obtained from the linearized approach for small values of the rarefaction parameter and complete evaporation-condensation at the interface. The solutions obtained from the DVM with and without the split method are presented for comparison.

Figure 4 shows the profile of the normalized distribution function f on the angle θ at r=1.975 and c p ≈ 2.3 obtained from the linearized and nonlinear approaches for small values of the rarefaction parameter and complete evaporation-condensation at the interface.

Appendix A. Split method to solve the linearized kinetic equation

The differential equations given by (54) for each thermodynamic force and subject to the corresponding boundary condition, given by ( 55) and (56), are solved analytically via the method of characteristics and the solutions are presented as follows

where

is the distance between a point in the gas flow domain with Cartesian coordinates (x, y, z) and a point on the spherical surface with Cartesian coordinates (x 0 , y 0 , z 0 ) directed towards -c. The angle θ 0 is given by

Thus, from ( 21)-( 24), the moments of the functions h (i) 0 are written as

Table 8: Comparison between linear and non-linear results for the bulk velocity and heat flux at the interface when T /T 0 =1.01 and n/n 0 =1.01.

u q δ linear non-linear linear non-linear σ=0.1 0.1 4.268×10 -4 4.284×10 -4 5.362×10 -3 5.364×10 -3 1 4.614×10 -4 4.617×10 -4 4.603×10 -3 4.679×10 -3 10 5.573×10 -4 5.626×10 -4 1.510×10 -3 1.549×10 -3 σ=0.5 0.1 2.152×10 -3 2.158×10 -3 4.867×10 -3 4.523×10 -3 1 2.454×10 -3 2.453×10 -3 3.795×10 -3 3.889×10 -3 10 3.242×10 -3 3.265×10 -3 1.277×10 -3 1.293×10 -3 σ=1. 0.1 4.360×10 -3 4.356×10 -3 3.403×10 -3 3.446×10 -3 1 5.339×10 -3 5.323×10 -3 2.679×10 -3 2.746×10 -3 10 8.157×10 -3 8.169×10 -3 8.500×10 -4 8.557×10 -4

Table 9: Comparison between linear and non-linear results for the bulk velocity and heat flux at the interface when T /T 0 =1.01 and n/n 0 =1.1.

Table 10: Comparison between linear and non-linear results for the bulk velocity and heat flux at the interface when T /T 0 =1.05 and n/n 0 =1.5.

1.803×10 -2 10 1.593×10 -2 1.669×10 -2 6.407×10 -3 6.394×10 -3 σ=0.5 0.1 7.502×10 -2 6.889×10 -2 -9.116×10 -3 -6.904×10 -3 1 8.086×10 -2 7.594×10 -2 -7.576×10 -3 -6.069 ×10 -3 10 9.270×10 -2 9.037×10 -2 -2.745×10 -4 1.414 ×10 -3 σ=1. 0.1 1.520×10 -1 1.125×10 -1 -4.690×10 -2 -3.977×10 -2 1 1.759×10 -1 1.472×10 -1 -4.437×10 -2 -3.376×10 -1 10 2.332×10 -1 1.972×10 -1 -1.250×10 -2 9.454×10 -3