CC BY
34
УЧЕНЫЕ ЗАПИСКИ КАЗАНСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА Том 148, кн. 3 Физико-математические пауки 2006
UDK 532.529.6
NUMERICAL STUDY OF GROWING DROPLETS DYNAMICS IN UNSTEADY THERMAL CONVECTION FLOW
, W. Hollander, S.K. Zaripov
Abstract
A mathematical model describing the unsteady thermal convection flow and deposition of growing droplets in an expansion-type Kelvin spectrometer is developed. The model includes the Navier Stokes equations to describe the gas dynamic processes and a droplet growth model accounting for vapor depletion due to condensation. The evolution of gas temperature and velocities and growing droplet dynamics inside the spectrometer chamber are studied using CFD (computational fluid dynamics) code Fluent. The relative droplet concentrations and mass densities of condensed vapor on all droplets as a time function are investigated. In so doing it has been taking into account the droplet growth by condensation, sedimentation, heat conduction, thermal convection, droplet evaporation. It is shown that the droplet life time in the measuring chamber is limited by two processes: deposition on the walls and evaporation. The relative contribution of evaporation and sedimentation depends on the initial saturation.
Introduction
Instruments capable of measuring droplet concentrations at given supersaturations are called Kelvin spectrometers. They are important for reliably characterizing atmospheric particles. An expansion-type Kelvin spectrometer has been designed in the laboratory of Aerosol Technology of Fraunhofer Institute Toxikologie und Experimentelle Medizin [2]. The design of this Kelvin spectrometer has a number of advantages compared with other measuring devices, especially, with diffusion chambers.
The Kelvin spectrometer consists of two cylindrical chambers. The measuring chamber is connected through a valve with a buffer chamber in which the gas is evacuated before the measuring cycle. Upon opening the valve between measuring and buffer chambers the pressures in both chambers equilibrate. The rapid pressure drop in the measuring chamber is accompanied by an adiabatic temperature decrease, which may-lead to sufficient supersaturation followed by droplet growth. The chamber wall always remains at the start temperature T10 while immediately after expansion the interior is at a temperature Tn. Heat conduction will increase the interior temperature and reduce the supersaturation ratio. Any chamber wall temperature nonuniformities may-lead to a natural convection flow that also can disturb the saturation ratio homogeneity. The droplets generated inside the measuring chamber due to nucleation process grow in vapor medium. They move inside the measuring chamber due to gravity and convection flow action and can deposit on the chamber walls. A compromise must be found since small L (L is the typical chamber dimension) are required for suppressing convection, while the large L are needed for obtaining large times of droplet growth undisturbed byheat conduction from the wall. To optimize the instrument it is required understanding the physics of the thermal convection processes activated by gas expansion and further temperature dynamics in measuring chamber. Such understanding can be obtained on
R.S. Gale.e.v
the base of the mathematical simulation of the natural convection flow and growing droplets dynamics in the chamber.
We developed a mathematical model describing the unsteady thermal convection flow and the deposition of growing droplets in the expansion-type Kelvin spectrometer. The model includes the Navier Stokes equations to describe the gas dynamic processes and droplet growth model accounting for vapor depletion due to condensation. Numerical calculations of gas flow by means of CFD code Fluent (Version 6.2. Fluent Inc.. NH) show that the evolution of processes inside the spectrometer chamber can be divided into two main stages. The first stage is the adiabatic gas expansion. During this stage the pressure in the chamber establishes very quickly. The gas temperature falls practically uniformly in the entire chamber volume. The next stage is characterized by a slow temperature increase beginning from regions near the walls. The temperature gradient causes the development of a gas flow upwards near the wall which is then transformed into a global vortex flow inside the chamber. To study numerically the droplet dynamics we neglect the expansion stage in view of its small duration and calculate the temperature and convection flow development starting with the uniform temperature distribution inside the chamber and a constant wall temperature larger than the interior gas temperature. Using the gas flow velocity and temperature distribution in the spectrometer chamber from CFD calculations we calculate the droplets paths taking into account their growth resulting from the local super-saturation.
1. Model of motion of growing droplets
A Lagrangian approach is used to calculate the droplet trajectories in a supersaturated vapor. The equations of motion of a single droplet neglecting all forces except the aerodynamic drag and gravity are written as
dvp u — vp — vs dr p
dt t ’ dt
= v
p-
(1)
where vp = vp(rp, t) is the droplet velocity, rp is the radius vector of droplet position. U = u(rp,t) is the gas velocity, vs = tq is the settling velocity, t = ppdp/\8^, pp is the droplet density, dp is the droplet diameter, ^ is the gas viscosity, g is the gravity acceleration vector. The relaxation time t relation holds for the continuum regime without Reynolds correction.
The droplet growth is described by the equation [5]
ddp 8A
<2'
where L is the latent heat of condensation, A is the carrier gas heat conductivity, Tx is the medium temperature.
For the droplet temperature T0 we use
LDVM fPs(Tq) ^
o D \ I /T1 ''T1
Rg A \ T o la
where Dv is the vapor diflusivity, M is the molecular weight, Rg is the universal gas constant, ps(T) is the saturation vapor pressure. The partial vapor pressure pTO with account made for vapor pressure depletion due to condensation is expressed through the formula
pTO = (S — !3)ps(T^),
where 3 = Nncdppp/6ps, N is the initial concentration of growing droplets, ps = = ps(T^)M/RgTTO. The saturation ratio S is defined as the ratio of the actual
ambient vapor pressure pTO0 after expansion to the saturation vapor pressure ps :
S = p^o/ps(T).
At the actual droplet location we solve the equation (3) to find T0 and include it into (2). The distribution Tx(rp, t) and the gas velocity u(rp, t) are found from the gas flow calculation inside the measuring chamber. The solution of equations (1) and (2) allows us to find the actual position and diameter of the droplet.
2. Deposition of droplets without convection flow influence
Droplets move in the measuring chamber due to gravity and thermal convection flow and can reach the chamber walls and deposit on it. The ratio of the number of droplets deposited on walls to their initial number characterizes losses during one measuring cycle. In the absence of convection flow and uniform initial distribution of temperature the part of deposited droplets can be easily theoretically estimated. For the horizontal orientation of the cylinder we will consider deposition of droplets in one cross section as a two-dimensional problem neglecting influence of the cylinder ends due to large ratio of cylinder height to its diameter. Let n0 be the droplets concentration that can be defined as the ratio of the total number of droplets N0 to the area of cross section of cylinder n0 = N0/nR2 (R is the cylinder radius). The droplets are uniformly distributed over a cross section of the cylinder and they fall down by gravity action with the settling velocity vs = rg. The relative number nd of droplets that deposits during the time period from t = 0 to t = td (deposition rate) can be expressed by the formula [3]
nd = ^ ^0 + cos , (4)
where 0 = arcsin (h/2R), h is the distance traveled by the droplet. Formula (4) can be easily obtained by analysing the Fig. 1 (shaded circle is the area of deposited droplets).
Fig. 1. Geometry of the area of deposited droplets
h
settling velocity and time h = vstd. If droplets grow the distance h is obtained byintegration of the variable settling velocity
td
h = j vs(t)dt.
o
3. Adiabatic relations of ideal gas
Tb calculate the initial conditions for measuring chamber after expansion the adiabatic relations of ideal gas can be used. Under adiabatic conditions in an ideal gas the following relations (x « 1.4 is the adiabatic coefficient) can be written (V is the gas volume)
p (1-k)/k t = const,
PVK = const, (5)
T 1/(-1+k) V = const.
Using relations (5) it is possible to calculate the gas parameters in the measuring chamber after expansion. Let T10 and P10 be ambient temperature and ambient pressure, the measuring chamber of volume V10 is filled with air. The buffer chamber has the volume V20 and also starts at ambient pressureP10. Next, it will be quickly (and, hence adi-abatically) evacuated to the measured pressure P200. New temperature T200 of buffer gas can be found from the relation
rp _____ p(1 — k)/k P — (1 — K)/Krp jp\
T 200 = P10 P200 T10. (0)
Obviously, T200 < T10 but heat conduction makes T200 reach T10 after some time
P200 P20
V20
P_ _ P200 _ 'P20_ T To 00 To o
P20
const,
P20 = Pl-1/K P200 . (8)
Next, the valve between the measuring and the buffer chambers is opened and pressure is equlibrated to P11 in both chambers. Immediately after pressure equilibration, the air mass originally contained in V10 of the measuring chamber at T10 is adiabat-ically cooled to T11 % expansion to Vu, while the air originally contained in V20 of the buffer chamber at T10 is adiabatically heated to T21 by compression to V21 where
V10 + V20 = V11 + V21.
Using the adiabatic relations we get
T/ V20(Vi0 + V20)
V21 ~ ~—717^-------------------------------------------------> (9)
V10 + V20
(10)
^ P20 J
Inserting V21 into the above adiabatic equation we finally get
/
V10 + V20
Uft)1'" vu + v™.
P11 = P10 ViKO (V10 + V20 - V2i)-K . (11)
Equations (10), (11) can be used to calculate the temperature and pressure in measuring chamber after expansion. To find the saturation immediately after expansion we can use the formula
S00 P11 psat(T10)
*0 - ------5-------(rp \ > (1^)
PlO psat(Tll)
S00
l
4. Results of calculations without convection (analytical solution for temperature distribution)
As was noted below for horizontally orientated chamber we can neglect the influence of the cylinder ends. Therefore the mathematical problem is formulated as two-dimensional (Fig. 2). The gas temperature is taken equal to T10 at t =0. The wall temperature at r = R is Tw (R is the cylinder radius, r = \Jx2 + y2). To simulate the temperature field development by conduction action in the absence of convection flow the analytical solution of heat conduction problem for cylindrical cavity can be used [1]. The corresponding solution for cross section of the chamber can be written as
T(r, i)=Tw-2jr (Tw - TW) Jo(^r) exp (-alt), (13)
^ anJi(anr)
n=1 v v
where r = r/R is the relative chamber radius, J0, J1 are the Bessel functions of the first kind of order zero and one respectively, the an are the positive zeros of J0, t = t/(R2X/pcp).
Fig. 2. Schematics of conduction problem in a cylindrical enclosure
Table 1. Physical data for system of H2 O-vapor
Ps lg ps(T) = A - B/T A = 21.18; B = 5367
Dv 2.58 • 1CT5 nr • s^1
M 0.018 kg • mol-1
Rg 8.3144 J • moP1 • KT1
L 2.3-10b J-kg
k 0.0261 W • m_i • KT1
Pp 1000 kg • in-"1
A' 1.34 • 10~b kg • m_i • s-1
^9 0.65 • 10“' m
Using solution (13) and physical data for considered system of H2 O-vapor (Table 1) we studied the temperature and saturation dynamics in the cylinder. The saturation in the chamber before expansion was taken equal to S00 = 0.99 that gives us the initial saturation after expansion S0 = 1.058. The temperature and saturation as a time function are given in Fig. 3 (T10 = 291.89 K, Tw = 293 K, R = 6 cm). The temperature increases with time until it reaches the value of wall temperature. The saturation decreases and tends to value S00 before expansion. It indicates the evaporation possibility in the chamber when the saturation will be smaller than unity. The depent
Fig. 4. The saturation decreases in time starting from the walls due to gas temperature establishment as a result of heat conduction.
The droplet growth in the center of the cylinder was investigated for temperature field given by the solution (13). The comparison of droplet diameters with experimental
Fig. 3. Temperature (a) and saturation (b) as a time function for cylinder center point
r r
ft) b)
Fig. 4. Temperature (a) and saturation (6) as a chamber radius function
30-
25-
20-
E 154
10-
5
0
N=59; S0=1.043
, ------- T=const
.... T=T(t) due to heat conduction
• experiment
10
20
t,s
30
40
a)
14-i
12-
10-
8-
6-
4-
2-
N=856; SQ=1.058 ------ T=const
■ ■ ■ ■ T=T(t) due to heat conduction
• experiment
10
t,s
15
20
Fig. 5. Comparison of droplet diameter with experiments for various initial saturation and concentrations cm-3
t,s
a)
1.06
1.05'
1.04'
1.03'
1.02'
1.01
1.00
0.99
------N=1
----- N=300
------N=1000
- - - -N=3000
80 100 120 140 16
t,s
Fig. 6. Droplet diameter (a) and saturation (b) for various initial droplet concentrations N. cm -3)
data is given in Fig. 5. Two figures correspond to various values of initial saturation and concentration of particles. The curves T const correspond to the case with constant gas temperature T = T10^d T = T(t) does so when heat conduction is taken into account using the solution (13). Experimental data were taken from the laboratory Aerosol Technology of Fraunhofer Institute Toxikologie und Experimentelle Medizin. We see satisfactory agreement between the theoretical and experimental data.
The droplet diameter and saturation as a time function for various initial droplet concentration N are shown in Fig. 6. The curves correspond to droplets starting at the point x = 0, y = 0 and are calculated without gravity. Increasing influence of pressure depletion is observed with the rise of particle concentration. In the absence of water pressure depletion (N =1) the droplet grows without limit. On the other hand, for N
of saturation after droplet evaporation coincide. (Fig. 6, b)
5. Results of calculations with convection (CFD solution)
In order to identify the influence of natural convection on temperature, supersaturation and droplet growth we have studied the flow using Fluent. The mathematical model includes the Navier Stokes equations [4]. The gas is initially steady and adheres to all solid walls (non-slip conditions). The gas flow in the measuring cylinder can be assumed to be laminar because of the small Rayleigh number. We use the same boundary conditions for the temperature as in the problem without convection.
The initial values of gas parameters after adiabatic expansion are calculated from the ideal gas relations (10), (11). The gas temperature T10 = 300 K after adiabatic expansion decreases to the value T11 = 295.7 K. We have two values of initial saturation S0 = 1.0^d S0 = 1.2 that correspond to the values of saturation S00 before expansion 0.86 and 0.975.
The average temperature within the cylinder cross-section calculated with and without convection as a time function is shown in Fig. 7. In the presence of convection flow the gas temperature faster approaches the wall temperature. This means that the saturation decreases faster (Fig. 8).
The time evolution of the temperature distribution and of contours of velocity magnitude is shown in Fig. 9 10. The temperature gradient causes the development of a gas flow upwards near the wall which in the course of time transforms into a global vortex flow inside the chamber (Fig. 9). The values of the velocity magnitude of the
1.20
298.5
1 - S0=1.06
2 - So=12
298.0-
297.5-
% 297.0-i-
296.5-
296.0-
0 2 4 6 8 10
t,s
Fig. 7. Average gas temperature with and without convection
295.5
0.90-1---------.------1------.------1------.------1------.------1 T—■
0 2 4 6 8 10
t,s
Fig. 8. The saturation based on the average temperature gas temperature with and without convection influence
vortex flow can achieve ~ 0.05 m/s. The convection flow violates the axially symmetric distribution of temperature and accelerates the increase of temperature in upper half of the chamber (Fig. 10).
The distributions of gas temperature Tx(rp,t) and velocity u(rp,t) found from CFD calculation and solution of equations of (1). (2) allow calculating the local droplet coordinates and diameters. The gas velocity and temperature distributions obtained from the CFD calculations are stored in files at several times. To provide accuracy of calculations the interpolation between the distributions TOT(fp^d u(rp) is applied for intermediate times.
To calculate the particle deposition we assume that the droplets are uniformly distributed initially over a cross section of the cylinder. Let the initial droplet diameter be dp0 = 1 ^m. We observe the path of each droplet. A few snapshots of the motion of the droplets are shown in Fig. 11. Legends show the droplet diameters. One can see that droplet motion in the chamber is influenced by both gravity and convection flow action. The time of droplet life in the measuring chamber is limited by two processes: deposition on the walls and evaporation. The thermal convection influences the growing droplet dynamics by two mechanisms. Firstly, the convection flow along with heat conduction increases the gas temperature and. hence, reduces saturation and may-lead to droplet evaporation. Secondly, the convection flow together with gravity affect the droplet path. In the vortex flow generated by thermal convection in the measuring chamber the droplet path becomes complex. The droplet can reach the regions with higher temperature and evaporate or can deposit on the walls. For the smaller initial saturation the evaporation dominates but the majority of the droplets reaches the cylinder wall for larger S0 ■ The example of two trajectories starting from one point but with various initial saturations is given in Fig. 12. The droplet with So = 1.06 will evaporate near the cylinder wall and droplet with S0 = 1.2 will deposit.
The relative droplet concentrations (n = N(t)/N(0)) obtained taking into account all processes (droplet growth by condensation, sedimentation, heat conduction, thermal convection, droplet evaporation) and without gravity influence (g = 0) are shown in Fig. 13. The corresponding liquid water content as a function of time is given in Fig. 14. Curves obtained without gravity allow estimating the concentration dynamics due to
6. Droplet dynamics
velocity-magnitude
n
0.05
0.045
0.04
0.035
0.03
0.025
0.02
0.015
0.01
5 s 10 s
Fig. 9. Time evolution of the contours of velocity magnitude in the cross section of the horizontal cylinder
5 s 10 s
Fig. 10. Time evolution of the temperature distribution in the cross section of the horizontal cylinder
Fig. 11. Droplet positions at various times at So = 1.06 and So = 1.2
Fig. 12. The droplet paths for two initial saturations
t,s
Fig. 13. The relative droplet concentration for N = 1000 cm-3 and various initial saturations
Fig. 14. The mass density of condensing vapor on all droplets for N = 1000 cm-3 and various initial saturations
the evaporation only. For the smaller initial saturation the losses of droplets generated by nucleation is determined mainly by evaporation. But for larger initial saturation the sedimentation becomes the main mechanism of droplet losses in the measuring chamber. The decrease of the slop of the curve n(t) with large saturation at t ~ 5 s (Fig. 13) is connected with motion of droplets in the vortex. The vortex flow modifies the relative concentration distribution by gravitational deposition of droplets. Therefore, thermal convection reduces the deposition rate.
Note that the results of calculations in Fig. 13 14 closely reproduce the experimental decay times under gravity and in parabolic flights where the microgravity conditions are realized.
Conclusions
The mathematical model developed here represents the gas dynamic processes and growing droplet dynamics in an expansion type Kelvin spectrometer. The difference between the temperature of the measuring chamber walls and temperature of interior gas after gas expansion leads to thermal convection flow inside the chamber. This flow
affects the temperature establishment and hence saturation ratio. In the presence of convection flow the droplet paths become very complex and natural thermal convection in the spectrometer chamber may increase as well as decrease droplet losses. The carrier gas temperature along the droplet trajectories changes and droplets can enter the regions where the saturation ratio falls below unity. It means that in these regions droplets can evaporate and the evaporation is one of the processes that decrease the cloud droplet concentration in the measuring chamber. The relative contribution of evaporation and sedimentation depends on the initial saturation.
The work was supported by RFBR (project 05-01-00794) and DLR (project 50WM0437).
Резюме
P. С. Галеев , В. Холла,пдер, Ш.Х. Зарипов. Численное исследование динамики рас-
тущих капель в условиях нестационарной тепловой конвекции.
Развита математическая модель движения растущих аэрозольных частиц в спектрометре Кельвина под воздействием силы тяжести в условиях нестационарной тепловой конвекции. Модель включает в себя систему уравнений Навье Стокса для описания газовой смеси при адиабатическом расширении и выравнивании температуры в камере за счет теплопроводности и конвекции и уравнения движения капель с учетом изменения их размера за счет конденсации или испарения. Уравнения движения капель интегрируются в поле температур и скоростей, пайдеппом из решения нестационарной газодинамической задачи с помощью CFD программы Fluent.. Исследованы нестационарные газодинамические и температурные поля, распределения насыщенности паровоздушной смеси и временная зависимость общей массы капель. Построены временные зависимости доли осажденных капель па стейках с учетом и без учета тепловой конвекции. Показано, что время жизни капли в камере ограничивается оседанием под действием силы тяжести и испарением. Относительный вклад двух процессов зависит от начального насыщения.
Literature
1. Carslaw H.S., Jaeger J.C. Conduction of Heat in Solids. - Clarendon, 1980. - 510 p.
2. Holiader W., Dunkhorst W., Lodding H., Windt H. Theoretical simulation and experimental characteristic of an expansion-type Kelvin Spectrometer with intrinsic calibration // J. of Atmospheric and Oceanic Technology. - 2002. - V. 19. - P. 1811-1825.
3. Kojic M., Tsuda A. A simple model for gravitational deposition of non-diffusing particles in oscillatory laminar pipe flow and its application to small airways // J. of Aerosol Science. - 2004. - V. 35. - P. 245-261.
4. Patankar S. V. Numerical Heat Transfer and Fluid Flow. - N. Y.: McGraw-Hill, 1980.
5. Seinfeld J.H., Pandis S.N. Atmospheric Chemistry and Physics. - Wiley, 1998. - 1326 p.
Поступила в редакцию 04.09.06
Холландер Вернер профессор лаборатории аэрозольных технологий Института токсикологии и экспериментальной медицины. Ганновер. Германия.
E-mail: hollaendereitem.fraunhofer.de
Зарипов Шамиль Хузеевич доктор физико-математических паук, профессор кафедры моделирования экосистем Казанского государственного университета.
E-mail: Shainil.ZaripovQksu.ru
