CFD Study of Thermal Comfort in Urban Area

The reduction of energy consumption for cooling should be addressed on a city level since the urban environment where the building operates increases its cooling needs and deteriorates the heat pumps operation. In the present work, a Computational Fluid Dynamics (CFD) model is used to study the improvement of microclimate conditions in an urban agglomeration of the city of Volos, through bioclimatic interventions. The Boussinesq approximation is used in order to take into account the thermal buoyancy while the Discrete Ordinate (DO) model is used for the radiation transport. The ground temperature is calculated by solving an energy balance model while the ground is also considered an infrared diffusive radiation emitter. The planted surfaces are considered finite thickness heat sinks and the roofs and building walls are also took as finite thickness isothermal walls and heat sources (the heat pumps and electrical device operation, lighting and, residences presence). Trees are regarded as porous volumes. The surfaces of water elements are considered isothermal heat sink walls. The local microclimate before and after the bioclimatic interventions is presented by comparing the temperature and wind speed distributions as well as the thermal comfort indices for a typical summer day.


Introduction
In the past, energy in buildings was consumed mainly for heating but today is also consumed to a high degree for cooling. This increase is due to: a) the development and widespread use of the cooling technology, b) the deterioration of the urban environment. Dense construction and lack of greenery are responsible for the development of the urban heat island. The reduction in energy consumption for heating can be dealt with on a building level, but the case of energy reduction for cooling should be dealt with on a city or neighborhood level. And this is because the urban environment, where a building operates, burdens considerably its cooling needs and the operation of the used heat pumps. The study of the microclimate developed in the urban core is a multivariate problem since the formation of the final pressure, air speed, temperature and humidity field depend upon factor such as: a) local climatic conditions, b) geometry and orientation of buildings, streets and openings, c) used materials in buildings, streets and public spaces, d) existence of greenery and water elements, e) the produced heat in buildings due to residents presence, device operation, lights and heat pumps etc.
One way to describe this complex problem is through field measurements. However those measurements give information about specific points and require a large number of assumptions to be adopted in order to draw conclusions for a whole region. Nevertheless this does not reduce their value and the information they provide. Furthermore, those measurements can be used for validation of produced numerical models for the description of the developed microclimate [1,2]. Measurement techniques are also used by numerical models for the preparation of the urban space computational grid [3].
The methods of Computational Fluid Dynamics (CFD), is another modern way to simulate the developed transport phenomena and to determine the urban microclimate behavior. Those methods can offer a full depiction of the developed microclimate in an area and even consider improvement scenarios by using appropriate materials, adding water elements and green areas, shaping public spaces etc. In this context, the flow field developed in an urban core was studied with CFD without taking into account the energy but focusing on the study of wind and turbulence modelling in paper [2]. A combination of measurements, CFD calculations for the flow aerodynamic and TRNSYS software for the thermodynamics are used in [4] in order to study bioclimatic measures aiming in the rehabilitation of specific open spaces in Tirana during winter and summer. The concept of different techniques coupling is also used in [5] where a CFD model is used for convective heat, air and vapor transfer with a Building Envelope Heat and Moisture (BE-HAM) model and a radiation model (RAD) for the study 9 of the microclimate of urban heat island effect in an urban street canyon. In the paper [6] also a CFD model was developed for the description of the urban microclimate considering the heat transfer through an energy balance model (EBM) for the definition of the boundary conditions but without direct simulation of the radiative transport equation (RTE). In another developed model [7], the buildings are described as sources of heat and radiation estimating the contribution of urban planting to the reduction of the air temperature in the urban environment. Another CFD model was also developed for the greenery and water elements simulation as sinks in the momentum and energy equations [8]. Since a number of individual models describing the effect of several factors on the microclimate were developed [9], CFD codes that combine some of these models were used for the study of real urban segments (eg Athens [10], Florina [11]) in order to study methods for mitigation of the urban heat island effect.
In the present work a CFD model was developed for the study of transport phenomena in the urban core of the city of Volos. An EBM was used for the description of buildings, streets, greenery and water elements as sources or sinks of heat. The greenery is simulated as porous material, while the radiation equation is solved directly for two wave length bands (solar radiation with λ=0 -0.76 μm and thermal radiation with λ=0.76 -100 μm). This model was used for the description of the current microclimate in the studied area and for the prediction of the microclimate evolution after an integrated intervention during a typical summer day when the urban heat island effect is most intensive.

Mathematical Model
The flow field and the transport phenomena developed in the computational domain are described by the Reynolds Averaged Navier-Stokes equations [12,13]. The continuity, momentum and energy conservation equations are solved in 3D (x,y,z), as well as the equations that correspond to the turbulent model. The flow is considered 3D, steady state, incompressible and turbulent. The turbulence effect in the developed flow and heat transfer is modeled through the standard k-ε high Reynolds model [14]. The effect of thermal buoyancy is modeled using the Boussinesq approximation. The radiation is modeled using the Discrete Ordinate (DO) model [15,16]. The radiation transfer equations are solved in two wave length bands, one corresponding to the solar radiation (short wave with λ=0-0.76 μm) and one corresponding to the thermal radiation (long wave λ=0.76 -100 μm). The greenery is considered porous materials that cause pressure drop and further more a heat sink. In the momentum conservation equations two source terms are added: one corresponds to thermal buoyancy and the other corresponds to the pressure drop inside the porous media.
Where, S Ui is the sink term in the i-direction of the momentum conservation equation, (C D =0.2) the drag coefficient, (LAD=7) the leaf area density, u i the local speed in the i-direction and, ͞ u the magnitude of local speed.
In the energy equation are added: a) a source term corresponding to the thermal radiation calculated by the DO model, b) a source term corresponding to the heat produced by buildings, c) a sink term that corresponds to the air temperature reduction due to evapotranspiration of greenery and d) a sink term that corresponds to temperature reduction due to evaporation of water elements.
The heat is produced by buildings and comes from: a) the resident's presence, b) the operation of the domestic equipment, c) the heat produced by lights operation d) the heat rejected by heat pumps.
Taking into account the above an average heat production is considered to be 28 W/m 2 . According to FAO instructions [18], the average evapotranspiration rate (ETc) for trees in the Mediterranean area with low humidity is 6-9 mm/day for temperatures higher than 30o C. If an average value is taken into account this leads to a heat consumption of 212 W/m 2 . The final heat consumption depends on the trees' leaf area index (LAI). This is considered to have an average value of 3.
The water element surface is considered an isothermal wall with a temperature of 25 o C, as well as a heat sink due to water evaporation. For a typical summer day with an average temperature of 33 o C and a relative humidity of 47% this leads to a heat sink of 200 W/m 2 . Finally to the turbulence models' equations are added the following source terms due to greenery elements [17].
Where, S k the sink term of the turbulence kinetic energy equation and S ε the sink term of the turbulent dissipation.

Numerical Model
The finite volume method is used for the resolve of the 3D Navier Stokes partial differential equations [13] using the OpenFOAM software. A hybrid multiblock grid of 6538192 cells was used, created by 22 blocks. The discretized continuity and momentum equations are coupled using the SIMPLE algorithm [12]. Due to the fact that in the solid boundaries wall functions are used, the grid is formed securing that the nondimensional distance y + is between 10 and 40. The diffusion terms are discretized using central differences while for the discretization of the convection terms 2 nd order upwind scheme are used. The absolute convergence criteria for the momentum, pressure and concentration equations are set to e=10 -5 and for the energy and radiation is set to e=10 -6 .

Geometry
In Figure 1, the geometry of the simulation domain is given. In order to avoid side effect in the simulation results, the computational domain boundaries are expanded by 8H from the interest area and 6H vertically over the buildings surface (H, is the buildings height) [14,17]. This way all the lateral sides can be used either as inlet or outlet surfaces given that during a day the wind direction varies.

Boundary Conditions
Depending on the time of the day, 3 of the lateral sides are considered inlet boundary conditions (known values of air velocity components, u i , turbulent kinetic energy, k, dissipation rate, ε, and air temperature T).The same boundary condition is set to the upper boundary. On the surface from which the air is expected to leave the simulation domain, a pressure outlet condition is set (known pressure equal to the atmospheric).
The streets, the sidewalks, the school courtyard as well as the under trees surfaces are modeled as isothermal walls. This corresponds to no-slip and no penetration boundaries for momentum equations while as far as the energy equation is concerned a predefined value of temperature is set. Concerning the radiation transport equation those surfaces are considered to emit thermal diffusive radiation according to their temperature, material and color.
The low planting surfaces are regarded as finite thickness walls and heat sinks due to evapotranspiration. The building's roofs and side walls are also considered isothermal walls with finite thickness. Inside this finite thickness it is assumed that the heat due to residents' activation, devices, lights and heat pumps operation is produced. The trees are considered porous volumes working as source and sink terms in the momentum, turbulence and energy equations. Finally, the water element surface is also considered to be an isothermal wall and heat sink. The wind speed in the boundary inlet is assumed to follow an exponential distribution profile up to 11 m and above that a logarithmic one. So the wind speed components are given according to the following relationships [19].
The values of turbulent kinetic energy, k, and of the dissipation rate, ε, are calculated for the inlet boundary conditions with the following relationships [17,8].
The temperature of the ground, which acts as an isothermal wall, for every considered cover material is calculated for each simulated hour using an energy balance model [6,8] taking into consideration the incident solar radiation, the air temperature, the underground temperature (in a depth of 15 cm), the thermal conductivity of the material, the absorption coefficient in the solar radiation band and the emission coefficient in the thermal radiation band. This way in the definition of the surface temperature the temperature rise due to incident solar radiation, the heat transfer to and from the ground through conduction, the heat exchange with the environment due to convection and the thermal radiation exchange with the sky are taken into account. The radiation exchange with the rest surfaces of the computational domain is modeled through the DO model. The incident solar radiation on every solid and porous surface of the computational domain is calculated according to the inclination and orientation of each surface taking shading into account [20].
The buildings' roofs are deemed to have a thickness of 15cm while the lateral walls have a thickness of 25cm. The buildings internal temperature used in the EBM was 28 o C. As far as the temperature of the lateral surface is concerned, two types are considered: a) glass surface, b) light color masonry and finally a weighted average value is calculated.

Material Properties
Surfaces corresponding to lateral walls according to their material are divided into masonry and glazing. The thermophysical properties of these materials are given in the following table. Assuming that 25% of the lateral walls area is covered by glazing, the weighted average values are calculated for all thermophysical properties. In Table 2, the thermophysical properties for the used materials are presented. For trees and greenery the assuming porosity is set to 40%. The air dynamic viscosity is μ=1.789x10 -5 (Pas) and the thermal expansion coefficient is β=0.00343. The optical properties of the materials are given in Table 3.

Parametric Study
The object of the study is an area in the urban core of Volos, a segment of the quarter named 'Palia' close to the harbor with stigma (φ=39 ο 21.65', L=22 o 56.15' and elevation of 3 m). In this area exist: i) public and private buildings of 2 or 3 stores, ii) a school courtyard covered by asphalt, iii) a few trees, iv) pedestrian walkways covered by dark color tiles, v) streets covered by common asphalt and vi) a square covered with dark tiles and cement, sparse trees and a small fountain.
The area was studied for a typical summer day. For this day steady state simulations were conducted every 1 hour, with the following climatic data.
The examined area is 54000 m 2 of which 25000 m 2 are covered by buildings. The interventions are summarized as following: a) replace 12% of streets and pavements with tiles characterized by favorable optical properties, b) replace 1.6% of the streets with cold asphalt, c) addition trees to 100 m 2 which also increase shading. In Table 5, the optical properties of the proposed materials are presented.

Results
In Figures 2 to 4 the temperature distribution before and after the proposed intervention at 1.8 m is given for the hours 6:00, 12:00 and 16:00 of the typical day.   The average temperature reduction at the level of 1.8m is 1.7 o C, the total cooling degree hours are decreased from 68.07 to 47.39 and the CPI is improved by 17%, but it remains at levels indicating an extremely hot situation. The basic reasons are the low air velocities caused by the trees located in the south boundary of the examined segment which prevents cooling from the sea wind. Nevertheless, the area is characterized by low humidity during the hot hours of the summer. For this reason, the Wet-bulb Globe Temperature (WBGT) is calculated before and after the intervention [21], taking into account the climatic data of the area for the relative humidity.
From the WBGT calculations, an improvement of 5% was observed, while there is a reduction in time interval in which precautions are required [22] as it is presented in Figure 9. Since the relative humidity remains low and the available shading is expanded, the discomfort feeling is not too intensive for the biggest part of the day.

Conclusions
The prediction of the microclimate in an urban segment is a multivariate problem that requires the modeling of various factors. In the present work, a numerical model was developed taking into account not only the geometry and local climatic conditions (wind speed, temperature and solar radiation) but also the materials used, the existence of greenery and water elements and the human activity in the buildings (human presence and operating devices).
According to the calculated thermal comfort indices and the developed temperature and flow fields, the interventions' necessity is fully justified and documented. The proposed interventions which are the replacement of the streets' and pavement cover materials with materials with favorable optical properties and the increase in the green coverage of public spaces demonstrate an improvement in the local microclimate in terms of average temperature and thermal comfort indices. Nevertheless, the arrangement of the trees in the south part of the domain prevents the sea breeze air current and its cooling convection impact, weakening somehow the total effect of the whole proposed intervention.