CFD Simulations AC3-01: Difference between revisions

From KBwiki
Jump to navigation Jump to search
m (Dave.Ellacott moved page Gold:CFD Simulations AC3-01 to CFD Simulations AC3-01 over redirect)
 
(48 intermediate revisions by 2 users not shown)
Line 1: Line 1:
{{AC|front=AC 3-01|description=Description_AC3-01|testdata=Test Data_AC3-01|cfdsimulations=CFD Simulations_AC3-01|evaluation=Evaluation_AC3-01|qualityreview=Quality Review_AC3-01|bestpractice=Best Practice Advice_AC3-01|relatedUFRs=Related UFRs_AC3-01}}
='''Buoyancy-opposed wall jet'''=
='''Buoyancy-opposed wall jet'''=


Line 20: Line 22:
Model 1: Eddy-Viscosity scheme with wall functions (KEWF)
Model 1: Eddy-Viscosity scheme with wall functions (KEWF)


This model applies the standard eddy-viscosity method, wherein transport equations are solved for the turbulence energy, κ and its dissipation rate ε and from the values of k and ε values of the turbulent viscosity μt are obtained. When the turbulent viscosity is multiplied by the mean strain rate, one obtains the Reynolds stresses. This model, however, cannot directly handle the low Reynolds number region near the wall and necessitates the use of a wall function to define the behaviour in the near-wall control volumes. In addition to the conventional wall-function strategies, that assume a logarithmic variation for the near-wall velocity, a novel approach is also tested here, in which the near-wall velocity variation is obtained from the analytical solution of the near-wall form of the momentum equation. Both these strategies are further discussed later on in this section.
This model applies the standard eddy-viscosity method, wherein transport equations are solved for the turbulence energy, &kappa; and its dissipation rate &epsilon; and from the values of &kappa; and &epsilon; values of the turbulent viscosity &mu;<sub>t</sub> are obtained. When the turbulent viscosity is multiplied by the mean strain rate, one obtains the Reynolds stresses. This model, however, cannot directly handle the low Reynolds number region near the wall and necessitates the use of a wall function to define the behaviour in the near-wall control volumes. In addition to the conventional wall-function strategies, that assume a logarithmic variation for the near-wall velocity, a novel approach is also tested here, in which the near-wall velocity variation is obtained from the analytical solution of the near-wall form of the momentum equation. Both these strategies are further discussed later on in this section.


   
   
Line 32: Line 34:
Model 3: Basic Second Moment Closure (Basic2mc)
Model 3: Basic Second Moment Closure (Basic2mc)


The second pair of turbulence closures solve transport equations for the four non-zero Reynolds stresses <uiuj>. In producing exact transport equations for these terms, one arrives at some terms which can be implemented exactly and others terms which require some form of turbulence model in order to close the equations. The first closure considered here is generally referred to as the ‘Basic Model’ adopting the Launder, Reece and Rodi (1975) proposal to handle the pressure strain jij, assumes isentropic dissipation εij, and used the generalised gradient diffusion hypothesis of Daly and Harlow (1970) for the diffusion dij. Two additional features are needed in order to handle the effects of the rigid wall correctly. Firstly, since this is essentially a high Reynolds number scheme, wall functions are required as described above for the KEWF model. Secondly, one needs to adapt the pressure-strain model jij, introducing a wall-reflection component ωwij. The form adopted here is that of Craft and Launder (1992), which is a development of the earlier and better known Gibson and Launder (1978) model, but which has been designed to handle flows impinging on, as well as parallel to the wall.
The second pair of turbulence closures solve transport equations for the four non-zero Reynolds stresses <&mu;<sub>i</sub>&mu;<sub>j</sub>>. In producing exact transport equations for these terms, one arrives at some terms which can be implemented exactly and others terms which require some form of turbulence model in order to close the equations. The first closure considered here is generally referred to as the ‘Basic Model’ adopting the Launder, Reece and Rodi (1975) proposal to handle the pressure strain j<sub>ij</sub>, assumes isentropic dissipation &epsilon;<sub>ij</sub>, and used the generalised gradient diffusion hypothesis of Daly and Harlow (1970) for the diffusion d<sub>ij</sub>. Two additional features are needed in order to handle the effects of the rigid wall correctly. Firstly, since this is essentially a high Reynolds number scheme, wall functions are required as described above for the KEWF model. Secondly, one needs to adapt the pressure-strain model j<sub>ij</sub>, introducing a wall-reflection component &omega;<sup>w</sup><sub>ij</sub>. The form adopted here is that of Craft and Launder (1992), which is a development of the earlier and better known Gibson and Launder (1978) model, but which has been designed to handle flows impinging on, as well as parallel to the wall.


   
   
Line 38: Line 40:
Model 4: Two-Component Limit Second Moment Closure (TCL2mc)
Model 4: Two-Component Limit Second Moment Closure (TCL2mc)


The fourth turbulence model considered here is a more advanced form of second-moment closure. Workers at UMIST over the last decade have evolved a model for the pressure-strain jij and dissipation εijprocesses (Kidger (2000)). This model exhibits two features: it is implicitly realisable, in that it cannot produce non-physical values for the Reynolds stresses, particularly when the turbulence approaches a two-component state as occurs near a wall or free-surface. Its second feature is that it eliminates the need for a wall-reflection process, or more importantly, the need to define a wall-normal vector y. In its first decade of testing this TCL (i.e. two-component limit) model has been found to yield superior results to the Basic model for both wall jets and flows which exhibit buoyant effects, both of which are very relevant to the present study (Craft et al, 1996; Craft and Launder, 2001).
The fourth turbulence model considered here is a more advanced form of second-moment closure. Workers at UMIST over the last decade have evolved a model for the pressure-strain j<sub>ij</sub> and dissipation &epsilon;<sub>ij</sub>processes (Kidger (2000)). This model exhibits two features: it is implicitly realisable, in that it cannot produce non-physical values for the Reynolds stresses, particularly when the turbulence approaches a two-component state as occurs near a wall or free-surface. Its second feature is that it eliminates the need for a wall-reflection process, or more importantly, the need to define a wall-normal vector y. In its first decade of testing this TCL (i.e. two-component limit) model has been found to yield superior results to the Basic model for both wall jets and flows which exhibit buoyant effects, both of which are very relevant to the present study (Craft et al, 1996; Craft and Launder, 2001).


   
   
Line 46: Line 48:
When the linear eddy-viscosity model is used for momentum transport, the scalar fluxes are consistently represented by
When the linear eddy-viscosity model is used for momentum transport, the scalar fluxes are consistently represented by


<math>
<center><math>
\overline{u_i\theta} = -{\nu_t \over {\sigma_\theta}}{{\partial\Theta}\over{\partial x_i}}
\overline{u_i\theta} = -{\nu_t \over {\sigma_\theta}}{{\partial\Theta}\over{\partial x_i}}
</math>
</math></center>


With a second-moment closure, however, rather than the differential transport equations for <math>\overline{u_i\theta}</math>, the so-called generalized gradient diffusion hypothesis (GGDH) is used instead:
With a second-moment closure, however, rather than the differential transport equations for <math>\overline{u_i\theta}</math>, the so-called generalized gradient diffusion hypothesis (GGDH) is used instead:


<math>
<center><math>
\overline{u_i\theta} = -c_\theta\ {\kappa \over \varepsilon}\ \overline{u_i u_j}\  
\overline{u_i\theta} = -c_\theta\ {k \over \varepsilon}\ \overline{u_i u_j}\  
{{\partial\Theta}\over{\partial x_j}}
{{\partial\Theta}\over{\partial x_j}}
</math>
</math></center>




Line 207: Line 209:
Computations have been carried out with different lengths of the computational domain, both below (Y<sub>L</sub>) and above (Y<sub>U</sub>) the splitter plate, in order to ensure that the boundary locations did not affect the computations. In the isothermal case the lower entry plane was located 1.3m below the splitter plate, so as not to influence the predicted penetration length of the jet. In non-isothermal flows the upper exit boundary was moved to 1.8 m above the splitter plate to prevent back flow across the exit boundary.
Computations have been carried out with different lengths of the computational domain, both below (Y<sub>L</sub>) and above (Y<sub>U</sub>) the splitter plate, in order to ensure that the boundary locations did not affect the computations. In the isothermal case the lower entry plane was located 1.3m below the splitter plate, so as not to influence the predicted penetration length of the jet. In non-isothermal flows the upper exit boundary was moved to 1.8 m above the splitter plate to prevent back flow across the exit boundary.


Different values for the turbulent kinetic energy and its dissipation rate have also been tried at the lower inlet boundary. The same case (isothermal, V<sub>chan</sub> /V<sub>jet</sub>=0.077) was computed with either k<sub>in</sub> = 0.01V<sub>in</sub><sup>2</sup> and εin = rcμkin2/(50μ), or kin = 0.025Vin2 and εin = kin3/2/lin, where lin is prescribed from DNS data for channel flow. The resulting predictions were identical.
Different values for the turbulent kinetic energy and its dissipation rate have also been tried at the lower inlet boundary. The same case (isothermal, V<sub>chan</sub> /V<sub>jet</sub>=0.077) was computed with either k<sub>in</sub> = 0.01V<sub>in</sub><sup>2</sup> and &epsilon;<sub>in</sub>&nbsp;=&nbsp;rc<sub>&mu;</sub>&kappa;<sub>in</sub><sup>2</sup>/(50&mu;), or &kappa;<sub>in</sub>&nbsp;=&nbsp;0.025V<sub>in</sub><sup>2</sup> and &epsilon;<sub>in</sub>&nbsp;=&nbsp;&kappa;<sub>in</sub><sup>3/2</sup>&nbsp;/&nbsp;l<sub>in</sub>, where l<sub>in</sub> is prescribed from DNS data for channel flow. The resulting predictions were identical.


With regard to the mesh spacing adjacent to the wall, with the low-Re model a maximum dimensionless near-wall mesh size (y*) of 1 at the jet exit, based on the near-wall turbulent kinetic energy, has been used. In the high-Re computations, mesh sizes have been used which result in y* values of between 65 and 160 at the jet exit, with no discernible differences between the predictions. On the wall opposite the jet, mesh sizes have been used which result in y* values of the order 5 to 10 at the bottom of the domain and 40 to 80 at the top, again with no discernible differences between the predictions.
With regard to the mesh spacing adjacent to the wall, with the low-Re model a maximum dimensionless near-wall mesh size (y*) of 1 at the jet exit, based on the near-wall turbulent kinetic energy, has been used. In the high-Re computations, mesh sizes have been used which result in y* values of between 65 and 160 at the jet exit, with no discernible differences between the predictions. On the wall opposite the jet, mesh sizes have been used which result in y* values of the order 5 to 10 at the bottom of the domain and 40 to 80 at the top, again with no discernible differences between the predictions.
Line 215: Line 217:
The cases computed are tabulated in Tables 2 and 3. In all cases the jet Reynolds number is 4,000.
The cases computed are tabulated in Tables 2 and 3. In all cases the jet Reynolds number is 4,000.


=='''Assessment of Computations'''==
'''Isothermal Flows'''
The assessment first focuses on comparisons between the LES and experimental data for the case of Vchan /Vjet=0.077, in the form of contours of the vertical, (V) velocity components, shown in Figure 13. The comparison shows that there is close agreement between LES and the experiment as far as the width and the penetration length of the downward jet is concerned. The jet reaches a distance of about 0.54m along the right side of the channel (from the splitter plate at y=1.3m to about 0.76m). The validity of the RANS computations can thus be assessed through comparisons with the more detailed LES predictions as well as through comparisons with the available experimental data.
In Figure 14, the LES contours of the vertical velocity are compared with those of the RANS computations. The low-Re k-ε model returns a jet penetration length of about 0.74m, while use of the high-Re k-ε model with the standard wall-function approach results in a penetration length of about 0.98m. These first two comparisons would suggest that neither the effective viscosity approximation nor the assumptions involved in the standard wall function are appropriate for opposed wall-jet flows. The introduction of the analytical wall function results in good agreement between the high-Re and low-Re k-ε predicted V contours, indicating that this novel wall-function strategy is general enough to be used in the prediction of opposed wall-jet flows.
The remaining four sets of contours in Figure 14 were produced with second-moment closures. They indicate that while the introduction of the basic second-moment closure results in only modest predictive improvements in comparison with the effective-viscosity model, the more elaborate TCL version produces more substantial improvements. Indeed the combination of the TCL second-moment closure with the analytical wall function results in a jet penetration length that is very close to that of the LES computations. Comparisons between measured and computed profiles of the vertical mean velocity, shown in Figures 15, 16 and 17, also confirm these conclusions. There is generally close agreement between the measured profiles and those computed using either LES, or the TCL second-moment closure with the analytical wall function. Comparisons between computed and measured profiles of the cross-channel component of the mean velocity, shown in Figures 18 to 20, also lead to similar conclusions, though for this component the LES profiles are closer to those measured than those of the TCL model with the analytical wall functions. Finally, for this case, Figures 21 to 23 show comparisons between the computed and measured profiles of the turbulent kinetic energy. The measurements show that k reaches its maximum levels along the right hand side of the channel, in the region where the wall jet turns upwards. The EVM computations manage to capture both the overall trend and the actual levels. The second-moment closures return similar distributions, but where the jet turns upwards somewhat lower levels. The LES computations produce k distributions that are different to those measured, possibly because the simulation has not advanced a sufficiently large number of time steps to produce statistically independent values of the turbulent kinetic energy.
While the opposed wall jet flow is geometrically very simple, the resulting flow structure is complex enough for both the effective viscosity approximation and also the assumption of a logarithmic near-wall velocity variation to break down. In reviewing earlier numerical studies of wall jets developing in stagnant surroundings, Launder and Rodi (1983) showed that the effective-viscosity model produced a more rapid spreading rate than second-moment closures. This was attributed to the prediction of stronger mixing between the jet and the surroundings by the EVM, which in the case of the opposed wall jet, in contrast to the findings of the present study, should have resulted in the computation of a shorter penetration length. The question that therefore needs to be addressed is why EVM models in this case return a longer penetration length than second-moment closures. One major difference between the case of a self-similar wall jet in stagnant surroundings and a confined, opposed wall jet in a channel, is the presence of an adverse pressure gradient. A closer look at the vertical velocity profiles of Figures 15 to 17 shows that the TCL-AWF modelling combination above and just below the splitter plate (yjet-y < 0.3 m) predicts stronger downward flow on the outer side of the jet. Further downstream (yjet-y = 0.4 m) this model starts to predict a lower downward velocity along the right wall and also a lower upward velocity along the left wall. The implication is that the loss of the wall-jet’s downward momentum arises principally from entrainment of fluid into the outer layer of the wall jet (which reduces the level of the downward velocity) accompanied by a strongly rising pressure in the downward direction. We note from the vector plots of the velocity field, Figure 24, that the TCL closure does lead to a more rapid thickening of the wall jet than the k-ε scheme. Indeed, for the TCL predictions, notice that external fluid is drawn from well above the origin of the wall jet, for entrainment into the shear flow. The extra entrainment means that the TCL wall jet is less able to withstand the strong adverse pressure gradient, Figure 25(a). As shown by the shear stress profiles of Figure 26, the second-moment closures produce higher levels of the turbulent shear stress above the splitter plate than the EVM, most likely due to the effects of convective transport. This is consistent with flow entrainment shown in Figure 24. A further contributing factor to the loss of the wall jet’s momentum is the wall friction. We note in Figure 25 that, for the first 0.4m downstream from discharge, the TCL closure gives a friction factor that, on average, is 25% higher than the k-ε scheme. It can also be seen, in the same figure, that the introduction of the analytical wall function also increases the wall shear stress at the initial stages of the jet, which is again consistent with the shorter penetration length produced by this wall-function strategy.
For the case of Vchan /Vjet=0.15, Figure 27 shows comparisons between the measured and the predicted contours of the vertical component of the mean velocity. As also seen in the earlier comparisons, the high-Re k-ε with the standard wall function produces a longer penetration length than the low-Re k-ε and with the introduction of the analytical wall function the high-Re k-ε predictions are brought to close agreement with those or the low-Re model. Replacing the EVM model with the TCL second-moment closure results in a further reduction of the predicted penetration length, but this reduction is not as substantial as in the case of Vchan/Vjet=0.077. Moreover, the measured contours show that the penetration length is more than a factor of two shorter. The measured contours, however also show that the upward inlet velocity at the lower boundary is higher than the nominal inlet velocity, by about a factor of two. Whether this is the result of three-dimensional variations, or of experimental uncertainties in the determination of the inlet flow rates it is hard to say. Comparisons between the measured and predicted mean velocity profiles are presented in Figures 28 to 33. They also raise the same questions about the experimental conditions. The corresponding k-profile comparisons are presented in Figures 34 to 36. For the first 0.2m after the splitter plate all RANS models return the correct k levels and distribution. From 0.2m to 0.5m after the splitter plate, where the measured jet has turned upwards while the predicted jets continue to proceed downwards, the turbulence models over predict k levels.
'''Non-Isothermal Flows'''
So far, for the non-isothermal cases, computations have only been carried out using EVM models.
For the case of Vchan/Vjet=0.077, the vector plot comparisons of Figure 37, show that the three EVM models employed (low-Re k-ε, high-Re k-ε with standard wall-functions, and high-Re k-ε with analytical wall functions) return similar penetration lengths that are in agreement with that implied by the measurements. The spreading of the jet across the channel is, however, under-estimated. While the measurements show that once it separates from the wall, the jet turns gradually upwards, thus spreading across the channel, in the EVM predictions, the jet turns sharply upwards, making a tight 180o turn, and attaches itself on the other side of the splitter plate. In contrast to the available data, above the jet entry point, the fluid on the left half of the channel moves downwards. Comparisons of contour plots for temperature, turbulent kinetic energy and of the vertical velocity, shown in Figures 38 to 40, also present a similar picture. The measurements show that as soon as the jet turns high turbulence and temperature levels spread across the channel. The EVM predictions on the other hand show that initially the high temperature fluid is confined to the right side of the channel and only crosses to the left side further up above the jet entry location, as a result of the recirculating motion. High turbulence levels are also initially confined to the left side of the channel and as the second, upward directed wall jet develops along the splitter plate, high turbulence levels gradually spread to the left side of the channel. More detailed comparisons, that include profiles of the two mean velocity components and also the turbulent kinetic energy, are presented in Figures 41 to 43. They also confirm that the three EVM models produce similar predictions, especially the low-Re model and the high-Re model with the analytical wall function, which return the correct penetration length, but under-estimates the mixing of the jet. The latter trend is most obvious in the k profile comparisons.
The differences between the EVM predictions and the measured behaviour appear to originate from the failure of the EVM model to return the rapid spread of high turbulence levels across the channel as the wall jet turns upwards. This in turn inhibits the transfer of thermal energy across the channel, and leaves the cooler, heavier fluid on the left side causing the weak down flow present in the predictions. A possible cause of the failure of the EVM model to return the rapid spread of turbulence levels is the way that the contribution of buoyancy to the generation rate of turbulence is taken into account. In buoyant flows, there is an additional contribution to the generation rate of turbulence Pg, the exact expression of which is
<math>\[P_{g} = -\overline{\rho'u_{i}g{i}} = \overline{\rho u_{i}\theta g_{i}}/\Theta\]</math>
<math>\[P_{g} = -\frac{g_{i}}{\Theta}\frac{\mu_{t}}{\sigma_\theta}\frac{\partial\Theta}{\partial x_{i}}\]</math>
Introduction of the effective diffusivity approximation for the turbulent heat fluxes leads to
<math>\[P_{g} = -\frac{g}{\Theta}\frac{\mu_{t}}{\sigma_\theta}\frac{\partial\Theta}{\partial {y}}\]</math>
In this specific case where the y direction is the vertical direction
What the above expression shows is that as a consequence of the effective diffusivity approximation, the buoyant generation rate term is only active in regions with a temperature gradient in the vertical direction. Such a region occurs at the interface between the downward directed hot jet fluid and the upward-directed, cold channel fluid. Since at that interface the temperature gradient ¶Θ/¶y is positive, Pgwould be negative, causing a reduction in turbulence levels.
<math>\[P_{g} = -\rho\frac{g_{i}}{\Theta}c_\theta\frac{k}{E}\overline{u_{i}u_{j}}\frac{\partial\Theta}{\partial x_{j}}\]</math>
A more realistic way of representing the turbulent heat fluxes would be through the generalized gradient diffusion hypothesis, GGDH. This is routinely used with second-moment closures for the Reynolds stresses and, as shown by Ince and Launder (1989), can also be used with EVM closures for the flow field, though with a different value for the constant cΘ. This approximation will result in
<math>\[P_{g} = -\rho\frac{g}{\Theta}c_\theta\frac{k}{E}\left({\overline{u\nu}\frac{\partial\Theta}{\partial x} + \overline{\nu^2}\frac{\partial\Theta}{\partial {y}}\right)\]</math>
For the wall-opposed jet geometry the above expression becomes
Consequently, use of GGDH will give rise to two contributions to the buoyant generation rate of turbulence, one due to the vertical temperature gradient, which is similar to that present when the effective diffusivity approximation is used, and one due to a horizontal temperature gradient. The latter will become active along the middle of the channel, after the jet has turned upwards, where ¶Θ/¶x will be positive and the turbulent shear stress <uv> will be negative. This additional contribution would thus tend to increase turbulence levels, which means it is likely to improve agreement with the experimental data.
The corresponding comparisons for the higher velocity ratio, reveal a different behaviour, as far as the prediction of the penetration length is concerned. The vector and contour plots of Figures 44 to 47 show that the low-Re k-ε and also the high-Re k-ε with the analytical wall functions produce a penetration length, which is similar to what was measured. In these predictions, and the measurements, the jet starts turning upwards as soon as it is injected into the channel. These two models (the low-Re k-ε model in particular) also predict a spreading rate across the channel closer to that shown in the experimental data. The high-Re k-ε with the standard wall function, on the other hand, produces a penetration length that is twice as long and a jet that, as in the previous case of the lower velocity ratio, does not spread across the channel. Arguably, because with the standard wall function the jet is predicted to remain attached for longer, the predicted turbulence levels are higher. Along the left wall the EVM models again predict that the upward flow eventually reverses its direction, though at this higher velocity ratio this happens further up the channel. More detailed profile comparisons are shown in Figures 48 to 50. They confirm that there are now substantial differences between the predictions produced by the high-Re model with the standard wall function and those returned by the either the low-Re model or the high-Re model with the analytical wall function. It is interesting to note that the double peaks in the measured cross-channel velocity over the initial stages of the jet development, (yjet-y<0.06 m in Figure 49) are faithfully reproduced by the two EVM computations that do not make use of the standard wall-function. The resulting abrupt changes in the gradient of the measured vertical velocity over the same region (Figure 48) are also equally well predicted by the same two models. As far as the k profiles are concerned, the low-Re EVM model and the high-Re model with the analytical wall-function, in contrast to the high-Re model with the standard wall function, under-estimate the k levels within the jet. They do, however, predict that the high turbulence region extends further across the channel in comparison with the predictions obtained with the standard wall function, though not as far as what is indicated by the measurements.
One question that arises here, is why at the higher velocity ratio the predictions based on the standard wall function result in a longer penetration length than those of the other two near-wall treatments. The main difference as far as the flow dynamics are concerned is that at the higher velocity ratio the wall jet is subjected to a stronger adverse pressure gradient. Not surprisingly, the assumptions of the standard wall function (log-law, constant shear stress, fixed dimensionless thickness of the viscous sub-layer etc) become more restrictive as the adverse pressure gradient becomes stronger. It thus appears that modelling of near-wall turbulence is critical to the computation of both isothermal and non-isothermal opposed wall-jet flows and that the analytical wall function provides an effective and inexpensive alternative to the use of low-Reynolds-number models.


A further question that is currently being addressed is whether the introduction of second-moment closures for the Reynolds stresses and the GGDH for the turbulent heat fluxes can result in further predictive improvements for non-isothermal flows.


© copyright ERCOFTAC 2004
© copyright ERCOFTAC 2004
Line 281: Line 227:
Site Design and Implementation: [[Atkins]] and [[UniS]]
Site Design and Implementation: [[Atkins]] and [[UniS]]
        [[#Top|Top]]              Next
{{AC|front=AC 3-01|description=Description_AC3-01|testdata=Test Data_AC3-01|cfdsimulations=CFD Simulations_AC3-01|evaluation=Evaluation_AC3-01|qualityreview=Quality Review_AC3-01|bestpractice=Best Practice Advice_AC3-01|relatedUFRs=Related UFRs_AC3-01}}

Latest revision as of 15:48, 11 February 2017

Front Page

Description

Test Data

CFD Simulations

Evaluation

Best Practice Advice

Buoyancy-opposed wall jet

Application Challenge 3-01 © copyright ERCOFTAC 2004


Overview of CFD Simulations

Introduction

A brief overview of the computational approach is given in paragraph 1.1 of this report. The preliminary computational studies described in issue 1 had considered only the isothermal case. This enabled an early evaluation of the different turbulence models considered, as well as providing verification that the implementation of the computational methods was correct. The more extensive studies, reported here, covered both isothermal and non-isothermal cases at two different velocity ratios (Vchan / Vjet 0.077 and 0.15). Moreover, two alternative wall functions have been used with all the high-Reynolds-number models. For the ratio Vchan / Vjet of 0.077, LES computations have also been carried out by Prof. Laurence’s group.


RANS Turbulence Modelling

The four turbulence models used in more detail are:

Model 1: Eddy-Viscosity scheme with wall functions (KEWF)

This model applies the standard eddy-viscosity method, wherein transport equations are solved for the turbulence energy, κ and its dissipation rate ε and from the values of κ and ε values of the turbulent viscosity μt are obtained. When the turbulent viscosity is multiplied by the mean strain rate, one obtains the Reynolds stresses. This model, however, cannot directly handle the low Reynolds number region near the wall and necessitates the use of a wall function to define the behaviour in the near-wall control volumes. In addition to the conventional wall-function strategies, that assume a logarithmic variation for the near-wall velocity, a novel approach is also tested here, in which the near-wall velocity variation is obtained from the analytical solution of the near-wall form of the momentum equation. Both these strategies are further discussed later on in this section.


Model 2: Low-Reynolds-number Eddy-Viscosity scheme (LRKE)

One need not, necessarily, be reliant on wall functions to compute the flow behaviour near to a rigid wall. If one adapts the turbulence model such that it can properly compute the flow in regions of low turbulent Reynolds numbers, one may solve the transport equations through the viscous sublayer, right up to the wall, without the need for a truncated grid and wall functions. The low Reynolds-number eddy-viscosity scheme used here is essentially that of Launder and Sharma (1974). There is a small modification to the transport equation for the dissipation ε, from that originally proposed. Further details are provided in Media:Appendix_A_issue2e.pdf


Model 3: Basic Second Moment Closure (Basic2mc)

The second pair of turbulence closures solve transport equations for the four non-zero Reynolds stresses <μiμj>. In producing exact transport equations for these terms, one arrives at some terms which can be implemented exactly and others terms which require some form of turbulence model in order to close the equations. The first closure considered here is generally referred to as the ‘Basic Model’ adopting the Launder, Reece and Rodi (1975) proposal to handle the pressure strain jij, assumes isentropic dissipation εij, and used the generalised gradient diffusion hypothesis of Daly and Harlow (1970) for the diffusion dij. Two additional features are needed in order to handle the effects of the rigid wall correctly. Firstly, since this is essentially a high Reynolds number scheme, wall functions are required as described above for the KEWF model. Secondly, one needs to adapt the pressure-strain model jij, introducing a wall-reflection component ωwij. The form adopted here is that of Craft and Launder (1992), which is a development of the earlier and better known Gibson and Launder (1978) model, but which has been designed to handle flows impinging on, as well as parallel to the wall.


Model 4: Two-Component Limit Second Moment Closure (TCL2mc)

The fourth turbulence model considered here is a more advanced form of second-moment closure. Workers at UMIST over the last decade have evolved a model for the pressure-strain jij and dissipation εijprocesses (Kidger (2000)). This model exhibits two features: it is implicitly realisable, in that it cannot produce non-physical values for the Reynolds stresses, particularly when the turbulence approaches a two-component state as occurs near a wall or free-surface. Its second feature is that it eliminates the need for a wall-reflection process, or more importantly, the need to define a wall-normal vector y. In its first decade of testing this TCL (i.e. two-component limit) model has been found to yield superior results to the Basic model for both wall jets and flows which exhibit buoyant effects, both of which are very relevant to the present study (Craft et al, 1996; Craft and Launder, 2001).


Modelling of the Turbulent Heat Fluxes

When the linear eddy-viscosity model is used for momentum transport, the scalar fluxes are consistently represented by

With a second-moment closure, however, rather than the differential transport equations for , the so-called generalized gradient diffusion hypothesis (GGDH) is used instead:


The experience at UMIST is that this form is usually adequate as long as the influential turbulent scalar and momentum fluxes are horizontal. It is emphasized, however, that in horizontally directed shear flows affected by buoyancy where the turbulent fluxes are vertical, it would be essential to use a complete second-moment closure, as described in Craft et al (1996). Indeed, in situations of extreme stable stratification third-moment closure is essential, Craft and Launder (2002), Kidger (2000).


Wall Function Strategies

In near-wall flows, the turbulence field undergoes rapid changes across a very thin near-wall sub-layer, known as the viscous sub-layer. Integration of the mean flow equations from the fully turbulent outer region to the wall is possible, but requires the use of specially extended models of turbulence, known as low-Reynolds number models, and also very fine near-wall meshes. The latter requirement makes this approach very expensive, especially in three-dimensional flows. A more economical alternative, the wall-function approximation, has consequently been developed and has become the most widely used practice for the mathematical modelling of the effects of near-wall turbulence. Its popularity lies in the fact that it does not need to resolve the viscous sub-layer. The near-wall grid node is located far enough from the wall to be in the fully turbulent region. The semi-empirical co-relations are then used to account for the wall effects on turbulence, resulting in estimates for the wall shear stress and also for the generation and dissipation rates of turbulence over the wall-adjacent control volumes. The main drawback of such strategies, at least until recently, has been the assumption that the near-wall region is one of uniform shear stress and with the production and dissipation rate of k in balance. This has limited the range of flows for which the use of wall functions can result in reliable predictions. A more refined alternative has been recently developed by the UMIST group, Craft et al (2002), which relies on assumptions made at a deeper, more generally valid level. In the present work this alternative has been used in combination with the three high-Reynolds-number models, in addition to the standard wall function. Both the standard and the Analytical (UMIST-A) wall functions are briefly outlined below.


Standard Wall-Function

Image209.gif


Figure 11. Prescribed semi-logarithmic velocity and temperature profiles


In this widely used approach, the variation in velocity between the near-wall node and the wall is assumed to follow the log-law. The total shear stress is assumed to remain constant and equal to the wall shear stress. The turbulent shear stress is assumed to be constant (also equal to the wall shear stress) up to the edge of the viscous sub-layer where it is assumed to abruptly fall to zero. Outside the viscous sub-layer, the dissipation rate of turbulence is assumed to be proportional to the inverse of the wall distance, and within the viscous sub-layer it is assumed to remain constant. Moreover the dimensionless thickness of the viscous sub-layer is also assumed to be independent of the mean flow field. Inevitably this rather large set of assumptions makes this approach unsuitable for many flows, including those affected by buoyancy.


Analytical Wall Function, Craft et al (2002).


Image210.gif


Figure 12. Prescribed variation of turbulent viscosity


The overall strategy is the same as in the standard wall function, namely the use of a large wall-adjacent control-volume and the wall shear stress and the generation and dissipation rates of turbulence estimated from certain assumptions. The approximations involved, however, are not as limiting as those of the standard wall-function. The wall shear stress and the appropriate thermal boundary condition are obtained through the analytical solution of simplified near-wall versions of the transport equations for the wall-parallel momentum and the thermal energy respectively. These analytical solutions are obtained by prescribing the variation of the turbulent viscosity and also by assuming that convective transport and the wall-parallel pressure gradient do not change across the near-wall control volumes. Over an inner viscosity-dominated sub-layer, determined by the value of the dimensionless wall distance y* (º ykP1/2/n), the turbulent viscosity is assumed to be zero. Outside the viscosity-dominated sub-layer, the turbulent viscosity is taken to be proportional to the distance from the edge of this viscosity-dominated sub-layer. The dimensionless thickness of the viscosity-dominated sub-layer is determined with reference to fully developed pipe flows. The analytical velocity variation and the prescribed turbulent viscosity lead to an expression for the local generation rate of turbulence, which can then be integrated over the wall-adjacent control volume. Further refinements have also been added to this approach. In heated flows the molecular viscosity is allowed to vary with temperature across the viscosity-dominated sub-layer. The analytical solution of the momentum equation also takes into account the buoyancy force, which is introduced through the analytical temperature variation. The dissipation rate of turbulence over the near-wall control volume is made sensitive to the mean flow development through the introduction of a parameter based on the ratio between the wall shear stress and the shear stress at the edge of the viscosity-dominated sub-layer. The wall-parallel momentum flux across faces of wall-adjacent control volumes is evaluated from the analytical velocity variation. In high Prandtl number fluids an effective molecular Prandtl number is introduced, to correct for the fact that the conduction sub-layer is thinner than the viscosity-dominated sublayer. These refinements make the analytical wall function applicable over a wider range of flows than the standard approach.

The equations involved together with the physical basis of the assumptions are included in Appendix B.


LES

LES computations have been carried out using two different codes. One is a research code, Saturne, developed at EDF and the other is the commercial code Star-CD. Both codes use second-order centered schemes in time and space and they employ pressure correction methods for the velocity-pressure coupling. Both the classical, Smagorinsky, and the dynamic, Germano et al sub-grid models have been tested.


Computational grid

A computational domain was chosen which covered the total flow measurement region (Figure 1). This necessitated having a cut-out in the flow to handle the region of the thick splitter plate. The inlet profiles of the two fluid streams were not reported in the experimental work, and the strategy adopted here has been to assume a ‘top-hat’ velocity profile both at entry to the splitter plate and at the bottom of the test section where the counter-flowing main stream enters. It is assumed that the solver will yield an acceptable boundary layer profile by the time the jet exits into the main cavity at the bottom of the splitter plate and where the channel flow reaches the mixing region.

The above strategy, however, has shortcomings, because of the presence of the splitter and the fact that the TEAM code uses a structured staggered grid. At the edges of the computational domain, the boundaries are arranged in such a way as to make the application of the appropriate boundary conditions easy since all dependent variables have a node on the boundary. However, this arrangement becomes somewhat complicated as a result of having the splitter plate within the computational domain since there will not be any nodes on the boundary on the splitter surface. For that reason, the strategy adopted was to use two distinct calculations, the first to compute the developing profile of the inlet jet between the splitter and the outer wall and the second to resolve the flow in the primary domain taking the values from the former to define the inlet jet.



Table 2. Isothermal Flow Computations
Vch/Vjet Tch °C Tjet °C Turbulence Model Wall Function Flow Domain Grids
YL YU
O.077 42 42 Low-Re κ-ε NA 1.3m 0.4m 130x280

103x165

0.077 42 42 High-Re κ-ε Standard 1.3m 0.4m 102x280

75x165

0.077 42 42 High-Re κ-ε Analytical 1.3m 0.4m 102x280

75x165

0.077 42 42 Low-Re κ-ε NA 1m 1m 130x300
0.077 42 42 High-Re κ-ε Standard 1m 1m 102x300
0.077 42 42 High-Re κ-ε Analytical 1m 1m 102x300
0.077 42 42 Basic 2mc Standard 1.3m 0.4m 102x280

75x165

0.077 42 42 Basic 2mc Analytical 1.3m 0.4m 102x280

75x165

0.077 42 42 TLC 2mc Standard 1.3m 0.4m 102x280

75x165

0.077 42 42 TLC 2mc Analytical 1.3m 0.4m 102x280

75x165

0.077 42 42 LES NA 1m 0.4m 58x142x27

59x141x40 601.6x103

0.15 42 42 Low-Re κ-ε NA 1.3m 0.4m 130x280

103x165

0.15 42 42 High-Re κ-ε Standard 1.3m 0.4m 102x280

75x165

0.15 42 42 High-Re κ-ε Analytical 1.3m 0.4m 102x280

75x165

0.15 42 42 Basic 2mc Standard 1.3m 0.4m 75x165
0.15 42 42 Basic 2mc Analytical 1.3m 0.4m 75x165
0.15 42 42 TLC 2mc Standard 1.3m 0.4m 75x165
0.15 42 42 TLC 2mc Analytical 1.3m 0.4m 75x165


YL and YU are the lengths of the computational domain below and above the splitter plates respectively

† Unstructured Mesh



Table 3. Non-Isothermal Flow Computations
Vch/Vjet Tch °C Tch °C Turbulence Model Wall Function Flow Domain Grids
YL YU
0.077 34 42 Low-Re κ-ε NA 0.6m 1.8m 130x280 130x150
0.077 34 42 High-Re κ-ε Standard 0.6m 1.8m 102x280 75x150
0.077 34 42 High-Re κ-ε Analytical 0.6m 1.8m 102x280 75x150
0.15 34 42 Low-Re κ-ε NA 0.6m 1.8m 130x280 103x150
0.15 34 42 High-Re κ-ε Standard 0.6m 1.8m 102x280 75x150
0.15 34 42 High-Re κ-ε Analytical 0.6m 1.8m 102x280 75x150
0.15 34 42 Low-Re κ-ε NA 0.6m 2.4m 130x300 103x150
0.15 34 42 High-Re κ-ε Standard 0.6m 2.4m 102x300 75x150
0.15 34 42 High-Re κ-ε Analytical 0.6m 2.4m 102x300

75x150


Computations have been carried out with different lengths of the computational domain, both below (YL) and above (YU) the splitter plate, in order to ensure that the boundary locations did not affect the computations. In the isothermal case the lower entry plane was located 1.3m below the splitter plate, so as not to influence the predicted penetration length of the jet. In non-isothermal flows the upper exit boundary was moved to 1.8 m above the splitter plate to prevent back flow across the exit boundary.

Different values for the turbulent kinetic energy and its dissipation rate have also been tried at the lower inlet boundary. The same case (isothermal, Vchan /Vjet=0.077) was computed with either kin = 0.01Vin2 and εin = rcμκin2/(50μ), or κin = 0.025Vin2 and εin = κin3/2 / lin, where lin is prescribed from DNS data for channel flow. The resulting predictions were identical.

With regard to the mesh spacing adjacent to the wall, with the low-Re model a maximum dimensionless near-wall mesh size (y*) of 1 at the jet exit, based on the near-wall turbulent kinetic energy, has been used. In the high-Re computations, mesh sizes have been used which result in y* values of between 65 and 160 at the jet exit, with no discernible differences between the predictions. On the wall opposite the jet, mesh sizes have been used which result in y* values of the order 5 to 10 at the bottom of the domain and 40 to 80 at the top, again with no discernible differences between the predictions.

In the LES computations three grids have been tested, two of which are fully structured and one that involves local grid refinement.

The cases computed are tabulated in Tables 2 and 3. In all cases the jet Reynolds number is 4,000.


© copyright ERCOFTAC 2004


Contributors: Jeremy Noyce - Magnox Electric

Site Design and Implementation: Atkins and UniS


Front Page

Description

Test Data

CFD Simulations

Evaluation

Best Practice Advice