Flow, Turbulence and Combustion 73: 277–305, 2004. © 2004 Kluwer Academic Publishers. Printed in the Netherlands.

277

Direct Numerical Simulation of a Turbulent Axisymmetric Jet with Buoyancy Induced Acceleration AMIT AGRAWAL1, B.J. BOERSMA2 and AJAY K. PRASAD3 1 Discipline of Mechanical Engineering, University of Newcastle, NSW 2308, Australia;

E-mail: [email protected] 2 Laboratory of Aero and Hydrodynamics, Delft University of Technology, 2628 CA Delft, The Netherlands; E-mail: [email protected], 3 Department of Mechanical Engineering, University of Delaware, Newark, DE 19716, U.S.A.; E-mail: [email protected] Received 6 November 2002; accepted in revised form 20 January 2004 Abstract. Direct numerical simulations of an axisymmetric jet with off-source volumetric heat addition are presented in this paper. The system solved here involves a three-way coupling between velocity, concentration and temperature. The computations are performed on a spherical coordinate system, and application of a traction free boundary condition at the lateral edges allows physical entrainment into the computational domain. The Reynolds and Richardson numbers based on local scales employed in the simulations are 1000 and 12 respectively. A strong effect of heat addition on the jet is apparent. Heating causes acceleration of the jet, and an increased dilution due to an increase in entrainment. Further, the streamwise velocity profile is distorted, and the cross-stream velocity is inward for all radial locations for the heated jet. Interestingly, the maximum temperature is realized off-axis and a short distance upstream of the exit of the heat injection zone (HIZ). The temperature width is intermediate between the scalar and velocity widths in the HIZ. Normalized rms of the concentration and temperature increases in the HIZ, whereas that of streamwise, crossstream and tangential velocities increases rapidly after decreasing. Both mass flux and entrainment are larger for the heated jet as compared to their unheated counterparts. The buoyancy flux increases monotonically in the HIZ, and subsequently remains constant. Key words: buoyancy flux, cumulus clouds, direct numerical simulation (DNS), entrainment, jets and plumes, mass flux, off-source volumetric heating.

1. Introduction An interesting and well known phenomenon is observed in cumulus clouds: as the heated column of air (plume) rises through the atmosphere and passes through the cloud base (the vertical location at which water vapor begins to condense to form droplets, thus rendering the cloud visible), its entrainment behavior changes drastically. Field observations of Colorado cumuli [19] clearly show that entrained air in a cumulus cloud originates either from the cloud base or from the cloud top, i.e., lateral entrainment of ambient air into the cloud virtually shuts down,

278

A. AGRAWAL ET AL.

in contrast to the normal (as expected for an ordinary plume) lateral entrainment at elevations below the cloud base. If one applies the standard plume value of the entrainment coefficient of 0.083 to model a cloud, one can either predict the vertical rise of the cloud, or its water vapor content correctly, but not both [24]. See [22] for more discussion on entrainment in geophysical flows, and [12] for an excellent review of entrainment in cumulus clouds. Bhat and Narasimha [6] suggested that this reduced entrainment is linked to the release of latent heat of condensation (and the concomitant addition of buoyancy) that commences at the cloud base and continues as the cloud ascends beyond it. They constructed a laboratory analogue consisting of an upward-pointed turbulent axisymmetric jet
to simulate the ascending cloud, and a device to inject volumetric heat in an off-source manner to simulate the latent heat release effect in cumulus clouds. The jet fluid (water) was selectively heated in the experiment in the following manner. A small amount of acid was added to the jet fluid to make it electrically conducting, while the ambient fluid was deionized (non-conducting) water. At some downstream distance, the jet was passed through a series of electrically energized wire grids, such that a current flowed only through the body of the jet, resulting in selective ohmic heating of the jet fluid. Agrawal et al. [2, 3] applied wholefield measurement techniques for velocity and temperature to the volumetrically heated jet problem within the heat injection zone (HIZ), and also extended the range of the non-dimensional heating parameter. The results in [2, 3] such as the reduction in scalar width and the disruption of coherent structures were in agreement with those in [6]. However, there was disagreement on the issues of rms and mass flux. See [3] for a detailed comparison between previous and our experimental results, and [10] for detailed discussion on normal jets and plumes. Our goal here is to undertake direct numerical simulations (DNS) of a turbulent jet with off-source buoyancy addition. Specifically, we wish to examine the behavior of the jet under a much stronger heating rate. The study would shed light on the effect of buoyancy addition on entrainment and other characteristics of free-shear flows. Experimental measurements on volumetrically heated jets are particularly difficult because the heating grids can obstruct the view and cast glare and shadows, and because experimental run times are restricted by the eventual recirculation from the top of the jet fluid into the ambient fluid which destroys their differential conductivities. Furthermore, it is difficult to experimentally attain large accelerations due to the excessively large input power required. In contrast, DNS studies permit the calculation of all the flow characteristics and allow an easier control of the parameters. However compared to experiments, fewer DNS studies of free shear flows have been performed in the past chiefly because it is difficult to define the computational domain due to an absence of clear flow boundaries. Second, these flows develop along the spatial coordinate, and the expanding physical boundaries of the flow must be suitably addressed. Boersma et
Clouds are perhaps better modeled as plumes rather than jets. However, jets are easier to create in the laboratory and are therefore somewhat better characterized than plumes in the literature.

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

279

al. [9] and Lubbers et al. [16] circumvented these difficulties by choosing a large computational domain and applying appropriate boundary conditions. Their code therefore correctly replicated the mean and turbulent characteristics of velocity and scalar concentration for an axisymmetric turbulent jet. The code of [9, 16] has been appropriately modified for the purpose of the present study. A literature survey immediately reveals that very few DNS studies have been undertaken on free-shear flows with off-source buoyancy addition. For example, Basu and Narasimha [5] performed DNS computations at Re = 1600 to investigate the effect of heating on jets, and the effect of different heating profiles on the jet. Their formulation consisted of a temporal analogue of the spatial problem. Therefore, they actually solved a temporally evolving cylindrical mixing layer; consequently, their characteristic local Reynolds number decreased with time. Further, their heating profile was not coupled with the instantaneous concentration field. Basu and Narasimha found a reduction in the normalized fluctuations with heat addition. They found that a narrow heating profile affects the jet much more strongly, and found an increase in vorticity with heat addition for all cases investigated. They also found large expulsive motions at certain transverse cross-sections, and attributed this as the reason for reduced entrainment in a heated jet. Their computations employed periodic boundary conditions along the axial coordinate, and therefore could not conclusively provide evidence for a change in entrainment rate with heat injection. As will be evident from the governing equations, the system solved here involves a three-way coupling between velocity, concentration and temperature. Such a system is both numerically challenging and interesting, and a numerical solution of this type of problem has not been previously reported to the best of our knowledge. We have not attempted to replicate moist, compressible convection that prevails in the atmosphere. In addition, the range of scales that exists here is much smaller than in clouds due to the much smaller value of Reynolds number employed in our simulations. Also, computational constraints described subsequently require us to employ a Richardson number which is about one order of magnitude higher than in clouds. Despite these shortcomings, our numerical simulation shows some features that are qualitatively similar to cumulus clouds. These results can further serve as a benchmark for modeling and computational efforts of buoyancy-added free-shear flows. 2. Computational Details We will now provide details of the simulations by examining the dimensional and non-dimensional versions of the governing equations, the resulting non-dimensional parameters and their prescribed values, the geometry of the computational domain, the computational scheme, and the initial and boundary conditions.

280

A. AGRAWAL ET AL.

2.1. GOVERNING EQUATIONS The governing equations for the problem are the equations for mass and energy conservation, and momentum and scalar transport equations. The fluid properties are assumed to be constant, and the Boussinesq approximation is used to restrict the effect of variation of density with temperature exclusively to the body force term. These equations are: ∇ · u = 0,

(1)

1 ∂u + (u · ∇)u = − ∇p + ν∇ 2 u + gβ(T − Ts ), ∂t ρ

(2)

∂C + (u · ∇)C = αD ∇ 2 C, ∂t

(3)

∂T + (u · ∇)T = α∇ 2 T + ∂t




q  (C) ρCp

0

in HIZ, otherwise.

(4)

In the above equations, u is the velocity vector, ρ is the fluid density at Ts , p is the pressure, ν is the kinematic viscosity, T is the temperature, Ts is the reference temperature which was chosen as the initial temperature, g is the gravitational vector, β is the coefficient of thermal expansion, C is the scalar concentration, α is the thermal diffusivity, αD is the scalar diffusivity, q  is the volumetric power, and Cp is the specific heat. Note the presence of a source term in both the momentum and energy equations. The magnitude of the local heat injection (or the strength of the energy source term) is a function of the local scalar concentration. This replicates the experimental situation in which the presence of acid lowers the local electrical resistance and increases the ohmic heating rate. Further, heating is confined to the jet fluid which corresponds to a non-zero value of scalar concentration, i.e., selective heating of the jet fluid occurs. The presence of the gβT term in the momentum equation which suggests a plume-like behavior has to be interpreted with care. A plume experiences buoyancy addition only at its source, and not off-source as in the present study. Consequently, T is zero in the pre-HIZ region, and a jet-like behavior applies in this region. Moreover, the temperature can increase with the axial coordinate in the HIZ, unlike for a plume which experiences an axial decay in temperature at the centerline. It is however not immediately obvious whether the behavior will be closer to a jet or a plume at the exit of the HIZ, although it is known that a jet with buoyancy addition will eventually evolve into a plume [20]. Advection of the scalar depends on the velocity components. As seen above, the velocities are dependent on the local temperature in the HIZ (and post-HIZ); the temperature rise in turn depends on the scalar concentration in the HIZ. Therefore, there is a three-way coupling between temperature, velocities, and scalar concentration in the HIZ. Such a three-way coupling of the momentum, scalar and energy

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

281

equations has not been attempted by previous DNS investigations of heated jets. For example, Basu and Narasimha [5] did not solve for the scalar field and instead used a time-averaged heat injection profile in their computations. However, in the post-HIZ, only the velocity and temperature are coupled. This difference in the amount of coupling in the HIZ and post-HIZ is apparent in the different behaviors for the two regions. 2.1.1. Non-Dimensional Governing Equations Orifice scales are employed for non-dimensionalization. Therefore, all distances are normalized by the orifice diameter, d, velocities by the velocity at the orifice exit, U0 , and concentration by the concentration at the orifice exit, C0 . Time is non-dimensionalized as t = t

U0 . d

Pressure is non-dimensionalized as p . p = ρU02 Using the non-dimensional scheme and dropping the primes for convenience, the governing equations may be written as follows. ∇ · u = 0,

(5)

1 2 gβθ0 d ∂u + (u · ∇)u = −∇p + ∇ u+ θ, ∂t Re U02

(6)

1 ∂C + (u · ∇)C = ∇ 2 C, ∂t Sc Re

(7)

1 ∂θ + (u · ∇)θ = ∇ 2θ + ∂t Pr Re


q  (C)

d ρCp U0 θ0

0

in HIZ, otherwise,

(8)

where θ = (T − Ts )/θ0 , Re is the Reynolds number (= U0 d/ν), Pr is the Prandtl number (= ν/α), and Sc is the Schmidt number (= ν/αD ). θ0 is a temperature scale that will be described shortly. We have set Re = 1000 to approximately match the value used in our experiments. On the other hand, we set Pr = 1, and Sc = 1 which do not match our experiments which have Pr = 6 and Sc ≈ 40. Our choice of Pr and Sc was essentially determined by computational constraints – the computations become very expensive for larger values of Pr and Sc.

282

A. AGRAWAL ET AL.

2.2. NON - DIMENSIONAL PARAMETERS Defining the Grashof number as gβθ0 d 3 , ν2 the momentum source term may be rewritten as Gr =

Gr gβθ0 d θ = 2 θ. 2 Re U0

(9)

Similarly, within the HIZ the non-dimensional energy source term is q  (C) d , ρCp U0 θ0 where q  (C) is the volumetric heat generation within the HIZ. We will assume that the local rate of heat generation is directly proportional to the local acid concentration (which is a reasonably close approximation to the experimental situation). In non-dimensional form, we can write q  (C) d = kC, ρCp U0 θ0

(10)

where k is a non-dimensional constant with an appropriate value. The value of k essentially sets the rate at which heat is injected into the body of the jet. As will be seen shortly, increasing k is equivalent to increasing the value of Richardson number. The total heat supplied to the HIZ can be obtained by integrating the volumetric heat generation over the volume of the HIZ (VHIZ ).  q  (C) dV . Q= VHIZ

We can now integrate both sides of Equation (10) over the volume of the HIZ to obtain  Q d =k C dV = kI, (11) ρCp U0 θ0 VHIZ

where



I=

C dV .

(12)

VHIZ

It is now possible to define a Richardson number as Ri =

gβ Q . ρCp dU03

(13)

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

283

Using Equations (9), (11) and (13) it is also easily shown that kI Ri Re2 = 3. Gr d In our simulations, we fixed the value of Gr/Re2 at 0.0468. The value of Ri cannot be set a priori because I can be determined only in post-processing. Note that I is a complicated function of k as well. We set k = 1 in our simulations. For our chosen geometry of the HIZ and the distribution of C therein, we obtain I ≈ 8d 3 from post-processing, which implies that Ri is about 0.4. It is now straightforward to compute the Richardson number based on local scales (Ri∗ ) as follows: Ri∗ = Ri

U03 d , U∗3 b∗

where U∗ and b∗ are the local velocity and length scales at the entry to the HIZ. We estimate Ri∗ to be about 12. Note that this value of Ri∗ is higher than that employed in the experiments which had Ri∗ = 0.3 [2, 3]. Such a high value was necessary because the jet takes much longer to achieve stationarity for lower values of Ri. It will be explained shortly that the initial thermal takes much longer to exit the computational domain, increasing the computational cost for small values of Ri. Despite the mismatch in Ri∗ , there is good general agreement between the computational and experimental results. Note that the temperature scale θ0 can be related to the temperature rise experienced by the jet fluid as it transits the HIZ, i.e., θ0 =

Q , µCp

(14)

where µ is a suitable average value of the changing mass flux of the jet within the HIZ. Let µ = f µ0 , where µ0 is the mass flux at the nozzle exit (µ0 = π d 2 U0 ρ/4) and f is a factor that depends on the simulation conditions. From Equations (11) and (14) one can show that 4 kI . (15) π d3 We will show in the results section that the numerical value of f for our chosen conditions is consistent with our expectations. f =

2.3. COMPUTATIONAL GEOMETRY AND NUMERICAL SCHEME We use DNS to provide a solution of Equations (6)–(8). These equations have been recast in spherical coordinates (for example, see [7]) for the purpose of these computations with, r, θ, φ denoting the radial, azimuthal, and tangential directions.
A
These directions correspond respectively to the streamwise, cross-stream, and out of plane components in the experimental setup of [2, 3].

284

A. AGRAWAL ET AL.

Figure 1. The computational domain (adapted from [9, 16]).

spherical geometry is superior to a cylindrical domain because it allows a more efficient utilization of the computational grid. The computational domain is a conical volume segment of a spherical shell (Figure 1). The computational domain spans between 50 and 93 nozzle diameters in the streamwise direction as measured from the origin of the spherical coordinate system. The lateral edge of the domain is angled at π/40 to the centerline. The chosen angle is close to the spreading rate of a jet, implying that the distance between the lateral edges of the jet and the lateral boundaries of the computational domain remains roughly the same for all axial locations. The governing equations are discretized in this spherical coordinate system on a three-dimensional staggered grid with help of a second order finite volume method. The singularity at the centerline (pole) of the system is removed by the finite volume method because all terms in the equations are multiplied by the Jacobian r 2 sin θ (for example, see [17]). The grid is non-uniform in the streamwise direction which allows for accurate calculation near the orifice without using too many grid points in the far-field. Therefore, this coordinate system allows for a well balanced numerical resolution both near the inflow as well as in the far field of the jet. The scalar field is discretized with a Total Variation Diminishing (TVD) scheme. Such a scheme has to be used because standard finite difference methods can cause an overshoot of the scalar concentration, i.e. values higher than the source concentration or even negative values, which is unphysical. The numerical scheme for the scalar is based on a limiter which interpolates the values of the concentration at the cell faces. Here we use the following limiter, see [15] 1 ci+1 − ci , ci+1/2 = ci + φ(ri+1/2 )(ci − ci−1 ), ri+1/2 = 2 ci − ci−1     1 2 , r, 2 , φ(r) = max 0, min 2r, min 3 3 where ci+1/2 is the concentration at the cell face. A standard second order scheme would give ci+1/2 = (ci +ci+1 )/2. The scheme is known as the κ −1/3 scheme and

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

285

has a global accuracy O(h2 ), where h is the gridspacing. Additional details may be found in [15, 25]. We use the second order Adams–Bashforth method for time integration. This method gives a prediction for the velocity at the new time level, here denoted by u∗   3 1 ∗ n n n−1 F (u ) − F (u ) , with F = −A − P + D, u = u + t 2 2 where u is the velocity vector and F denotes the sum of advection (A), diffusion (D) and pressure (P ) operators and t is the time step of the numerical scheme. The velocity u∗ is in general not divergence free and therefore has to be corrected to obtain un+1 . This is done with help of the pressure correction method in which an additional step is taken to obtain un+1 un+1 = u∗ − t∇pc , where pc is the pressure correction. Taking the divergence of the equation above and assuming that ∇ ·(un+1 ) = 0 gives a Poisson equation ∇ 2 pc = ∇ ·u. This Poisson equation is solved with a fast Poisson solver [9]. Once the pressure correction is known, the calculation of un+1 is straightforward. The computational grid consists of 270 × 80 × 48 points in the r, θ, φ directions respectively. As mentioned earlier, this discretization scheme is able to follow the streamwise spreading of the jet and to allow a well-balanced resolution of the flow field with a reasonable number of grid points. In fact, Boersma et al. [9] have used a similar grid resolution as in the present computations, and showed that the resolution is adequate for resolving all the scales in the flow. Using a similar analysis, we will show here that we are resolving length scales which are comparable to the Kolmogorov length scale, η [18]. The Kolmogorov length scale, η = (ν 3 / )1/4 where the dissipation energy, can be estimated from = 0.5Uc3 /(z−z0), where Uc is the local centerline velocity [14]. From the above equation, η/d is estimated to vary as 0.13, 0.09, 0.09 at z/d = 20, 30 and 40 in our computation. The grid resolution (R/d, Rθ/d, R sin θφ/d) near the centerline at these axial locations are, respectively, (0.17, 0.07, 0.015), (0.21, 0.08, 0.017), and (0.26, 0.09, 0.02). So all the grid spacings are of the same order as the Kolmogorov length scale, implying that we are indeed resolving sufficient energy of the energy spectrum. Further, because Sc = Pr = 1, the concentration and temperature scales are also sufficiently resolved. The time step for the computation is calculated from the following criterion: 0.15  . t =  u   u   uφ   1 2  1 2  1 2 r θ  + + + ν r + rθ + r sin θφ r rθ r sin θφ  The non-dimensional time step used in the simulations was approximately 10−3 . This time step is also sufficiently smaller than the advection time scale of the smallest eddies, implying sufficient temporal resulotion as well.

286

A. AGRAWAL ET AL.

For the results presented in this paper, the statistical average is performed over approximately 40 and 12 eddy turnover times near the orifice exit and end of the HIZ respectively (the time scale of the flow becomes large with increasing distance to the jet nozzle). As shown in Figure 1, the jet orifice is confined to a central area of radius equal to 10 grid units in the cross-stream direction on the lower spherical surface of the computational domain, i.e., the orifice area is 1/64 of the total lower surface area. The domain extends from 0 ≤ z/d ≤ 43 in the streamwise direction, and the HIZ is located between 23 ≤ z/d ≤ 31. The ratio of the streamwise length of the HIZ to the nominal jet width is ≈ 3 and matches that of the experiments [2, 6]. The computations extend beyond the HIZ to capture the jet behavior upon exiting the HIZ. However the results for the last few d’s are not reliable due to the influence of the outflow boundary condition. 2.4. INITIAL AND BOUNDARY CONDITIONS The boundary conditions for the inflow, outflow and lateral boundaries of the jet have been discussed in detail in [9, 16]. Therefore, these are only briefly discussed here. 2.4.1. Inflow Boundary The axial velocity component and concentration are equal to U0 and C0 at the jet orifice, and zero everywhere else on the boundary, i.e., we have a top-hat profile at the orifice for these two quantities. A small sinusoidal perturbational with an amplitude equal to 2% of the velocity at the orifice is also added to enable an easier transition to turbulent flow. The frequency of the perturbation expressed in non-dimensional terms corresponds to a Strouhal number of 0.3 which is close to the natural frequency of the jet. The pressure on the inflow plane is left free. 2.4.2. Outflow Boundary At the outflow, we apply the boundary conditions proposed by Akselvoll and Moin [4]. These are ∂ur ∂ur = −U¯ , ∂t ∂r ∂uθ ∂uθ = −U¯ , ∂t ∂r ∂uφ ∂uφ = −U¯ , ∂t ∂r ∂C ∂C = −U¯ , ∂t ∂r

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

287

∂T ∂T = −U¯ . ∂t ∂r Here, U¯ is the mean streamwise velocity over the outflow boundary. Note that as the underlying coordinate system is spherical, we have written the velocity components in terms of radial, azimuthal and tangential velocities. As stated earlier, these velocity components correspond respectively to the streamwise (U ), cross-stream (V ), and out of plane (W ) components in the experimental setup of [2, 3]. 2.4.3. Lateral Boundary Traction free boundary conditions are applied on the lateral boundaries [13]. Mathematically,   ∂ui ∂uj + , σij .nj = 0 with σij = −pδij + ν ∂xj ∂xi where σij is the stress tensor, and nj is the unit normal on the boundary, and δij is the Kronecker delta. This boundary condition has been evaluated extensively (see [8, 9] for details) for a turbulent jet. The main advantage of this traction-free condition over a free-slip or no-slip boundary condition is that a velocity across the boundary is allowed, thereby allowing physical entrainment into the computational domain. The gradient of concentration along the normal of the lateral boundary is set to zero. It should be noted that the lateral boundary condition is stable for fluid flowing into the computational domain only if the cell Reynolds number at the lateral boundary, defined as Rec =

uθ rθ ν

is smaller than 2. However, this is not a serious constraint because usually the cross-stream velocity is small, and therefore grid refinement near the lateral edges is not required. 2.4.4. Initial Condition A self-similar jet as in [9] was first obtained. This was achieved by allowing the jet to develop in a quiescent environment for a long time. Heat injection began only after the jet had achieved complete self-similarity. Our results confirm that even the second order moments have reached self-similarity before the jet entered the HIZ. 3. Results and Discussion We will now present results from the simulations and discuss them in detail. The results have been organized as instantaneous and time-averaged behavior, variation

288

A. AGRAWAL ET AL.

of the fluctuating quantities at the centerline, and streamwise variation of the integral quantities. Under each subheading, the results are further separated for the three zones: pre-HIZ, HIZ, and post-HIZ.

3.1. INSTANTANEOUS AND TIME - AVERAGED CONTOURS Figure 2a shows a contour plot for instantaneous scalar concentration in the jet. (We have arbitrarily chosen one of the 24 axial planes for illustration.) Dyed fluid can be seen ejecting from the orifice with a top-hat profile. The potential core is apparent uptil about z/d = 5. The shear layer at the edge of the jet is seen to roll up beginning around z/d = 2 to form ring vortices. The vortices appear to be break up beyond z/d = 7. The scalar concentration drops continuously with both the axial and radial coordinates. In fact, the drop is quite monotonic, and is due to the dilution of jet fluid caused by entrainment of ambient fluid as the jet proceeds downstream. The drop in scalar concentration appears somewhat sharper in the HIZ and beyond it, and suggests an increased amount of entrainment for these regions. (The HIZ corresponds to 23 ≤ z/d ≤ 31.) Note, that the contour levels have been deliberately chosen to be non-uniform to highlight the concentration values in the HIZ and beyond. The jet appears to ascend like a column in the HIZ, i.e., with little lateral spread [2, 6]. We have included an instantaneous laser induced fluorescence (LIF) image of the heated jet in Figure 3 for comparison with the simulated instantaneous scalar concentration field of Figure 2a. Although there is a difference between the parameters (as indicated in the caption of Figure 3), the general agreement between the experimental and numerical results is very good. Figure 2b shows a contour plot for instantaneous streamwise velocity. The velocity width of the jet is seen to continuously increase downstream, although there is some indication of the velocity width being constant in the HIZ. The velocity at the jet centerline however shows a non-monotonic variation – decreasing first and then increasing, contrary to the centerline velocity of a normal jet which decreases with the axial coordinate. Heat injection first arrests the decay of the centerline velocity and then causes the jet to accelerate. A similar observation has been made experimentally [2, 6]. (We refer the reader to [2, 3] for a comparison with experimental data. The mismatch in the flow parameters between the experiments and simulation should be kept in mind while comparing the results.) The contour plot for temperature in Figure 2c indicates that the temperature is unaffected before the HIZ is encountered. This is expected because initially the temperatures of the jet and ambient fluids are identical. Volumetric heating is encountered upon entering the HIZ, and therefore a temperature rise starts to develop from the beginning of the HIZ. Due to a direct coupling between the temperature rise and concentration, the temperature field in the HIZ appears to be well correlated with the scalar field (compare Figure 2a with Figure 2c). Note that Figures 2a–2c correspond to the same simulation time. The temperature field,

Figure 2. Instantaneous (a) scalar concentration, (b) streamwise velocity, and (c) temperature contours in a volumetrically heated jet. The HIZ is located between 23 ≤ z/d ≤ 31.

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

289

290

A. AGRAWAL ET AL.

Figure 3. An LIF image representing instantaneous scalar concentration for comparison with Figure 2a. The experiment corresponds to Re = 2000, Ri∗ = 0.5, with the HIZ located between 195 ≤ z/d ≤ 260. The dark band in the image is part of the support for the six wire grids which are visible further downstream. Note that the illuminating laser sheet is directed downward from the top of the tank and its focal plane is situated within the HIZ, resulting in a peak intensity in this region.

although highly correlated with the concentration field, appears uncorrelated with the velocity field. The maximum temperature appears to occur within the HIZ, and temperature decreases as expected in the post-HIZ. A thermal-like structure formed for the initial time steps after heat addition began. A vortex-ring was apparent at the head of this thermal. The lateral spread rate of the thermal as it traveled downstream appeared to be larger than for normal jets. The formation and travel of the thermal corresponds to the transient part of this problem. However, the jet soon achieved stationarity once the thermal exited

Figure 4. Time- and azimuthally-averaged (a) scalar concentration, (b) streamwise velocity, and (c) temperature contours in a volumetrically heated jet. The HIZ is located between 23 ≤ z/d ≤ 31.

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

291

292

A. AGRAWAL ET AL.

the computational domain. The plots for average quantities presented here do not include this transitory flow condition. Figure 4a shows a contour plot for time-averaged scalar concentration in the jet. This plot was obtained after averaging in the azimuthal direction as well, therefore only one half-plane needs to be shown. The ring vortices in Figure 2a are not evident after time-averaging. The potential core is still apparent uptil about z/d = 5. The large dilution in the HIZ is also apparent from the low value of concentration in this region, and beyond it. The scalar flux is directly related to both the mass flux and scalar concentration. As seen in a subsequent figure, the mass flux is rather large in this region. This implies that the volume averaged concentration in the HIZ has to drop to conserve the scalar flux. Figure 4b shows a contour plot for time-averaged streamwise velocity. For a decelerating velocity field the contours will tend to converge towards the jet centerline, and this is evident in the pre-HIZ. Within the HIZ and beyond, the contours clearly diverge away from the centerline and this is indicative of an accelerating flow in these regions. A roughly constant velocity width in the HIZ is more evident after averaging. An increasing streamwise velocity is seen in the post-HIZ. The acceleration of the jet is due to the added buoyancy (as confirmed later the buoyancy flux attains its maximum value towards the end of the HIZ). Figure 4c shows a contour plot for time-averaged temperature. The temperature width seems nearly uniform (also evident in Figure 2c). There is some evidence of the maximum temperature occurring slightly off-axis. We will highlight this off-axis temperature peak by presenting radial profiles (see Figure 6 below). Figure 5 shows a plot for time-averaged radial velocity. The velocity is indeed small in the pre-HIZ. The time-averaged radial velocity profile in the pre-HIZ is zero at the centerline and again experiences a zero-crossing at r/z ≈ 0.1, i.e., there is a region of outflow in the neighborhood of the axis, and a region of inflow for larger radial locations in the pre-HIZ. The outflow near the jet axis is due to the decay of the centerline velocity with the axial coordinate [1, 2]. However, heat addition causes a dramatic change in the radial velocity profile. Figure 5 shows that the radial velocities are large in magnitude and negative in the HIZ and post-HIZ. In other words, the neighborhood of outflow surrounding the centerline vanishes and the radial velocities turn inwards for all radial locations with the addition of heat. A similar behavior was observed experimentally for a much smaller heating rate [2]. This inflow is related to the streamwise acceleration of the jet with heat addition. The absence of an outflow region and large magnitude of inflow velocity suggests that the mass flux should increase at a faster rate than for normal jets for the simulated conditions.

3.2. TIME - AVERAGED PROFILES The profiles for time-averaged concentration, velocity and temperature at different downstream locations for the jet reveal several interesting features. Note that once

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

293

Figure 5. Time- and azimuthally-averaged profile for radial velocity. The HIZ is located between 23 ≤ z/d ≤ 31.

again, the profiles are also averaged in the azimuthal plane to obtain smoother curves. Although not shown, the time-averaged concentration profile shows a Gaussian behavior in the pre-HIZ. However heat injection causes the concentration profile to change in shape from a Gaussian to a double-peaked Gaussian in the HIZ (Figure 6). (Note that we only see an off-axis peak in Figure 6, however if the data were reflected on to the opposite side of the jet-axis, a second peak would appear there. Therefore, the term “double-peaked Gaussian” introduced in [2] is retained here.) The reason for this change in profile shape is provided below. Also note the rather small value of the concentration in Figure 6. Similar to the concentration profile, it was verified that the streamwise velocity is well described by a Gaussian in the pre-HIZ. Heat injection causes a change in

294

A. AGRAWAL ET AL.

Figure 6. Time- and azimuthally-averaged profile at z/d = 32 for streamwise velocity, temperature, and scalar concentration. The HIZ is located between 23 ≤ z/d ≤ 31.

the shape of the streamwise velocity profile. A double-peaked Gaussian can again be discerned (Figure 6). Similarly, the time-averaged temperature profile shows an interesting double-peak in the HIZ, i.e., the maximum temperature rise occurs some distance away from the jet centerline as also noted above. A strong correlation between the average concentration and temperature profiles is evident. The reason for the change to a double-peaked shape of the temperature profile is because of the two competing effects of heat injection and residence time in the HIZ [2]. For a Gaussian velocity and scalar distribution, the fluid near the centerline exhibits the highest scalar concentration but has the least residence time in the HIZ; similarly, the fluid at the jet edge has the largest residence time but the smallest scalar concentration. Now, the temperature rise increases with residence time, and

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

295

Figure 7. Axial variation in the centerline streamwise velocity, temperature rise and concentration. The vertical lines indicate the location of the HIZ.

the rate of heat injection (which is proportional to scalar concentration). These two competing effects result in a temperature rise which is maximized away from the jet centerline. The centerline variation in average velocity, concentration and temperature for the heated jet is shown in Figure 7. The curves exhibit a small amount of jitter and some non-monotonic variations because of inadequate averaging. Values at the centerline do not benefit from azimuthal averaging because azimuthal averaging contributes nothing right at the centerline; therefore, these quantities are essentially only time-averaged. Yet, it is apparent that the centerline variation of these quantities is different in the pre-HIZ, HIZ and post-HIZ. Essentially, the jet

296

A. AGRAWAL ET AL.

behaves like a normal jet in the pre-HIZ, implying that Uc ∼ z−1 , Cc ∼ z−1 , and θ = 0. The constants Bu and Bc characterizing the axial decay of the centerline velocity and concentration are respectively 5.6 and 4.2 in excellent agreement with previous researchers [23]. The decay of the centerline velocity for the jet is reversed after the jet enters the HIZ – the jet actually accelerates in the HIZ. The centerline value of concentration continuously drops in the HIZ to a value substantially smaller than for the unheated jet. The temperature starts to rise at the beginning of the HIZ. In fact, the rise is very rapid towards the beginning of the HIZ. The temperature rise is expected to reach a maximum at the exit of the HIZ due to a continuous and cumulative addition of heat during the passage of fluid through the HIZ. However it is found that the location of maximum temperature rise unexpectedly occurs a short distance upstream of the exit of the HIZ. So what causes this rather anomalous temperature behavior inside the HIZ? There are several dynamic factors which determine the overall behavior of the jet as it passes through the HIZ. The dominant factor is the increased amount of mixing of the ambient fluid with the jet fluid. The enhanced mixing can lower the temperature distribution in the HIZ in three ways: (i) The ambient fluid is always cooler than the jet fluid, and therefore mixing directly lowers the temperature of the jet. (ii) Mixing lowers the scalar concentration in the jet which in turn reduces the strength of the source term in the temperature equation. (iii) Temperature being an active scalar, increases the velocities in the jet, thereby reducing the residence time of the fluid in the HIZ, and therefore affect the temperature. It so happens that the combined effect of these factors results in a temperature drop which more than compensates for the temperature rise due to heat addition. The maximum temperature is therefore attained before the fluid exits the HIZ. The decay in the centerline concentration and temperature is apparent in the post-HIZ. This is caused by mixing of the jet fluid with the non-acidified (i.e., no scalar), cold ambient fluid. The velocity however seems to continuously increase due to the presence of heated, buoyant fluid. The scalar, velocity and temperature widths for the jet are plotted in Figure 8. The width is defined as the radial distance at which the value of the variable drops to 1/e of the centerline value. Both the scalar and the velocity widths increase linearly in the pre-HIZ. The ratio of scalar to velocity width is approximately 1.2 in agreement with previous measurements [22, 23]. This behavior for scalar and velocity widths is expected for a normal jet. A large increase in the scalar width can be noted just upstream of the HIZ. A bulge in the scalar width has been observed experimentally at the beginning of the HIZ for lower heating rates [3, 6]. Interestingly, the scalar width is less than the velocity width in the HIZ, and becomes larger again in the post-HIZ. The difference arises because concentration changes from an active scalar in the HIZ to a passive scalar in the post-HIZ, and therefore begins to display a behavior in post-HIZ consistent with normal jets and plumes.

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

297

Figure 8. Axial variation in the velocity, scalar and temperature widths for the jet. The vertical lines indicate the location of the HIZ.

The jet remains isothermal in the pre-HIZ, therefore, the temperature width is undefined in this region. The temperature width is intermediate between the scalar and velocity widths in the HIZ. As stated earlier, this is because temperature is dependent on the rate of heat injection (concentration) and the residence time (velocity) in the HIZ. In the post-HIZ the temperature width is smaller than the velocity width which in turn is smaller than the concentration width, i.e., this observation is consistent with plumes. However, this observation itself does not imply a plume-like behavior for the jet in the post-HIZ. The difference in the temperature width for the two regions arises because the temperature is decoupled from the concentration in the post-HIZ.

298

A. AGRAWAL ET AL.

Figure 9. Axial variation in the normalized fluctuations at the jet centerline for the streamwise velocity and concentration. The vertical lines indicate the location of the HIZ.

3.3. PROFILES OF FLUCTUATING QUANTITIES Figure 9 shows the normalized fluctuations for concentration, and the streamwise velocity at the jet centerline, while Figure 10 depicts it for the radial and azimuthal velocities and temperature. (We have divided these results into two plots for clarity.) Normalization is performed by dividing by the mean centerline concentration, streamwise velocity, or temperature as appropriate. It is immediately apparent that the amount of jitter is much larger for these cases, owing to three causes: (i) second-order moments take longer to converge than first-order moments, (ii) we are normalizing by mean centerline values which themselves suffer from a small amount of jitter, and (iii) as in the case of the mean centerline values, the rms

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

299

Figure 10. Axial variation in the normalized fluctuations at the jet centerline for the radial and azimuthal velocities, and temperature. The vertical lines indicate the location of the HIZ.

values of the centerline cannot benefit from azimuthal averaging. Nevertheless, it is apparent that the rms values for concentration, radial and azimuthal velocities are very small for the first five diameters of the jet, and correspond to the region of the potential core. Beyond the potential core, the rms values gradually increase to their fully developed values. The fully developed normalized rms for concentration and for the three velocity components from our computations is approximately 0.2, in good agreement with previous measurements [14]. As noted earlier, the data at the centerline exhibits some jitter due to inadequate sampling. The zero rms for temperature in the pre-HIZ reflects the isothermal conditions in this region. The normalized rms for concentration, streamwise velocity and temperature increases rapidly inside the HIZ (note that the streamwise rms velocity shows

300

A. AGRAWAL ET AL.

an initial decrease). Further, the temperature fluctuations are well correlated with those of concentration similar to their mean values. This strong correlation indicates that the temperature at a point in the HIZ is strongly dependent on the local value of concentration (through the source term in the energy equation) at least for the coefficients employed here. The large value of the normalized rms for streamwise velocity indicates that relaminarization of the jet does not occur with heat injection. The experiments of Bhat and Narasimha [6] had revealed a complex behavior of rms of streamwise velocity fluctuations. Specifically, [5, 6] had reported a decrease in their normalized fluctuations, whereas the experimental findings in [3] have indicated that the normalized rms for streamwise velocity increases. The present computations, albeit at a higher heating rate, suggest a complex behavior of the rms fluctuations for the streamwise velocity within the HIZ that possibly contributed to the contradictory results reported earlier. To explain the rather dramatic increase in normalized fluctuations for streamwise velocity in the latter part of the HIZ, we propose that the jet fluid entering the HIZ could experience periods of higher velocity followed by periods of slow moving fluid. Note that the same explanation applies for the large fluctuations in scalar concentration and temperature. The normalized rms for the cross-stream velocity component remains virtually unchanged with heat addition. The rms of azimuthal velocity shows a decrease in the HIZ and then rises to the original value. The initial decrease is because the velocity fluctuations do not increase to the extent that the mean streamwise velocity increases resulting in a drop of the normalized value. However, subsequently, the fluctuations also increase leading to a nearly constant value for the normalized rms. The normalized rms for concentration and temperature decreases somewhat in the post-HIZ. However, the value is much higher than for a normal jet. This indicates that although the mean values for concentration and temperature are reducing, the fluctuations are reducing by a larger extent. Again note the strong correlation between concentration and temperature. As explained above, these two fields are very similar in the HIZ, and due to identical governing equations they will be advected in a similar manner in the post-HIZ. The fluctuations for streamwise velocity remain at almost the same value as in the HIZ. Similarly, the normalized rms for both cross-stream and azimuthal velocity components remains virtually unchanged from the value attained at the end of the HIZ.

3.4. MASS , MOMENTUM AND BUOYANCY FLUX The contour plot for the concentration field, and the radial velocity profile has suggested that the mass flux for the heated jet is larger than the corresponding unheated case. The actual calculations for the mass flux for the heated jet are shown in Figure 11. The mass and momentum fluxes were obtained by integrating the streamwise velocity. The variation for the corresponding normal jet is also shown in the figure. The mass flux for the heated jet increases linearly in the pre-HIZ as

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

301

Figure 11. Variation of mass and momentum flux for a heated and a normal jet. The vertical lines indicate the location of the HIZ.

expected for a turbulent axisymmetric jet. Further, the variation of mass flux in the pre-HIZ compares well with previous findings. The small bulge in the mass flux between 6 ≤ z/d ≤ 14 probably arises from the transition to turbulence. The mass flux for the heated jet increases rapidly once the jet enters the HIZ. This suggests an increased amount of entrainment for the heated jet in the HIZ. Further, the dilution rate for the heated jet in and even beyond the HIZ is much larger than a normal jet, and in fact exceeds the dilution rate in the HIZ. The large mass flux is related to acceleration of the heated jet. As the heated jet begins to accelerate due to the added buoyancy, continuity dictates that a larger amount of ambient fluid is drawn inwards to the centerline.

302

A. AGRAWAL ET AL.

We will now compare the mass flux of heated and unheated jets. The value of I can be estimated from Equation (12), where the mean value of C can be obtained from Figure 7. Inserting I ≈ 8d 3 into Equation (15) we can estimate the value of f as about 10, implying that the average mass flux through the HIZ is about 10µ0 . This matches well with the mid-HIZ mass flux value in Figure 11. Note that the corresponding value for an unheated jet would be about 7. It is also possible to compare the value of the entrainment coefficient α for the heated jet relative to its unheated counterpart by using the following equation [22] dµ = 2π bUc α, dz and inserting the values of relative mass-flux gradients, velocity widths, and centerline velocities from the respective plots. Our simulations indicate that the entrainment coefficient increases substantially in the HIZ and beyond it. For example, the entrainment coefficient for the heated jet approximately doubles at the end of the HIZ relative to an unheated jet. Thus, for the value of Ri used in this simulation (about 12) we find that both mass flux and entrainment increase with heat addition. A direct comparison with experiments is not possible owing to the smaller Ri (about 0.3) used therein. Nevertheless, it is relevant to note that Bhat and Narasimha [6] reported a decrease in mass flux and entrainment with heat addition, whereas the measurements of [3, 11] indicated a larger mass flux for the heated jet. A larger mass flux with heating is also predicted by the empirical model of [2]. Figure 11 also shows the momentum flux (calculated using the mean velocity) variation for the heated jet. The momentum flux is nearly constant in the pre-HIZ, as expected for a normal jet. Buoyancy addition in the HIZ results in a continuous increase in the momentum flux for the jet in the HIZ. The added buoyancy drives a continuous increase in the momentum flux in the post-HIZ (somewhat similar to that in a plume). (A slight drop in the mass and momentum flux beyond z/d = 40 is due to the effect of the outflow boundary condition.) Figure 11 also shows the axial buoyancy flux, B obtained by integrating the product of streamwise velocity and temperature over the cross-sectional area of the jet, i.e., 2π ∞ B(z) =

U θr dr dφ, φ=0 r=0

where dφ is the incremental change in the azimuthal angle. The buoyancy flux is zero due to isothermal conditions in the pre-HIZ. However the buoyancy flux increases monotonically in the HIZ due to continuous heat addition into the jet. Because buoyancy is neither added or removed in the post-HIZ, one expects the buoyancy flux to be conserved in this region. Our computations indicate that this is indeed true, i.e., a decrease in the product of velocity and temperature is exactly

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

303

compensated by an increase in the cross-sectional area of the jet. Further, it can be shown that the energy balance dictates that the buoyancy flux at the end of the HIZ must equal Q/ρCp . This reduces to kI Bend HIZ = 3 . d As stated earlier, we found I ≈ 8d 3 , therefore Bend HIZ ≈ 8 (because k = 1). This result is consistent with Figure 11. 4. Conclusions A three-way coupled system involving momentum, temperature and scalar concentration is solved using direct numerical simulations in a flow field consisting of an off-source volumetrically heated jet. A strong jet acceleration is noted for the value of the parameters chosen here, and this strong acceleration has a marked effect on the jet behavior. These simulations complement previously reported experiments which used laboratory analogues to explore cumulus cloud entrainment, where such large accelerations were difficult to achieve. The profiles for concentration, streamwise velocity and temperature are found to be double-peaked Gaussians. Further, there is a strong correlation between concentration and temperature. The appearance of a double peaked Gaussian profile is related to the competing effects of residence time and heat injection rate, and elucidates the strong coupling between the three fields. The centerline streamwise velocity decreases in the pre-HIZ, increases due to the cumulative effect of heat addition in the HIZ, and then shows a non-monotonic trend in the post-HIZ. The centerline concentration decreases monotonically throughout. The centerline temperature attains its maximum value inside the HIZ, while subsequent mixing of cold ambient fluid with the jet fluid results in a drop in temperature. The radial velocity turns inward for all radial locations in the HIZ due to the acceleration of the jet. Heat injection causes an eventual increase in the normalized fluctuation of the concentration, streamwise velocity and temperature, while that of the radial and azimuthal velocities remains virtually unchanged. The mass and momentum fluxes increase with heat addition. The entrainment rate is also larger for the heated jet as compared to its unheated counterpart. The buoyancy flux behaves as theoretically predicted. Combining the insights gained experimentally and numerically at both low [2, 3] and high heating rates, it is evident that there is a change in the shape of U , V and T profiles and an interesting evolution of these quantities along the centerline. There is an eventual increase in the normalized rms of streamwise velocity. The mass flux of the heated jet is larger than its unheated counterpart, and at least for high heating rates, heating leads to increased entrainment. The insights gained by these experiments and simulations have helped in the development of a conceptual model for cumulus cloud entrainment [1].

304

A. AGRAWAL ET AL.

Acknowledgements This work was supported by the National Science Foundation under grant NSFATM-0095122. The first author gratefully acknowledges the travel grant by the Burgers Center, Delft University of Technology, and the use of the computing facilities at Delft. References 1. 2.

3. 4. 5. 6. 7. 8. 9.

10. 11. 12. 13. 14. 15.

16. 17. 18. 19. 20. 21.

Agrawal, A., Effect of off-source volumetric heat addition on entrainment in a turbulent jet with application to cumulus clouds. Ph.D. Dissertation, University of Delaware (2002). Agrawal, A., Sreenivas, K.R. and Prasad, A.K., Velocity and temperature measurements in an axisymmetric jet with cloud-like off-source heating. Internat. J. Heat Mass Transfer 47 (2004) 1433–1444. Agrawal, A. and Prasad, A.K., Evolution of a turbulent jet subjected to volumetric heating. J. Fluid Mech. (in review). Akselvoll, K. and Moin, P., Large-eddy simulation of turbulent confined coannular jets. J. Fluid Mech. 315 (1996) 387–411. Basu, A.J. and Narasimha, R., Direct numerical simulation of turbulent flows with cloud-like off-source heating. J. Fluid Mech. 385 (1999) 199–228. Bhat, G.S. and Narasimha, R., A volumetrically heated jet: Large eddy structure and entrainment characteristics. J. Fluid Mech. 325 (1996) 303–330. Bird, R.B., Stewart, W.E. and Lightfoot, E.N., Transport Phenomena. Wiley, New York (1960). Boersma, B.J., Entrainment boundary conditions for free shear flows. Internat. J. Comput. Fluid Dynam. 13 (2000) 357–363. Boersma, B.J., Brethouwer, G. and Nieuwstadt, F.T.M., A numerical investigation on the effect of the inflow conditions on the self-similar region of a round jet. Phys. Fluids 10 (1998) 899– 909. Chen, C.J. and Rodi, W., Vertical Turbulent Buoyant Jets – A Review of Experimental Data. Pergamon Press, Oxford (1980) pp. 11–12. Elavarasan, R., Bhat, G.S., Narasimha, R. and Prabhu, A., An experimental study of a jet with local buoyancy enhancement. Fluid Dynam. Res. 16 (1995) 189–202. Emanuel, K.A., Atmospheric Convection. Oxford University Press, Oxford (1994). Gresho, P.M., Incompressible fluid dynamics: Some fundamental formulations issues. Annu. Rev. Fluid Mech. 23 (1991) 413–453. Hussein, H.J., Capp, S.P. and George, W.K., Velocity measurements in a high Reynolds number, momentum-conserving axisymmetric turbulent jet. J. Fluid Mech. 258 (1994) 31–75. Koren, B., A robust upwind discretization method for advection, diffusion and source terms. In: Vreugdenhill, C.B. and Koren, B. (eds), Numerical Methods for Advection-Diffusion Problems, Notes on Numerical Fluid Mechanics, Vol. 45. Vieweg, Braunschweig (1993) pp. 117–138. Lubbers, C.L., Brethouwer, G. and Boersma, B.J., Simulation of the mixing of a passive scalar in a round turbulent jet. Fluid Dynam. Res. 28 (2001) 189–208. Mohensi, K. and Colonius, T., Numerical treatment of polar coordinate singularities. J. Comput. Phys. 157 (2000) 787–795. Moin, P. and Mahesh, K., Direct numerical simulation: A tool in turbulence research. Annu. Rev. Fluid Mech. 30 (1998) 539–578. Paluch, I.R., The entrainment mechanism in Colorado cumuli. J. Atmos. Sci. 36 (1979) 2467– 2478. Rodi, W., Private communication (2003). Schlichting, H., Boundary-Layer Theory, seventh edition. McGraw-Hill, New York (1979).

DNS OF A TURBULENT AXISYMMETRIC JET WITH BUOYANCY

22. 23. 24. 25.

305

Turner, J.S., Turbulent entrainment: The development of the entrainment assumption, and its application to geophysical flows. J. Fluid Mech. 173 (1986) 431–471. Wang, H. and Law, A.W.-K., Second order integral model for a round turbulent buoyant jet. J. Fluid Mech. 459 (2002) 397–428. Warner, J., On steady-state one-dimensional models of cumulus convections. J. Atmos. Sci. 27 (1970) 1035–1040. Zijlema, M. and Wesseling, P., Higher-order flux-limiting schemes for the finite volume computation of incompressible flow. Internat. J. Comput. Fluid Dynam. 9 (1998) 89–109.