Albert F. Kurbatskii

1. Modeling and Simulation of Turbulent Transport of Active and Passive Scalars above Urban Heat Island in Stably Stratified Environment Albert F. Kurbatskii Institute of Theoretical and Applied Mechanics of Russian Academy of Sciences, Siberian Branch Russia, Novosibirsk Shallow Flows Symposium, TU Delft, 2003

2. Introduction For stratified atmospheric flows the LES models and third-order closure models should be considered as fundamental research tools because of their large computer demands. A growing need for detailed simulations of turbulent structures of stably stratified flows motivates the development and verification of computationally less expensive closure models for applied research in order to reduce computational demands to a minimum.

3. Objectives The principal aim of this investigation is the development of turbulent transport model for the simulation of the urban-heat-island structure and pollutant dispersion in the stably stratified environment.

4. Objectives The algebraic modeling techniques can be used in order to obtain for buoyant flows the fully explicit algebraic model for turbulent fluxes of the momentum, heat and mass.  The principal object of this work is the development of three-four-parametric turbulence model minimizes difficulties in simulating the turbulent transport in stably stratified environment and reduces efforts needed for the numerical implementation of model.

5. The mathematical model of the urban heat-island The penetrative turbulent convection is induced by the constant heat flux H0 from the surface of a plate with diameter D. It simulates a prototype of an urban heat island with the low-aspect-ratio plume (zi / D « 1) under near calm conditions and stably stratified atmosphere.

6. The mathematical model of the urban heat-island The mixing height, zi , is defined as a height where the maximum negative difference between the temperature in the center of the plume and the ambient temperature Ta is achieved.

7. The mathematical model of the urban heat-island The problem of the development and evolution turbulent circulation above a heat island is assumed to be axisymmetric.  At the initial moment the medium is at rest and it is stably stratified.

8. Limitations of Laboratory Measurements for Full-scale Simulation There are important limitations utilized in the laboratory experiment and simulation of the real urban heat-island in the nighttime atmosphere: Very large heat fluxes from the heater surfaces Very strong temperature gradients that required to obtain the low aspect ratios (zi/D) and small Froude numbers.

9. Governing Equations Governing equations describing the turbulent circulation above a heat island can be written in the hydrostatic approximation at absence of the Coriolis force and radiation with use a Boussinesq approximation.

10. Governing Equations in RANS-approach

11. Transport Equations for heat and mass fluxes

12. Explicit Algebraic Expressions for Turbulent Fluxes The explicit algebraic models for the turbulent heat flux vector and turbulent mass vector were derived by truncation of the closed transport equations for turbulent fluxes of heat and concentration by assuming weak equilibrium, but retaining all major flux production terms. For turbulent stresses we applied eddy viscosity expression.

13. CLOSURE:full explicit turbulent fluxes models foractive (heat) and passive (mass) scalars

14. Three-Equation ModelE - - 2 The closure of expressions for the turbulent stresses and heat flux vector is achieved by solving the equations for turbulent kinetic energy, its dissipation rate and temperature variance, resulting in three-equation model E--2

15. CLOSURE:three-equation model for active (heat) scalar field

16. Four-equation modelE - - 2 - c For the closure of expression for turbulent flux vector of a passive scalar is involved the equation for covariance concentration – temperature.    Thus, for the description of a concentration field is formulated the four-equation model E--2-c

17. CLOSURE : four-equation model for passive scalar field

18. Boundary Conditions The domain of integration is a cylinder of a given height H. The heated circular disc of diameter D is located at the center of the cylinder bottom.

19. At the plume axis and at its outer boundary symmetric conditions (/r) = (/r) = (/r) (2/r) = 0 areprescribed. (Ur=0 at r = 0 and at r =1.8R) At the top boundary the zero-flux condition V/z = /z = /z = = 2/z = 0 is enforced. Domain of integration is a cylinder HEAT TRANSFER BOUNDARY CONDITIONS z Top boundary Heat Flux, H0 r s o u r c e

20. The surface heat source is placed on the bottom (z = 0) has the size 0  r / D  0.5. Boundary conditions at thebottom are specified as heat flux H0 is prescribed values of E,  and 2 at the first level above surface are chosenaccording to Kurbatskii (JAM, 2001, vol.40, No.10) Domain of integration is a cylinder HEAT TRANSFER BOUNDARY CONDITIONS Z Top boundary Heat Flux, H0 0 r s o u r c e

21. Atthe plume axis and at outer boundary, (C/r) = (c/r) = 0. At the top, Constantflux of mass, is prescribedinside a source. At the bottomand outside of a source Domain of integration is a cylinder MASS TRANSFER BOUNDARY CONDITIONS Z Top boundary L = 0.5 D r mass source

22. MASS TRANSFER BOUNDARY CONDITIONS Z Domain of integration is a cylinder Top boundary The same boundary conditions are used for source of small length located at the center of a heat island. L=0.1D r mass source

23. MASS TRANSFER BOUNDARY CONDITIONS Z Domain of integration is a cylinder Top boundary and the same boundary conditions are used for source of small length located at the periphery of a heat island. L= 0.1D r mass source

24. Numerical Method Fr , Fz – turbulent fluxes of momentum, heat and mass Semi-implicit alternating direction scheme

25. The numerical method uses a staggered mesh. The difference equations are solved by the three-diagonal-matrix algorithm. Staggered mesh          Numerical Procedure z r r/2 z z/2 r 0 Ur Uz E, , T, <2>, C, <c>

26. Main Results of Simulation The results of simulation correspond to a quasi-steady state of circulation over an area heat source in stable stratified environment. Figure (c):shadowgraph picture at t = 240 sec when the full circulation is established.

28. Calculation of Normal Turbulent Stresses In this problem a simple gradient transport model preserves certain anisotropy of the normal turbulent stresses is turbulent viscosity.

31. RESULTS: Temperature profiles Calculated temperature profiles inside the plume have characteristic “swelling”: the temperature inside the plume is lower than the temperature outside at the same height creating an area of negative buoyancy due to the overshooting of the plume at the center. This behavior indicates that the plume has a dome-shaped upper part in the form of a “hat”.  

34. Pollutant Dispersion Investigation of buoyancy effects on distribution of mean concentration in mixing and inversion layers of urban heat island was the main goal in modeling and simulation of pollutant dispersion from a continuous surface source.

35. Pollutant Dispersion Experimental measurements were not available for the quantitative validation of simulation data. Instead, we present some preliminary results that illustrate interesting properties of pollutant dispersion from a continuous source located inside the urban heat island and on its periphery.

39. Pollutant Dispersion One characteristic feature is observed in all three cases. The contaminant penetrates not only into the inversion layer but even higher beyond its boundary. This behavior was recently reproduced in the laboratory experiment of Snyder et al. (BLM,2002, vol.102, 335-413.). These experimental data clearly show penetration of the continuous buoyant plumes into inversion above the convective boundary layer.

40. Difference in Turbulent Diffusion Between Active and Passive Scalars 1 Sometimes assumed that in the stratified atmospheric boundary layer the eddy diffusivity of heat (KH) is equal to the eddy diffusivity of contaminant (KC). However, the stratificationcauses a larger difference in the eddy diffusivities between active heat and passive mass. Indeed, for the ratio of the vertical eddy diffusivities of heat and mass can be written the following expression 

41. Difference in Turbulent Diffusion Between Active and Passive Scalars 2 [1-D1bg <cq>(C/z)-1] [1-D2bg <q2>(Q/z)-1] If both mass and heat are passive additives (the buoyancy terms in this expression are negligible) then it is evident from this expression that KC=KH.  It appears that cases for which the largest deflection of KC/ KH from unit will occur are when either T or C is acting as a passive additive.  In our case the mass is acting as a passive additive.

42. Turbulent Fluxes of Active and Passive Scalars

43. Ratio of Eddy Diffusivities of Passive Mass to that Active Heat

44. CONCLUSIONS 1 The three-equation model of turbulent transport of heat reproduces structural features of the penetrative turbulent convection over the heat island in a stably stratified environment. This model minimizes difficulties in describing the non-homogeneous turbulence in a stably stratified environment and reduces computational resources required for the numerical simulation.

45. CONCLUSIONS 2 The four-equation model for the description of pollutant dispersion in the stable stratified atmospheric boundary layer is formulated. Favorable comparison the numerical results of pollution dispersion from the continuous surface source above the urban heat island with laboratory measurements in the convective boundary layer showing penetration of the continuous buoyant plumes into inversion above the convective boundary layer is found.

46. The friction velocity u(r)/ wD The friction velocity on the underlying surface can be obtained on the calculatedmean velocity asu (r) =  ( Ur / z ).

47. Turbulent Velocity Scale Turbulent velocity scale Uf was estimated as 1/30 of “a mean wind velocity ” – velocity scale wD of the mean inflow velocity: Uf 1/30wD. This value was used as characteristic scale of the turbulent velocity for boundary conditions for E1 and e1 at the first level of a grid above an underlying surface.

48. Numerical Procedure It took about 2.8·104 time steps to drive the numerical solution to a quasi-steady state.  Computations were performed on a mesh with 25 (and 50) points in radial direction.  In vertical direction 116 (and 232) mesh points were used.  The time step was chosen so that the numerical accuracy was preserved.