Unsteady Near-Field Plumes

Underlying Flow Regime 1-07

Best Practice Advice

Best Practice Advice for the UFR

Key Physics

The key physics of this UFR is the transient, unsteady behaviour in the near-field of a turbulent buoyant helium-air plume. The flow features two key instabilities. Firstly, the Rayleigh-Taylor instability related to the presence of dense fluid above less-dense fluid, which gives rise to fingers or spikes of dense fluid separated by rising bubbles of lighter fluid. Secondly, the Kelvin-Helmholtz instability related to the shear-layer interface between the rising plume and the ambient fluid, which produces roll-up vortex sheets on the boundary between the two layers of fluid travelling at different velocities. The flow is very challenging to predict using CFD, due to the sharp density gradients at the plume exit which produce flow conditions where small scales of turbulent motion feed into the larger scales.

Numerical Modelling

For LES, the flow cannot be treated as two-dimensional or axisymmetric. Full three-dimensional time-dependent simulations must be performed.

For simulation of the selected UFR test case, open boundaries should be used on all sides of the flow domain except for the floor. Constant pressure boundaries may be used, although if a fully-compressible code is used, care will need to be taken to ensure that the boundaries are non-reflective.

For simulation of the selected UFR test case, the domain should extend at least 4 metres radially and vertically to minimize any effects of the open boundaries on the development of the plume. Ideally, tests should be performed to ensure that the location of the open boundaries has no significant effect on the results.

The finest mesh should be used given the available computing resources. The results discussed above suggest that a mesh of around 4 million nodes should give good agreement with the experiments in terms of mean flow quantities, but may still be insufficient for good predictions of fluctuations or RMS values. Tieszen^[1] noted that at least 75 cells across the base diameter of the plume are necessary to avoid significant differences in the vertical centreline velocity compared to the measured values. Ideally, a grid-dependence study should be undertaken to investigate the magnitude of these effects.

Physical Modelling

Either the fully-compressible or the low-Mach-number form of the Navier-Stokes equations can be used. The fully-compressible N-S equations require more careful treatment to avoid acoustic waves reflecting back into the domain from open boundaries. Furthermore, they will require a very short time-step, based on the speed of sound instead of the local flow speed, unless special treatments are used. For details of a fully-compressible N-S treatment, see DesJardin et al. [1].

The baroclinic torque is non-zero and therefore should not be neglected.

The Boussinesq approximation, where flow properties are assumed not to vary as a function of temperature or composition, and where buoyancy is only included as an additional body-force term in the momentum equations, should not be used. The Boussinesq approximation is only appropriate for modelling small density differences, equivalent to a temperature difference less than around 15°C in air [80].

If an LES approach is used, the effect of the unresolved small-scale turbulence on the resolved motion can either be accounted for by an explicit sub-grid-scale model, such as the dynamic Smagorinsky model, or by numerical damping in an implicit LES (a “no-model” approach). If an explicit approach is taken, central differencing should ideally be used for convection in the momentum equation but bounded upwind-biased schemes will probably be needed for the scalar equations to prevent unbounded under/overshoots. DesJardin et al. [1] obtained slightly better results with the implicit approach but this is likely to depend on grid resolution. If a coarse grid is used, an explicit LES model should probably be avoided. Both implicit and explicit approaches should ideally be tested to examine the sensitivity of results to the turbulence treatment. Recent work by Blanquart & Pitsch [40] has shown very good predictions for both mean momentum and concentration using the Lagrangian dynamic SGS model of Meneveau et al. [51] for turbulent diffusion terms in both the momentum and helium mass-fraction transport equations.

It is difficult to provide definitive guidance on use of RANS models, since to date it appears that there have only been two relevant studies for this flow, and they produced somewhat contradictory results. Chung & Devaud [39] found in the helium plume experiments of O‘Hern et al. [4] that steady flow behaviour was obtained using k – ε models with SGDH or GGDH and values of the model constant C_ε3 varying between 0 and 1. They also found that assuming the flow to be axisymmetric or using a fully three-dimensional approach gave practically identical results and an axisymmetric approach with 22,882 cells gave a grid-independent solution. In contrast, Nicolette et al. [38] found a “standard” k – ε model produced unsteady flow behaviour using all but the very coarsest of meshes, which had only 56,000 cells for the three-dimensional geometry. For meshes containing 500,000 to 2 million cells, the predicted flow behaviour was unsteady with the finer meshes resolving an increasing proportion of the unsteady flow structures. These differences in resolving steady or unsteady flow behaviour could in part be due to the former study using a steady solution method whilst the latter used a transient time-stepping approach. Nevertheless, Chung & Devaud [39] reported that the residuals in their steady simulations could be reduced to low levels (maximum residuals of 10^-5), whilst usually in flows where there is a tendency towards transient behaviour it is difficult to obtain such good convergence. Putting these differences to one side, the Chung & Devaud [39] study showed that good predictions of the steady flow behaviour in the near-field of buoyant plumes could be achieved provided that special care was taken over the choice of the model constant, C_ε3. Different optimum values of C_ε3 were found when using either SGDH or GGDH, and varying the value of C_ε3 produced very marked changes in the flow predictions. It is recommended therefore if studying the near-field flow behaviour of plumes using similar models to examine the sensitivity of the results to this parameter. The Nicolette et al. [38] study mainly focused on testing their newly-developed Buoyant Vorticity Generation (BVG) extension to the k – ε model. The model showed promising results in comparison to the helium plume experiments from NIST [31][32] and the Sandia FLAME facility with a low helium inlet velocity of 0.13 m/s. However, the model gave less encouraging predictions when compared to the O‘Hern et al. [4] experiments where the inlet velocity was higher, due to the delayed predicted onset of laminar to turbulent transition. There were also additional complications with mesh-dependent transient behaviour, as mentioned above. Overall, the two studies indicate that further work is needed before definitive best-practice advise can be provided on the use of RANS models in the near-field of buoyant plumes.

The above comments address the use of LES and RANS models in the near-field region, up to around five diameters downstream from the source. LES predictions of the fully-developed plume in the far-field are presented by Zhou et al. [45][46] and Pham et al. [50]. For information on best-practice modelling of the steady far-field behaviour of buoyant plumes, see the companion UFR.

Application Uncertainties

The mean velocity of the helium gas mixture flowing through the 1-metre-diameter inlet was different according to whether the value was Favre- or Reynolds-averaged. In the former case it was 0.339 m/s and in the latter case it was 0.325 m/s.

The mean inlet helium and air temperatures in the experiments were slightly different, 11 °C and 13 °C, respectively. However the change in density associated with the 2 °C temperature difference is very small in comparison to that associated with the difference in the molecular weight of the helium mixture and air. The flow can therefore be treated as isothermal, at approximately 12 °C.

The plume experiments involved the release of a helium, acetone and oxygen gas mixture with a molecular weight of 5.45 g/mol, rather than a pure helium with a molecular weight of 4.0 g/mol.

Turbulence levels in the flow issuing from the 1-metre-diameter inlet and in the surrounding entrained air flow were not directly measured in the experiments. However, flow visualization strongly suggested the conditions in the helium inlet were laminar, and there was only weak turbulence in the entrained air. O‘Hern et al. [4] considered that transition to turbulence in the plume was not driven by residual vorticity from the boundary layers in the inlet flows but instead came from gravitational and baroclinic torque in the plume. This was confirmed by DesJardin et al. [1], who found that superimposing turbulent fluctuations on the inlet velocity in their simulations did not affect the resulting predicted flow behaviour.

The experiments were carefully designed to mimic conditions where the plume is unconfined and surrounded at the source by an infinite flat ground plane [4]. Measurements were undertaken to ensure the uniformity of the flow through the annular air inlet into the chamber and the facility was designed following extensive CFD modelling studies to eliminate any disturbances to the plume. Without models reproducing exactly the same geometry as used in the experiments this remains a potential source of uncertainty. Whilst the majority of the published studies have chosen to simulate an unconfined plume on a flat plane, Chung & Devaud [39] and Blanquart & Pitsch [40] both modelled the complete geometry.

There are some discrepancies in the boundary conditions used by DesJardin et al. [1]. The plume experiments involved the release of a helium, acetone and oxygen gas mixture with a molecular weight of 5.45 g/mol, whereas DesJardin et al. instead modelled a pure helium plume with a molecular weight of 4.0 g/mol. They also used an inlet velocity of 0.351 m/s whereas the inlet Favre-averaged velocity measured in the experiments was 0.339 m/s [4]. Furthermore, DesJardin et al. [1] used a small co-flow velocity of 0.01 m/s around the plume source whereas the velocity in the experiments was zero, due to the presence of the 0.51 m wide steel plate. The simulations presented by DesJardin et al. [1] were carried out at the same time as the experiments based on an initial design, and there was insufficient time to repeat the CFD simulations once the experimental conditions were fully established^[2]. Due to these differences, there remains some uncertainty in the model predictions.

It is considerably more challenging to establish grid-independence for LES than it is for traditional RANS simulations. Grid independence can only ever be achieved for statistical quantities, such as mean velocities and Reynolds stresses, and not for the instantaneous field, which will change with grid resolution [83]. Unlike RANS, it is usually not possible to separate the effects of the turbulence model and the discretization errors. Most commercial CFD codes use implicit filtering, whereby discretization of the governing equations is assumed to act as a filtering process. The filter width, Δ, in the SGS model is then usually taken as the cube-root of the cell volume (or twice the cube-root in some codes). The turbulence model is therefore intrinsically linked to the grid resolution: refining the grid affects both the discretization errors and the turbulence model itself^[3]. In theory, a truly grid-independent solution could be obtained by refining the grid to the point where the solution becomes effectively a DNS, or by using an alternative approach based on explicit filtering, which separates the discretization errors from the SGS modelling effects (see for example Gullbrand & Chow [84]). However, these approaches are costly and are rarely used in practice. Consequently, most LES solutions use implicit filtering and the solutions often involve a complex mixture of numerical and modelling errors. The interaction between these errors has been studied by Geurts et al. [85][86][81], who found that in some flows their interaction can actually lead to results becoming worse as the grid is refined. They suggested that it may be more effective to run a number of simulations with relatively coarse meshes to help understand the grid-dependence issues and optimize the SGS model coefficients, rather than run just one or two simulations using the finest mesh possible (see also Klein et al. [87][88]). Their studies have so far been limited to relatively simple flows, such as homogeneous isotropic turbulence. In the plume studies presented here, there appears to be a consistent trend for results to improve as the grid is made progressively finer. This appears to be due to key physical processes controlling the development of turbulent structures in plumes being present only in the very small scales. The need for a very fine grid to resolve the Rayleigh-Taylor instabilities in plumes is discussed in some detail in DesJardin et al. [1] and Tieszen et al. [2], }and is summarized above. This matter introduces uncertainties for industrial LES practitioners seeking to reproduce the results from the above plume studies using their own commercial or in-house CFD codes that rely on different numerical schemes and meshing practices. Although various measures have been proposed to indicate whether the LES grid is sufficiently fine, such as the LES Index of Quality of Celik et al. [89][90], these are still the subject of ongoing research and have yet to be proven widely applicable in practice [91].

The usefulness of the QNET UFRs is to provide best-practice guidance for CFD simulations of various different generic flows. It is important that the guidance is general rather than being specific to just one particular set of experiments. With this in mind for the present UFR, it is perhaps worthwhile reflecting on the differences between non-reacting helium plumes and fire plumes to help identify possible uncertainties that may arise in extrapolating the advice from this UFR to the modelling of fire plumes. In helium plumes the strongest gradients are at the base of the plume and the mean driving force decreases with height. In fire plumes, however, due to combustion being dependent upon fuel and air mixing, the driving force first increases with height as the combustion rate increases and only further downstream starts decreasing. Most fuel vapours are denser than air and produce a fuel core or vapour dome around the source [42]. They may not therefore be subject to the same Rayleigh-Taylor instabilities that occur in helium plumes. This suggests that to model fires, it may not be necessary to model the bubble and spike structures that are important in helium plumes and grid dependence may be less of an issue. Further work is necessary to confirm this.

Recommendations for further work

Whilst the studies discussed in this UFR provide some guidance on appropriate grid resolution for the O‘Hern et al. [4] helium plume, it is difficult to formulate generic rules from this for other plumes involving for example different density ratios, Richardson and Reynolds numbers. Such guidance would be extremely welcome in the industrial CFD community where studies are often made of flows where there is no experimental data available and where CFD models are used as truly predictive tools.

The discussion of grid-dependence issues in the works of Tieszen et al. [2][42] suggests two possible criteria that could be investigated in future work to help formulate guidance on grid-resolution for LES simulations of plumes. In [2], the characteristic length scale of the Rayleigh-Taylor instability is estimated as:

\lambda =\left({\frac {D^{2}}{Ag}}\right)^{1/3}

\left(20\right)

where D is the diffusion coefficient between plume and air, A is the Atwood number (the density difference divided by the sum of the densities) and g is the gravitational acceleration. A possible rule for grid refinement could be based on the local grid size as a function of λ in much the same way that Baggett et al. [92] suggest that the grid size should be a certain fraction of the integral length scale (see also [93]).

Secondly, Tieszen et al. [42] noted that simulations with fine meshes resolved the air-spike structures characteristic of the Rayleigh-Taylor instability but coarse-grid simulations did not. A second measure could therefore be based on the appearance (or not) of air spike structures in the simulations. Clearly it would be necessary beforehand to have some physical criteria for the actual appearance of air spikes in real plumes, perhaps based on the Reynolds, Froude or Richardson number, or the width of the plume source (see discussion here).

It would be useful to investigate further the performance of steady and unsteady RANS models in predicting the near-field behaviour of turbulent plumes, in particular algebraic and differential stress models. The only relevant studies to date appears to be the work the works of Nicolette et al. [38] and Chung & Devaud [39], who tested variants of the k – ε model.

Finally, further well-instrumented and carefully controlled experiments would be welcome for industrially-relevant turbulent plumes. This would preferably be combined with a cohort of CFD studies where different CFD models and numerical methods are compared to experimental data.

Footnotes

↑ S. Tieszen, Private Communication, March 2010.
↑ DesJardin, Personal Communication, 2010.
↑ This should not be confused with the MILES approach [82], where the numerical scheme is devised to account for the effects of the unresolved turbulence and there is no explicit SGS model, although the same considerations apply.