Large-Eddy Simulation on turbulent flow and plume dispersion over a 2-dimensional hill

Abstract. The dispersion analysis of airborne contaminants including radioactive substances from industrial or nuclear facilities is an important issue for air quality maintenance and safety assessment. In Japan, many nuclear power plants are located at complex coastal terrains. In these cases, terrain e ffects on the turbulent flow and plume dispersion should be investigated. In this study, we perform Large-Eddy Simulation (LES) of turbulent flow and plume dispersion over a 2-dimensional hill flow and investigate the characteristics of mean and fluctuating concentrations.


Introduction
The dispersion analysis of airborne contaminants including radioactive substances from industrial or nuclear facilities is an important issue in air quality maintenance and safety assessment.In a flat terrain, the mean concentration of a plume can be easily predicted by using a Gaussian plume model.However, many nuclear power plants in Japan are located in complex coastal terrains.In this case, the effects of the terrain should be taken into consideration when predicting plume dispersion.Terrain effects on turbulent flow and/or plume dispersion have been investigated in many studies.For example, Jackson and Hunt (1975) proposed a theoretical model for the analysis of turbulent structures over a gentle hill.In one of the earliest wind tunnel experimental studies, Khurshudyan et al. (1981) conducted wind tunnel experiments of turbulent flow and plume dispersion over 2-dimensional hills with different slopes and investigated the influence of hill slope and source location on the flow and concentration fields.Castro and Snyder (1982) conducted wind tunnel experiments on plume dispersion over 3dimensional hills with various ratios of spanwise hill breadth to height and investigated the effect of hill shape on mean concentration distributions.Arya and Gadiyaram (1986) also conducted wind tunnel experiments on flow and dispersion in the wakes of 3-dimensional low hills and investigated the difference between the ground-level mean concentration of a plume released from a point source with and without the hill.
Correspondence to: H. Nakayama (nakayama.hiromasa@jaea.go.jp)Sada (1991) investigated the streamwise variation of vertical profiles of mean wind velocity, turbulence kinetic energy, and mean concentration over a 2-dimensional gentle hill using wind tunnel experiments.In one of the earliest numerical studies, Hino (1968) conducted numerical experiments for plume dispersion over a complex terrain and showed the contaminant distribution patterns.Recently, Castro and Apsley (1997) performed numerical simulations of plume dispersion over 2-dimensional hills using the standard k-ε turbulence model and compared this with the experimental data.
These previous studies have focused on the terrain effects on the characteristic of turbulent flow and/or mean concentration.Another issue related to plume dispersion is the potential problem of the accidental release of hazardous and flammable materials.During the assessment of the hazard to human health from toxic substances, the existence of high concentration peaks in a plume should be considered.For safety analysis of flammable gases, certain critical threshold levels should be evaluated.In such a situation, not only mean but also fluctuating concentrations should be estimated.In this study, we perform a numerical simulation of unsteady behaviors of turbulent flow and plume dispersion over a 2dimensional hill by using an LES that can give more detailed information on the flow and concentration fields in comparison with the wind tunnel experiments of Sada (1991) and investigate the characteristics of mean and fluctuating concentrations.
Published by Copernicus Publications.

Numerical model and computational settings
The basic equations for the LES model are the spatially filtered continuity equation, Navier-Stokes equation and the transport equation for concentration.The subgrid-scale (SGS) Reynolds stress is parameterized by using the standard Smagorinsky model (Smagorinsky, 1963) with a Van Driest damping function (Van Driest, 1956), where the Smagorinsky constant is set to 0.12 (Iizuka and Kondo, 2004) for estimating the eddy viscosity.The subgrid-scale scalar flux is also parameterized by an eddy viscosity model and the turbulent Schmidt number is set to 0.5.Various SGS models for LES have been proposed besides the standard Smagorinsky model.For example, dynamic Smagorinsky models have been proposed by Germano et al. (1991), Lilly (1992), and Meneveau et al. (1996).However, Iizuka and Kondo (2004) examined the influence of various SGS models on the prediction accuracy of LESs of turbulent flow over a hilly terrain and showed that the prediction accuracy of LESs with the standard Smagorinsky model is better than that of the LESs with the dynamic Smagorinsky type models.This indicates that the dynamic Smagorinsky type models are not always effective for determining model constant.As a static type SGS model, Nicoud and Ducros (1999) proposed the walladapting local eddy-viscosity (WALE) model.This model can capture the effects of both the strain and the rotation rate of the small-scale turbulent motions without a damping function from the wall.According to Temmerman et al. (2003), the WALE model shows better performance when compared with the dynamic/standard Smagorinsky model.However, the conventional Smagorinsky model that has the advantage of simplicity and low computational costs is adopted in our LES model because the focus of our research is not the small order effects of turbulent flow.
The coupling algorithm of the velocity and pressure fields is based on the Simplified Marker and Cell (SMAC) method with the second-order Adams-Bashforth scheme for time integration.The SMAC method is the algorithm for solving numerically the Navier-Stokes equation and was proposed by Amsden and Harlow (1970).The Poisson equation is solved by the Successive Over-Relaxation (SOR) method which is an iterative method for solving a Poisson equation for pressure.For the spatial discretization in the governing equation of the flow field, a second-order accurate central difference is used.For the dispersion field, Cubic Interpolated Pseudoparticle (CIP) method proposed by Takewaki et al. (1985) and a second-order accurate central difference method are used for the advection and diffusion terms, respectively.The time step interval ∆tU ∞ /H is 0.005 (∆t: time step).
To perform LES of plume dispersion, a thick turbulent boundary layer (TBL) flow with strong velocity fluctuations should be simulated in order to mimic a plume released from an elevated point source.Various kinds of methods for producing a realistic TBL flow have been proposed.For example, Mochida et al. (1992) performed a preliminary LES of a channel flow for the inlet boundary condition.This method can be easily applied, but the flow is driven by a pressure gradient under the periodic boundary condition in the streamwise direction.Because it is known that the influence of a pressure gradient on TBL characteristics is large (Kline, 1967), it is desirable to reproduce an approaching flow without a pressure gradient.
Figure 1 shows a schematic illustration of the computational regions for plume dispersion over a hill immersed in a fully-developed TBL flow.First, a spatially-developing TBL flow with strong velocity fluctuations in the driver region is generated by incorporating the inflow turbulence generation method proposed by Kataoka and Mizuno (2002) into an upstream small fraction of the driver region with a 2-dimensional roughness bar placed at the ground surface.Next, the inflow turbulence data obtained near the exit of the driver region is imposed at the inlet of the main region at each time step, and calculation of the turbulent flow and plume dispersion over a hill is performed.The total size and number of grid points for the computational regions are 81.5H×12.5H ×25.0H (H: hill height) with a Cartesian grid system and 810 × 120 × 100 in x-, y-, and z-directions, respectively.At the exit of the driver and main regions, a Sommerfeld radiation condition (Gresho, 1992) is imposed.At the top, a free-slip condition for streamwise and spanwise velocity components is imposed and the vertical velocity component is 0. At the side, a periodic condition is imposed and at the flat ground surface, a non-slip condition for each velocity component is imposed.Here, in the flow field, Immersed Boundary Method proposed by Fadlun et al. (2000) is used in order to consider the terrain effects.In the concentration field, a zero gradient is imposed at all the boundaries and a fractional volume term is introduced into the scalar conservation equation to account for the space occupied by terrain.
A 2-dimensional hill with a mean slope angle of 20 • is placed in the main region.The position of the center of the hill is located at 12.5 H downstream of the inlet in the main region.
The origin of the coordinates is the ground surface at the center of the hill.A release point of a tracer gas is located at a distance of 9.1 H from the center of the hill at an elevation of 0.45 H.The Reynolds number based on hill height and free-stream velocity (U ∞ ) at the hill height is almost 5000.

Previous wind tunnel experiments for evaluation of the model performance
Many wind tunnel experimental studies of the dispersion characteristics of a plume over a hilly terrain have been conducted.For example, Sada et al. (1991) investigated turbulence structures and dispersion characteristics of a plume over a 2-dimensional hill.The experiments were conducted under conditions of a neutrally stratified TBL flow in the wind tunnel of Central Research Institute of Electric Power Industry.The test section of the wind tunnel is 20-m long, 3-m wide and 1.5-m high.A TBL flow with strong velocity fluctuations from the ground surface to the upper level was generated using roughness elements with an L-shaped cross sections placed on the floor at the entrance of the wind tunnel section.The 2-dimensional hill was 0.11-m high, 0.6-m long, with a hill aspect ratio (ratio of hill length to height) of about 5.5.A release point was located at a distance of 1.0 m upstream of the center of a hill at a height of 0.05 m.
The thickness of the TBL at the downstream position of the release point is 0.3 m.Using those parameters, the vertical profiles of mean wind velocity, turbulence kinetic energy and mean concentration over a 2-dimensional hill were obtained.
In this study, in order to evaluate the model performance, we compare the LES data of turbulent flow and plume dispersion over a 2-dimensional hill with the wind tunnel experimental data by Sada et al. (1991), described above.

Turbulence characteristics of approaching flow
Figure 2 shows a comparison of LES results with the experimental data obtained by Sada (1991) for vertical profiles of mean wind velocity, turbulence intensities and Reynolds stress in the driver region.Each turbulence statistics obtained by an LES is found to be almost consistent with the experimental data.Therefore, it is considered that the TBL flow is successfully simulated.

Turbulence characteristics over hill
Figures 3 and 4 show a comparison of LES results with the experimental data (Sada, 1991) for vertical profiles of mean wind velocity (U/U ∞ ) and turbulence kinetic energy (k/U 2 ∞ ) at the positions of x/H = −2.7,0.0, 1.8, 3.6, 8.2, 11.8, 14.5 and 18.2.There are slight discrepancies of mean wind velocity between the LES results and the experimental data at x/H = 1.8 and 3.6.The values for the reattachment length behind a hill for the previous experiment (Sada, 1991) and the LES are L/H = 5.9 and 6.9 (L: reattachment length), respectively.Iizuka and Kondo (2006) also performed an LES with the standard Smagorinsky model for a 2-dimensional hill with a hill aspect ratio of 5.0 and compared it with the reattachment length of the experimental data of Ishihara et al. (2001).The results show that the computed reattachment           the wake region.Figure 6 shows a comparison of LES results with the experimental data (Sada, 1991) for streamwise variations of the mean concentration near ground level.The mean concentration (C) is normalized by the free-stream velocity, hill height and source strength (Q).It is found that the basic characteristics such as the increase of mean concentration towards the crest, the rapid decrease behind a hill and the gradual decrease with downwind distance are almost similar to the experimental data.Figure 7 shows a comparison of the LES results with the previous experimental data for mean concentration at x/H = −4.5, −2.7, 0.0, 1.8, 4.5, 8.2, 11.8, 14.5 and 18.2.The mean concentration is normalized by the maximum of mean concentration (C max ) at each downstream position.The peak locations of the mean concentration obtained from the LES are lower than those obtained from the experiment, and the mean concentration values are overestimated in the wake region.However, the basic characteristics of plume dispersion such as the rapid vertical spread of a plume behind a hill and the formation of uniform mean concentration profiles with downwind distance are similar to the experimental data.Figure 8 shows streamwise variation of vertical profiles of r.m.s.concentration (c r.m.s. ).The r.m.s.concentration is normalized by the maximum of the r.m.s.concentration (c r.m.s.max ) at downstream position.
Here, we show the r.m.s.concentration profiles obtained by the LES, since the fluctuating characteristic of concentration was not discussed in the previous experiment (Sada, 1991).We found that its values become small inside the wake region due to the smoothing effects of recirculating flows, compared with the mean concentration values.

Conclusions
In this study, we performed a numerical simulation of unsteady behaviors of turbulent flow and plume dispersion over a 2-dimensional hill immersed in a fully-developed TBL by using an LES and investigated the characteristics of mean and fluctuating concentrations.First, a spatially-developing TBL flow with strong velocity fluctuations in the driver region was generated by incorporating Kataoka's method (2002) into an upstream small fraction of the driver region with a 2-dimensional roughness bar placed at the ground surface.Then, we obtained the inflow turbulence data at the inlet of the main region, and calculated the flow and dispersion over a hill.As compared to the experimental data (Sada, 1991), the main characteristics such as complex behaviors of turbulent flow and plume dispersion behind a hill were successfully simulated.Furthermore, the r.m.s.concentration is found to be smaller than the mean concentration inside the wake region due to the smoothing effects of recirculating flows.

Fig. 1 .Fig. 3 .
Fig. 1.Schematic illustration of computational model.(a)Driver region for generation of a spatially-developing TBL flow and (b)main region for plume dispersion over hill.

Figure 1 .
Figure 1.Schematic illustration of computational model.(a) Driver region for generation of a spatially-developing TBL flow and (b) main region for plume dispersion over hill.

Fig. 1 .Fig. 4 .Figure 2 .
Fig. 1.Schematic illustration of computational model.(a)Driver region for generation of a spatially-developing TBL flow and (b)main region for plume dispersion over hill.

Fig.Figure 3 .
Fig. Streamwise variation of vertical profiles of mean wind velocity.

Fig. 5 .
Fig. 5. Instantaneous plume dispersion field.The yellow areas on the isosurface indicate 0.01% of initial concentration.

Fig. 5 .
Fig. 5. Instantaneous plume dispersion field.The yellow areas on the isosur indicate 0.01% of initial concentration.

4. 3 Fig. 7 .
Figure5the shows instantaneous plume dispersion field at times t * (= tU ∞ /H) = 10.8, 21.6 and 43.2 after the plume release.It shows that the plume moves upward above the hill at first, then a portion of the plume is entrained into the wake region.At t * = 43.2, the plume is entirely dispersed in