**3. Mathematical model and computational procedure**

#### **3.1 System considered and assumptions**

4 Advanced Fluid Dynamics

The PDIA system (VisiSizer N60V, Oxford Laser Ltd. Didcot, United Kingdom) schematically illustrated in Fig. 1 is a spatial multiple counting apparatus that captures instantaneously (i.e., in 4 ns) shadow images of the droplets moving through a thin (~400 m) sampling volume and analyzes them in real time (Minchaca et al., 2011). A dual head Nd:YAG laser sends light pulses (15 mJ at 532 nm) through a fluorescent diffuser to illuminate the region of interest from behind while a high resolution camera placed in front captures the shadow images of the objects passing in between. The disposition of these elements is illustrated in the figure. Operating in dual pulse mode the laser and camera are triggered to capture image pairs separated by a time interval of 1.7 s, the figure displays a single pair extracted from superposed frames. The analysis of single and superposed frames allows, respectively, the simultaneous determination of the size and velocity of the drops appearing. The criteria employed for the consideration of single drops and drop pairs have been described elsewhere and were validated by off-line analysis of single and superposed frames (Minchaca, 2011; Minchaca et al., 2011). Lenses with two magnifications (2×, 4×) were employed to resolve the whole spectrum of drop sizes. With each magnification 1000 frames were captured to obtain samples with over 5500 drops that ensured statistical confidence limits of 95 % (Bowen & Davies, 1951). The magnifications allowed resolving drops with sizes ranging from 5 m to 366 m and velocities of up to 185 m/s. The field of view with both magnifications was 2.561×2.561 mm2, which allowed combination of the samples obtained from both to carry out statistical analysis of the data. The calibration (i.e., m/pixel) provided for the camera, lens and magnifications used was validated measuring standard circles in a reticule and standard line spacings in a grating and the agreement was better than 0.5 %. The traversing rail shown in the figure moved the diffuser and camera to 7 prescribed x-positions (0.0013, 0.030, 0.059, 0.088, 0.116, 0.145 and 0.174 m), while the y and zs positions where maintained constant at 0 m and 0.175 m, respectively. Differently from the measurements with the patternator the measurements with the PDIA system were carried out with the nozzle oriented vertically downward, but it was experimentally verified that the distributions of dd and u obtained with

The water for the pneumatic nozzle was supplied from a reservoir using an immersion pump instrumented with a digital turbine flow-meter, a valve and a digital pressure gauge. A compressor provided the air and this line was instrumented with an automatically

Fig. 1. Schematic of experimental setup

both orientations were not significantly different.

The 3-D system domain considered in the model is shown in Fig. 3(a), it includes the twophase free-jet issuing from a pneumatic nozzle and the surrounding environment; the mixing chamber is excluded from the analysis. Since the visualization of the jets and the

Fig. 3. (a) Schematic of system considered and computational domain, (b) quadruple exposure PIV image of drops in the neighborhood of the nozzle orifice, (c) schematic of assumed air-velocity profiles at nozzle exit and (d) schematic of assumed drop velocity profiles and water flux distribution at nozzle exit

An Experimental and Computational Study of the Fluid Dynamics of Dense Cooling Air-Mists 7

The meaning of the symbols appearing in the equations is given in Section 7. The turbulence of the air was treated by the k- model for low Reynolds flows of Lam-Bremhorst modified by Yap, 1987. The constants and functions appearing in Eqs. (4) and (5) are listed in Table 1. Also, under the considerations done in Sec. 3.1 the equation of motion for individual drops

iiii

i <sup>i</sup> u dt

(7)

3 2 2

'z t <sup>k</sup> Re; 'zk Re

> <sup>d</sup> ii

uUd

0 x lx' Ux = 0, Uy= 0, Uz= Uz,max

Ux= Uz,max (tan x), Uy= 0 Uz= Uz,max(lx - x)/(lx - lx') x = (x - lx')/(lx - lx')

Uz,max= (A+W)/ [2ly(lx + lx')] \* <sup>0</sup>xlx yxyx )ll4/()ll(k2;U01.0k 5.1 oo <sup>2</sup> <sup>o</sup> ,maxz

<sup>i</sup> g)1()uU(uU <sup>d</sup>

The drag coefficient CD was assumed to vary with the particle Reynolds number, Red, according to the expressions given in Table 1. The trajectory of the drops was computed

> Constants and functions involved in the turbulence model C1= 1.44; C2= 1.92; Cdk= 1.0; = 1.3

'12 20.5 0.05 (1 exp( 0.0165Re )) 1 ; 1 ; 1 exp( Re ) Re *z t*

> Drag coefficient expressions Stokes law region, Red 2 CD = 24/Red Intermediate region, 2 ≤ Red ≤ 500 CD = Re/10 <sup>d</sup> Newton's law region, 500 ≤ Red ≤ 2105 CD = 0.44

> > d

Boundary conditions at nozzle exit

Physical properties = 1.02 kg/m3; = 1.8×10-5 Pa s; d= 998 kg/m3 Nozzle dimensions and parameters lx= 0.01 m; ly= 0.00325 m; lx'= 0.00585 m zs= 0.175 m; = 45°; = 10°

Re

*t*

<sup>k</sup> f ; C

where

Drop-phase Air-phase

lx'x lx

2

*f f <sup>f</sup> <sup>f</sup>*

t d

i d

dx (8)

2

where,

in the mist and under the effects of aerodynamic drag and gravity is expressed as,

from the variation with time of the components of the position-vectors, according to,

dd

D

C 4 3

dt du

ux,k = uz,t sin(x/lx) uy,k = uz,t sin(y/ly) 2/1 <sup>2</sup>

2 k,x

\*A is computed at local conditions: 25°C, 86 kPa.

2 k,z t,z uuuu

k,y

Table 1. Auxiliary equations, properties and dimensions

measurements of their impact footprints indicated double symmetry over the x*–*z and y–z planes the computational domain involved just one quarter of the physical domain, as seen in the figure. Additionally, since the air- and water-flow rates were stable it was assumed that on a time average basis, the flow characteristics of the two-phase jet could be simulated in steady-state conditions. For their treatment the continuous air-phase was considered in an Eulerian frame of reference and the discrete droplets were regarded in a Lagrangian frame.

The assumptions for the model were: (a) the liquid emerges from the nozzle as drops. This is supported by PIV observations done close of the nozzle orifice, as that displayed in Figure 3(b). This figure shows a quadruple-exposure photograph with trails of 4-images of droplets. Additionally and in agreement with PDIA observations the drops are assumed spherical; (b) the size distribution of the drops exiting the orifice is equal to the distribution measured at a distance z= zs (i.e., at the typical working distance of a given nozzle). This is reasonable since drop coalescence and break-up are rare events. The low volume fraction, d of the drops prevents coalescence and the PDIA images, taken at different positions in the free mists, rarely show droplet break-up; (c) drops of all the specified sizes leave the orifice at the terminal velocity reached in the mixing chamber while dragged by the air. Calculations indicate that this would be the case for individual drops and since d is low the assumption would seem reasonable for the dilute multi-drop system moving within the chamber; in the mixing chamber d< 0.08; (d) the droplets in the jet do not interact with each other and only interact with the air through interfacial drag, and (e) the air and the droplets are at room temperature and condensation and vaporization are negligible.

#### **3.2 Governing equations**

Under the considerations just described, the governing equations for the motion of the air are: the continuity equation (2), the Navier-Stokes equations (3) and the turbulence transport equations (4) and (5), which are expressed as follows,

$$\frac{\partial \mathbf{U}\_{\mathrm{i}}}{\partial \mathbf{x}\_{\mathrm{i}}} = \mathbf{0} \tag{2}$$

$$\rho \mathbf{U}\_{\mathbf{i}} \frac{\partial \mathbf{U}\_{\mathbf{j}}}{\partial \mathbf{x}\_{\mathbf{i}}} = -\frac{\partial \mathbf{p}}{\partial \mathbf{x}\_{\mathbf{j}}} + \frac{\partial}{\partial \mathbf{x}\_{\mathbf{j}}} \left[ (\mu + \mu\_{\mathbf{t}}) \left( \frac{\partial \mathbf{U}\_{\mathbf{i}}}{\partial \mathbf{x}\_{\mathbf{j}}} + \frac{\partial \mathbf{U}\_{\mathbf{j}}}{\partial \mathbf{x}\_{\mathbf{i}}} \right) \right] + \mathbf{S}\_{\mathbf{i}} \tag{3}$$

$$\rho \mathbf{U}\_{\mathbf{i}} \frac{\partial \mathbf{k}}{\partial \mathbf{x}\_{\mathbf{i}}} = \frac{\partial}{\partial \mathbf{x}\_{\mathbf{i}}} \left[ \left( \mu + \frac{\mu\_{\mathbf{t}}}{\sigma\_{\mathbf{k}}} \right) \frac{\partial \mathbf{k}}{\partial \mathbf{x}\_{\mathbf{i}}} \right] + \mu\_{\mathbf{t}} \left( \frac{\partial \mathbf{U}\_{\mathbf{i}}}{\partial \mathbf{x}\_{\mathbf{j}}} + \frac{\partial \mathbf{U}\_{\mathbf{j}}}{\partial \mathbf{x}\_{\mathbf{i}}} \right) \frac{\partial \mathbf{U}\_{\mathbf{i}}}{\partial \mathbf{x}\_{\mathbf{j}}} - \rho \varepsilon \tag{4}$$

$$\rho \mathbf{U}\_{\mathbf{i}} \frac{\partial \mathbf{z}}{\partial \mathbf{x}\_{\mathbf{i}}} = \frac{\partial}{\partial \mathbf{x}\_{\mathbf{i}}} \left[ \left( \mu + \frac{\mu\_{\mathbf{t}}}{\sigma\_{\mathbf{c}}} \right) \frac{\partial \mathbf{z}}{\partial \mathbf{x}\_{\mathbf{i}}} \right] + \mathbf{f}\_{\mathbf{l}} \mathbf{C}\_{1} \mu\_{\mathbf{t}} \frac{\mathbf{c}}{\mathbf{k}} \left( \frac{\partial \mathbf{U}\_{\mathbf{i}}}{\partial \mathbf{x}\_{\mathbf{j}}} + \frac{\partial \mathbf{U}\_{\mathbf{j}}}{\partial \mathbf{x}\_{\mathbf{i}}} \right) \frac{\partial \mathbf{U}\_{\mathbf{i}}}{\partial \mathbf{x}\_{\mathbf{j}}} - \rho \mathbf{f}\_{2} \mathbf{C}\_{2} \frac{\mathbf{c}^{2}}{\mathbf{k}} \tag{5}$$

Si in Eq. (3) is the source term expressing the i-direction momentum transferred between the air and the drops in a given cell of the fixed Eulerian grid over a Lagrangian time step. It is equal to the change in the momentum (only due to interfacial drag) of the drops following all the trajectories traversing a cell over that time step (Crowe et al., 1977) and it is given as,

$$\mathbf{S}\_{\rm i} = \frac{\pi \rho\_{\rm d}}{6 \,\mathrm{v}\_{\rm cell}} \sum\_{k=1}^{n\_{\rm cyl}} (\mathbf{N}\_{\rm d} \,\mathrm{d}\_{\rm d}^3 (\mathbf{u}\_{\rm i,out} - \mathbf{u}\_{\rm i,in}))\_{\rm k} \tag{6}$$

6 Advanced Fluid Dynamics

measurements of their impact footprints indicated double symmetry over the x*–*z and y–z planes the computational domain involved just one quarter of the physical domain, as seen in the figure. Additionally, since the air- and water-flow rates were stable it was assumed that on a time average basis, the flow characteristics of the two-phase jet could be simulated in steady-state conditions. For their treatment the continuous air-phase was considered in an Eulerian frame of reference and the discrete droplets were regarded in a Lagrangian frame. The assumptions for the model were: (a) the liquid emerges from the nozzle as drops. This is supported by PIV observations done close of the nozzle orifice, as that displayed in Figure 3(b). This figure shows a quadruple-exposure photograph with trails of 4-images of droplets. Additionally and in agreement with PDIA observations the drops are assumed spherical; (b) the size distribution of the drops exiting the orifice is equal to the distribution measured at a distance z= zs (i.e., at the typical working distance of a given nozzle). This is reasonable since drop coalescence and break-up are rare events. The low volume fraction, d of the drops prevents coalescence and the PDIA images, taken at different positions in the free mists, rarely show droplet break-up; (c) drops of all the specified sizes leave the orifice at the terminal velocity reached in the mixing chamber while dragged by the air. Calculations indicate that this would be the case for individual drops and since d is low the assumption would seem reasonable for the dilute multi-drop system moving within the chamber; in the mixing chamber d< 0.08; (d) the droplets in the jet do not interact with each other and only interact with the air through interfacial drag, and (e) the air and the droplets

are at room temperature and condensation and vaporization are negligible.

jji

 

 

t

xx p

i

n

cell d

<sup>S</sup> cell

1k

 

 

 

equations (4) and (5), which are expressed as follows,

j

x U

 

xx

ii

xx x

 

 

U

i

i

 

ii

Under the considerations just described, the governing equations for the motion of the air are: the continuity equation (2), the Navier-Stokes equations (3) and the turbulence transport

> x U i <sup>i</sup>

0

t

 

 

t

k

 

 

 

x U

x U

j i

in,iout,i

<sup>k</sup> U (4)

j i

<sup>U</sup> )(

<sup>i</sup> S

U

 

 

 

x k

Cf

t11

Si in Eq. (3) is the source term expressing the i-direction momentum transferred between the air and the drops in a given cell of the fixed Eulerian grid over a Lagrangian time step. It is equal to the change in the momentum (only due to interfacial drag) of the drops following all the trajectories traversing a cell over that time step (Crowe et al., 1977) and it is given as,

> 3 dd

<sup>i</sup> ))uu(dN( v6

k i t

(2)

i j

 

 

x U

 

 

> x U

 

i j

x U

i j

k

(6)

j i

(3)

x

 

i

 

> k Cf

22

2

(5)

j i

x U

x U

j i

 

**3.2 Governing equations** 

The meaning of the symbols appearing in the equations is given in Section 7. The turbulence of the air was treated by the k- model for low Reynolds flows of Lam-Bremhorst modified by Yap, 1987. The constants and functions appearing in Eqs. (4) and (5) are listed in Table 1. Also, under the considerations done in Sec. 3.1 the equation of motion for individual drops in the mist and under the effects of aerodynamic drag and gravity is expressed as,

$$\frac{d\mathbf{u}\_{\rm i}}{dt} = \frac{3}{4} \mathbf{C}\_{\rm D} \frac{\rho}{\rho\_{\rm d} \mathbf{d}\_{\rm d}} |\mathbf{U}\_{\rm i} - \mathbf{u}\_{\rm i}| (\mathbf{U}\_{\rm i} - \mathbf{u}\_{\rm i}) + (1 - \frac{\rho}{\rho\_{\rm d}}) \mathbf{g}\_{\rm i} \tag{7}$$

The drag coefficient CD was assumed to vary with the particle Reynolds number, Red, according to the expressions given in Table 1. The trajectory of the drops was computed from the variation with time of the components of the position-vectors, according to,

$$\frac{d\mathbf{x\_i}}{dt} = \mathbf{u\_i} \tag{8}$$


\*A is computed at local conditions: 25°C, 86 kPa.

Table 1. Auxiliary equations, properties and dimensions

An Experimental and Computational Study of the Fluid Dynamics of Dense Cooling Air-Mists 9

1

Furthermore, based on the form of the mist footprint obtained with a patternator the water flow was assumed to be distributed in the orifice according to an obelisk-shaped distribution (Camporredondo et al., 2004), as that shown schematically in Figure 3(d). With the criteria given, the initial conditions for Eqs. (7) and (8) were assigned to each port of the nozzle orifice, such that at t = 0 the position and velocity of the drop are specified as

As mentioned in Sec. 3.1, the initial velocity, (ui)k, of the drops was prescribed by assuming that the drops exit with the terminal velocity that they reach in the mixing chamber (Minchaca et al., 2010). As suggested by the observed drop trajectories (Hernández et al., 2008), according to the position assigned to the drop the angle of the velocity varied from 0 deg to deg in the x-direction and from 0 deg to deg in the y-direction. A schematic representation of the velocity vectors of the drops is displayed in Fig. 3(d), and the expressions for the ux,k, uy,k and uz,k velocity components are given in Table 1. Also, the physical properties of the fluid, the dimensions of the orifice and the angles of expansion of

Fig. 4. Numeric and volume frequency drop size distributions measured for W= 0.58 L/s and pa= 279 kPa. The log-normal and Nukiyama-Tanasawa, NT, distributions fitted to the

The Eulerian and Lagrangian equations of the model were solved using the control volume method and the particle tracking facility implemented in Phoenics. The mesh used had 128×25×99 control volumes to achieve mesh independent results and the number of ports was 100×12 in the x-y directions. The convergence criterion specified a total residual for all

*nT*

*k W Nv* 

follows,

the jet are given in the table.

respective data are included

the dependent variables ≤ 10-3.

**3.4 Solution procedure** 

, ,

(12)

i i ii k k u u and x x (13)

*dk dk*

#### **3.3 Boundary and initial conditions**

At the boundaries of the calculation domain shown in Fig. 3(a), Eqs. (2) through (5) describing the turbulent motion for the air-phase were solved imposing the following conditions:


$$\mathbf{P} = \mathbf{P}\_{\text{amb}} \; ; \; \mathbf{k} \; = \mathbf{e} = \; \mathbf{0} \tag{9}$$

To approach these conditions, the boundaries were located far away from the jet.


$$\mathbf{U}\_{\mathbf{j}} = \frac{\partial \mathbf{U}\_{\mathbf{i}}}{\partial \mathbf{x}\_{\mathbf{j}}} = \frac{\partial \mathbf{k}}{\partial \mathbf{x}\_{\mathbf{j}}} = \frac{\partial \mathbf{c}}{\partial \mathbf{x}\_{\mathbf{j}}} = \mathbf{0} \tag{10}$$

where j represents the index for the coordinate normal to the respective symmetry plane.


$$\mathbf{U}\_{i} = \begin{array}{c} \mathbf{k} \end{array} = \boldsymbol{\varepsilon} = \begin{array}{c} \mathbf{0} \end{array} \tag{11}$$


$$\mathcal{W} = \sum\_{k=1}^{n\_{\overline{l}}} N\_{d,k} \upsilon\_{d,k} \tag{12}$$

Furthermore, based on the form of the mist footprint obtained with a patternator the water flow was assumed to be distributed in the orifice according to an obelisk-shaped distribution (Camporredondo et al., 2004), as that shown schematically in Figure 3(d).

With the criteria given, the initial conditions for Eqs. (7) and (8) were assigned to each port of the nozzle orifice, such that at t = 0 the position and velocity of the drop are specified as follows,

$$\mathbf{u}\_{\mathbf{i}} = \begin{pmatrix} \mathbf{u}\_{\mathbf{i}} \end{pmatrix}\_{\mathbf{k}} \text{ and } \mathbf{x}\_{\mathbf{i}} = \begin{pmatrix} \mathbf{x}\_{\mathbf{i}} \end{pmatrix}\_{\mathbf{k}} \tag{13}$$

As mentioned in Sec. 3.1, the initial velocity, (ui)k, of the drops was prescribed by assuming that the drops exit with the terminal velocity that they reach in the mixing chamber (Minchaca et al., 2010). As suggested by the observed drop trajectories (Hernández et al., 2008), according to the position assigned to the drop the angle of the velocity varied from 0 deg to deg in the x-direction and from 0 deg to deg in the y-direction. A schematic representation of the velocity vectors of the drops is displayed in Fig. 3(d), and the expressions for the ux,k, uy,k and uz,k velocity components are given in Table 1. Also, the physical properties of the fluid, the dimensions of the orifice and the angles of expansion of the jet are given in the table.

Fig. 4. Numeric and volume frequency drop size distributions measured for W= 0.58 L/s and pa= 279 kPa. The log-normal and Nukiyama-Tanasawa, NT, distributions fitted to the respective data are included

#### **3.4 Solution procedure**

8 Advanced Fluid Dynamics

At the boundaries of the calculation domain shown in Fig. 3(a), Eqs. (2) through (5) describing the turbulent motion for the air-phase were solved imposing the following

P Pamb ; k e 0 (9)

(10)

U k i ε 0 (11)

0

of the flat hollow portion of the flanged orifice, in the rest of the orifice

xx k

jjj


To approach these conditions, the boundaries were located far away from the jet.

x <sup>U</sup> <sup>U</sup>

i <sup>j</sup> 

where j represents the index for the coordinate normal to the respective symmetry plane. - Non-penetration and non-slip conditions were specified at the external wall of the nozzle,


the profiles decreased to zero varying in angle from 0 deg to deg at the edge; the distributions were the same throughout the whole thickness (y-)direction of the orifice. These velocity profiles were suggested by the geometry of the flanged orifice and are supported by the results presented in Section 4. The expressions describing the profiles are listed in Table 1, together with the expressions for the turbulence kinetic energy and the dissipation rate of turbulence kinetic energy at this boundary, ko and o, respectively. - Positions and velocities were specified to the droplets as initial conditions for the solution of their motion (7) and trajectory (8) equations. For doing this, a first step was to decide a series of criteria to distribute throughout the orifice drops of different size and velocity in a random fashion that reflected that the water flux profile decreases from its center to its edges. To do this, the orifice was simulated as a grid of ports, k, releasing drops satisfying the diameter distribution measured at z= zs (according with the assumption indicated in Section 3.1) for the particular set of nozzle operating conditions under consideration. For deciding the number of ports assigned to each drop size category it is important to establish what type of distribution to use, number or volume frequency? Figure 4 shows both size distributions measured for a representative set of operating conditions W and pa. The number frequency distribution shows that droplets smaller than 25 m account for a large number percentage of the drops (82.85 % of them), but that they represent only 6.82 % of the volume of the drops in the sample. Since the number of ports that can be used cannot be excessively large the assignment of the ports according to the number frequency distribution would leave many sizes unrepresented. The volume frequency distribution does not present this disadvantage and was chosen to designate the number of ports for each size category. The size assigned to each port was done through a random number generator to simulate the stochastic emergence of drops with different characteristics from distinct sites of the orifice. The drops with volume vd,k= (dd3)k/6 exiting the ports k with a number frequency, Nd,k, had to satisfy the water-flow rate, W,


**3.3 Boundary and initial conditions** 

conditions:

over the length lx'

according to the following expression,

The Eulerian and Lagrangian equations of the model were solved using the control volume method and the particle tracking facility implemented in Phoenics. The mesh used had 128×25×99 control volumes to achieve mesh independent results and the number of ports was 100×12 in the x-y directions. The convergence criterion specified a total residual for all the dependent variables ≤ 10-3.

An Experimental and Computational Study of the Fluid Dynamics of Dense Cooling Air-Mists 11

located at the different x-positions. The regression curves fitted to the experimental and computational results show an excellent agreement. The dispersion exhibited by the computed results displayed in Figs. 5 and 6 is smaller than the experimental due to the very different number of drop trajectories traversing the sampling volumes, thousands in the experimental case versus a few tens in the computational case. Despite of this the model is

(a) (b)

**4.2 Effect of nozzle operating conditions on the velocity of the drops** 

**4.2.1 Effect of water flow rate at constant air inlet pressure** 

are for W= 0.58 L/s and pa= 279 kPa

their behavior in the free jet.

Fig. 6. Measured and computed velocity components as a function of x-position for drops of all sizes: (a) tangential velocity component and (b) normal velocity component. The results

In the application of air-mist nozzles for the cooling of surfaces at high temperature it is common to vary the water flow rate maintaining constant the air inlet pressure. This procedure would be equivalent to move along the curves of constant pa appearing in the operating diagram of Figure 2. The reason behind this is that the spray cooling intensity is commonly associated only to the flux of water impinging upon the hot surface, when actually there is another mist parameter that plays an important role and this is the velocity of the drops (Hernández et al., 2011). Experimentally, it has been found that the velocity of the drops increases with W up to a certain value, but once this value is exceeded the opposite effect takes place and the drop velocity decreases markedly (Minchaca et al., 2011). With the increase in W at constant pa the drops generated by the nozzle become larger (Minchaca et al., 2011) and the air flow rate, A, gets smaller as indicated by Figure 2. Both factors will alter the terminal velocities that the drops will reach at the nozzle exit and also

The multivariate effects that the droplet velocity experiences when changing W at constant pa are complex. Therefore, it was important to examine the predictions of the CFD model in

able to represent very well the trend in the behavior of the actual system.
