Numerical calculation of the condensational growth of liquid particles in non-dilute and non-ideal media

International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

Contents lists available at ScienceDirect

International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt

Numerical calculation of the condensational growth of liquid particles in non-dilute and non-ideal media F.J. Fernández a,⇑,1, M.M. Prieto a,1, J.M. Fernández-Díaz b a b

Department of Energy, University of Oviedo, Campus de Viesques, 33204 Gijón, Spain Department of Physics, University of Oviedo, C/Calvo Sotelo, 33007 Oviedo, Spain

a r t i c l e

i n f o

Article history: Received 21 January 2014 Received in revised form 9 June 2014 Accepted 12 June 2014

Keywords: Particle Condensational Growth Simultaneous heat and mass transfer Numerical solution Finite volume method

a b s t r a c t Condensational growth in aerosols involves simultaneous heat and mass transfer between particles and their surrounding medium. There are several ways of modelling this phenomenon: dilute medium in the vapour species, ideal gas behaviour, ideal mixture, constant thermophysical properties, etc. Most of the solutions proposed for this problem are analytical, at the cost of simplifying many aspects of the phenomenon. This study develops a numerical method for the solution of the steady-state condensational growth of aerosol liquid particles submerged in a gaseous mixture of condensable vapour and inert (noncondensing) gas. Following the principles of the finite volume method, the mass and energy conservation and the droplet heat balance equations are discretized in a spherical mesh around the particle, obtaining their solution by means of a matrix procedure. Heat and mass are transferred satisfying the first-order phenomenological equations. Transport and thermophysical properties of the mixture can be calculated independently, avoiding other assumptions required in analytical procedures. In order to test the method, it was compared with two analytical solutions for the non-dilute condensational growth problem. One of these methods was applied as formulated in the literature, but the other was significantly improved, applying new factors to separate heat and mass transfer as functions of temperature and composition, respectively. The comparison was performed for four condensing substances (H2 O, R-134a, n-pentane and n-octane) with air as the inert gas. The results show good agreement under conditions for the ideal gas law, with differences in the case that deviates slightly from this behaviour. Ó 2014 Elsevier Ltd. All rights reserved.

1. Introduction The interaction between liquid particles and their surrounding gaseous medium is a key aspect in the study of aerosols. Consisting of a continuous gas phase with a discrete liquid phase suspended within it, aerosols of this kind, also known as fogs, are very important in atmospheric physics, but can also be present in other systems. The term ‘condensational growth’ describes the process by means of which a droplet grows or shrinks due to the condensation or evaporation on its surface of the species that forms it. The study of the condensational growth of aerosol particles can be addressed in various ways. Some authors use techniques based on the method of moments [1,2], while others use molecular kinetics [3,4]. If changes in droplet size or composition and/or in the bulk medium are expected to be very fast, non-steady state ⇑ Corresponding author. Tel.: +34 985182112; fax: +34 985182143. E-mail addresses: [email protected] (F.J. Fernández), manuelap@uniovi. es (M.M. Prieto), [email protected] (J.M. Fernández-Díaz). 1 Group associated with the Spanish Scientific Research Council (CSIC). http://dx.doi.org/10.1016/j.ijheatmasstransfer.2014.06.034 0017-9310/Ó 2014 Elsevier Ltd. All rights reserved.

versions of the energy balance are proposed, with a simpler formulation of heat and mass transfer than in first-order transport equations [5–7]. Among the latter, the relative velocity between the medium and the particle is sometimes considered and the problem is studied in terms of convective heat and mass transfer, like studies in drug administration [8,9], particles in combustion engines [10] and other more general studies [11]. In other cases, purely diffusive heat and mass transfer are considered [12,13]. Non-steady particle growth caused by diffusive heat and mass transport between phases in gas–liquid contact devices with multiple species involved was studied by Schaber et al. [14–16]. Condensational growth of aerosols in atmospheric physics is usually based on the assumption of low vapour concentration (dilute medium). It is formulated in this way in benchmark books in this field [17,18] and is a frequent assumption in particles surrounded by moist air. There are many other situations in which the composition of the gaseous mixture can have a broad range, so the dilute assumption is not applicable. The first theoretical basis of droplet growth under the dilute medium assumption was established by Maxwell [19], with heat

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

1103

Nomenclature c cp Dg DT dc E FT ; FM G g^ h IT T IT Y i j ^| K1; K2 Kn Kt k M m _ m N NC nm p Q_ q Rp RU r S T

Mole concentration Specific heat capacity at constant p Binary diffusion coefficient Thermal diffusion coefficient Collision diameter Energy Factors of dependence in eqs. with separated variables Geometrical factor in mesh definition Molar Gibbs free energy Specific enthalpy Energy conservation indep. terms Continuity independent terms Cell counter Mass flux Mole flux Factors for linear law in Knudsen number Ratio of diffusivities Thermal conductivity Molecular Mass Mass Mass flow rate Number of moles Number of mesh cells Molecular concentration Pressure Heat transfer rate eat flux Particle radius Ideal gas universal const. Radial position Saturation ratio Temperature

and mass transport described by Fourier’s law of heat conduction and Fick’s law of diffusion, respectively, considering a motionless droplet within the uniform medium following a steady state process. The more general first-order transport equations in Wagner [20] formulate condensational growth considering simultaneous heat and mass transfer, Stefan-flow, thermal diffusion and the Dufour effect around the droplet. Researchers have applied different hypotheses to obtain analytical solutions of these first-order transport equations over time. Wagner [20] assumed a zeroth-order (Fourier and Fick’s laws) diffusion problem to obtain a solution. Kulmala and Vesala [21] and Heidenreich [22] assumed a relationship between concentration and temperature fields based on a pure diffusion conception of the problem (although Heidenreich [22] also provided a more complex option). The well-known solution of Barrett and Clement [23] starts from the transport equations without thermal diffusion or the Dufour effect and includes radiative heat exchange on the droplet surface. An important simplification in analytical methods arises from the way they address the formulation of thermophysical properties, together with the equation of state defining the behaviour of the state variables. The problem is defined by the fields of temperature and composition of the gas mixture around the droplet, so the value of thermophysical properties should be obtained as a function of these variables. Difficult to formulate heat and mass transport properties are also needed. For example, Wagner [20] and Kulmala and Vesala [21] consider the exponential temperature dependence of the diffusion coefficient, while Heidenreich [22]

t Ur ~ ur V v v^ y

Time Energy change by surface tension Unitary radial vector Volume Velocity Specific molar volume Gas phase mole fraction

Greek Symbols DHfg Specific latent heat of condensation CT ; HT ; XT Energy conservation matrix terms CY ; HY ; XY Continuity matrix terms at Thermal diffusion factor aT ; aM Mass and heat accommodation coeffs. bT ; bM Factors for non-continuum regime hX ; hT Factors for heat transfer eq. k Mean Free Path l Dynamic viscosity nX ; nT Factors for mass transfer eq. q Density rlg Liquid–gas surface tension sT Factor for heat transfer eq. u Fugacity coefficient Subscripts b g in l p sat v

bulk gas mixture inert gas liquid particle saturation vapour

calculates thermal conductivity as a function of temperature and composition. In these cases, the ideal gas equation of state is implicit in the solutions. When the composition and temperature profiles around the droplet are considered to remain stationary for the solution of the heat and mass transfer, the process can also be considered as stationary. This is a common approach in condensational growth that can be applied when the relaxation time for temperature and composition profiles [18] is much shorter than the time taken by significant changes in the boundary conditions [24,25]. Known as the quasi-steady state approach, it is assumed in classical studies like Wagner [20], Heidenreich [22] and Kulmala and Vesala [21]. It is also assumed in studies dealing with binary condensational growth like in Kulmala et al. [26], whose analytical solution was recently applied by Davies et al. [25], and in others related to simulation of the behaviour of cloud condensation nuclei counters [27,28] in which condensational growth is calculated via the theory of diffusional growth [18]. During a study of aerosols in condensers of pure substances mixed with inert gas, the need arose for an accurate tool to calculate the heat and mass exchange between continuum and discrete phases. Inside a condenser, the assumption of the vapour being the minority species is not applicable; it could even be the main species. Obtaining one solution of the first-order phenomenological equations in Wagner [20] as precise as possible became the objective that has given rise to this paper. A numerical method based on the finite volume method (FVM) has been developed to achieve this aim, resulting in a highly flexible and precise procedure. The

1104

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

FVM proposes the discretization of equations over finite control volumes, ensuring conservation of mass and energy over the full calculation domain [29]. This makes the method suitable for solving problems of a convective–diffusive nature. The use of a numerical method allows the thermophysical properties to be calculated separately from the solution of the transport equations, assigning local values to each point of a mesh scheme around the particle. Ideal gas and ideal mixture assumptions can be avoided. These simplifications are not acceptable under high pressure conditions, high molecular weight vapours, important differences between condensing and inert species, etc. The method was validated by comparing its results with those obtained using analytical solutions of the problem. The first method used for this validation is based on the solution given by Wagner [20]. The original method separates fields of composition and temperature in the transport equations by means of simplifications that have been avoided here, obtaining improved factors for the equations. The second analytical method can be found in two important benchmark studies in condensational droplet growth in non-dilute medium: Kulmala and Vesala [21] and Heidenreich [22]. The comparison is carried out with results corresponding to four condensing species (H2 O, R-134a, n-octane, n-pentane) mixed with air. This paper proposes the numerical solution of the non-dilute first-order transport equations for unary condensational growth, using FVM techniques. Isothermal conditions will be assumed inside the droplet. Particles will move floating with the gas stream with no significant relative velocity. The conservative equations and the heat and mass balances will be formulated in the stationary form. No equations of state or simplification of fluid properties will be inserted in the equations. The appropriate corrections derived from the non-continuum regime and the effect of curvature in the equilibrium of partial pressure in the droplet interface will be considered. The aim is to contribute to the study of aerosol behaviour when attention is focused on special conditions of the gas mixture: high pressures, precise formulation of mixture properties, wide range of compositions, etc. 2. Theory 2.1. Heat and mass transfer around a liquid particle The physical problem is defined by a spherical liquid particle with radius Rp surrounded by a gas mixture of vapour and inert gas. At a distance from the droplet, the bulk temperature and mole fraction of the vapour are T b and yv ;b , respectively. Heat and mass transfer rates and fluxes are defined as positive when exiting the particle. The mass flux density of vapour can be expressed in terms of the velocity of the vapour relative to the center of the droplet: ~ |v ¼ qv ~ vv . The inert gas condition for the other species results in zero velocity relative to the droplet. According to the formulation of the first-order phenomenological equation for the mass transfer in the absence of pressure gradients and external forces [20], we have:

~ vv ¼ 

Dg yv yin



ryv þ

 Kt rT : T

Kt ¼

DT : Dg

ð2Þ

This ratio can be calculated by means of the thermal diffusion factor, at , and both are related through the product of the mole fractions of the species [30]:

K t ¼ at yv yin :

ð3Þ

Within the temperature gradient, the two species move in opposite directions, so the sign of at is not arbitrary. Heavier molecules generally move to cooler areas. Hence, with the sign defined in Eq. (1), the thermal diffusion factor, at , is positive when the vapour is heavier than the inert gas, and negative otherwise. The first-order phenomenological equation for energy transfer in the binary mixture, also formulated by Wagner [20], once the inert gas condition has been applied, is:

~ q ¼ kg rT þ pg K t ~ v v þ q v hv ~ vv :

ð4Þ

The first term on the right-hand side of Eq. (4) is conductive. The last term is convective, the energy being transported in form of the enthalpy of the vapour moving with velocity ~ vv . The term in the middle of the right-hand side of Eq. (4) is the heat flux density originated by the Dufour effect: reciprocal to the Soret effect, within the gradient of composition, a thermal gradient will appear in opposition to the original. Both effects are determined by the same ratio, K t . The heat and the mass transfer depend on the gradients of both temperature and composition. Eqs. (1) and (4) can be modified so as to be directly expressed in terms of these gradients:

~ |v ¼ ðnY ryv þ nT rT Þ; ~ q ¼ ðhY ryv þ hT rT Þ: Factors nY ; nT ; hY ; hT and

ð5Þ ð6Þ

sT have been defined for this purpose:

Dg qv ; yv yin Dg at qv nT ¼ ; T hY ¼ pg Dg at sT ;

nY ¼

pg a2t yv yin Dg sT hT ¼ kg þ ; T q v hv sT ¼ 1 þ : pg at yv yin

ð7Þ ð8Þ ð9Þ ð10Þ ð11Þ

2.2. Energy balance in the droplet The heat and mass transfer must achieve an energy balance at the droplet surface. The droplet energy experiences a change resulting from these fluxes. Its mass and temperature and the surface area of the liquid–gas interface vary. Following the signs in Fig. 1, the heat transfer and mass flow rate on the droplet surface are:

ð1Þ

The first term on the right-hand side of Eq. (1) describes the diffusion of vapour through the non-condensing gas induced by the concentration gradient, following Stefan flow. The second term corresponds to the flux originated by thermal diffusion. This phenomenon, also known as the Soret effect, expresses the movement of the species in the gas mixture when a thermal gradient exists. The transport property involved is the thermal diffusivity, DT , though Eq. (1) is expressed in terms of the ratio:

Fig. 1. Heat flux density and mass flux density in a radial position around the particle.

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

Q_ p ¼ Q_ ðRp Þ ¼ 4 p R2p qðRp Þ; _ v ðRp Þ ¼ 4 p _ v ;p ¼ m m

R2p jv ðRp Þ :

ð12Þ ð13Þ

The sum of rates of change in the particle energy due to the change in mass, temperature and surface tension equals the net energy transfer rate to the droplet. The energy change by unit mass due to the change in the liquid–gas interface surface area is:

Ur ¼

2 rlg : q l Rp

ð14Þ

Q_ p includes the convective term, so the energy balance in the particle is:

mp

 dhl ðT p Þ  _ v ;p ¼ Q_ p : þ hl ðT p Þ þ U r m dt

ð15Þ

The first term on the left-hand side of Eq. (15) is the variation in energy by temperature change in the droplet. It is not considered in steady state studies. Ignoring this term and the interfacial energy change, the balance given in Wagner [31] is obtained. The form of the energy balance in the droplet surface as used in the numerical analysis is:

_ v ;p ¼ Q_ p : ðhl ðT p Þ þ U r Þ m

ð16Þ

2.3. Conservative equations The equations for the conservation of mass and energy can be formulated in the space around the droplet:

dqg |v ¼ 0; þ r ~ dt  d q ¼ 0: q hv þ qin hin þ r  ~ dt v

ð17Þ

r ~|v ¼ 0; q ¼ 0: r ~

p0v ;sat ðT p Þ ¼ pv ;sat ðT p Þ exp

ð19Þ ð20Þ

2.4. Steady and dilute condensational growth: the Maxwell solution The Maxwell solution of condensational growth assumes a dilute medium (yv ! 0) and the heat and mass transfer equations are simplified to Fourier’s law of heat conduction: ~ q ¼ kg rT and ^|v ¼ Dg rcv . The steady state conFick’s law of binary diffusion: ~ tinuity and energy conservation equations can be integrated from the particle surface to an ‘infinite distance’ (bulk). Vapour mole concentrations and temperatures are ðcv ;p ; T p Þ and ðcv ;b ; T b Þ at the droplet surface and bulk positions, respectively. Considering a constant binary diffusion coefficient (Dg ) and thermal conductivity of the gas mixture (kg ), the vapour mole concentration and temperature result in hyperbolic profiles:



 2 rlg Mv : Rp RU T p q l

ð23Þ

Frequently referred to as Kelvin’s equation, it is obtained from the state of metastable equilibrium during the homogeneous nucleation of new droplets within a supersaturated medium [32,18]. It is deduced from the Gibbs free energy change when a liquid droplet is formed:

DG ¼ Np ðg^l  g^v Þ þ 4 p R2p rlg :

ð24Þ

^v Þ is integrated at constant temperature from a flat If ðg^l  g interface with vapour pressure pv ;sat to a droplet with a curved interface and vapour pressure pv , given that ðg^l  g^v Þpv ;sat ¼ 0 and ^ dp: dg^ ¼ v

Z

pv

pv ;sat

dðg^l  g^v Þ ¼ g^l  g^v ¼

Z

pv

pv ;sat

ðv^ l  v^ v Þdp:

ð25Þ

In the deduction of the Kelvin equation, the specific molar volume of liquid is neglected and ideal gas behaviour is assumed at this point. For the sake of consistency, these assumptions should be avoided. The integral on the right-hand side of Eq. (25) must be calculated. An iterative process and external calculation of the thermodynamic properties of fluids may be advantageously employed to solve it. If particles already exist, saturation rates:

S¼

_ v, Eqs. (19) and (20) mean that heat and mass transfer, Q_ and m remain constant at any radial position outside the droplet.

cv ðrÞ  cv ;b Rp ¼ ; cv ;p  cv ;b r TðrÞ  T b Rp ¼ : Tp  Tb r

equilibrium with a convex liquid surface is higher than the partial pressure on a flat surface at the same temperature [18]. This effect is known as the Kelvin effect and is usually considered by means of a factor that corrects the saturation vapour pressure on a flat surface, pv ;sat , to a curved interface, p0v ;sat :

ð18Þ

The density of the binary mixture is the sum of the densities of its components: qg ¼ qv þ qin . The first terms on the left-hand side of Eqs. (17) and (18) are time derivatives. For the steady state, the conservative equations are:

1105

pv ; pv ;sat ðTÞ

ð26Þ

as high as in homogeneous nucleation processes are not expected, letting pv ;sat and pv be very close. Under these conditions the integral can be approximated by ðv^ l  v^ v Þðpv þpv ;sat Þ=2 ðpv  pv ;sat Þ. Obtaining the maximum of DG in a similar way to the deduction of the Kelvin equation, a new expression can be obtained for the saturation vapour pressure of curved interfaces for non-ideal gas behaviour:

2 rlg qv p0v ;sat ’ pv ;sat þ ; Rp ðql  qv Þ

ð27Þ

where qv and ql are evaluated at ðpv ;sat þ p0v ;sat Þ=2 and the interface temperature, T p . Due to the simplification in the calculation of the integral, the expression in Eq. (27) is valid for droplets much bigger than the particles in processes of homogeneous nucleation surrounded by a medium whose saturation rate is close to unity. A real approach to the gas mixture would force the formulation of the equilibrium in the interface in terms of the fugacity coefficient of the vapour, uv [33]:

yv ;p ¼

p0v ;sat ðT p ; Rp Þ : pg uv ðyv ; T p ; pg Þ

ð28Þ

ð21Þ

2.6. Corrections for non-continuum effects

ð22Þ

The phenomenological Eqs. (1) and (4) are valid in the continuum regime. The effects of kinetic or transition regimes must be considered since the droplet size could be comparable with the mean free paths of the vapour or the gas mixture. The nondimensional Knudsen number, Kn, evaluates the dynamic regime of particles through the comparison of the mean free path, k, and the particle radius: Kn ¼ k=Rp .

2.5. Kelvin effect in a non-ideal gas A boundary condition for vapour partial pressure is needed at the vapour–liquid interface. The vapour partial pressure in

1106

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

Following the criterion proposed by Jacobson [17], the corrections to mass and heat transfer can be implemented through correction factors bM and bT by multiplying the parameters of molecular diffusion and heat conduction, Dg and kg :

in the gas–liquid interface is the energy balance equation in the droplet, Eq. (16). The other boundary condition at this position is Eq. (28), that relates temperature and composition at this point.

D0g ¼ bM Dg ;

ð29Þ

3.2. Spatial discretization

0 kg

ð30Þ

¼ bT kg :

The corrections in [17] are due to collision geometry and sticking probability. Expressions for these two corrections are [34]:

1 þ Knv  ; 1 þ 3 aM þ 0:377 Knv þ 3 a4M Kn2v 1 þ Kng   : bT ¼ 1 þ 3 4aT þ 0:377 Kng þ 3 4aT Kn2g

bM ¼



4

ð31Þ ð32Þ

After experimental verification, Laaksonen et al. [35] recommend unitary values of the mass and heat accommodation coefficients (aM ; aT ) for condensation theory in cloud modelling. Two different Knudsen numbers are given for mass and heat transfer. In mass transfer, the vapour mean free path should be considered, resulting in:

Knv ¼

kv ;in : Rp

ð33Þ

The appropriate mean free path in heat transfer is that of the gas mixture:

Kng ¼

kg : Rp

ð34Þ

kv ;in can be calculated using the expression given in Davis [36], kg can be calculated with the expression from Seinfeld and Pandis [18]:

1 kv ;in ¼ pffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffi  2 ; d þd 2 2p nm;v dc;v þ p MvMþMin nm;in c;v 2 c;in in sffiffiffiffiffiffiffiffiffiffiffiffiffi 2 l g p RU T kg ¼ : 8 Mg p

ð35Þ

ð36Þ

3. Numerical method 3.1. General description of the numerical method The phenomenological Eqs. (5) and (6) for heat and mass transfer can be included in the steady state conservative Eqs. (19) and (20):

r  ðnY ryv Þ þ r  ðnT rT Þ ¼ 0; r  ðhY ryv Þ þ r  ðhT rT Þ ¼ 0:

ð37Þ ð38Þ

The conservative equations are then formulated as functions of gradients of temperature and composition. Authors seeking analytical solutions [20–22] have used different strategies to separate variables. One numerical method based on the FVM can be applied to these equations, thereby avoiding further assumptions. The vapour mole fraction gradient can be considered as the only variable in the continuity Eq. (37), leaving the terms in rT as independent. The opposite is done in Eq. (38), i.e. the energy conservation equation. The yv field obtained from Eq. (37) is transferred to Eq. (38). Eq. (37) is updated with the temperature field obtained from Eq. (38). Successive applications of this procedure will converge in temperature and composition fields that meet both equations. The bulk temperature and composition (T b ; yv ;b ) are boundary conditions of the problem. One of the boundary conditions needed

The space around the droplet is divided into control volumes. Exploiting the spherical symmetry of the problem, the division is performed by means of concentric spherical surfaces. The control volumes are the regions delimited by two successive spherical surfaces, as shown in Fig. 2. A geometric law with factor Gwas selected for the definition of the radial position of the nodes. N C is the total number of cells or nodes. The lower and upper limits of a control volume i are at radial positions ri and r iþ , situated on the geometrical mean of the radial positions of adjacent nodes.

8 r ¼ G ri > < iþ1 pffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffi ffi r iþ ¼ ri riþ1 ¼ r i G > :  pffiffiffiffiffiffiffiffiffiffiffiffi priffiffiffi r i ¼ ri1 r i ¼ G

ði ¼ 1; 2; . . . N C  1Þ; ði ¼ 1; 2; . . . N C Þ;

ð39Þ

ði ¼ 2; 3; . . . N C Þ:

The radius of the particle, RP ð¼ r 1 Þ, and a distant radial position, r b ð¼ r NC þ Þ, far enough away from the particle, constitute the lower and upper limits, where the boundary conditions are to be applied. 3.3. Discretization of conservative equations Following the normal procedure in the numerical FVM, the conservative equations are integrated over each control volume and the resulting equations are discretized and rearranged in a tridiagonal matrix system. In the integration of the continuity Eq. (37), the terms multiplying mole fractions at the nodes are the factors in the tridiagonal matrix, the terms multiplying temperatures at the nodes, together with the mole fraction boundary conditions (yv ;p ; yv ;b ) are incorporated in the column matrix of independent terms, IT Y;i . The unknown variables are the mole fractions at the nodes, yv ;i :

2

CY;1 HY;1 6 XY;2 CY;2 6 6 6 0 XY;3 6 6 . 6 .. 6 6 4 0

. ...

0

...

..

0

...

HY;2 CY;3

...

..

.

... .. .

XY;NC 1 CY;NC 1 0 XY;NC

3 2 32 y 3 IT Y;1 v ;1 7 6 yv ;2 7 6 0 7 IT Y;2 7 76 7 7 6 76 7 6 7 6 0 76 yv ;3 7 6 IT Y;3 7 76 7 7¼6 7 7: 6 .. 7 6 .. .. 76 7 7 6 . . . 76 7 6 7 76 7 6 7 5 4 HY;NC 1 4 yv ;NC 1 5 IT Y;NC 1 5 CY;NC yv ;NC IT Y;NC 0

ð40Þ Similarly, the energy balance, Eq. (38), integrated over the control volumes of the mesh scheme yields another tridiagonal system. In this case, the temperature terms are the factors in the

Fig. 2. Scheme of spatial discretization around the particle.

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

tridiagonal matrix, the temperatures at the nodes, T i , are the unknown variables and the independent terms IT T;i include the terms multiplying mole fractions and the temperature boundary conditions (T p ; T b ) at the mesh borders:

2

CT;1 HT;1 CT;2 XT;3

6 XT;2 6 6 6 0 6 6 . 6 .. 6 6 4 0 0

..

.

... ...

0

...

HT;2 CT;3

...

..

.

... .. .

XT;NC 1 CT;NC 1 0 XT;NC

3 2 3 IT T;1 T1 7 7 6 6 0 76 T 2 7 6 IT T;2 7 7 76 7 6 7 0 76 T 3 7 6 IT T;3 7 76 7 6 7 76 . 7 ¼ 6 7: .. .. . 7 7 7 6 6 . . 76 . 7 6 7 76 7 6 7 HT;NC 1 54 T NC 1 5 4 IT T;NC 1 5 0

CT;NC

32

T NC

The integration of the conservative equations to obtain the factors in the tridiagonal matrix systems is detailed in Appendix A. 3.4. Application of the energy balance The droplet energy balance under steady state condition, Eq. (16), must be satisfied at the liquid–gas interface as a boundary condition, but it can also be applied at the nodes and surfaces of the mesh scheme around the droplet. Total heat and mass fluxes, _ v , are constant with distance outside the droplet, so the Q_ and m equation may be transformed into:

ð42Þ

Substituting heat and mass fluxes as they are expressed in Eqs. (5) and (6), Eq. (42) can be easily discretized into one equation that relates mole fraction and temperature changes at the nodes:

DT ¼ 

ðhl ðT p Þ þ U r Þ nY  hY Dyv ; ðhl ðT p Þ þ U r Þ nT  hT

where the size of the first cells adjacent to the particle depend on G. The main temperature and composition changes take place in the region close to the particle, hence the mesh scheme must be dense here. Fig. 3(c) shows the convergence with G for a fixed mesh size (rb =Rp ¼ 500). Convenience of mesh sizes above 50 times the particle radius can be deduced from Fig. 3(b) . Values of G below 1.05 are also convenient to obtain precise solutions in accordance with Fig. 3(c). The results presented in this work were calculated with a grid formed by N c ¼ 160 cells and a geometrical ratio G ¼ 1:04, ensuring in all the cases that results are grid independent.

IT T;NC ð41Þ

_ v ðrÞ ¼ Q_ ðrÞ: ðhl ðT p Þ þ U r Þ m

1107

ð43Þ

being DT ¼ T i  T i1 and Dyv ¼ yv ;i  yv ;i1 . The terms nY ; hY ; nT ; hT are evaluated at the border of one control volume. Actually, it is not strictly necessary to use this equation for every node in the mesh. It could be applied as a boundary condition in the droplet surface and the application of conservative Eqs. (37) and (38) would ensure its accomplishment over the entire domain. However, the application of the balance at the nodes is convenient to achieve a faster numerical convergence of the solution.

4. Analytical solutions Two analytical solutions of the phenomenological equations are used to compare the results and validate the correctness of the numerical method. The non-linearity of these solutions mean they have to be solved by means of iterative techniques. Molecular and transition regimes are considered in both methods, correcting with bM and bT in the calculation of Dg and kg , in accordance with Eqs. (29) and (30). 4.1. Analytical solution 1 (AS1) The first analytical solution is based on the method devised by Wagner [20]. This solution is obtained starting from expressions of the heat and mass transfer in terms of separate variables and introducing two factors of dependence, F T and F M :

~ q ¼ kg F T rT; ~ |v ¼ Dg F M rqv :

ð44Þ ð45Þ

In order to separate the variables, Wagner includes a number of assumptions avoided here. Alternative expressions for F T and F M were obtained. The detailed procedure to their obtaining, where ideal gas law is assumed, is the subject of a separate study.

FT ¼

1 

1  hhv

ð46Þ

;

l

FM ¼



pin þ ðhv  

pg ;  Dg M v pv ðpg þpin at Þ hl Þ kg RU T 2



where hl and hv are: 

3.5. Mesh study and boundary conditions A study of the spatial discretization around the particle is performed to achieve a compromise between accurateness and computing time. The mesh is defined from its geometrical ratio, G, and the total number of cells, N C , while its total size can be expressed in terms of the particle radius by the ratio rb =Rp , which is fixed by the geometrical scheme as rb =Rp ¼ GNc . Fig. 3(a) shows the results of the heat flux on the particle surface for several mesh schemes. Large sizes are necessary to obtain a certain degree of stability, so a large number of cells and high computation times are needed. This is mainly due to the position of the boundary condition, since the temperature and gas composition are fixed at a finite distance in their bulk (infinite) values. This effect easily disappears using boundary conditions of a different kind, based on the lower order, also called the dilute or Maxwell solution of heat and mass diffusion, Eqs. (22) and (21). The equivalent results using this technique are shown in Fig. 3(b). The solution converges to the same value, but with much smaller mesh sizes. The solution converges to different values depending on the geometrical ratio, G. This is due to the scheme of discretization

ð47Þ

hl ¼ hl ðT p Þ þ U r ; p p  hv ¼ hv ðTÞ þ at v in : pg qv

ð48Þ ð49Þ

Heat and mass transfer from the droplet to the medium are calculated by means of [20]: Q_ p ¼ 2 p Rp ðkg;b F T;b þ kg;p F T;p ÞðT p  T b Þ; _ v ;p ¼ 2 p Rp Dg;b m

ð50Þ

! T lp D 1 Tb  Tp  ðF M;b þ F M;p Þðqv ;p  qv ;b Þ: ð51Þ ð l  1Þ T b T lD 1  T plD 1 D b

The term lD is the exponent of the temperature in a law of dependence of the diffusion coefficient, Dg , in terms of temperature. F T;p and F T;b are the values of F T evaluated on the particle surface and in the bulk gas, respectively. The solution forces the fulfilment of the energy balance on the droplet surface, including the term of the change in energy with surface tension, _ v ;p ¼ Q_ p . ðhl ðT p Þ þ U r Þ m 4.2. Analytical solution 2 (AS2) Another alternative solution was formulated by Heidenreich [22] very similar to the solution given by Kulmala and Vesala

1108

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

(a)

(b)

1890

(c)

1880

1870

G=1.02 G=1.04 G=1.06 G=1.08 G=1.1

1860

G=1.02 G=1.04 G=1.06 G=1.08 G=1.1

1850 10

100

1000

10000 10

100

rb / Rp

1000

10000 1

rb / Rp

1.05

1.1

1.15

1.2

G

Fig. 3. Heat flux density on the particle surface (H2 O þ air; Dp ¼ 5 lm;T b ¼ 35 C;S ¼ 1:2; yv ;b ¼ 0:4). Solutions for different mesh schemes: (a) T b and yv ;b as boundary conditions in r b ; (b) first-order solution (dilute) as boundary conditions in r b ; (c) convergence of the solution with the geometrical ratio, G, for r b =Rp ¼ 500.

[21]. Heidenreich considers the gas thermal conductivity as a function of temperature and composition, while Kulmala and Vesala consider the variation of the diffusion coefficient with temperature. Both studies assume the simplified relationship between fields of temperature and composition shown below, based on the profiles obtained from Fourier’s and Fick’s laws for heat and mass transfer, respectively.

@T T p  T b @pv ¼ : @r pv ;p  pv ;b @r

ð52Þ

A formulation of the heat balance in the particle was considered that neglects the effect of changes in energy with surface tension, not including the term U r in Eq. (42). The variation in vapour enthalpy is calculated with a constant value for cp : hv ðTÞ ¼ hv ðT p Þ þ cp;v ðT  T p Þ. The thermal conductivity of the gas mixture obeys a linear law (kg ðTÞ ¼ K 1 þ K 2 T) with factors K 1 and K 2 according to Lindsay and Bromley [37]. Ideal gas behaviour is also assumed. Using the formulation proposed by Heidenreich [22], the temperature in the particle can be calculated by means of:

  _ v ;p ðcp;v þ at RU AH Þ DHfg ðT p Þ  cp;v T p K 1 m þ  4 p Rp K 2 cp;v þ at RU AH K2   DHfg ðT p Þ þ at RU AH T p ;  log hv ðT b Þ  hl ðT p Þ þ at RU AH T p

Tp ¼ Tb 

ð53Þ

where:

AH ¼

yv yin yv M v þ yin M in



1þ

 yin Min : yv M v

ð54Þ

To calculate the mass transfer Heidenreich use:

_ v ;p ¼  m

4pRp M v Dg pg RU T b

log

! 1  yv ;p at ðT p  T b Þðyv ;p þ yv ;b Þ  : 2 Tb 1  yv ;b

5. Fluid properties Calculation of the fluid thermophysical properties must be included in the models. The thermal conductivity, heat capacity, viscosity, density and specific enthalpy of the gas phase are needed. The specific enthalpy and density of the liquid in the particle are also needed. Since there is phase coexistence, surface tension, latent heat and saturation vapour pressure are also required. The gas phase is a mixture of vapour and inert gas, so the mass binary and thermal diffusion coefficients (Dg and DT ) are also very important properties. 5.1. Pure substances Air was considered the inert gas in all cases. The results of condensational growth were obtained for several condensable species: H2 O, 1-1-1-2-tetrafluoroethane (R-134a), n-pentane and n-octane. The aim was to test substances with different saturation pressure and molecular mass. The properties of pure substances were taken from the National Institute of Standards and Technology [38] site and stored in text files as a function of pressure and temperature. In the case of n-pentane and n-octane, the properties were obtained from the Engineering Equation Solver [39] software. The properties of pure substances can be easily deduced at local pressure and temperature by programming the appropriate subroutine included in the calculation scheme. The properties of the supersaturated vapour are extrapolated from the values corresponding to the saturated and superheated vapour at the same partial pressure. The saturation pressure of water is obtained from the formulation of Preining et al. [40]. 5.2. Gas mixture

ð55Þ Then the heat transfer can be easily calculated using the aforementioned heat balance in the particle.

The properties of the gas mixture are calculated from the properties of the pure components and its composition. To calculate partial pressures of components and to obtain the density of the gas mixture, unitary fugacity coefficients were assumed for each

1109

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

pair of substances. The specific enthalpy and heat capacity of the mixture are obtained from its mass composition. Viscosity is obtained using the formulation of Wilke [41], applying a similar law [33] for the thermal conductivity of the mixture. Some more detailed explanations relative to the calculation of the binary and thermal diffusion coefficients are needed.

6

5

Rp (µm)

4

5.3. Binary diffusion coefficient The empirical formulation proposed by Marrero and Mason [42] is a good option for the pair of substances H2 O-air. There are equations based on molecular interaction models which are more generally applicable to other pairs. Avoiding the ideal gas law, the expression by Chapman and Enskog in Poling et al. [33] may be written as:

Dg ¼ 3:203 106



T M v ;in

1=2

1 Mv

2 þ M1in

2

1

0 0

Mg

qg r2v ;in XD

20

40

:

60

80

100

120

Growth time(ms)

ð56Þ

;

where M g is the apparent molecular mass of the gas mixture, and:

M v ;in ¼ 

3

Fig. 5. Condensational growth: evolution of a droplet radius. T b ¼ 25 C; pg ¼ 1atm; S ¼ 2:0. — numerical method (continuum regime),- - - numerical method (transition regime), —  —  — Heidenreich [22], case at ¼ 0:034.

ð57Þ

rv ;in ¼ ðrv þ rin Þ=2, incorporated in Eq. (56) in Å. The formulation

A law of intermolecular forces is needed to evaluate rv ;in and XD . Following the Lennard-Jones 12–6 potential, the characteristic length parameter for molecular interaction, rv ;in , may be taken as the arithmetic mean of the species length parameters,

Table 1 Energy and length parameters for the Lennard-Jones 12–6 law of intermolecular forces.

Air H2 O R-134a n-pentane n-octane

αt

kB

e ðKÞ

r(Å)

78.6 809.1 301.1 308.3 357.7

3.711 2.641 5.016 5.916 7.024

(a)

by Neufeld et al. [43] was used to calculate the dimensionless collision integral for diffusion, XD . 5.4. Thermal diffusion coefficient The thermal diffusion coefficient, DT , is locally calculated for each pair of species considering the temperature and composition of the mixture. The thermal diffusion factor, at , is obtained using an expression from the molecular kinetics theory:

at ¼

Kt yv S v  yin S in ¼ 5C 2 : yv yin yv Qv þ y2in Qin þ yv yin Qv ;in

Further information about terms C; S and Q is available in Chapman and Cowling [30]. The Lennard-Jones 12–6 potential law of intermolecular forces is used to calculate them and for the

(b)

0.30 n-octane+Air R134a+Air n-pentane+Air Water+Air

0.20

bar

,

(c)

100

10

0.25

ð58Þ

10 1 n-octane+Air R134a+Air n-pentane+Air Water+Air

0.15

0.1

1

0.10

0.1

T = 40ºC S =1.0

0.05

0.01

n-octane R134a n-pentane Water

0.01

0.00

T = 40ºC S =1.0

-0.05

0.001

0.2

0.4

0.6

yv

0.8

1

0.001

0.2

0.4

0.6 yv

0.8

1

0

20

40

60

80

100

T(ºC)

Fig. 4. (a) Calculated thermal diffusion factor, at and (b) binary diffusion coefficient, Dg , at T ¼ 40 C for a saturated gas mixture (S ¼ 1:0); (c) saturation vapour pressures.

1110

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

Results of the thermal diffusion factor, at , binary diffusion coefficient, Dg , together with the saturation vapour pressures, are shown in Fig. 4 for the substances under study.

Table 2 Mass flux on the droplet surface and droplet temperature calculated by the numerical method (NUM) and the analytical solutions (AS1, AS2). Time (ms) 0.0 3.5 17.5 40.0 71.5

NUM Rp ðlmÞ 0.5000 1.0102 2.0255 3.0142 4.0109   jv ;p s kg m2

AS1

AS2

0.5000 1.0103 2.0253 3.0139 4.0105

0.5000 1.0092 2.0233 3.0112 4.0071

0.0 3.5 17.5 40.0 71.5

0.201129 0.105135 0.053828 0.036476 0.027529

0.201166 0.105125 0.053820 0.036471 0.027524

0.200527 0.104993 0.053779 0.036446 0.027507

0.0 3.5 17.5 40.0 71.5

T p ð Þ C 35.377 35.424 35.446 35.453 35.457

35.378 35.422 35.443 35.450 35.453

35.348 35.391 35.412 35.418 35.421

6. Results 6.1. Growth of a water droplet under atmospheric conditions The results are first compared with others presented in Heidenreich [22]. The evolution of the size of a water droplet under atmospheric conditions (T b ¼ 25 C;pg ¼ 1 atm) in a supersaturated gas mixture (S ¼ 2:0) is calculated. The conservative equations were solved and the droplet size was calculated starting from a 1lm diameter and discretizing in time steps Dt ¼ 0:5 ms. Assuming constant vapour mass flux density throughout the time step, the particle size changes with:

j Dt DRp ¼  v :

ð59Þ

ql

The evolution of the particle size over a growth period of 100 ms is presented in Fig. 5. The similarity with the result obtained by Heidenreich [22] can be appreciated in this figure. The case chosen for this comparison was the one in which the author considered a thermal diffusion factor at ¼ 0:034; i.e. the closest one to the at ¼ 0:0123 calculated with our procedure.

ð1Þ

collision integrals (XD and X12 ð1Þ). The energy and length parameters in this law (keB and r) were taken from Poling et al. [33] for H2 O and air; from Mourits and Rummens [44] for n-alkanes and from Shibasaki-Kitakawa et al. [45] for R-134a (Table 1).

H2O+air

0.0020

n-octane+air 0.008

0.0015

0.006

0.0010

0.004

0.0005

0.002

0.0000

0.000

yyv=0.2 v,b=0.2 yyv=0.3 v,b=0.3 yyv=0.4 v,b=0.4 yyv=0.5 v,b=0.5 yyv=0.6 v,b=0.6 yyv=0.7 v,b=0.7 yyv=0.8 v,b=0.8 yyv=0.9 v,b=0.9

-0.0005 -0.0010 -0.0015

yyv=0.2 v,b=0.2 yyv=0.3 v,b=0.3 yyv=0.4 v,b=0.4 yyv=0.5 v,b=0.5 yyv=0.6 v,b=0.6 yyv=0.7 v,b=0.7 yyv=0.8 v,b=0.8 yyv=0.9 v,b=0.9

-0.002 -0.004 -0.006

-0.0020

-0.008

0.8

0.9

1 S

1.1

1.2

0.8

0.9

0.012

0.04

0.008

0.02

0.004

yyv=0.5 v,b=0.5

1.1

1.2

yyv=0.4 v,b=0.4 yyv=0.5 v,b=0.5 yyv=0.6 v,b=0.6 yyv=0.7 v,b=0.7 yyv=0.8 v,b=0.8 yyv=0.9 v,b=0.9

-0.004

yyv=0.6 v,b=0.6 yyv=0.7 v,b=0.7

-0.04

1.2

0.000

yyv=0.4 v,b=0.4

-0.02

1.1

n-pentane+air

R134a+air 0.06

0.00

1 S

-0.008

yyv=0.8 v,b=0.8 yyv=0.9 v,b=0.9

-0.012

-0.06

0.8

0.9

1 S

1.1

1.2

0.8

0.9

1 S

Fig. 6. Results of the numerical method. Mass flux on the surface of a liquid particle (50 lm diameter); liquid-vapour species: H2 O, n-octane, R-134a and n-pentane; air as inert gas; T b ¼ 40 C.

1111

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

n-octane+air

H2O+air

0.0020

0.008

0.0015

0.006

0.0010

0.004

0.0005

0.002

0.0000

0.000

-0.0005

-0.002

-0.0010

-0.004

-0.0015

-0.006 -0.008

-0.0020 0.8

0.9

1

1.1

S

1.2

0.04

0.8

0.9

1

1.1

S

1.2

n-pentane+air

R134a+air

0.012

0.03 0.008 0.02 0.004

0.01

0.00 0.000 -0.01 -0.02

-0.004

-0.03 -0.008

-0.04 -0.05

-0.012 0.8

0.9

1 S

1.1

1.2

0.8

0.9

1 S

1.1

1.2

Fig. 7. Mass flux on the surface of a liquid particle (50 lm diameter). Comparison among methods: — numerical, - - - AS1 and   AS2. Liquid–vapour species: H2 O, n-octane, R134a and n-pentane; air as inert gas; yv ;b ¼ 0:5 T b ¼ 40 C.

These calculations were also performed by means of AS1 and AS2. As the results obtained by each method were very similar, the curves corresponding to the analytical solutions were not plotted in Fig. 5, since differences would not be appreciated on this scale. Table 2 explains the closeness of the results obtained using the three methods. This is illustrated by the mass flux density on the droplet surface and by the droplet temperature, obtained by the numerical method, AS1 and AS2. The initial situation (Rp ¼ 0:5 lm) and instants corresponding to a radius very close to 1, 2, 3 and 4 lm were selected to be shown here.

46

44

Tp ( ºC)

42

40

38

36

6.2. Condensational growth for different condensing substances

34

32

0.8

0.9

1

1.1

1.2

S Fig. 8. Results of the numerical method: Temperature of the liquid particle (50 lm diameter); liquid–vapour species: — H2 O, - - - n-octane,    R-134a and      npentane; air as inert gas; yv ;b ¼ 0:5 T b ¼ 40 C.

The curve labeled as ‘Continuum Regime’ does not include the corrections of the non-continuum effects (bM ¼ 1:0 and bT ¼ 1:0), as they were not included in Heidenreich’s study either.

Particles of other condensing substances were tested, keeping air as the inert gas. The conditions were also varied in terms of the saturation ratio and gas composition. The coincidence in results obtained with numerical and analytical solutions was tested. Fig. 6 shows the vapour mass flux on the surface of a 50 lm diameter particle, calculated by means of the numerical method. The proportion of vapour to inert gas in the gas mixture and the saturation ratio were varied, while the bulk temperature was kept at 40° C in all cases. The variation in the saturation ratio is hence obtained by changing partial and total pressures.

1112

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

It can be seen that droplets grow under supersaturation conditions (S > 1) and shrink for superheated vapour (S < 1), since the outgoing fluxes are positive. Mass fluxes decrease as the vapour molar fraction increases, mainly due to the change in binary diffusion coefficients with pressure (Fig. 4). Fig. 8 shows the results of droplet temperature calculated with the numerical method for constant bulk gas composition and temperature (yv ;b ¼ 0:5; T b ¼ 40 C). As expected, growing droplets are hotter than the surroundings while shrinking droplets are colder. Higher temperature differences are found in substances in which higher mass fluxes are found, coinciding with those having higher saturation vapour pressures. Under these conditions, the results of the numerical method are very similar to those obtained using the analytical solutions. This is likewise illustrated in Fig. 7 by the mass flux density on the droplet surface; for constant bulk gas composition and temperature (yv ;b ¼ 0:5; T b ¼ 40 C). The plots of different calculation methods for n-pentane and n-octane are indistinguishable. Analytical solution 2 shows a slight difference with respect to the other methods when the pair is H2 O þ air. Another appreciable discrepancy is present in the case of the pair R134a þ air. In this case, there is even a small difference between the numerical method and AS1. This difference could be due to the departure of these gases from ideal gas behaviour, since the pressures for the phase change of R-134a are higher (the saturation vapour pressure of R-134a at 40 C is 10.16 bar). The differences observed in the results of the numerical and analytical methods when the gases depart from ideal behaviour is worth highlighting.

The calculation model was originally conceived as a part of a model for the study of condensation processes in the presence of an aerosol and an inert gas, but it may be applied in a wide variety of aerosol systems. Condensers or evaporators, containment heat exchangers employed in nuclear safety, aerosol generators and nuclei counters in aerosol experimental science are only some examples, although other systems exist where this approach to aerosol condensational growth will surely be of interest. Conflict of interest None declared. Appendix A. Factors in the discretized equations The factors multiplying ryv and rT make up the terms in the tridiagonal matrix systems in Eqs. (40) and (41): CY;i ; HY;i ; XY;i ; IT Y;i ; CT;i ; HT;i ; XT;i and IT T;i . The terms come from the integration of the conservative Eqs. (37) and (38) in each control volume.

Z CV i

Z

CV i

r  ðnY ryv ÞdV þ r  ðhY ryv ÞdV þ

Z

r  ðnT rTÞdV ¼ 0;

ðA:1aÞ

r  ðhT rTÞdV ¼ 0:

ðA:1bÞ

CV i

Z

CV i

After applying the divergence theorem and substituting the spatial derivatives at the borders of the control volume (see Fig. A.9) by a Taylor-series expansion, the expressions of the discretized equations at a central node are obtained:

ðr 2i nY;i GÞyv ;i1 þ ðr 2i nY;i G  r2iþ nY;iþ Þyv ;i þ ðr2iþ nY;iþ Þyv ;iþ1 7. Conclusions

þ ðr 2i nT;i GÞT i1 þ ðr 2i nT;i G  r 2iþ nT;iþ ÞT i þ ðr 2iþ nT;iþ ÞT iþ1 ¼ 0; ðA:2aÞ

A numerical method has been developed to calculate the condensational growth of particles in a non-dilute medium in the presence of an inert gas. Steady state and no significant droplet– gas relative velocity were assumed. The FVM was applied to solve conservative equations in a discretized domain around the particle. The first-order phenomenological equations were included in the conservative equations considering species diffusion by Stefan flow, heat conduction, mass convection, and Soret and Dufour effects, together with a detailed energy balance in the particle. Special attention was paid to the calculation of fluid properties, avoiding the inclusion of parameters from ideal gas law or another equation of state in the equations to be solved. An alternative formulation of the partial pressure equilibrium in the curved liquid– gas interface was formulated to avoid the ideal gas law in Kelvin’s equation. Classical correction factors were included in order to consider non-continuum kinetic and transition regimes. The result is a calculation procedure that can be applied to any condensinginert pair of substances in any range of physical states where their properties are known. It can be easily included in one more general simulation model of systems containing a discrete phase. A first comparison with a classical study of droplet growth has confirmed the validity of the method. Heat and mass transfer results were compared with those obtained by two classical analytical solutions of the first-order phenomenological equations for several condensing substances in a wide range of compositions and saturation rates. In the case of the analytical solution with separate variables, two new factors of dependence were formulated, finding a high level of coincidence with the numerical method. The second analytical solution produces similar results, although with a slightly poorer coincidence. Future research would comprise a more in-depth comparison, testing substances and states that diverge further from ideal gas behaviour, for which the differences would be expected to be larger.

ðr 2i hY;i GÞyv ;i1

þ

ðr 2i hY;i G

þ ðr 2i hT;i GÞT i1 þ

 r2iþ hY;iþ Þyv ;i þ 2 ðr i hT;i G  r 2iþ hT;iþ ÞT i

ðr2iþ hY;iþ Þyv ;iþ1 þ ðr 2iþ hT;iþ ÞT iþ1

¼ 0; ðA:2bÞ

where nY;i is nY in Eq. (7) evaluated at a distance r i , corresponding to the inner border of the control volume, and likewise for the remaining similar terms. The integrations in the first (i ¼ 1) and in the last (i ¼ N C ) control volumes constitute special cases. In the control volume closest to the particle (i ¼ 1):

pffiffiffiffi pffiffiffiffi GÞÞ yv ;p þ ðR2p nY;p ðG þ GÞ  r 21þ nY;1þ Þ yv ;1 pffiffiffiffi þ ðr 21þ nY;1þ Þ yv ;2 þ ðR2p nT;p ðG þ GÞÞ T p pffiffiffiffi þ ðR2p nT;p ðG þ GÞ  r 21þ nT;1þ Þ T 1 þ ðr21þ nT;1þ Þ T 2 ¼ 0; pffiffiffiffi pffiffiffiffi ðR2p hY;p ðG þ GÞÞ yv ;p þ ðR2p hY;p ðG þ GÞ  r 21þ hY;1þ Þ yv ;1 pffiffiffiffi þ ðr 21þ hY;1þ Þ yv ;2 þ ðR2p hT;p ðG þ GÞÞ T p p ffiffiffi ffi þ ðR2p hT;p ðG þ GÞ  r 21þ hT;1þ Þ T 1 þ ðr21þ hT;1þ Þ T 2 ¼ 0:

ðR2p nY;p ðG þ

Fig. A.9. Control volume corresponding to node i : CV i .

ðA:3aÞ

ðA:3bÞ

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

1113

Atpffiffiffithe last node (i ¼ N C ), where the geometric relation ffi r b ¼ G rNC was maintained:

References

pffiffiffiffi ðr2N nY;NC GÞ yv ;NC 1 þ ðr2N nY;NC G  r2b nY;b ð G þ 1ÞÞ yv ;NC C C pffiffiffiffi þ ðr2b nY;b ð G þ 1ÞÞ yv ;b þ ðr2N nT;NC GÞ T NC 1 pffiffiffiffi C pffiffiffiffi þ ðr 2N nT;NC G  r 2b nT;b ð G þ 1ÞÞ T NC þ ðr2b nT;b ð G þ 1ÞÞ T b ¼ 0;

[1] C. Yuan, F. Laurent, R.O. Fox, An extended quadrature method of moments for population balance equations, J. Aerosol Sci. 51 (2012) 1–23. [2] M. Yu, J. Lin, Binary homogeneous nucleation and growth of water-sulfuric acid nanoparticles using a {TEMOM} model, Int. J. Heat Mass Transfer 53 (2010) 635–644. [3] D. Liang, Q. Guo, S. Gong, A new splitting wavelet method for solving the general aerosol dynamics equation, J. Aerosol Sci. 39 (2008) 467–487. [4] V. Levdansky, J. Smolik, P. Moravec, Trapping of impurity molecules in condensation from mixtures of gases, Int. J. Heat Mass Transfer 43 (2000) 629–637. [5] G. Hubbard, V. Denny, A. Mills, Droplet evaporation: effects of transients and variable properties, Int. J. Heat Mass Transfer 18 (1975) 1003–1008. [6] B. Abramzon, W. Sirignano, Droplet vaporization model for spray combustion calculations, Int. J. Heat Mass Transfer 32 (1989) 1605–1618. [7] M.D. Protheroe, A. Al-Jumaily, R.J. Nates, Prediction of droplet evaporation characteristics of nebuliser based humidification and drug delivery devices, Int. J. Heat Mass Transfer 60 (2013) 772–780. [8] P.W. Longest, J. Xi, Condensational growth may contribute to the enhanced deposition of cigarette smoke particles in the upper respiratory tract, Aerosol Sci. Technol. 42 (2008) 579–602. [9] P.W. Longest, M. Hindle, Cfd simulations of enhanced condensational growth (ecg) applied to respiratory drug delivery with comparisons to in vitro data, J. Aerosol Sci. 41 (2010) 805–820. [10] A. Warey, S. Balestrino, P. Szymkowicz, M.R. Malayeri, A one-dimensional model for particulate deposition and hydrocarbon condensation in exhaust gas recirculation coolers, Aerosol Sci. Technol. 46 (2012) 198–213. [11] E.W.C. Lim, S.H. Koh, L.K. Lim, S.H. Ore, B.K. Tay, Y. Ma, C.-H. Wang, Experimental and computational studies of liquid aerosol evaporation, J. Aerosol Sci. 39 (2008) 618–634. [12] S.P. Fisenko, W.-N. Wang, M. Shimada, K. Okuyama, Vapor condensation on nanoparticles in the mixer of a particle size magnifier, Int. J. Heat Mass Transfer 50 (2007) 2333–2338. [13] R. Saleh, A. Shihadeh, Hygroscopic growth and evaporation in an aerosol with boundary heat and mass transfer, J. Aerosol Sci. 38 (2007) 1–16. [14] R. Ehrig, O. Ofenloch, K. Schaber, P. Deuflhard, Modelling and simulation of aerosol formation by heterogeneous nucleation in gas–liquid contact devices, Chem. Eng. Sci. 57 (2002) 1151–1163. [15] K. Schaber, J. Körber, O. Ofenloch, R. Ehrig, P. Deuflhard, Aerosol formation in gas–liquid contact devices—nucleation, growth and particle dynamics, Chem. Eng. Sci. 57 (2002) 4345–4356. [16] A. Wix, K. Schaber, O. Ofenloch, R. Ehrig, P. Deuflhard, Simulation of aerosol formation in gas–liquid contact devices, Chem. Eng. Commun. 194 (2007) 565–577. [17] M.Z. Jacobson, Fundamentals of Atmospheric Modeling, Cambridge University Press, 2005. [18] J.H. Seinfeld, S.N. Pandis, Atmospheric Chemistry and Physics: From Air Pollution to Climate Change, A Wiley-Intersciencie publications, Wiley, 2006. [19] J.C. Maxwell, in: W.D. Niven (Ed.), The Scientific Papers of James Clerk Maxwell, vol. II, Cambridge University Press, 1890, pp. 636–640. [20] P.E. Wagner, Aerosol growth by condensation, in: W.H. Marlow (Ed.), Aerosol Microphysics II: Chemical Physics of Microparticles, Topics in Current Physics, Springer-Verlag, 1982, pp. 129–178. [21] M. Kulmala, T. Vesala, Condensation in the continuum regime, J. Aerosol Sci. 22 (1991) 337–346. [22] S. Heidenreich, Condensational droplet growth in the continuum regime–a critical review for the system air–water, J. Aerosol Sci. 25 (1994) 49–59. [23] J.C. Barrett, C.F. Clement, Growth rates for liquid drops, J. Aerosol Sci. 19 (1988) 223–242. [24] T. Vesala, M. Kulmala, R. Rudolf, A. Vrtala, P.E. Wagner, Models for condensational growth and evaporation of binary aerosol particles, J. Aerosol Sci. 28 (1997) 565–598. [25] J.F. Davies, A.E. Haddrell, J.P. Reid, Time-resolved measurements of the evaporation of volatile components from single aerosol droplets, Aerosol Science and Technology 46 (2012) 666–677. [26] M. Kulmala, T. Vesala, P.E. Wagner, An analytical expression for the rate of binary condensational particle growth, Proc.: Math. Phys. Sci. 441 (1993) 589– 605. [27] G.C. Roberts, A. Nenes, A continuous-flow streamwise thermal-gradient ccn chamber for atmospheric measurements, Aerosol Sci. Technol. 39 (2005) 206– 221. [28] T.L. Lathem, A. Nenes, Water vapor depletion in the dmt continuous-flow ccn chamber: effects on supersaturation and droplet growth, Aerosol Sci. Technol. 45 (2011) 604–615. [29] S. Patankar, Numerical heat transfer and fluid flow, Series in Computational Methods in Mechanics and Thermal Sciences, Hemisphere Pub. Corp., 1980. [30] S. Chapman, T.G. Cowling, Viscosity, thermal conduction and diffusion: general expressions, in: The Mathematical Theory of Non-Uniform Gases: An Account of The Kinetic Theory of Viscosity, Thermal Conduction and Diffusion in Gases, Cambridge mathematical library, Cambridge University Press, 1970, pp. 149– 166. [31] P.E. Wagner, Energy balance during condensational drop growth, J. Aerosol Sci. 22 (1991) 789–791.

C

ðr2N C

pffiffiffiffi G þ 1ÞÞ yv ;NC

ðr2N C

ðA:4aÞ

r2b hY;b ð

hY;NC GÞ yv ;NC 1 þ hY;NC G  pffiffiffiffi 2 þ ðrb hY;b ð G þ 1ÞÞ yv ;b þ ðr2N hT;NC GÞ T NC 1 pffiffiffiffi C pffiffiffiffi þ ðr 2N hT;NC G  r 2b hT;b ð G þ 1ÞÞ T NC þ ðr 2b hT;b ð G þ 1ÞÞ T b ¼ 0: C

ðA:4bÞ The parameters for the equations in the tridiagonal matrix systems are obtained identifying Eqs. (A.2)–(A.4) with the matrix system Eqs. (40) and (41): At central nodes (i ¼ 2; . . . N C  1):

XY;i ¼ r 2i nY;i G;

ðA:5aÞ

CY;i ¼ r 2i nY;i G  r2iþ nY;iþ ;

ðA:5bÞ

HY;i ¼ r 2iþ nY;iþ ;

ðA:5cÞ

IT Y;i ¼ ðr 2i nT;i GÞ T i1  ðr2i nT;i G  r 2iþ nT;iþ Þ T i  ðr 2iþ nT;iþ Þ T iþ1 ;

ðA:5dÞ

XT;i ¼ r 2i hT;i G;

CT;i ¼

r 2i

hT;i G  

ðA:5eÞ r2iþ

ðA:5fÞ

hT;iþ ;

HT;i ¼ r 2iþ hT;iþ ;

ðA:5gÞ

IT T;i ¼ ðr 2i hY;i GÞ yv ;i1  ðr 2i hY;i G  r 2iþ hY;iþ Þ yv ;i  ðr 2iþ hY;iþ Þ yv ;iþ1 :

ðA:5hÞ

At the first node (i ¼ 1):

pffiffiffiffi

CY;1 ¼ R2p nY;p ðG þ GÞ  r 21þ nY;1þ ; r 21þ

HY;1 ¼

nY;1þ ;

pffiffiffiffi pffiffiffiffi GÞÞ yv ;p  ðR2p nT;p ðG þ GÞÞ T p p ffiffiffi ffi  ðR2p nT;p ðG þ GÞ  r 21þ nT;1þ Þ T 1  ðr21þ nT;1þ Þ T 2 ; pffiffiffiffi ¼ R2p hT;p ðG þ GÞ  r 21þ hT;1þ ;

ðA:6aÞ ðA:6bÞ

IT Y;1 ¼ ðR2p nY;p ðG þ

CT;1

HT;1 ¼

r 21þ

hT;1þ ;

pffiffiffiffi pffiffiffiffi GÞÞ T p  ðR2p hY;p ðG þ GÞÞ yv ;p p ffiffiffi ffi  ðR2p hY;p ðG þ GÞ  r 21þ hY;1þ Þ yv ;1  ðr 21þ hY;1þ Þ yv ;2 :

ðA:6cÞ ðA:6dÞ ðA:6eÞ

IT T;1 ¼ ðR2p hT;p ðG þ

ðA:6fÞ

At the last node (i ¼ NC ):

XY;NC ¼ r2N nY;NC G; C

CY;NC ¼

r 2N C

IT Y;NC

ðr 2b nY;b ð

pffiffiffiffi G þ 1Þ;

r2b nY;b ð

nY;NC G  pffiffiffiffi ¼ G þ 1ÞÞ yv ;b  ðr 2N nT;NC GÞ T NC 1 pffiffiffiffi C  ðr 2N nT;NC G  r 2b nT;b ð G þ 1ÞÞ T NC C pffiffiffiffi  ðr 2b nT;b ð G þ 1ÞÞ T b ;

XT;NC ¼

r2N C

hT;NC G;

CT;NC ¼

r 2N C

IT T;NC

ðr 2b hT;b ð

pffiffiffiffi G þ 1Þ;

r2b hT;b ð

hT;NC G  pffiffiffiffi ¼ G þ 1ÞÞ T b  ðr 2N hY;NC GÞ yv ;NC 1 pffiffiffiffiC  ðr 2N hY;NC G  r 2b hY;b ð G þ 1ÞÞ yv ;NC C pffiffiffiffi  ðr 2b hY;b ð G þ 1ÞÞ yv ;b :

ðA:7aÞ ðA:7bÞ

ðA:7cÞ ðA:7dÞ ðA:7eÞ

ðA:7fÞ

1114

F.J. Fernández et al. / International Journal of Heat and Mass Transfer 77 (2014) 1102–1114

[32] P. Reist, Aerosol Sci. Technol., McGraw-Hill, 1993. [33] B. Poling, J. Prausnitz, J. O’Connell, The properties of gases and liquids, McGraw-Hill Professional Publishing, 2001. [34] N. Fuchs, A. Sutugin, Highly dispersed aerosols, Ann Arbor Science Publishers, 1970. [35] A. Laaksonen, T. Vesala, M. Kulmala, P.M. Winkler, P.E. Wagner, Commentary on cloud modelling and the mass accommodation coefficient of water, Atmospheric Chemistry and Physics 5 (2005) 461–464. [36] E. Davis, Transport phenomena with single aerosol particles, Aerosol Sci. Technol. 2 (1982) 121–144. [37] A.L. Lindsay, L.A. Bromley, Thermal conductivity of gas mixtures, Ind. Eng. Chem. 42 (1950) 1508–1511. [38] NIST, Natinal institute of standards and technology, 2013, . [39] EES, Engineering equation solver, 2013, . [40] O. Preining, P. Wagner, F. Pohl, W. Szymanski, Heterogeneous Nucleation and Droplet Growth, Institute for Experimental Physics of the University of Vienna, 1981.

[41] C.R. Wilke, A viscosity equation for gas mixtures, J. Chem. Phys. 18 (1950) 517– 521. [42] T.R. Marrero, E.A. Mason, Gaseous diffusion coefficients, J. Phys. Chem. Ref. Data 1 (1972) 3–118. [43] P.D. Neufeld, A.R. Janzen, R.A. Aziz, Empirical equations to calculate 16 of the transport collision integrals for the Lennard-Jones (12–6) potential, J. Chem. Phys. 57 (1972) 1100–1102. [44] F.M. Mourits, F.H.A. Rummens, A critical evaluation of Lennard-Jones and Stockmayer potential parameters and of some correlation methods, Can. J. Chem. 55 (1977) 3007–3020. [45] N. Shibasaki-Kitakawa, M. Takahashi, C. Yokoyama, Viscosity of gaseous hfc134a (1,1,1,2-tetrafluoroethane) under high pressures, Int. J. Thermophys. 19 (1998) 1285–1295.