J2.2 DEVELOPMENT AND APPLICATIONS OF CFD SIMULATIONS IN SUPPORT OF AIR QUALITY STUDIES INVOLVING BUILDINGS EPA/600/A-04/078 *1Alan Huber Atmospheric Sciences Modeling Division, ARL/NOAA, RTP, NC 1 on assignment to National Exposure Research Laboratory, US EPA, RTP, NC Wei Tang National Research Council Research Associate at National Exposure Research Laboratory, US EPA, RTP, NC Anita Flowe Independent Contractor, Roxboro, NC and Brian Bell, Karl Kuehlert, and Walter Schwarz Fluent Inc, Lebanon, NH 1. INTRODUCTION There is a need to properly develop the application of Computational Fluid Dynamics (CFD) methods in support of air quality studies involving pollution sources near buildings at industrial sites. CFD models are emerging as a promising technology for such assessments, in part due to the advancing power of computational hardware and software. CFD simulations have the potential to yield more accurate solutions than other methodologies because they are a solution of the fundamental physics equations and includes the effects of detailed three-dimensional geometry and local environmental conditions. However, the tools are not well validated for environmental flows and best-practice methodologies have not been established. Fluent, Inc and the US EPA National Exposure Research Laboratory are working cooperatively to demonstrate CFD model simulations as a proven and applied tool in support of environmental assessment studies. See also Huber et al 2000a, 2000b, and 2001 for additional perspectives related to this project. The results of CFD simulations can both be directly used to better understand specific case studies as well as be used to support the development of better-simplified algorithms that may be generally applied. Unlike most currently used regulatory air quality models, CFD simulations are able to include specific details of building structures as well as a range of physical processes that affect atmospheric turbulent boundary layers. Dispersion in absence of buildings is demonstrated in this paper to be comparable with standard plume dispersion models for point and line source pollutant emissions. Boundary layer turbulence is being simulated as characterized by surface roughness (characterized by * Corresponding author address: Alan H Huber, NOAA/ERL/ARL/AMD, Mail-Code E243-03, US EPA NERL, RTP, NC 27711 e-mail: huber.alan@epa.com u*) and surface heat flux (characterized by the Obukhov length L). This paper discusses ongoing development and application of CFD simulations through case studies using CFD software for simulating air pollutant concentrations from sources near buildings. Comparisons of CFD simulations to reference wind tunnel data and field measurement studies are being studied to provide model evaluation/validation. CFD simulations should be shown to be comparable with simple proven air dispersion models being reliably applied today in routine air quality studies. This is critical to demonstrate that the complex numerical techniques that are part of CFD software are well behaved under simple conditions. We do not need CFD software to support studies where simple analytical solutions are possible. We do want to extend CFD applications for complex conditions where we know simple analytical solutions are not appropriate. While we have already explored many of the basic elements of CFD software there is much ongoing that needs to be completed before we make recommendations. An overview on progress with our evaluation and development of CFD applications appropriate in supporting of air quality studies involving buildings is presented herein. 2. CFD SOFTWARE A brief overview of the numerical methods is provided here. The discussion is meant only to present an introduction for someone that may be new to CFD software. CFD software involves many layers of coding with complex interactions. For this reason any CFD software should be carefully examined and have a history of quality assurance testing before one begins to apply it to support air quality studies. Those interested in additional introductory reading on CFD issues can find many good reference books (for example: Ferziger and Peric, 1997; Wesseling, 2000; and Wilcox, 1998) ------- The FLUENT software (Fluent Inc, 2003) solves the governing equations for the conservation of mass, momentum, energy, and scalars such as a pollutant. The study domain is divided into discrete control volume cells using a computational grid mesh. Unstructured meshing supports variable volume cell sizes throughout the domain. This allows for better computational efficiencies by being able to concentrate the grid mesh in areas where finer mesh is most critical in resolving complex flows. Algebraic equations for discrete dependent variables such as velocities and pollutants are constructed and solved. There are options for both a coupled equation solver using either an implicit and explicit discretization, or a segregated equation solver having implicit discretization. For atmospheric flows the segregated solver using implicit discretization is appropriate and is being used for our studies. The momentum equations are solved, and then a pressure-correction is applied to update the pressure field to support calculation of mass fluxes to ensure conservation of mass. The solutions for energy, turbulence and other scalar equations (i.e., pollutants) follow separately. In the implicit discretization for a given variable the unknown value in each cell represented at the cell center is calculated using both existing and unknown values from neighboring cells. Overall the software uses an algebraic multigrid method to solve the resultant system of equations for the dependant variable in each cell. The calculations continue and update all the cell properties until selected criteria for a converged solution is reached. There are options for obtaining volume face values by applying first-order, second-order, power-law, and for quadrilateral/ hexahedral grid mesh the QUICK (Quadratic Upstream Interpolation for Convective Kinematics) scheme. There are specific options for pressure interpolation including linear, second-order, body- force-weighted, and PRESTO (PREssure Staggering Option). For pressure-velocity coupling the options are SIMPLE (Semi-Implicit Method for Pressure- Linked Equations), SIMPLEC, and PISO (Pressure- Implicit with Splitting of Operators). We have not noticed a significant effect among these different choices for our simulations to date. The software has options for either steady or unsteady (time-varying) solutions. There are options for a first order and higher order implicit schemes for temporal discretization of the time derivative. To date we have not been examining unsteady flow solutions. We have started with the simplest applicable CFD models for supporting air quality studies involving buildings. We have been evaluating solutions for the Reynolds-Averaged Navier-Stokes (RANS) governing equations for momentum. Solutions require a selection of boundary conditions and a model for turbulence. The software has options for the wall (ground surface) boundary conditions and several turbulence models. We have been evaluating the performance of standard k-e (turbulent kinetic energy: k; turbulent energy dissipation rate: e, epsilon) turbulence model. This is our base case. In the future we plan to examine higher order turbulence closure models including Reynolds Stress Models (RSM) and Large Eddy Simulation (LES) along within the framework of unsteady solutions. The computational requirements of these higher order solutions may not be practical for support of routine air quality studies but may be useful for special cases studies where detailed analyses of pollutant dispersion is important in human exposures. Also, higher order simulations may support the development of reliable simplified models of human exposures to pollution. 3. AMOSPHERHIC BOUNDARY LAYER AND PLUME DISPERSION Simulation of the atmospheric boundary layer is critical to modeling plume dispersion. While the primary interest in application of CFD methods is to simulate flow around buildings the CFD code should first be demonstrated to correctly model a plume in absence of any buildings. The CFD simulated flow should be simple and well defined in absence of building influences. If there are problems within the CFD code they can be more easily identified. Flat plate and atmospheric boundary layer theory provides a basis for testing the sensitivity of CFD code parameters over a range of boundary conditions. Monin-Obukhov similarity theory is applicable to atmospheric boundary layers. For CFD simulations a surface heat flux (H, W/m2) is fixed as the bottom boundary condition to simulate non-neutral stability conditions. The CFD boundary layer flow is set up by using a finely resolved grid with application of the "law of the wall" near the bottom. The surface friction velocity (u*) is estimated from the resulting wind profile. Figure 1 presents example profiles of mean streamwise velocity and temperature with and without added heat flux. Figure 2 presents a summary of simulated Obukhov length (L, m) versus surface friction velocity that result from a range of simulations. These results are found to compare well with Monin- Obukhov theory (see Figure 11.1, Arya 2001). Simulation of pollutant plume dispersion requires good models for both the bulk transport and the turbulent dispersion of the pollutant. Good simulation of the bulk transport is expected if the mean flow field is correctly modeled. Good turbulent dispersion requires that the turbulent flow be correctly modeled. For this study we are only evaluating RANS simulations and using several k-e turbulence models, which produces turbulent kinetic energy (TKE) driving turbulent pollutant dispersion. Dispersion from line and point sources are being studied to evaluate the performance of the turbulence models. The standard turbulence model has been generally working well using standard code default parameters for simple ------- 14 12 10 U) E 8 > o 6 o 4 2 0 -velocity, heat flux = 200W/mA2 -velocity, no heat flux ^temperature, heat flux = 200W/mA2 292 291 290 289 288 3 TO c o _l > o 3 •fi O 10 1 —400W/mA2 —300W/mA2 -*-200W/mA2 100W/mA2 —40W/mA2 t / 0.01 10 m/s 0.1 1 Friction Velocity, u* Figure 2. Monin-Obukhov theory applied to a range of case studies. flows. Some additional evaluations and refinements are ongoing to improve performance, especially for highly buoyant cases. No work has yet been started to evaluate strongly stable stratified flow. Figure 3 presents a comparison with a line source assuming P-G Gaussian plume urban dispersion parameters (Stability C). The concentrations are normalized by the wind speed at 10 m (u, m/s) and the source strength (q, gm/s). Matching both the full lateral and full vertical profiles of the plume beyond the centerline present challenges which are being examined in more detail than can be covered in this presentation. The CFD simulations include effects of wind shear and characteristics of the turbulence model that must be more carefully evaluated. While plume centerline peak concentrations are of primary interest in air quality studies, when the study includes a series of ranging wind directions the location of overall peak concentrations can be altered significantly by off- centerline concentrations. Also, peak concentrations near the ground can be significantly altered by off centerline plume concentrations for elevated sources. 1.0E+00 - heat flux = 200W/mA2, U100m = 3.0m/s Ri = -0.94 heat flux = 400W/mA2, U100m = 7.9m/s Ri = -0.19 urban line source — 1.0E-01 o 1-0E-02 w 1.0E-03 E 1.0E-04 1.0E+01 1.0E+02 1.0E+03 Distance from source, m 1.0E+04 Figure 3. Comparison of normalized concentrations (Cu/q, m"2) with urban P-G line source. 4. PROJECT PRAIRIE GRASS SIMULATION Field measurements in the atmospheric boundary layer contain inherent variability due to unsteady winds that are case specific. This factor is not fully captured by the simple Gaussian straight-line plume formulations. Development of Gaussian straight-line plume models in regulatory use has been formed, in part, on an early field study, Project Prairie Grass (Barad, 1958). Field measurements from this study are being used to evaluate methods for direct application of CFD simulations of specific cases for simple atmospheric flows. For these case studies there are ranges of wind directions that are not part of the steady-state RANS CFD simulations. We are working with these field measurements to evaluate methods for the best CFD simulations of all measurements for each case. Methods include accounting for variation in wind direction by smoothing the steady solution over the wind distribution or by enhancing the lateral dispersion internally. Figure 4 presents an example simulation for one case having minimal wind variation (Case 55, Barad 1958). In the figure for each of 4 distances downstream from the near ground source the CFD simulation is compared with the field measurements, routine P-G (stability D) straight-line plume estimates, and the CFD solution weighted by ------- a) Arc distance = 100 m 5 4 !¦ — CFD-weighted • measurements — CFD plume model • • 1 L ¦»1 W*. 320 330 340 Direction, degrees b) Arc distance = 200 m i CFD-weighted • measurements CFD —» plume model • 1 • I • A A IV -A v f/T N • 320 330 340 Direction, degrees c) Arc distance = 400 m 5 4 — CFD-weighted • measurements — CFD plume model J 1 * I • | If ¦J A- -y\ lv\ yV. 310 320 330 340 Direction, degrees d) Arc distance = 800 m 4.5 4 3.5 i 3 £ S 25 — CFD-weighted • measurements CFD .... plume model ml m • f • Jr A ft JJr -.v\ ¦ .rii 330 340 Direction, degrees Figure 4, Example Prairie Grass case. the distribution in wind direction as a smoothing function. The CFD simulation was matched to the measured vertical profile of wind speed from 2 m to 16 m with a friction velocity u*=0.44 m/s and roughness height zO = 0.009m. Similar results are being observed for other cases. It appears that the measurements lie between the default steady-state RANS simulation and these simulations smoothed over a function of the wind direction. Best methods will be developed based on the whole database. While the Project Prairie Grass is an especially good database for examining the horizontal plume it has only a few vertical plume profiles. Additional field measurements including more vertical profiles are desirable. There are few vertical profiles because they are naturally more challenging to collect. Additional databases for better examining the vertical plume profiles are being searched. The preliminary CFD simulations of flow and dispersion for an atmosphere-like boundary layer has been determined sufficient to begin applications to areas with buildings. 5. BUILDING SIMULATIONS Fortunately there are many databases on flow near buildings from scaled physical model studies in wind and water tunnels. The boundary layers are simple without many of the chaotic and complicating factors in actual field situations. Simple idealized buildings can be systematically studied. These data are ideal for evaluating the performance of CFD models because boundary conditions are well controlled. CFD models should be demonstrated to simulate the scaled model studies before moving forward to full-scale field situations. These simulations should help identify potential errors in model coding or identify limitations of physical models. Data from studies conducted in the US EPA's Meteorological wind tunnel are being used to initially evaluate CFD code for our project. This EPA wind tunnel study (Lawson et al, 2000) was conducted in collaboration with the Los Alamos National Laboratory and the Lawrence Livermore National Laboratory. The study included measurements of velocity and tracer concentrations within arrays of two-dimensional (2D) and three-dimensional (3D) buildings. These data are being used by others to evaluate other model performance (Brown et al.. 2000; Chan et al., 2000; Kastner-Klein et al, 2000). The buildings for the 2D study were set with seven square cross sections spanning the width of the wind tunnel. The separation between the buildings was equal to the cross-section building scale (0.15 m). The buildings for the 3D study covered the same area as that for the 2D study but consisted of cubes separated by open space of the same volume as each cube. ------- For these studies the wind tunnel ceiling was adjusted to reduce the horizontal pressure gradient. Three triangular fins (spires) at the inlet and 0.19 mm blocks on the floor were used to develop a 165 m deep wind tunnel boundary layer. The roughness blocks are absent on the floor in the study area with the model buildings. The CFD simulations begin by developing a set up trying to match the wind tunnel boundary layer at the inlet. The wind tunnel measurements of mean velocity can be matched at the inlet well as long as the velocity profile is set correctly. Matching the TKE field requires a good estimate of the dissipation rate, for which there are no direct measurements. Figures 5 and 6 present comparisons of the measured (green) with the CFD (black) velocity vectors near the leading two rows of buildings. The coordinate origin (x=0, y=o, z=0) is located at the base of the leading building (cross-stream center of the leading building. The region in front of the leading building and over its roof has the greatest difference in the flow field between 2D versus 3D and is the most challenging to CFD simulation because of the small regions of recirculation flow. The CFD simulations resulted from using standard default set up. The grid resolution was with 30 cells per building face. Studies of grid resolution demonstrated grid independence at this fine scale. The measured flows were very similar and well matched for the flow over the roof and street canyons between rows 3 to last row 7. Figure 7 presents comparisons for a vertical profile in the first building street canyon (between building 1 and 2). Included are simulation profiles at the same location without the buildings. The CFD simulations resulted from default set up parameters. The mean velocities both with and without buildings are well matched with the measurements. The TKE without buildings is well matched except for the near wall zone. The TKE profiles with buildings are similar to measures but too low in a near wall zone and too high in the upper zone. The wind tunnel study with the 3D buildings includes an examination of tracer dispersion from a source placed at the leeward base of the first building (x=0.15, y=0). Figures 8 and 9 present example comparisons between the measurements and the CFD simulations using the standard default set up. Concentrations gradients in the first (source containing) street canyon are greatest and are well matched by the CFD simulation. The comparisons are likewise good in the second street canyon (as well as the other street canyons not shown) where the gradients are greatly reduced due to the uniform mixing imposed by the flow through the streets. These comparisons demonstrate that with good simulations of the mean flow even without matching TKE, the transport and dispersion of tracer concentration in street canyons can be matched. Profiles of velocity and concentrations are more complex and have more complex gradients in the along-stream street canyons. Further examinations are ongoing for more complex building street canyon studies. Modifications to the standard default CFD set up are being examined to identify how best to improve the simulations of TKE for this study. These include assessing performance of different turbulence models, boundary conditions (especially inlet and top), surface wall models, and grid resolution. Some preliminary comparisons are presented in Figures 10 and 11 near the leading 2-D building where differences have been most noticeable. Figure 10 shows how blockage effects that were observed for the 2-D case study can significantly affect the flow in front of the leading 2-D building. The roof in the wind tunnel is adjusted to minimize horizontal pressure gradients in the free-stream flow. Therefore having the correct ceiling boundary condition is critical. Also, further improvements in the CFD simulations are possible when dissipation rate (e - epsilon) is better estimated. Refined simulations of the wind tunnel boundary layer development are providing significant improvement as presented in Figure 10. In Figure 11 the present simulations show that TKE over the roof of the leading 2-D building is not significantly affected by the ceiling boundary conditions but is significantly affected by inlet dissipation rate, Black=FLUENT Green=Wind tunnel Reference, m/s J— CD "S £ c *-» .c o> "a> I 0.3-r 0.2- 0.1- Iff! fi! if if! $3 J I -0.1 0 0.1 0.2 0.3 0.4 0.5 Streamwise distance in meters Figure 5. Velocity near two-dimensional buildings. 0.6 Black=FLUENT Green=Wind tunnel Reference, m/s o" 3H 0 0.1 0.2 0.3 0.4 0.5 Streamwise distance in meters Figure 6. Velocity near three-dimensional buildings. ------- a) Mean Velocity (m/s) 2 — — — — FLUENT 2-0 flat plate FLUENT 2-D with buildings FLUENT 3-0 flat plate FLUENT 3-0 with buildings • # • Experimental 3-D with buildings • + • Experimental 2-0 with buildings 9 • # Experimental flat plate "i ¦—i 1 1 1 1 1 1 -2 0 2 4 6 Streamwise velocity in meters 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0.55 O.f Streamwise distance in meters Figure 8. Concentrations within 3-D building canyons at y=0 (CFD color contour, Wind tunnel data numbers). III; Wind tunnel x=O S25 FLUENT x=0.525 Wind tunnel x=0.225 FLUENT x =0.225 -0-2 O 0.2 Cross-stream, meters c) Turbulent Dissipation (m Is — FLUENT 2-0 flat plate FLUENT 2-0 with buildings FLUENT 3-0 flat plate FLUENT 3-0 with buildings —I 1 I II II11' I ¦ I ¦ I l 11| 1 1 I l 1 II l| 0.01 0,1 1 10 Turbulent dissipation mA3/sA3 Figure 9. Cross-stream concentration profiles with 1s building canyon (x=0.225) and 2"" building canyon (x=0.525). s < £ UJ Figure 7. Profiles within the first building canyon. Figure 10. Turbulent kinetic energy in front of leading 2-D building face (x=~0.038 m, y=0). b) Turbulent Kinetic Energy (m Is 2.5 -i FLUENT 2-D flat plate FLUENT 2-D with bu i Idings FLUENT 3-D flat ptale FLUENT 3-D with buildings • • $ Experimental 3-D with buildings • 9 9 Experimental 2-D with buildings 9 • Experimental flat plate r 0.4 0.8 TKE in m*2/sA2 Ceiling BC, estimated dissipation Ceiling BC. spires simulation 9 # t Wind tunnel data - Symmetry BC, estimated dissipation 0.1 0.2 0.3 0.4 Vertical position in meters ------- 2.5 -i 2 - (M JS 1.5 (M < £ LU * 1 " 0.5 - —T 1 1 1 1 T" 0.1 0.2 0.3 0.4 Vertical distance in meters 0.5 Figure 11. Turbulent kinetic energy on leading 2-D building roof edge (x=0.045 m, y=0). 6. OVERVIEW Much is being learned about how best to set up CFD simulations to support environmental simulations and the issues that most affect comparability with both physical model studies and field measurement studies. The choice of boundary conditions, grid resolution and structure, and turbulence models affect the outcome of a solution significantly. Transport and dispersion can be well simulated for flat plate boundary layers as used in physical model studies. Transport and dispersion simulations are more complicated for atmospheric flows due to the complex temporal-spatial wind fluctuations. To date the project has focused on RANS steady- state solutions and the standard k-e turbulence models. This is being extended to include unsteady solutions and higher order turbulence models. Detailed technical papers will be prepared as this project reaches significant conclusions. REFERENCES: Arya, S. Pal, 2001: Introduction to micrometeorology. 2nd Edition. Academic Press, San Diego CA. Barad, M.L. (Editor), 1958: Project Prairie Grass, A Field Program in Diffusion. Geophysical Research Paper, No. 59, Vol I and II, Report AFCRF-TR-58- 235, Air Force Cambridge Research Center, Bedford, MA. Brown, M., R. Lawson, D. Descroix, and R. Lee, 2000: Mean flow and turbulence measurements around an array of buildings in a wind tunnel. 11th AMS Conference on Applications of Air Pollution Meteorology. January, Long Beach, CA. Chan, S.T., D. Stevens, and R. Lee, 2000: A model for flow and dispersion around buildings and its validation using laboratory measurements. Proc. 3rd AMS Symposium on the Urban Environment, 14-18 August, Davis, CA., 56-57. Ferziger, J. and M. Peric 1997: Computational Methods for Fluid Dynamics. ISBN 3540594345. Springer-Verlag, New York. 364 p Fluent, Inc 2003: FLUENT 6.1 User's Guide. Fluent Inc, Lebanon, NH. Huber, A., S. Rida, E. Bish, and K. Kuehlert, 2000a: Addressing Environmental Engineering Challenges with Computational Fluid Dynamics. Proceedings of the 93rd Annual Meeting of the Air & Waste Management Association, June 18-22, 2000, Salt Lake City, UT. Air & Waste Management Association, Pittsburgh, Paper No. 1073. Huber, A., Bolstad, M. Freeman, S. Rida, E. Bish, and K. Kuehlert, 2000b: Addressing human exposures to air pollutants around buildings in urban areas with computational fluid dynamics (CFD) models, 3rd Symposium on the Urban Environment, 14-18 August, Davis, CA., 62-63. Huber, A, M. Freeman, S. Rida, K. Kuehlert, and E. Bish, 2001: Development and Applications of CFD in Support of Air Quality Studies of Roadway and Building Microenvironments. Proceedings of the 94th Annual Meeting of the Air & Waste Management Association, Orlando, Florida, June 24-28,2001. (CD-ROM). Air & Waste Management Association, Pittsburgh, Paper No. 1035. Kastner-Klein, P., M. Rotach, M. Brown, E. Fedorovich, and R. Lawson, 2000: Spatial variability of mean flow and turbulence fields in street canyons. Proc. 3rd AMS Symposium on the Urban Environment, 14-18 August, Davis, CA, 13- 14. Lawson, R., S. Perry, and R. Thompson 2000: Measurement of Velocity and Concentration Fields in Arrays of 2-dimensional and 3-dimensional Buildings in a Simulated Neutrally-Buoyant Atmospheric Boundary Layer, U.S. EPA, RTP, NC, 55 pp. Wesseling, P., 2000: Principals of Computational Fluid Dynamics. ISBN 3540678530. Springer- Verlag, New York. 642 p. Wilcox, D., 1998: Turbulence Modeling for CFD. ISBN 0963605151. DCW Industries, La Canada, CA. 540 p. ACKOWLEGMENT: Appreciation is extended to Dr Steven Perry and Dr David Heist of the US EPA Fluid Modeling Facility for their advice and assistance in understanding the wind tunnel model studies. Ceiling BC, estimated dissipation - Ceiling BC, spires dissipation 9 # • Wind tunnel data Symmetry BC, estimated dissipation ------- Disclaimer: The research presented here was performed under the Memorandum of Understanding between the U.S. Environmental Protection Agency (EPA) and the U.S. Department of Commerce's National Oceanic and Atmospheric Administration (NOAA) and under agreement number DW13921548. Although it has been reviewed by EPA and NOAA and approved for publication, it does not necessarily reflect their policies or views. ------- |