CFD Simulations AC3-11
Downward flow in a heated annulus
Application Challenge 3-11 © copyright ERCOFTAC 2004
Overview of CFD Simulations
Introduction
Earlier three-dimensional calculations of similar buoyancy-affected annulus flows have been reported by Guy et al (1999). One of the main objectives of this study was to establish whether buoyancy might cause the flow to become unstable, and consequently three-dimensional. The circumferential grid employed thus covered the entire 360o of the annulus cross-section. The linear k-ε scheme was employed, with three different near-wall treatments, namely standard wall functions; a low-Reynolds-number 1-equation near-wall model, and the low-Reynolds-number Launder-Sharma k-ε scheme.
At Reynolds numbers greater than 3000, the use of symmetric boundary conditions was found to always result in a symmetric solution, confirming the experimental findings that the flow does not exhibit large-scale three-dimensional asymmetries. The only exception to this behaviour was for one case, at a very low Reynolds number of 2000 (for which 2-dimensional calculations suggested the flow was laminar before the heated section). In this case it proved impossible to obtain stable 3-dimensional solutions when using the low-Reynolds-number k-ε scheme. Coupled with experimental evidence, it was argued that this might indicate flow instability in this case.
Three different wall treatments were tested by Guy et al. Their conclusion was that the zonal approach in which a 1-equation model was employed near the wall, or the low-Reynolds-number k-ε scheme of Launder-Sharma produced reasonable agreement with the available experimental data, whilst the wall function approach resulted in an underprediction of the Nusselt number of around 25%. It should be noted that in their low-Reynolds-number k-ε calculations, they reported that the inclusion of a lengthscale correction in the ε equation, based on that proposed by Yap (1987) or similar, was also beneficial in terms of predictive accuracy.
A further important finding of the above study was the effect that the inclusion of temperature-dependent molecular viscosity and Prandtl number had on the resulting levels of heat transfer. For the one case reported, inclusion of property variations lead to an increase in Nusselt number of almost 30%, Figure 2. The conclusion was that for water flows at the heating rates studied the variation of fluid properties with temperature should not be neglected.
In view of the above findings, the more recent calculations reported here have treated the problem as a two-dimensional axisymmetric flow, which obviously requires significantly less computing resources.
Four of the experimental cases have been examined, covering values of the buoyancy parameter Bo ranging from 0.22 to 2.89, and details of these are summarized in Table 4. As indicated in Section 1, particular attention has been given to the modelling of the near-wall, viscosity-affected, sublayer. In principle, as noted above, the most reliable method of treating this region is with a low-Reynolds-number model, with grids fine enough to fully resolve the rapid variation of mean and turbulence quantities across the layer. However, this is not always feasible in complex industrial flows, because of the associated computational cost. Instead, wall functions are often used which are based on the assumption that across this sublayer the turbulence is in a state of simple shear where the production and dissipation rates of turbulence energy are equal. However, these idealised conditions are not generally found in the flows arising in industry, where reliable CFD predictions are sought. There, the flow physics will be highly complex, whether due to a complicated 3-dimensional strain field, or to the action of force-fields such as buoyant effects considered here. Hence the standard wall functions have a rather limited width of applicability.
Figure 2: Effect of property variation on the predicted Nusselt number, using wall functions.
Over the past three years the Turbulence Mechanics Group at UMIST has been developing a new wall function strategy that is designed to remove some of the weaknesses of the standard approach. It has primarily been developed and tuned by reference to buoyancy aided and opposed flows in vertical pipes. Subsequently, it has also proved to be successful in predicting an opposed wall jet, and the natural convection vertical boundary layer. This approach has therefore been further tested in the present cases.
| CFD Tests | Re | Bo | Expt. Test No. (Table 1) |
|---|---|---|---|
| 1 | 6000 | 0.22 | 4 |
| 2 | 6000 | 0.78 | 9 |
| 3 | 4000 | 0.83 | 2 |
| 4 | 4000 | 2.89 | 7 |
Table 4: Summary of cases studied with CFD.
Numerical and Modelling Procedure
The computations have been carried out using the TEAM computer code (Huang & Leschziner, 1983) which is a finite-volume solver employing a staggered storage arrangement on a rectangular plane or axisymmetric Cartesian grid. The SIMPLE pressure correction scheme of Patankar (1980) is employed, whilst convection is approximated via the QUICK scheme of Leonard (1979) in equations for mean quantities, but with the PLDS scheme (Patankar, 1980) in the turbulence equations. Because the streamlines are nearly parallel with the constant-radius grid lines, false diffusion effects from using PLDS will be negligible.
As noted above, it is known to be important to include property variations with temperature in these buoyancy-affected water flows, and consequently, both the molecular viscosity and Prandtl number (as well as density) were allowed to depend on temperatrure.
The turbulence modelling approaches tested are:
Model 1: Low-Reynolds-Number Eddy-Viscosity Scheme (LRN)
The low-Reynolds-number scheme of Launder & Sharma (1974) was employed, on a grid fine enough to allow the near-wall flow to be resolved accurately. Although linear EVM’s are known to be unsuitable for many complex flows, this particular scheme has been originally developed (Jones & Launder, 1972) to account for changes in the dimensionless viscous sublayer thickness. Moreover, because the flow does closely approximate a simple shear flow, the model does actually return predictions in good agreement with experimental data for the flows considered here. Hence it provides a good modelling framework within which to test and compare the different wall function treatments. Although the earlier work of Guy et al (1999) highlighted the beneficial effects of a lengthscale correction term in the ε equation, this has not been included in the present calculations. Some initial tests were carried out which suggested that in these cases (with lower heat loadings than those simulated by Guy et al) its influence was negligible. Details of the model are given in Appendix A.
Model 2: Standard Wall Function (StWF)
In this case, Model 1 is applied in the main flow region. However a considerably coarser near-wall grid is employed, with the wall-adjacent node sufficiently far from the wall for it to be located in fully turbulent fluid. The widely employed wall function form of Spalding (1967) is used to obtain the wall shear stress and temperature, which assumes the near-wall flow to be a simple shear in local equilibrium. The form of equations employed in the scheme are given in Appendix B.
Model 3: Algebraic Wall Function (AWF)
This approach uses the same grid distribution and main flow region model as Model 2, but employs the new wall functions recently developed at UMIST to cover the near-wall cell. These wall functions are based on the analytical integration of simplified forms of the near-wall momentum and temperature equations, with a suitable assumption being made for the turbulent viscosity variation across the near wall cell. This allows buoyancy forces and fluid property variations to be taken into account. The model also includes a mechanism to account for the thickening or thinning of the viscous sublayer in non-equilibrium flows. Details of the model are summarized in Appendix C.
Computational Grid and Boundary Conditions
The domain employed for the calculations is shown in Figure 3, which was chosen to match the test section dimensions of the experiment. For the low-Reynolds-number model a grid of 252 (streamwise) by 60 (radial) nodes was employed, which was sufficient (at these low bulk Reynolds numbers) to ensure that the value of y+ at the near-wall node was less than unity, and to produce grid-independent results. With the wall function approaches, the number of radial grids could be reduced to 12 whilst still obtaining grid-independent results. Because of the low bulk Reynolds numbers of these flows, a fairly large near-wall cell of around 10% of the annular gap width was used. This ensured that, when the algebraic wall function was employed, the first grid node lay in the fully turbulent flow region. However, with the standard wall function, even these large near-wall cells did, in some cases, result in the near-wall node still lying within the viscosity-affected region. It is worth noting at this point that the analytical wall function is, in fact, designed to be fairly insensitive to the size of the near-wall cell (and can even be applied when the near-wall node does not lie in the fully turbulent region). Predictions obtained with the standard wall function, on the other hand, are known to show a dependence on the size of the near-wall cell.
Figure 3: Computational domain.
Fully developed pipe inlet conditions were applied, although tests with uniform inlet conditions returned indistinguishable results – because the 23 diameters over which the flow can develop before the heated section of pipe ensured that fully developed flow conditions were achieved by the beginning of the heated section. Zero streamwise gradients were applied at the outlet plane. Regarding the thermal boundary conditions, adiabatic conditions were applied on the outer wall, together with the initial and end sections of the inner wall, whilst a constant heat flux was applied to the middle section of the inner wall, to match the experimental conditions.
CFD Results and Assessment of Calculations
As indicated in Section 1, the parameter of principal interest is the variation of Nusselt number along the annulus, and for all the cases studied, the variation of this quantity, and the corresponding wall temperature distribution, will be compared with the experimental measurements. In the figures, the predictions are labelled as LRN, StWF and AWF, corresponding to the models described above. For the low-Reynolds-number approach a set of computations were also performed where the gravitational acceleration was set to zero, allowing the effect of (principally) the buoyant term in the vertical momentum equation to be seen.
Figure 4 presents results for the case where buoyant effects are weak (Bo=0.22 at a mean pipe Reynolds number of 6000, corresponding to experimental Test 4). According to the low-Reynolds-number model the buoyant terms have the effect of increasing the Nusselt number by 30% compared with the case where buoyant forces are suppressed. In this case both wall function approaches produce results which accord well with the low-Reynolds-number computations and the experimental data.
However, a different picture begins to emerge when the buoyancy parameter is increased. Figure 5 shows results of simulating experimental Test 9 (Bo=0.78, Re=6000). At this level of heating, the low-Reynolds-number model still returns levels of Nu in close agreement with the measured values, and indicates that buoyant effects raise the level of Nusselt number above the level of purely forced convection by around 70%. The standard wall function calculations, however, only capture about 40% of the increase in Nusselt number relative to the non-buoyant case. The AWF predictions do much better, and are again in close agreement with the low-Reynolds-number predictions and experimental data.
If the Reynolds number of the flow is reduced to 4000, whilst keeping the buoyancy parameter roughly the same, then Figure 6 shows that the AWF approach still reproduces the data with good accuracy over most of the development length, whereas the standard wall function approach now returns levels below the forced convection low-Reynolds-number scheme results.
The final case computed is for the situation where Bo is increased to 2.89 (experimental Test 7), which represents a very strong buoyant effect (although still not quite enough to cause near-wall flow reversal). In this case the low-Reynolds-number results suggest that the buoyancy raises the Nusselt number by 150% relative to the non-buoyant situation (Figure 7). Even under these strong buoyant influences the AWF scheme still gives a good account of the heat transfer coefficient variation, returning levels within 10% of the measurements over most of the flow. The standard wall function approach, on the other hand, leads to predicted levels of Nusselt number in this case of less than 40% of the measured levels. This last case studied corresponds to one of the configurations for which measured velocity profiles are also available. Figure 8 shows the predicted and measured vertical velocity profiles across the annular gap for both the non-buoyant (unheated) and buoyant cases. Despite the large differences in heat transfer results seen in Figure 7, both wall function approaches appear to produce similar mean velocity profiles, which do reproduce the asymmetry of the experimental data.
One reason for the large differences in heat transfer predicted by the two wall function approaches at high values of the buoyancy parameter can be traced to the importance of correctly representing the influence of buoyancy in the near-wall region of these flows. Although the velocity profiles in Figure 8 are similar, a different picture emerges if the very near-wall velocity profiles are examined. Figure 9 shows predicted profiles of the mean velocity and temperature very close to the heated wall (covering almost the width of the first two grid cells). In addition to the results returned by both wall function approaches at the main grid nodes, the figures also show the analytical near-wall profiles which the new wall function returns. As can be seen, the near-wall temperature profile returned by the analytical wall function is in reasonable agreement with the low-Reynolds-number model results. The analytical integration of the mean momentum equation in the new wall function does explicitly include buoyancy affects arising from this temperature profile, in addition to molecular property variations and other effects outlined above. As a result, the shape of the near-wall profile is again seen to mimic that returned by the full low-Reynolds-number model, with a region of significantly low velocity adjacent to the wall where the buoyancy effects act to slow down the near-wall fluid. Since the analytical wall function subsequently provides wall shear stress and temperature values from these near-wall analytical profiles, the increase in heat transfer due to buoyancy is also captured by this scheme.
© copyright ERCOFTAC 2004
Contributors: Mike Rabbitt - British Energy
Site Design and Implementation: Atkins and UniS