UFR 1-06 Test Case: Difference between revisions
No edit summary |
m (Dave.Ellacott moved page SilverP:UFR 1-06 Test Case to UFR 1-06 Test Case) |
||
(199 intermediate revisions by the same user not shown) | |||
Line 4: | Line 4: | ||
}} | }} | ||
= Axisymmetric buoyant far-field plume | = Axisymmetric buoyant far-field plume = | ||
Underlying Flow Regime 1-06 | Underlying Flow Regime 1-06 | ||
Line 20: | Line 20: | ||
* The hot air is discharged at a velocity of U0 = 67 cm/s with a approximately a top-hat profile. | * The hot air is discharged at a velocity of U0 = 67 cm/s with a approximately a top-hat profile. | ||
* Temperature and velocity fluctuations at the inlet are less than 0.1%. | * Temperature and velocity fluctuations at the inlet are less than 0.1%. | ||
* George ''et al.'' [[UFR_1-06_References|[3]]] present experimentally measured profiles | * George ''et al.'' [[UFR_1-06_References|[3]]] present experimentally measured profiles of both mean and fluctuating components of the temperature and axial velocity in the self-similar region at x/D = 8, 12 and 16 above the source. | ||
of both mean and fluctuating components of the temperature and axial velocity in the self-similar region | |||
at x/D = 8, 12 and 16 above the source. | |||
== Test Case Experiments == | == Test Case Experiments == | ||
The experiments used in this UFR are those of George ''et al''. [[UFR_1-06_References|[3]]] | |||
which were conducted in 1974 at the Factory Mutual Research Corporation and were subsequently repeated by | |||
Shabbir & George [[UFR_1-06_References|[34]]] at the University of Buffalo. | |||
The general arrangement is shown in Figure 4. | |||
Compressed air is passed through a set of heaters and porous mesh screens before exiting through a nozzle into | |||
the enclosure. The nozzle is stated as a 15:1 contraction in [[UFR_1-06_References|[3]]], | |||
a 12:1 contraction in [sg92] | |||
and appears to be different again in a drawing of the arrangement in [[UFR_1-06_References|[3]]] | |||
(see Figure 5). It resulted in a velocity | |||
profile through the exit which was uniform to within 2% outside the wall boundary layer. The velocity and | |||
temperature fluctuations at the exit were measured to be very low, less than 0.1% | |||
in [[UFR_1-06_References|[3]]] | |||
and 0.5% in [[UFR_1-06_References|[34]]]. | |||
The temperature of the source was 300°C and the ambient environment | |||
29°C. Both were controlled to an accuracy of within 1°C. The discharge velocity was 67 cm/s, as | |||
calculated from the measured heat flux. These source conditions corresponded to Reynolds number, ''Re''<sub>0</sub> = 870, and densimetric Froude number, ''Fr''<sub>0</sub> = 1.23 | |||
<ref name="ftn3"> | |||
The densimetric Froude number is calculated here from the source and | |||
ambient temperatures, the exit velocity and source diameter given by | |||
George ''et al''. [[UFR_1-06_References|[3]]], using Equation (1). | |||
However George ''et al''. [[UFR_1-06_References|[3]]] stated that the densimetric Froude number | |||
was 1.4. It is unclear how they determined this value. Using the approach taken by | |||
Chen & Rodi [[UFR_1-06_References|[1]]] in which the source density instead of the | |||
ambient density is used to make the density difference dimensionless, and Froude number is defined using | |||
the square of the expression given in Equation (1), this gives a Froude number of 0.80. | |||
</ref> | |||
There was no evidence of laminar flow behaviour at a position two inlet | |||
diameters downstream from the source. The effective origin of the plume, ''x''<sub>0</sub>, was | |||
found to be at the same location as the exit | |||
(see [[UFR_1-06_References|[3]]] for details of how this was determined). | |||
The screen enclosure around the plume exit was 2.44 × 2.44 metres in cross-section and 2.44 metres high | |||
(there is, presumably, an error in [[UFR_1-06_References|[3]]] | |||
which suggests that the enclosure is 2.44 × 2.44 × 2.44 mm). | |||
In the later Shabbir & George experiments, a 2 × 2 × 5 metre enclosure was used. The purpose | |||
of the screens was to minimize the effect of cross-draughts and other disturbances affecting the flow. | |||
Two-wire probes were used by George ''et al''. [[UFR_1-06_References|[3]]] | |||
to record velocities and temperature. | |||
[[Image:UFR1-06_figure5.gif|center|400px]] | |||
<center>'''Figure 5''' Schematic of the George ''et al.'' [[UFR_1-06_References|[3]]] experiments, | |||
from Shabir & George [[UFR_1-06_References|[11]]]</center> | |||
[[Image:UFR1-06_figure6.gif|center|400px]] | |||
<center>'''Figure 6''' Schematic of the plume generator used in the experiments, | |||
from George ''et al.'' [[UFR_1-06_References|[3]]]</center> | |||
George ''et al''. [[UFR_1-06_References|[3]]] reported that measurement errors, stemming from | |||
directional ambiguity of the hot wire and its thermal inertia, were around 3% for the velocity and lower for | |||
other mean and RMS values. The frequency response of the hot wires was estimated to be around 300 Hz compared | |||
to the frequency of the energy-containing eddies at around 50 Hz and the Kolmogorov microscale at 1 Khz. | |||
It was noted that measurement errors were likely to be higher on the outer edge of the plume where | |||
the velocity fluctuations were higher. | |||
In their review of plume experiments, Chen & Rodi [[UFR_1-06_References|[1]]] | |||
noted that the data from George ''et al.'' differed significantly from earlier measurements by | |||
Rouse ''et al''. [[UFR_1-06_References|[64]]]. | |||
However, they considered it to be more reliable due to its use of more | |||
sophisticated instrumentation. | |||
George [[UFR_1-06_References|[40]]], describes an experimental program at the University of | |||
Buffalo that was set up following publication of the original | |||
George ''et al.'' [[UFR_1-06_References|[3]]] paper to investigate | |||
possible causes of differences in experimental plume results. Possible sources of errors discussed | |||
included: | |||
* ambient thermal stratification | |||
* the size of the enclosure | |||
* the use of porous screens used to minimise disturbances from the far-field affected the plume source. | |||
* hot wire measurement errors | |||
The most significant concern was ambient thermal stratification. | |||
One of the features of buoyant plumes in neutral environments is that the integral of the | |||
buoyancy across the whole cross-section of the plume, ''F'', should remain constant and equal to the | |||
buoyancy added at the source, ''F''<sub>0</sub>. | |||
George [[UFR_1-06_References|[40]]] discussed how thermal stratification involving | |||
small temperature differences of the order of 1°C across a 3 metre vertical span would be sufficient | |||
to cause ''F'' to decrease to only 50% of the source value. | |||
This would be likely to cause differences in measured temperature and velocity plume profiles. | |||
In the initial experiments of George ''et al''. [[UFR_1-06_References|[3]]], | |||
the thermal stratification was not strictly controlled. | |||
However, results from later experiments published in the PhD thesis of | |||
Shabbir [[UFR_1-06_References|[32]]] | |||
(reproduced in [[UFR_1-06_References|[34]]] and [[UFR_1-06_References|[40]]]), | |||
which conserved buoyancy to within 10%, are in good agreement with the earlier results from | |||
George ''et al''. [[UFR_1-06_References|[3]]]. | |||
This suggests that, perhaps fortunately, ambient thermal stratification did not contaminate the | |||
George ''et al''. [[UFR_1-06_References|[3]]] results significantly. | |||
A summary of the original results from George ''et al''. [[UFR_1-06_References|[3]]] | |||
and those reproduced later by Shabbir & George [[UFR_1-06_References|[11]]] | |||
is presented in Table 3. | |||
Also shown are the recommended values from | |||
Chen & Rodi's review [[UFR_1-06_References|[1]]] and other studies. | |||
The parameters given in Table 3 relate to the following empirical formulae for the mean vertical velocity: | |||
{| | |||
|- | |||
|width="750" align="center"|<math>W=F_0^{-{1/3}}z^{-{1/3}}f_w</math> | |||
|width="40" align="right"|<math>\left(9\right)</math> | |||
|} | |||
and effective buoyancy acceleration: | |||
{| | |||
|- | |||
|width="750" align="center"|<math>g\frac{\Delta\rho}{\rho}=F_0^{2/3}z^{-{5/3}}f_T</math> | |||
|width="40" align="right"|<math>\left(10\right)</math> | |||
|} | |||
where <math>f_w</math> and <math>f_T</math> are Gaussian functions: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
f_w=A_w\exp\left[-B_w\left(\frac{r}{x}\right)^2\right]\qquad | |||
f_T=A_T\exp\left[-B_T\left(\frac{r}{x}\right)^2\right] | |||
</math> | |||
|width="40" align="right"|<math>\left(11\right)</math> | |||
|} | |||
The parameters, <math>l_{\Delta T/2}</math> and <math>l_{w/2}</math> are the dimensionless half-widths | |||
of the plume, as defined by the location where the normalized buoyancy or mean velocity falls to half | |||
its centreline value. | |||
The RMS temperature and axial velocity fluctuations normalized by their centreline mean values are | |||
denoted, <math>\left(\overline{t^2}\right)^{1/2}/\Delta T_c</math> | |||
and <math>\left(\overline{w^2}\right)^{1/2}/W_c</math>, respectively. | |||
As noted earlier, Dai ''et al.'' [[UFR_1-06_References|[10]]][[UFR_1-06_References|[37]]] | |||
[[UFR_1-06_References|[38]]][[UFR_1-06_References|[39]]][[UFR_1-06_References|[41]]] disputed the | |||
accuracy of the George ''et al.'' [[UFR_1-06_References|[3]]] experiments and suggested | |||
that they had made measurements too near the source, before the plume had reached a fully-developed state. | |||
Their arguments are disregarded by Shabbir & George [[UFR_1-06_References|[11]]] | |||
[[UFR_1-06_References|[34]]]. | |||
<center>'''Table 3''' Summary of mean flow parameters and turbulence intensities, from Shabbir & George [[UFR_1-06_References|[11]]]</center> | |||
[[Image:UFR1-06_table3.gif|center|650px]] | |||
== CFD Methods == | == CFD Methods == | ||
<math>\frac{\partial }{\partial x_{j}}\left(\overline{{\rho | == Van Maele & Merci: Description of CFD Work == | ||
=== Numerical Methods === | |||
Van Maele & Merci [[UFR_1-06_References|[2]]] | |||
used the finite-volume-based commercial CFD code, Fluent | |||
<ref name="ftn4">[http://www.fluent.com http://www.fluent.com]</ref>, | |||
to simulate the plume experiments of George ''et al.'' [[UFR_1-06_References|[3]]]. | |||
For the discretization of the convective terms in the momentum, turbulence and energy equations a second-order | |||
upwind scheme was used. Diffusion terms were discretized using second-order central differences and the SIMPLE | |||
algorithm was used for pressure-velocity coupling. The flow was treated as axisymmetric and elliptic | |||
calculations were performed used a Cartesian grid arrangement. | |||
The low-Mach-number form of the Favre-averaged Navier-Stokes equations were used. In this weakly-compressible | |||
approach, the density is treated as only a function of temperature and not pressure. Pressure only affects the | |||
flow field through the pressure-gradient term in the momentum equations. The ideal gas law is used to link the | |||
mean density, <math>\overline{\rho}</math>, to mean temperature, ''T'' as follows: | |||
{| | |||
|- | |||
|width="750" align="center"|<math>p_*=\overline{\rho}RT</math> | |||
|width="40" align="right"|<math>\left(12\right)</math> | |||
|} | |||
where ''p<sub>*</sub>'' is taken as constant and equal to the atmospheric pressure. The | |||
low-Mach-number approximation implies that the effect of the mean kinetic energy and the work done by viscous | |||
stresses and pressure are negligible in the energy equation. | |||
=== Turbulence Modelling === | |||
Two turbulence models were used by Van Maele & Merci [[UFR_1-06_References|[2]]]: | |||
the standard ''k – ε'' model of | |||
Jones & Launder [[UFR_1-06_References|[65]]] | |||
and the realizable ''k – ε'' model | |||
of Shih ''et al.'' [[UFR_1-06_References|[66]]]. | |||
In the former model, the eddy viscosity is given by: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
\mu=\overline{\rho}c_\mu\frac{k^2}{\varepsilon} | |||
</math> | |||
|width="40" align="right"|<math>\left(13\right)</math> | |||
|} | |||
where ''c<sub>μ</sub>'' is a constant equal to 0.09 and the standard ''k'' and ''ε'' | |||
equations are written: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
\frac{\partial }{\partial x_{j}}\left(\overline{{\rho | |||
}}U_{j}k\right)=\frac{\partial }{\partial x_{j}}\left[\left(\mu | |||
+\frac{\mu }{\sigma _{k}}\right)\frac{\partial k}{\partial | |||
x_{j}}\right]+P_{k}+G-\overline{{\rho }}\varepsilon | |||
</math> | |||
|width="40" align="right"|<math>\left(14\right)</math> | |||
|} | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
\frac{\partial }{\partial x_{j}}\left(\overline{{\rho | |||
}}U_{j}\varepsilon \right)=\frac{\partial }{\partial | |||
x_{j}}\left[\left(\mu +\frac{\mu }{\sigma _{\varepsilon | |||
}}\right)\frac{\partial \varepsilon }{\partial | |||
x_{j}}\right]+c_{\varepsilon 1}P_{k}\frac{\varepsilon | |||
}{k}-c_{\varepsilon 2}\overline{{\rho }}\frac{\varepsilon | |||
^{2}}{k}+S_{\varepsilon B} | |||
</math> | |||
|width="40" align="right"|<math>\left(15\right)</math> | |||
|} | |||
where ''c<sub>ε1</sub>'' = 1.44, | |||
''c<sub>ε2</sub>'' = 1.92, | |||
''σ<sub>k</sub>'' = 1.0, | |||
''σ<sub>ε</sub>'' = 1.3 and ''P<sub>k</sub>'' is the production term due to mean shear. | |||
The terms ''G'' and ''S<sub>εB</sub>'' are source terms related to the influence of buoyancy on | |||
the ''k'' and ''ε'' equations. The treatment of these terms is discussed below. | |||
The Shih ''et al.'' [[UFR_1-06_References|[66]]] model involves two changes to the standard | |||
''k – ε'' model. Firstly, ''c<sub>μ</sub>'' is made a function of strain and | |||
vorticity invariants to ensure that the model always returns positive normal Reynolds stresses and satisfies | |||
the Schwarz inequality for the turbulent shear stresses. The function form of ''c<sub>μ</sub>'' | |||
is given by: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
c_{\mu }=\left(A_{0}+A_{s}U^{(\ast {})}\frac{k}{\varepsilon | |||
}\right)^{-1} | |||
</math> | |||
|width="40" align="right"|<math>\left(16\right)</math> | |||
|} | |||
where: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
U^{(\ast {})}=\sqrt{S_{\mathit{ij}}S_{\mathit{ij}}+\Omega | |||
_{\mathit{ij}}\Omega _{\mathit{ij}}} \ \ \ | |||
S_{\mathit{ij}}=\frac{1}{2}\left(\frac{\partial U_{i}}{\partial | |||
x_{j}}+\frac{\partial U_{j}}{\partial x_{i}}\right)\ \ \Omega | |||
_{\mathit{ij}}=\frac{1}{2}\left(\frac{\partial U_{i}}{\partial | |||
x_{j}}-\frac{\partial U_{j}}{\partial x_{i}}\right) | |||
</math> | |||
|width="40" align="right"|<math>\left(17\right)</math> | |||
|} | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
A_{s}=\sqrt{6}\cos \phi \ \ \phi =\frac{1}{3}\arccos | |||
\left(\sqrt{6}W\right) \ \ \ | |||
W=2^{3/2}\frac{S_{\mathit{ij}}S_{\mathit{jk}}S_{\mathit{ki}}}{S^{3}} | |||
\ \ S=\sqrt{2S_{\mathit{ij}}S_{\mathit{ij}}} | |||
</math> | |||
|width="40" align="right"|<math>\left(18\right)</math> | |||
|} | |||
and ''A<sub>0</sub>'' is a constant equal to 4.04. | |||
Secondly, a different ''ε''-equation is used to resolve the problem of the round-jet/plane-jet | |||
anomaly (see Pope [[UFR_1-06_References|[67]]]): | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
\frac{\partial }{\partial x_{j}}\left(\overline{{\rho | |||
}}U_{j}\varepsilon \right)=\frac{\partial }{\partial | }}U_{j}\varepsilon \right)=\frac{\partial }{\partial | ||
x_{j}}\left[\left(\mu +\frac{\mu }{\sigma _{\varepsilon | x_{j}}\left[\left(\mu +\frac{\mu }{\sigma _{\varepsilon | ||
Line 34: | Line 330: | ||
x_{j}}\right]+c_{\varepsilon 1}\overline{{\rho }}S\frac{\varepsilon | x_{j}}\right]+c_{\varepsilon 1}\overline{{\rho }}S\frac{\varepsilon | ||
}{k}-c_{\varepsilon 2}\overline{{\rho }}\frac{\varepsilon | }{k}-c_{\varepsilon 2}\overline{{\rho }}\frac{\varepsilon | ||
^{2}}{k+\sqrt{\nu \varepsilon }}+S_{\varepsilon B}</math> | ^{2}}{k+\sqrt{\nu \varepsilon }}+S_{\varepsilon B} | ||
</math> | |||
|width="40" align="right"|<math>\left(19\right)</math> | |||
|} | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
c_{\varepsilon 1}=\mathit{max}\left[0.43,\eta /(\eta +5)\right] | |||
\ \ \ \ \eta =Sk/\varepsilon | |||
</math> | |||
|width="40" align="right"|<math>\left(20\right)</math> | |||
|} | |||
where ''S'' is the strain-rate invariant as before, ''c<sub>ε2</sub>'' = 1.9, | |||
''σ<sub>k</sub>'' = 1.0 and ''σ<sub>ε</sub>'' = 1.3. | |||
The Shih ''et al.'' [[UFR_1-06_References|[66]]] model was developed for high | |||
Reynolds number turbulent flows and therefore a zonal or wall-function | |||
approach must be used to bridge the viscous sub-layer near walls. | |||
Compared to the standard ''k–ε'' model, | |||
it has been shown to produce improved behaviour in a number of free | |||
shear flows, boundary-layer flows and a backward-facing step flow [[UFR_1-06_References|[66]]]. | |||
One of the major weaknesses of the standard ''k–ε'' model | |||
is that it produces too much turbulent | |||
kinetic energy at stagnation points [[UFR_1-06_References|[68]]]. | |||
The Shih ''et al.'' | |||
model should in principle suffer less from this weakness since the | |||
functional form of ''c<sub>μ</sub> should | |||
reduce the over-production of ''k''. However, its overall | |||
performance in stagnating flows will depend on the type of wall model | |||
used. | |||
==== Production due to Buoyancy, ''G'' ==== | |||
The term ''G'' in the ''k''-equation relates to the | |||
influence of buoyancy on the turbulent kinetic energy, and is given by: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
G=\overline{{\rho u_{j}}}g_{j} | |||
</math> | |||
|width="40" align="right"|<math>\left(21\right)</math> | |||
|} | |||
where ''g<sub>j</sub>'' is the gravitational | |||
acceleration vector. In stably stratified flows, where the temperature | |||
increases with height, ''G'' is negative. Conversely, in unstably | |||
stratified flows, where temperature decreases with height, ''G'' | |||
is positive and acts to increase ''k''. The unknown | |||
density-velocity correlation, <math>\overline{{\rho u_{j}}}</math>, must be | |||
modelled. The most common approximation of this term is the so-called | |||
Boussinesq Simple Gradient Diffusion Hypothesis (SGDH): | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
\overline{{\rho u_{j}}}=-{\frac{\mu _{t}}{\sigma | |||
_{t}}\frac{1}{\overline{{\rho }}}\frac{\partial \overline{{\rho | |||
}}}{\partial x_{j}}} | |||
</math> | |||
|width="40" align="right"|<math>\left(22\right)</math> | |||
|} | |||
The production due to buoyancy using SGDH is then as follows: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
G=-{\frac{\mu _{t}}{\sigma _{t}}\frac{1}{{\overline{{\rho | |||
}}}^{2}}\frac{\partial \overline{{\rho }}}{\partial x_{j}}\rho _{\infty | |||
}g_{j}} | |||
</math> | |||
|width="40" align="right"|<math>\left(23\right)</math> | |||
|} | |||
In their paper, Van Maele & Merci [[UFR_1-06_References|[2]]] | |||
erroneously included | |||
an additional pressure-gradient term <math>\left(\partial P/\partial | |||
x_{j}\right)</math> in Equation (23) related to the pressure-work rather than | |||
buoyancy (see Wilcox [[UFR_1-06_References|[69]]]). Since the term is negligible in | |||
incompressible flows, such as the buoyant plumes considered here, it | |||
has therefore been ignored. The ratio of the reference density to the | |||
mean density, <math>\rho _{\infty }/\overline{\rho }</math>, appears in | |||
Equation (23) due to the use of a non-Boussinesq approach and Favre-averaging, | |||
which are discussed later. Van Maele & Merci [[UFR_1-06_References|[2]]] assumed | |||
that ''σ<sub>t</sub>'' was constant and | |||
equal to 0.85. | |||
Instead of writing the buoyancy production in terms of the | |||
density-velocity correlation, <math>\overline{{\rho u_{j}}}</math>, the equation | |||
can be written in terms of the heat flux, <math>\overline{{u_{j}t^\prime}}</math>: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
\overline{{u_{j}t^\prime}}=-{\frac{\mu _{t}}{\sigma _{t}}\frac{\partial | |||
T}{\partial x_{j}}} | |||
</math> | |||
|width="40" align="right"|<math>\left(24\right)</math> | |||
|} | |||
and the ''G'' term is then written: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
G=-{\beta \frac{\mu _{t}}{\sigma _{t}}\frac{\partial T}{\partial | |||
x_{j}}g_{j}} | |||
</math> | |||
|width="40" align="right"|<math>\left(25\right)</math> | |||
|} | |||
where ''t′'' is the temperature fluctuation, ''T'' is the mean | |||
temperature and ''β'' is the volumetric expansion | |||
coefficient, <math>\beta =-\left(1/\overline{{\rho }}\right)\partial | |||
\overline{{\rho }}/\partial T</math>. Other equivalent expressions can also | |||
be formulated using the ideal gas law and the assumption that density | |||
is only a function of temperature, not pressure (the low-Mach-number | |||
approximation). The conversion from mean density to temperature | |||
gradients is then as follows: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
\frac{\partial \overline{{\rho }}}{\partial x_{j}}=\frac{\partial | |||
\overline{{\rho }}}{\partial T}\frac{\partial T}{\partial | |||
x_{j}}=\frac{-P_{\ast {}}}{RT^{2}}\frac{\partial T}{\partial | |||
x_{j}}=\frac{\overline{{\rho }}}{T}\frac{\partial T}{\partial | |||
x_{j}} | |||
</math> | |||
|width="40" align="right"|<math>\left(26\right)</math> | |||
|} | |||
The SGDH model predicts zero density-velocity correlation or heat flux | |||
components ( <math>\overline{{\rho u_{j}}}=0</math> or <math>\overline{{u_{j}t^\prime}}=0</math>) | |||
in situations were the density or temperature gradients are zero in | |||
that direction. However, as Ince & Launder [[UFR_1-06_References|[70]]] noted, in a | |||
simple shear flow in which there are only cross-stream temperature | |||
gradients, the heat flux in the streamwise direction actually exceeds | |||
that in the cross-stream direction. This shortcoming of the SGDH model | |||
was confirmed by the analysis of Shabbir & Taulbee [[UFR_1-06_References|[33]]], who | |||
showed that the model significantly underestimates the magnitude of the | |||
heat flux in vertical buoyant plumes. The underprediction of | |||
<math>\overline{{\rho u_{j}}}</math> or <math>\overline{{u_{j}t^\prime}}</math> by the SGDH model | |||
leads to an overly-small production term, ''G'', and hence a | |||
turbulent kinetic energy, ''k'', which is too small, producing too | |||
little mixing in the modelled plume. The study by | |||
Yan & Holmstedt [[UFR_1-06_References|[53]]] | |||
provides a clear example of how the ''k – ε'' | |||
model with SGDH produces buoyant plumes which | |||
are too narrow and with overly high temperatures and velocities in the | |||
core of the flow. | |||
Van Maele & Merci [[UFR_1-06_References|[2]]] | |||
examined a different model for | |||
''G'' based on the the Generalized-Gradient Diffusion Hypothesis | |||
(GGDH) of Daly & Harlow [[UFR_1-06_References|[52]]]. | |||
This was first used in the | |||
context of practical CFD calculations with the ''k – ε'' | |||
model by Ince & Launder [[UFR_1-06_References|[70]]], and is | |||
written as follows: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
\overline{{\rho u_{j}}}=-{\frac{3}{2}\frac{c_{\mu }}{\sigma | |||
_{t}}\frac{k}{\varepsilon }\left(\overline{{u_{j}u_{k}}}\frac{\partial | |||
\overline{{\rho }}}{\partial x_{k}}\right)} | |||
</math> | |||
|width="40" align="right"|<math>\left(27\right)</math> | |||
|} | |||
which Van Maele & Merci [[UFR_1-06_References|[2]]] expressed as follows: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
G=-{\frac{3}{2}\frac{\mu _{t}}{\sigma _{t}{\overline{\rho | |||
}}^{2}k}\left(\overline{{u_{j}u_{k}}}\frac{\partial \overline{{\rho | |||
}}}{\partial x_{k}}\right)}\rho _{\infty }g_{j} | |||
</math> | |||
|width="40" align="right"|<math>\left(28\right)</math> | |||
|} | |||
again with ''σ<sub>t</sub>'' = 0.85. As | |||
previously, Van Maele & Merci [[UFR_1-06_References|[2]]] included a | |||
pressure-gradient term in the above equation related to pressure-work | |||
but this can effectively be ignored in the present application. The | |||
advantage of the GGDH approach is that transverse density gradients | |||
affect the production term. | |||
In their plume simulations, | |||
Van Maele & Merci [[UFR_1-06_References|[2]]] used a | |||
slightly modified form of the above relation. They replaced the normal | |||
stress in the streamwise (vertical) direction, | |||
<math>\overline{{\mathit{ww}}}</math>, with the turbulent kinetic energy, ''k''. | |||
They justified this on the basis that the ''k – ε'' | |||
model gives poor predictions of normal | |||
stresses in plumes. Experimental measurements indicate that the | |||
streamwise normal stress is approximately twice the magnitude of the | |||
transverse components (i.e. <math>\overline{{\mathit{ww}}}\approx | |||
2\overline{{\mathit{uu}}}</math>) whereas the ''k – ε'' | |||
model predicts them to be roughly equal. Since | |||
''k'' can be approximated from <math>k\approx | |||
\frac{1}{2}\left(\overline{{\mathit{ww}}}+2\overline{{\mathit{uu}}}\right)</math> | |||
and the ''k – ε'' model predicts | |||
<math>\overline{{\mathit{ww}}}\approx \overline{{\mathit{uu}}}</math>, they | |||
suggest that it is more appropriate to use ''k'' rather than | |||
<math>\overline{{\mathit{ww}}}</math>, to artificially increase the stress to a | |||
more realistic value. This ad-hoc correction may not be appropriate in | |||
more complex flows where the gravitational vector is not aligned to the | |||
Cartesian axes. | |||
A simplification frequently made to the buoyancy treatments described | |||
above is to assume that the mean density is equal to the reference | |||
density, <math>\overline{{\rho }}\approx \rho _{\infty }</math>, an approach | |||
known as the Boussinesq approximation | |||
<ref> | |||
A more complete definition of the Boussinesq approximation is that the variation in fluid properties | |||
due to changes in temperature or pressure are assumed to be zero (i.e. not only constant density, | |||
but constant modular viscosity, thermal conductivity etc.). In Van Maele & Merci's | |||
work, although they stated that they tested the Boussinesq approximation, in fact their simplifications | |||
only consisted of assuming <math>\overline{\rho}\approx\rho_\infty</math> in the production term, ''G '' | |||
(Van Maele, private communication, 2007). It is unclear whether variation of the density in the other | |||
remaining terms, or variation of the viscosity, thermal conductivity etc., affected the results significantly. | |||
</ref>. | |||
For the SGDH written in | |||
terms of temperature gradients, this gives: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
G=-{\beta \frac{\mu _{t}}{\sigma _{t}}\frac{\partial T}{\partial | |||
x_{j}}g_{j}} | |||
</math> | |||
|width="40" align="right"|<math>\left(29\right)</math> | |||
|} | |||
Van Maele & Merci [[UFR_1-06_References|[2]]] examined the effect of this | |||
simplification on the prediction of the | |||
George ''et al.'' [[UFR_1-06_References|[3]]] buoyant plume experiments. | |||
==== Buoyancy Source Term in the ''ε–''Equation, ''S<sub>εB</sub>'' ==== | |||
The buoyancy source term, ''S<sub>εB</sub>'', in the ''ε–''equation is given by: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
S_{\varepsilon B}=c_{\varepsilon 1}\left(1-c_{\varepsilon | |||
3}\right)\frac{\varepsilon }{k}G | |||
</math> | |||
|width="40" align="right"|<math>\left(30\right)</math> | |||
|} | |||
Unlike other model constants in the ''k – ε'' | |||
model, there is still some controversy over the best value or | |||
formula for ''c<sub>ε3</sub>''. | |||
Different approaches have been proposed by different researchers, | |||
partly depending on whether the flows are horizontal or vertical and | |||
whether there is stable or unstable stratification. For a review of the | |||
performance of various models, see Rodi [[UFR_1-06_References|[71]]], | |||
Markatos ''et al.'' [[UFR_1-06_References|[72]]] | |||
or Worthy ''et al.'' [[UFR_1-06_References|[73]]]. | |||
In their paper, Van Maele & Merci [[UFR_1-06_References|[2]]] | |||
provided a summary of the | |||
values proposed previously in 20 published papers and, based on | |||
analysis of these studies, used a constant value for ''c<sub>ε3</sub>'' of 0.8. | |||
==== Favre-Averaging ==== | |||
Throughout their paper, Van Maele & Merci [[UFR_1-06_References|[2]]] refer to | |||
Favre-averaged mean velocity, enthalpy and temperature | |||
(<math>\widetilde{U}</math>, <math>\widetilde{h}</math> and <math>\widetilde{T}</math>) instead of the | |||
perhaps more familiar Reynolds-averaged values, (''U'', ''H'' and ''T''). | |||
The Favre average of a variable | |||
<math>\left.\phi\right.</math>, denoted <math>\widetilde{\phi }</math> is calculated from: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
\tilde {{\phi }}=\frac{\overline{{\rho \phi }}}{\overline{\rho | |||
}} | |||
</math> | |||
|width="40" align="right"|<math>\left(31\right)</math> | |||
|} | |||
where overbars represent long time or ensemble averages in the | |||
traditional Reynolds-averaged sense. The turbulent stresses appearing | |||
in the Favre-averaged Navier-Stokes equations are: | |||
<!-- Had to replace \widetilde by \tilde in eqn (32) - (34) until I resolve issue --> | |||
{| | |||
|- | |||
|width="750" align="center"|<math> - \widetilde{{u_{i}u_{j}}}=2\nu_{t} \tilde{{S_{\mathit{ij}}}}-\frac{2}{3}k\delta_{\mathit{ij}}</math> | |||
|width="40" align="right"|<math>\left(32\right)</math> | |||
|} | |||
This is the same as the usual Reynolds-averaged stress except that the | |||
strain-rate tensor, <math>\widetilde{{S_{\mathit{ij}}}}</math> is now | |||
Favre-averaged. The Favre-averaged strain-rate is calculated from: | |||
{| | |||
|- | |||
|width="750" align="center"|<math> | |||
\widetilde{{S_{\mathit{ij}}}}=\frac{1}{2}\left(\frac{\partial | |||
\widetilde{{U}}_{i}}{\partial x_{j}}+\frac{\partial | |||
\widetilde{{U}}_{j}}{\partial x_{i}}\right) | |||
</math> | |||
|width="40" align="right"|<math>\left(33\right)</math> | |||
|} | |||
where the Favre-averaged mean velocity, <math>\widetilde{{U_{i}}}</math>, is the | |||
parameter solved for in the momentum equations. | |||
Even though they solved Favre-averaged transport equations, Van Maele & Merci | |||
modelled the buoyancy production term, ''G'', using | |||
the Reynolds-averaged stress, not the Favre-averaged stress. Formally, | |||
one can expand the Reynolds-averaged stress as follows [[UFR_1-06_References|[74]]]: | |||
{| | |||
|- | |||
|width="750" align="center"|<math>{\overline{{u_{i}u_{j}}}= \widetilde{{u_{i}u_{j}}}-\frac{\overline{{\rho u_{i}u_{j}}}}{\overline{\rho }}-\frac{\overline{{\rho u_{i}}}\ \overline{{\rho u_{j}}}}{{\overline{{\rho}}}^{2}}}</math> | |||
|width="40" align="right"|<math>\left(34\right)</math> | |||
|} | |||
However, it is often assumed that the last two terms in this expansion | |||
can be neglected, in which case <math>\overline{{u_{i}u_{j}}}\approx \widetilde{{u_{i}u_{j}}}</math>. | |||
This was the approach adopted by | |||
Van Maele & Merci [[UFR_1-06_References|[2]]]. | |||
If one replaces all the Favre-averaged quantities in | |||
Van Maele & Merci’s equations with Reynolds-averaged quantities | |||
(i.e. substitute ~ symbols with ‾ ), their equations appear | |||
identical to the incompressible Reynolds-averaged Navier-Stokes | |||
equations. In terms of the actual coding of the model, the mean | |||
quantities in the flow equations can therefore be interpreted as either | |||
Reynolds-averaged or Favre-averaged variables. For this reason, the | |||
transport equations have been written in this UFR without the | |||
Favre-averaged (~) symbols. | |||
In the majority of the other papers reviewed for this UFR, the transport | |||
equations are stated as being Reynolds averaged. | |||
Hossain & Rodi [[UFR_1-06_References|[8]]] | |||
noted that correlations between fluctuating velocities | |||
and fluctuating density are important in combustion applications but | |||
are small in comparison to correlations between velocity fluctuations | |||
in simple buoyant plumes. | |||
Brescianini & Delichatsios [[UFR_1-06_References|[42]]] also | |||
commented that “depending on certain assumptions, the | |||
mean quantities ... can be interpreted as either time-averaged or | |||
Favre-averaged variables. The differences between these two types of | |||
averages is small when compared to the experimental uncertainties for | |||
the plumes examined in this study, and as a result, no large | |||
distinction is made between the two forms”. In the | |||
experiments of O'Hern ''et al.'' [[UFR_1-06_References|[75]]], | |||
reviewed for the companion UFR on unsteady plumes, simultaneous | |||
measurements were made of velocities and mass fraction. This enabled | |||
both Favre-averaged and \ Reynolds-averaged quantities to be derived. | |||
O'Hern ''et al.'' [[UFR_1-06_References|[75]]] found that the | |||
difference between Favre- and Reynolds-averaged velocities and | |||
second-order turbulent statistics was less than the uncertainty in the | |||
data throughout the flow field. Further details on Favre-averaging can | |||
be found in Chassaing ''et al.'' [[UFR_1-06_References|[74]]]. | |||
=== Boundary Conditions === | |||
Van Maele & Merci [[UFR_1-06_References|[2]]] | |||
modelled the plume source using a | |||
diameter, ''D<sub>0</sub>'' , of 6.35 cm, an inlet | |||
temperature of 573 K, an axial velocity of 0.67 m/s, and a turbulence | |||
intensity of 0.5\% . To define a value for ''ε'' at | |||
the inlet, they assumed a turbulence length scale | |||
<math>{(k^{3/2}/\varepsilon )}</math> of 4.5 mm, equivalent to | |||
''D<sub>0</sub>''/15. Around the source, on the radial | |||
plane, a wall boundary was specified. The domain was axisymmetric, | |||
extending 1 metre in the radial direction and 3 metres in the axial | |||
direction. On the side and top boundaries, static pressure conditions | |||
were specified to allow the flow into or out of the computational | |||
domain. Ambient air entrained through the open boundaries was assumed | |||
to have a temperature of 302 K, a turbulent kinetic energy of | |||
10<sup>-6</sup>m<sup>2</sup>/s<sup>2</sup> and | |||
a dissipation rate of 10<sup>-9</sup>m<sup>2</sup>/s<sup>3</sup>. | |||
The atmospheric pressure was | |||
taken as 101 325 Pa. No tests were undertaken to examine whether the | |||
location or conditions of the entrainment boundaries affected the flow | |||
solution. | |||
=== Grid Used === | |||
Van Maele & Merci [[UFR_1-06_References|[2]]] | |||
used a 40 × 100 node | |||
rectangular grid in the (radial × axial) directions. | |||
There were 10 equispaced cells across the source and 30 across the | |||
adjacent wall region. The grid was stretched in both the axial and | |||
radial directions from the source. | |||
To ensure that results were grid-independent, Van Maele & Merci | |||
performed simulations using 80 × 200 and 160 × 400 node grids. | |||
The predicted spreading rates and centreline values on | |||
the coarsest and finest grids differed by less than 2%. The 40 × 100 grid | |||
was therefore considered to provide adequately grid-independent results. | |||
=== Discussion === | |||
The CFD methodology employed by Van Maele & Merci appears to have been | |||
performed to a high standard. Details of the modelling and numerical | |||
techniques used in their work were recorded clearly in their paper. | |||
They did not show a plot of the grid that was used, however this may | |||
have been because it would not have been well reproduced in print and, | |||
in any case, the mesh was a simple Cartesian grid arrangement and is | |||
adequately described in the text. | |||
The sensitivity of results to the turbulence length scale at the inlet, | |||
and to the level of turbulence in the entrained air was not explored, | |||
although reasonable approximations appear to have been used for these | |||
values. Likewise, the sensitivity to the size of the domain and the | |||
entrainment boundaries was not explored. | |||
== Footnotes == | |||
<references/> | |||
Line 44: | Line 803: | ||
{{ACContribs | {{ACContribs | ||
|authors=Simon Gant | |authors=Simon Gant | ||
|organisation= | |organisation=UK Health & Safety Laboratory | ||
}} | }} | ||
© copyright ERCOFTAC 2010 | © copyright ERCOFTAC 2010 |
Latest revision as of 19:12, 11 February 2017
Axisymmetric buoyant far-field plume
Underlying Flow Regime 1-06
Test Case
Brief Description of the Study Test Case
The experiments used in this UFR are those of George et al. [3] which were conducted in 1974 at the Factory Mutual Research Corporation and were subsequently repeated by Shabbir & George [34] at the University of Buffalo.
- Heated air is discharged through a circular orifice into ambient air that is at rest.
- The plume source temperature is 300°C and the ambient air is 29°C.
- The source has diameter, D = 6.35 cm.
- The hot air is discharged at a velocity of U0 = 67 cm/s with a approximately a top-hat profile.
- Temperature and velocity fluctuations at the inlet are less than 0.1%.
- George et al. [3] present experimentally measured profiles of both mean and fluctuating components of the temperature and axial velocity in the self-similar region at x/D = 8, 12 and 16 above the source.
Test Case Experiments
The experiments used in this UFR are those of George et al. [3] which were conducted in 1974 at the Factory Mutual Research Corporation and were subsequently repeated by Shabbir & George [34] at the University of Buffalo.
The general arrangement is shown in Figure 4.
Compressed air is passed through a set of heaters and porous mesh screens before exiting through a nozzle into
the enclosure. The nozzle is stated as a 15:1 contraction in [3],
a 12:1 contraction in [sg92]
and appears to be different again in a drawing of the arrangement in [3]
(see Figure 5). It resulted in a velocity
profile through the exit which was uniform to within 2% outside the wall boundary layer. The velocity and
temperature fluctuations at the exit were measured to be very low, less than 0.1%
in [3]
and 0.5% in [34].
The temperature of the source was 300°C and the ambient environment
29°C. Both were controlled to an accuracy of within 1°C. The discharge velocity was 67 cm/s, as
calculated from the measured heat flux. These source conditions corresponded to Reynolds number, Re0 = 870, and densimetric Froude number, Fr0 = 1.23
[1]
There was no evidence of laminar flow behaviour at a position two inlet
diameters downstream from the source. The effective origin of the plume, x0, was
found to be at the same location as the exit
(see [3] for details of how this was determined).
The screen enclosure around the plume exit was 2.44 × 2.44 metres in cross-section and 2.44 metres high
(there is, presumably, an error in [3]
which suggests that the enclosure is 2.44 × 2.44 × 2.44 mm).
In the later Shabbir & George experiments, a 2 × 2 × 5 metre enclosure was used. The purpose
of the screens was to minimize the effect of cross-draughts and other disturbances affecting the flow.
Two-wire probes were used by George et al. [3]
to record velocities and temperature.
George et al. [3] reported that measurement errors, stemming from
directional ambiguity of the hot wire and its thermal inertia, were around 3% for the velocity and lower for
other mean and RMS values. The frequency response of the hot wires was estimated to be around 300 Hz compared
to the frequency of the energy-containing eddies at around 50 Hz and the Kolmogorov microscale at 1 Khz.
It was noted that measurement errors were likely to be higher on the outer edge of the plume where
the velocity fluctuations were higher.
In their review of plume experiments, Chen & Rodi [1]
noted that the data from George et al. differed significantly from earlier measurements by
Rouse et al. [64].
However, they considered it to be more reliable due to its use of more
sophisticated instrumentation.
George [40], describes an experimental program at the University of
Buffalo that was set up following publication of the original
George et al. [3] paper to investigate
possible causes of differences in experimental plume results. Possible sources of errors discussed
included:
- ambient thermal stratification
- the size of the enclosure
- the use of porous screens used to minimise disturbances from the far-field affected the plume source.
- hot wire measurement errors
The most significant concern was ambient thermal stratification.
One of the features of buoyant plumes in neutral environments is that the integral of the
buoyancy across the whole cross-section of the plume, F, should remain constant and equal to the
buoyancy added at the source, F0.
George [40] discussed how thermal stratification involving
small temperature differences of the order of 1°C across a 3 metre vertical span would be sufficient
to cause F to decrease to only 50% of the source value.
This would be likely to cause differences in measured temperature and velocity plume profiles.
In the initial experiments of George et al. [3],
the thermal stratification was not strictly controlled.
However, results from later experiments published in the PhD thesis of
Shabbir [32]
(reproduced in [34] and [40]),
which conserved buoyancy to within 10%, are in good agreement with the earlier results from
George et al. [3].
This suggests that, perhaps fortunately, ambient thermal stratification did not contaminate the
George et al. [3] results significantly.
A summary of the original results from George et al. [3]
and those reproduced later by Shabbir & George [11]
is presented in Table 3.
Also shown are the recommended values from
Chen & Rodi's review [1] and other studies.
The parameters given in Table 3 relate to the following empirical formulae for the mean vertical velocity:
and effective buoyancy acceleration:
where and are Gaussian functions:
The parameters, and are the dimensionless half-widths
of the plume, as defined by the location where the normalized buoyancy or mean velocity falls to half
its centreline value.
The RMS temperature and axial velocity fluctuations normalized by their centreline mean values are
denoted,
and , respectively.
As noted earlier, Dai et al. [10][37]
[38][39][41] disputed the
accuracy of the George et al. [3] experiments and suggested
that they had made measurements too near the source, before the plume had reached a fully-developed state.
Their arguments are disregarded by Shabbir & George [11]
[34].
CFD Methods
Van Maele & Merci: Description of CFD Work
Numerical Methods
Van Maele & Merci [2] used the finite-volume-based commercial CFD code, Fluent [2], to simulate the plume experiments of George et al. [3]. For the discretization of the convective terms in the momentum, turbulence and energy equations a second-order upwind scheme was used. Diffusion terms were discretized using second-order central differences and the SIMPLE algorithm was used for pressure-velocity coupling. The flow was treated as axisymmetric and elliptic calculations were performed used a Cartesian grid arrangement.
The low-Mach-number form of the Favre-averaged Navier-Stokes equations were used. In this weakly-compressible
approach, the density is treated as only a function of temperature and not pressure. Pressure only affects the
flow field through the pressure-gradient term in the momentum equations. The ideal gas law is used to link the
mean density, , to mean temperature, T as follows:
where p* is taken as constant and equal to the atmospheric pressure. The
low-Mach-number approximation implies that the effect of the mean kinetic energy and the work done by viscous
stresses and pressure are negligible in the energy equation.
Turbulence Modelling
Two turbulence models were used by Van Maele & Merci [2]: the standard k – ε model of Jones & Launder [65] and the realizable k – ε model of Shih et al. [66]. In the former model, the eddy viscosity is given by:
where cμ is a constant equal to 0.09 and the standard k and ε
equations are written:
where cε1 = 1.44,
cε2 = 1.92,
σk = 1.0,
σε = 1.3 and Pk is the production term due to mean shear.
The terms G and SεB are source terms related to the influence of buoyancy on
the k and ε equations. The treatment of these terms is discussed below.
The Shih et al. [66] model involves two changes to the standard
k – ε model. Firstly, cμ is made a function of strain and
vorticity invariants to ensure that the model always returns positive normal Reynolds stresses and satisfies
the Schwarz inequality for the turbulent shear stresses. The function form of cμ
is given by:
where:
and A0 is a constant equal to 4.04.
Secondly, a different ε-equation is used to resolve the problem of the round-jet/plane-jet
anomaly (see Pope [67]):
where S is the strain-rate invariant as before, cε2 = 1.9,
σk = 1.0 and σε = 1.3.
The Shih et al. [66] model was developed for high
Reynolds number turbulent flows and therefore a zonal or wall-function
approach must be used to bridge the viscous sub-layer near walls.
Compared to the standard k–ε model,
it has been shown to produce improved behaviour in a number of free
shear flows, boundary-layer flows and a backward-facing step flow [66].
One of the major weaknesses of the standard k–ε model
is that it produces too much turbulent
kinetic energy at stagnation points [68].
The Shih et al.
model should in principle suffer less from this weakness since the
functional form of cμ should
reduce the over-production of k. However, its overall
performance in stagnating flows will depend on the type of wall model
used.
Production due to Buoyancy, G
The term G in the k-equation relates to the influence of buoyancy on the turbulent kinetic energy, and is given by:
where gj is the gravitational
acceleration vector. In stably stratified flows, where the temperature
increases with height, G is negative. Conversely, in unstably
stratified flows, where temperature decreases with height, G
is positive and acts to increase k. The unknown
density-velocity correlation, , must be
modelled. The most common approximation of this term is the so-called
Boussinesq Simple Gradient Diffusion Hypothesis (SGDH):
The production due to buoyancy using SGDH is then as follows:
In their paper, Van Maele & Merci [2]
erroneously included
an additional pressure-gradient term in Equation (23) related to the pressure-work rather than
buoyancy (see Wilcox [69]). Since the term is negligible in
incompressible flows, such as the buoyant plumes considered here, it
has therefore been ignored. The ratio of the reference density to the
mean density, , appears in
Equation (23) due to the use of a non-Boussinesq approach and Favre-averaging,
which are discussed later. Van Maele & Merci [2] assumed
that σt was constant and
equal to 0.85.
Instead of writing the buoyancy production in terms of the
density-velocity correlation, , the equation
can be written in terms of the heat flux, :
and the G term is then written:
where t′ is the temperature fluctuation, T is the mean
temperature and β is the volumetric expansion
coefficient, . Other equivalent expressions can also
be formulated using the ideal gas law and the assumption that density
is only a function of temperature, not pressure (the low-Mach-number
approximation). The conversion from mean density to temperature
gradients is then as follows:
The SGDH model predicts zero density-velocity correlation or heat flux
components ( or )
in situations were the density or temperature gradients are zero in
that direction. However, as Ince & Launder [70] noted, in a
simple shear flow in which there are only cross-stream temperature
gradients, the heat flux in the streamwise direction actually exceeds
that in the cross-stream direction. This shortcoming of the SGDH model
was confirmed by the analysis of Shabbir & Taulbee [33], who
showed that the model significantly underestimates the magnitude of the
heat flux in vertical buoyant plumes. The underprediction of
or by the SGDH model
leads to an overly-small production term, G, and hence a
turbulent kinetic energy, k, which is too small, producing too
little mixing in the modelled plume. The study by
Yan & Holmstedt [53]
provides a clear example of how the k – ε
model with SGDH produces buoyant plumes which
are too narrow and with overly high temperatures and velocities in the
core of the flow.
Van Maele & Merci [2]
examined a different model for
G based on the the Generalized-Gradient Diffusion Hypothesis
(GGDH) of Daly & Harlow [52].
This was first used in the
context of practical CFD calculations with the k – ε
model by Ince & Launder [70], and is
written as follows:
which Van Maele & Merci [2] expressed as follows:
again with σt = 0.85. As
previously, Van Maele & Merci [2] included a
pressure-gradient term in the above equation related to pressure-work
but this can effectively be ignored in the present application. The
advantage of the GGDH approach is that transverse density gradients
affect the production term.
In their plume simulations,
Van Maele & Merci [2] used a
slightly modified form of the above relation. They replaced the normal
stress in the streamwise (vertical) direction,
, with the turbulent kinetic energy, k.
They justified this on the basis that the k – ε
model gives poor predictions of normal
stresses in plumes. Experimental measurements indicate that the
streamwise normal stress is approximately twice the magnitude of the
transverse components (i.e. ) whereas the k – ε
model predicts them to be roughly equal. Since
k can be approximated from
and the k – ε model predicts
, they
suggest that it is more appropriate to use k rather than
, to artificially increase the stress to a
more realistic value. This ad-hoc correction may not be appropriate in
more complex flows where the gravitational vector is not aligned to the
Cartesian axes.
A simplification frequently made to the buoyancy treatments described
above is to assume that the mean density is equal to the reference
density, , an approach
known as the Boussinesq approximation
[3].
For the SGDH written in
terms of temperature gradients, this gives:
Van Maele & Merci [2] examined the effect of this
simplification on the prediction of the
George et al. [3] buoyant plume experiments.
Buoyancy Source Term in the ε–Equation, SεB
The buoyancy source term, SεB, in the ε–equation is given by:
Unlike other model constants in the k – ε
model, there is still some controversy over the best value or
formula for cε3.
Different approaches have been proposed by different researchers,
partly depending on whether the flows are horizontal or vertical and
whether there is stable or unstable stratification. For a review of the
performance of various models, see Rodi [71],
Markatos et al. [72]
or Worthy et al. [73].
In their paper, Van Maele & Merci [2]
provided a summary of the
values proposed previously in 20 published papers and, based on
analysis of these studies, used a constant value for cε3 of 0.8.
Favre-Averaging
Throughout their paper, Van Maele & Merci [2] refer to Favre-averaged mean velocity, enthalpy and temperature (, and ) instead of the perhaps more familiar Reynolds-averaged values, (U, H and T). The Favre average of a variable , denoted is calculated from:
where overbars represent long time or ensemble averages in the
traditional Reynolds-averaged sense. The turbulent stresses appearing
in the Favre-averaged Navier-Stokes equations are:
This is the same as the usual Reynolds-averaged stress except that the
strain-rate tensor, is now
Favre-averaged. The Favre-averaged strain-rate is calculated from:
where the Favre-averaged mean velocity, , is the
parameter solved for in the momentum equations.
Even though they solved Favre-averaged transport equations, Van Maele & Merci
modelled the buoyancy production term, G, using
the Reynolds-averaged stress, not the Favre-averaged stress. Formally,
one can expand the Reynolds-averaged stress as follows [74]:
However, it is often assumed that the last two terms in this expansion
can be neglected, in which case .
This was the approach adopted by
Van Maele & Merci [2].
If one replaces all the Favre-averaged quantities in
Van Maele & Merci’s equations with Reynolds-averaged quantities
(i.e. substitute ~ symbols with ‾ ), their equations appear
identical to the incompressible Reynolds-averaged Navier-Stokes
equations. In terms of the actual coding of the model, the mean
quantities in the flow equations can therefore be interpreted as either
Reynolds-averaged or Favre-averaged variables. For this reason, the
transport equations have been written in this UFR without the
Favre-averaged (~) symbols.
In the majority of the other papers reviewed for this UFR, the transport
equations are stated as being Reynolds averaged.
Hossain & Rodi [8]
noted that correlations between fluctuating velocities
and fluctuating density are important in combustion applications but
are small in comparison to correlations between velocity fluctuations
in simple buoyant plumes.
Brescianini & Delichatsios [42] also
commented that “depending on certain assumptions, the
mean quantities ... can be interpreted as either time-averaged or
Favre-averaged variables. The differences between these two types of
averages is small when compared to the experimental uncertainties for
the plumes examined in this study, and as a result, no large
distinction is made between the two forms”. In the
experiments of O'Hern et al. [75],
reviewed for the companion UFR on unsteady plumes, simultaneous
measurements were made of velocities and mass fraction. This enabled
both Favre-averaged and \ Reynolds-averaged quantities to be derived.
O'Hern et al. [75] found that the
difference between Favre- and Reynolds-averaged velocities and
second-order turbulent statistics was less than the uncertainty in the
data throughout the flow field. Further details on Favre-averaging can
be found in Chassaing et al. [74].
Boundary Conditions
Van Maele & Merci [2] modelled the plume source using a diameter, D0 , of 6.35 cm, an inlet temperature of 573 K, an axial velocity of 0.67 m/s, and a turbulence intensity of 0.5\% . To define a value for ε at the inlet, they assumed a turbulence length scale of 4.5 mm, equivalent to D0/15. Around the source, on the radial plane, a wall boundary was specified. The domain was axisymmetric, extending 1 metre in the radial direction and 3 metres in the axial direction. On the side and top boundaries, static pressure conditions were specified to allow the flow into or out of the computational domain. Ambient air entrained through the open boundaries was assumed to have a temperature of 302 K, a turbulent kinetic energy of 10-6m2/s2 and a dissipation rate of 10-9m2/s3. The atmospheric pressure was taken as 101 325 Pa. No tests were undertaken to examine whether the location or conditions of the entrainment boundaries affected the flow solution.
Grid Used
Van Maele & Merci [2] used a 40 × 100 node rectangular grid in the (radial × axial) directions. There were 10 equispaced cells across the source and 30 across the adjacent wall region. The grid was stretched in both the axial and radial directions from the source.
To ensure that results were grid-independent, Van Maele & Merci performed simulations using 80 × 200 and 160 × 400 node grids. The predicted spreading rates and centreline values on the coarsest and finest grids differed by less than 2%. The 40 × 100 grid was therefore considered to provide adequately grid-independent results.
Discussion
The CFD methodology employed by Van Maele & Merci appears to have been performed to a high standard. Details of the modelling and numerical techniques used in their work were recorded clearly in their paper.
They did not show a plot of the grid that was used, however this may have been because it would not have been well reproduced in print and, in any case, the mesh was a simple Cartesian grid arrangement and is adequately described in the text.
The sensitivity of results to the turbulence length scale at the inlet, and to the level of turbulence in the entrained air was not explored, although reasonable approximations appear to have been used for these values. Likewise, the sensitivity to the size of the domain and the entrainment boundaries was not explored.
Footnotes
- ↑ The densimetric Froude number is calculated here from the source and ambient temperatures, the exit velocity and source diameter given by George et al. [3], using Equation (1). However George et al. [3] stated that the densimetric Froude number was 1.4. It is unclear how they determined this value. Using the approach taken by Chen & Rodi [1] in which the source density instead of the ambient density is used to make the density difference dimensionless, and Froude number is defined using the square of the expression given in Equation (1), this gives a Froude number of 0.80.
- ↑ http://www.fluent.com
- ↑ A more complete definition of the Boussinesq approximation is that the variation in fluid properties due to changes in temperature or pressure are assumed to be zero (i.e. not only constant density, but constant modular viscosity, thermal conductivity etc.). In Van Maele & Merci's work, although they stated that they tested the Boussinesq approximation, in fact their simplifications only consisted of assuming in the production term, G (Van Maele, private communication, 2007). It is unclear whether variation of the density in the other remaining terms, or variation of the viscosity, thermal conductivity etc., affected the results significantly.
Contributed by: Simon Gant — UK Health & Safety Laboratory
© copyright ERCOFTAC 2010