Air Flow Study around Isolated Cubical Building in the City of Athens under Various Climate Conditions

Abstract

This study focuses on the airflow and pollutant dispersion around an isolated cubical building located in a warm Mediterranean climate, taking into account the local microclimate conditions (of airflow, albedo of building and soil, and air humidity) using a large-eddy simulation (LES) numerical approach. To test the reliability of computations, comparisons are made against the SILSOE cube experimental data. Three different scenarios are examined: (a) Scenario A with adiabatic walls, (b) Scenario B with the same constant temperature on all the surfaces of the building, and (c) Scenario C using convective and radiative conditions imposed by the local microclimate. For the first two cases the velocity and temperature fields resulting are almost identical. In the third case, the resulting temperature on the surfaces of the building is increased by 19.5%, the center (eye) of the wake zone is raised from the ground and the maximum pollutant concentration is drastically reduced (89%).

1. Introduction

The excessive population concentration in megacities has resulted in environmental pollution problems and overcrowded living conditions. The air quality in these cities is determined by the airflow in the complex urban terrain, the temporal and spatial conditions of pollutant sources, and the local urban microclimate parameters. The airflow and pollutant dispersion in an urban environment can be studied in different geometrical scales, such as blocks of buildings, street canyons, and isolated buildings. Several wind-tunnel, real-scale experiments and numerical studies examine the airflow and pollutant dispersion around blocks of buildings [1,2,3,4,5] and street canyons [6,7,8,9,10,11,12].

The city’s local microclimate parameters are air and surface temperatures, humidity, direct and diffuse solar irradiation, and wind speed at different directions [13]. The knowledge of the urban microclimate is necessary to prevent and reduce pollution problems. The local microclimate conditions are playing an important role in the airflow and pollution distribution in an urban environment [14]. For this reason, complex models can be used to calculate mandatory variables such as the heat capacity of a building, the shortwave and longwave radiation according to the solar position, the emissivity, and the air humidity. According to these variables, the increment or decrement of the thermal comfort rate can be found to define the human comfort conditions. Another important issue is the values of the heating rate and of the thermal absorption of the surfaces of the building. Using complicated microclimate models, the distribution of the environmental temperature and their effect on pedestrians can be defined [15,16,17,18,19]. Small-scale meteorological models are widely used to define the comfort conditions inside the urban environment for winter and summer seasons [20,21,22,23].

A phenomenon by which warm layers of air are accumulated in densely populated centers is defined as an Urban Heat Island (UHI) creating uncomfortable thermal conditions for urban living [24]. A high air-temperature increment leads to the need for higher energy consumption for cooling in the summer, with a consequence of increasing air pollutants in the atmosphere [25].

It is also identified that the rapid augmentation of the air temperature depends on the airflow velocity. At calm wind conditions, where the airflow velocities are lower than 3 m/s, the density of hot gas masses is much higher than at heavy wind gusts [26]. As a result, wind gusts with high airflow velocities can remove the air hot masses from the urban environments, which leads to the decrement of the UHI phenomenon [27].

The effect of the microclimate in densely populated urban centers is also an important topic for study. Gronemeier and Raasch [28] studied the urban flow for the city of Hong Kong with PALM software. The numerical results showed a rapid increment in wind velocity at the pedestrian crossing areas and high changes of the airflow direction [28,29]. Another study by Pfafferott and Rißmann [30] with the PALM code examined the interaction between buildings and the urban microclimate using an energy balance solver. They found that the temperature increment of the urban environment was caused by the heated building, as the outcome of energy from the indoor model [30]. Resler and Eben [31] examined the airflow over the city of Prague in the Czech Republic with PALM code. Their study showed that the temperature and the concentration of NOx are in a good agreement with real field measurements during the summer season. In the winter season, the deviations between numerical results and experimental data pointed to inaccuracies in modeling the atmospheric boundary layer [31]. Thus, real field urban measurements can give important information for both real-scale meteorological phenomena and their impact on an urban boundary layer. To better understand the complex phenomena of the airflow in an urban environment, it is important to combine real field measurements and wind tunnel experiments so that these two methods can be interrelated [32,33]. Numerical methods can be also used to simulate the turbulence modelling of the airflow around both small- and large-scale cubical obstacles [34,35,36,37,38].

Several wind-tunnel experiments exist that describe the airflow characteristics around a cubical building. The wind-tunnel experiments are suitable for collecting detailed spatial experimental data [32,33,39,40,41] under stable conditions. However, experimental wind tunnel experiments are not efficient for the prediction of real-scale meteorological conditions.

The present study is focusing on the simulation of the flow and pollution distribution around an isolated building [42] and, more importantly, on the influence of a meteorological model on these results, i.e., what differences exist between the use of real boundary conditions and those employed for reasons of convenience. In addition, the use of an isolated building instead of a large urban area ensures good reliability of the numerical results for the available resources. In this respect, such studies do not seem to exist in the scientific literature.

The isolated building of this study is in the climatological area of Athens, Greece. Three different scenarios are examined. Scenario A studies the neutral flow and pollutant dispersion in the wake area of an isolated cubical building with adiabatic walls, without any buoyancy forces. Scenario B studies the same building with specified constant temperature of 50 °C on its surfaces. Scenario C studies the same building with convective and radiative conditions imposed by the local microclimate.

2. Configuration and Pollutant Dispersion Modelling

2.1. Problem Configuration

Figure 1 presents the geometry of the computational domain around the isolated cubical building of this study, with exactly the same dimensions as the SILSOE cube. The SILSOE cube field experiment for the study of the airflow around it [43], with a 6 m height cube placed in an open area, has been used as a test case in several other studies [44]. This geometry remains unchanged in all the different simulated scenarios of this work.

The global coordinate system origin is the frontal left edge of the computational domain and the distance between the frontal area of the computational domain and the building is 5 H. The distance between the lateral boundaries of the computational domain and the building’s surface is also 5 H [45]. The distance from the rear surface of the building up to the outlet boundary of the computational domain is set at 12 H and the total height of the computational domain is 5 H. The resulting blockage effect value is approximately 1.81%. Under the recommendations of the German Association of Engineers (VDI), the blockage effect should be maintained below 10% during the simulations [46] and according to [46,47] should be even lower than 3%, which is satisfied in the present study.

The Reynolds number is kept constant at 4.03 × 10⁶, based on the height of the cube and the free stream velocity. The specified flow field at the inlet of the computational domain is in the form of a logarithmic profile with the no-slip condition at the ground. A passive source of pollutants is located on the floor of the computational domain at a distance H from the rear surface of the building, using as a pollutant methane (CH₄) gas with a concentration and release of 600 ppm and 18.5 L/h, respectively. The dimensionless concentration coefficient K of the passive pollutant is defined as [48,49]:

$$K=\frac{\left(\raisebox{1ex}{$c_{Measured}$}\!\left/ \!\raisebox{-1ex}{$c_{Source}$}\right.\right)u_H{H}^2}{Q_{Source}}$$

(1)

where $c_{Measured}$ is the tracer concentration at any position, $c_{Source}$ is the tracer concentration on the pollutant source, $Q_{Source}$ is the release rate of the pollutant, and $u_H$ is the velocity of the flow at the height of the building.

In Scenarios B and C, the temperature is affecting the flow field through the usual Boussinesq approximation. In these cases, the Grashof number, i.e., the ratio of thermal buoyancy forces over viscous forces is defined as [50]:

$$G{r}_H=\frac{g\beta \left(T_{Surface}-T_{\infty }\right)H^3}{v^2}$$

(2)

where $T_{Surface}$ is the temperature on the surfaces of the building, $T_{\infty }$ is the bulk temperature, $\beta$ is the coefficient of thermal expansion, and $v$ is the kinematic viscosity of the fluid. In Scenarios Β and C, the Grashof number based on the height of the building is 2.43 × 10¹¹ and 7.68 × 10¹¹, respectively. Thus, the buoyancy forces are intensive and the boundary layers close to the heated surfaces are highly turbulent.

2.2. Governing Equations

The PALM model system, an open-source software that simulates atmospheric and oceanic boundary layers, is used for the present simulations. The Navier-Stokes equations in a non-hydrostatic, filtered, and incompressible form under the Boussinesq approximation are solved. The basic equations for the conservation of mass, momentum, energy, and moisture are defined as:

$$\frac{\partial u_i}{\partial x_j}=0$$

(3)

$$\frac{\partial u_i}{\partial t}=-\frac{\partial u_i{u}_j}{\partial x_j}-{\epsilon }_{ijk}{f}_j{u}_k+{\epsilon }_{i3j}{f}_3{u}_{g,j}-\frac{1}{{\rho }_0}\frac{\partial {\pi }^{\ast }}{\partial x_i}+g\frac{{\theta }_v-\langle {\theta }_v\rangle }{\langle {\theta }_v\rangle }{\delta }_{i3}-\frac{\partial }{\partial x_j}\left(\overline{u_i^{″}{u}_j^{″}}-\frac{2}{3}e{\delta }_{ij}\right)$$

(4)

$$\frac{\partial \theta }{\partial t}=-\frac{\partial u_j\theta }{\partial x_j}-\frac{\partial }{\partial x_j}\left(\overline{u_j^{″}{\theta }^{″}}\right)-\frac{L_v}{C_p\Pi }{{\Psi }_q}_v$$

(5)

$$\frac{\partial q_v}{\partial t}=-\frac{\partial u_j{q}_v}{\partial x_j}-\frac{\partial }{\partial x_j}\left(\overline{u_j^{″}{q}_v^{″}}\right)+{{\Psi }_q}_v$$

(6)

$$\frac{\partial s}{\partial t}=-\frac{\partial u_{j}s}{\partial x_j}-\frac{\partial }{\partial x_j}\left(\overline{u_j^{″}{s}^{″}}\right)+{\Psi }_s$$

(7)

where $i, j, k\in \left\{1, 2, 3\right\}$, the components of the velocity $\left(u_1=u, u_2=v \mathrm{and} u_3=w\right)$ are defined by the $u_i$ variable at specific positions in the flow field $x_i$, where $x_1=x, x_2=y$ and $x_3=z$, and $t$ is time. $q_v$ and $s$ are the latent heat transfer and moisture, respectively. The gravitational acceleration is denoted by $g$, the density of the dry air is defined by ${\rho }_0$, and $L_v$ is the latent heat of vaporization. The source term of the variable $q_v$, is defined as ${{\Psi }_q}_v$ and the sink term of the variable $s$ is defined as ${\Psi }_s$. The angle brackets correspond to the horizontal averages of the flow field and the double prime indicated the subgrid-scale variables.

The subgrid-scale turbulence kinetic energy is defined as:

$$e=\frac{1}{2}\left(\overline{u_i^{″}{u}_i^{″}}\right)$$

(8)

The modified perturbation pressure based on the perturbation pressure $p^{\ast }$ is defined as:

$${\pi }^{\ast }=p^{\ast }+\frac{2}{3}{\rho }_{0}e$$

(9)

The potential temperature is defined from the equation:

$$\theta =\frac{T}{\Pi }$$

(10)

where $T$ is the instantaneous absolute temperature and the Exner function is calculated by the equation:

$$\Pi ={\left(\frac{p}{p_0}\right)}^{\frac{R_d}{C_p}}$$

(11)

where $p$ is the hydrostatic pressure of the air, $p_0=1000 \mathrm{hPa}$ is the reference pressure, $R_d$ is the specific gas constant for dry air, and $C_p$ is the specific heat capacity under constant pressure.

The virtual potential temperature is defined by the equation:

$${\theta }_v=\theta \left[1+q_v\left(\frac{R_v}{R_d}-1\right)-q_l\right]$$

(12)

where $R_v$ is the specific gas constant for water vapors, and $q_l$ is the liquid water mixing ratio.

2.3. Urban Surface Model (USM) of the PALM Model System

The urban surface model (USM) of the PALM model system is used for the energy balance in the computational domain. The energy balance solver is driven by three main procedures: (a) the solver predicts the temperature of the outer surfaces, (b) the turbulence dispersion of the sensible temperature is calculated at the surfaces near the walls, and (c) the calculation of the sublayer heat flux caused by convection. The first two procedures of the energy balance solver are executed simultaneously. The third procedure is executed under the calculations of a subgrid-scale model that predicts the thermal diffusion from the bluff body.

The energy balance solver is defined as:

$$C_0\frac{d{T}_0}{dt}=R_n-H-LE-G$$

(13)

where $C_0$ is the specific heat capacity, $T_0$ is the radiative temperature of the surface skin layer, $R_n$ is the net radiation, $H$ is the sensible heat, $LE$ is the latent heat, and $G$ is the heat flux at the surface of the ground. On the USM, the energy balance is calculated separately for each surface, and the three different types of radiations (sensible heat, latent heat, and ground heat flux) from the surface heat are combined.

The calculation of the convective heat transfer between the air and the outer surfaces is defined by the equation:

$$H=h\left({\theta }_1-{\theta }_0\right)$$

(14)

where ${\theta }_0$ is the temperature of the outer surface of the building, ${\theta }_1$ is the temperature of the air mass in contact with the outer surface of the building, and $h$ is the coefficient of convective heat transfer and is parameterized only for the vertical surfaces inside the computational domain [51]. For the horizontal surfaces inside the computational domain, the coefficient of convective heat transfer is parameterized under the Monin-Obukhov similarity, which involves the calculation of the local friction velocity [52,53]. The friction velocity is used for the calculation of flow momentum near the surface, for each surface separately.

3. Initial and Boundary Conditions

The inlet velocity profile is described by the logarithmic law of the wall with a free-stream velocity value of 10.13 m/s. Figure 2 presents the normalized velocity profiles at both the inlet of the computational domain and in the X/H = 2 position from the inlet for Scenarios A–C, where A and B should be identical.

As observed in Figure 2a, an important deviation appears in Scenario C that corresponds to the application of the meteorological model in comparison with the other two scenarios. This difference appears also in the free-stream flow region and is caused by the strong buoyancy effects caused by the increased thermal heating of the air.

As shown in Figure 3, the present numerical results for the pressure coefficient are in good agreement with the corresponding experimental data of Scenarios A and B.

The pressure coefficient at the surfaces of the building for all scenarios is defined as:

$$C_P=\frac{p-p_{REF}}{\frac{1}{2}{\rho }_{REF}{u}_{\infty }^2}$$

(15)

where $p$ is the static pressure of the fluid at a point, $p_{REF}$ and ${\rho }_{REF}$ are the static pressure and density of the fluid at the free-stream at the inlet of the computational domain, respectively, and $u_{\infty }$ is the free-stream velocity at the inlet.

As observed in Figure 3, important differences exist between Scenarios A–C because of the heating of the surfaces of the building and the microclimate model application. The pressure coefficient distribution on the frontal surface of the building for Scenarios B and C is higher than the corresponding for Scenario A because of the influence of the mixed convection heat transfer, which causes reduction on the density of the streamlines in this region. In Scenario C, this influence is expected to be larger than for Scenario B. In Scenario B, the pressure distribution at the rear surface of the building significantly approaches the corresponding for Scenario A. In Scenario C, the pressure distribution at the rear surface of the building is the highest compared with other scenarios because of the domination of the buoyant heating forces there.

Figure 4 presents the vertical profile of the turbulent kinetic energy at the inlet of the computational domain for Scenarios A–C, as compared with the corresponding numerical data on the same normalized upstream position from the cube [18].

As described in Figure 4, the turbulent kinetic energy distribution at the inlet of the computational domain is almost identical for Scenarios A and B and different from Scenario C because of the effect of the buoyancy forces from the microclimate conditions.

On the lateral boundaries of the computational domain, periodic boundary conditions are employed because of the temporal and spatial periodicity of the flow. Neumann boundary conditions are used for the turbulent kinetic energy $\left(e\right)$, temperature $\left(\theta \right),$ and the perturbation pressure $\left(p^{\ast }\right)$ to calculate the following equations concerning the building’s height:

$$\overline{e}\left(-\frac{\Delta y}{2}\right)=\overline{e}\left(+\frac{\Delta y}{2}\right)$$

(16)

$$\overline{\theta }\left(-\frac{\Delta y}{2}\right)=\overline{\theta }\left(+\frac{\Delta y}{2}\right)$$

(17)

$$\overline{p^{\ast }}\left(-\frac{\Delta y}{2}\right)=\overline{p^{\ast }}\left(+\frac{\Delta y}{2}\right)$$

(18)

The Dirichlet boundary condition is used for the velocity at the ground of the computational domain for the no-slip condition. Therefore, the velocity components at the ground are:

$$\overline{u}\left(z=0\right)=0$$

(19)

$$\overline{v}\left(z=0\right)=0$$

(20)

$$\overline{w}\left(z=0\right)=0$$

(21)

On the staggered computational grid, the velocity components in the X and Y directions, respectively, are defined at the specific height of $z=\pm \frac{\Delta z}{2}$. Thus, the symmetry boundary condition is used as:

$$\overline{u}\left(-\frac{\Delta z}{2}\right)=-\overline{u}\left(+\frac{\Delta z}{2}\right)$$

(22)

$$\overline{v}\left(-\frac{\Delta z}{2}\right)=-\overline{v}\left(+\frac{\Delta z}{2}\right)$$

(23)

On the upper surface of the computational domain, the perturbation pressure $\left(p^{\ast }\right)$ obeys the Neumann boundary condition. For the velocity to maintain the free-stream regime over the atmospheric boundary layer, Neumann boundary conditions are also applied to the velocity components of the velocity in the X and Y directions:

$${{\partial }_z\overline{u}|}_{top}=constant$$

(24)

$${{\partial }_z\overline{v}|}_{top}=constant$$

(25)

Additionally, on the upper surface of the computational domain, Dirichlet boundary conditions are used for both the turbulence kinetic energy $\left(e\right)$ and the temperature $\left(\theta \right)$.

4. Nested Computational Grid

PALM provides the capability of computational grid self-nesting. In the self-nesting mode, there is a root/parent computational domain with the ability of nesting up to 63 levels of nested/child computational subdomains. Each subdomain may be root/parent to another nested subdomain. As a result, the second subdomain is simultaneously nested referring to the first parent, and also root, to the third nested. The nested computational subdomain receives all the appropriate and mandatory information for the three components of the velocity and all the prognostic variables of the vectors from its boundaries with the root computational domain. The flow field data are interpolated from the coarse to the finer computational grid. By the end of each time step, the corrected solution from the solver is reversely interpolated to the root/parent domain.

Figure 5 shows the grid arrangement used in the simulations. The root computational grid (coarse grid) is shown at full height of the computational domain in contrast with the nested (finer grid-darker parts), which is shown at two building heights.

For the simulations, three different computational grids are used with an incremental resolution from the coarse to the finer: the coarse computational grid with 924,400 cells, the medium grid with 1,248,400 cells, and the fine grid with 1,842,400 cells. The grid errors using two computational grids are estimated from the following equation [4]:

$$GCI=\frac{f_2-f_1}{1-r^p}$$

(26)

where $f_2$ is the numerical solution that results from the medium computational grid and the corresponding from the finer is $f_1$. $r$ is the refinement factor between the two computational grids and $p$ is the accuracy of the algorithm, which is 3 for the present study. The refinement factor between the medium and finer computational grid is approximately 1.47.

Figure 6a presents the profile of the error bars of the normalized velocity, at the axial position 0.5 H from the rear surface of the building for scenario A. Additionally, Figure 6b shows the profile of the normalized concentration of the pollutant at the same axial position and scenario.

The mean value of the errors of the normalized velocity between the coarse and the medium computational grids is approximately 3.77% and the corresponding value between the medium and the fine grids is about 1.54%. The mean value of the errors of the pollutant between the coarse and the medium computational grids is approximately 2.77% and the corresponding value between the medium and the fine grids is about 1.24%. For both of the aforementioned variables, the decrement of the errors going from the coarse to the finer computational grids is clearly observed.

5. Numerical Details

The global discretization of the computational grid is accomplished using finite differences on a staggered Cartesian Arakawa-C grid [55]. The equations are spatially discretized on a fifth-order differential upwind scheme [56], while a third-order Runge-Kutta scheme is used on the temporal discretization [57]. The equations are implicitly filtered by the discretization of the computational grid and the subgrid-scale processes are calculated by the 1.5-order Deardorff numerical scheme [58]. Thus, it is assumed that the energy that is transported by the subgrid-scale vortices is proportional to the local gradients of the mean quantities [59,60].

The convergence criterion maintained on each simulation is kept below the 10⁻⁴ value for each computed variable based on the error. The time step is automatically adjusted for the CFL constant value of 0.9.

Every simulation is carried out until flow stationarity has been achieved. Figure 7 presents the normalized velocity components with respect to the free-stream velocity (u_x/u_∞, u_y/u_∞, u_z/u_∞), at the point with Cartesian coordinates of X: 10 H, Y: 5.5 H, and Z: 0.5 H and for the period from 600 to 1400 s for Scenario A.

As it is shown in Figure 7, the flow seems to have attained stationarity conditions from the time point of 700 (s) onward. The attainment of stationarity conditions is confirmed by finding out that the mean value and higher order statistics are independent of the time of initiation of the measurements. If this is valid for the mean value and the autocorrelation function, the process is said to be weakly stationary. In this respect, the mean, the variance, and the autocovariance with a time delay of 2 (s) values are computed and given in tabulated form below (

Table 1. Mean Value, Variance and Autocovariance for Specific Time Periods.

Averaging Period (s)	Mean Value	Variance	Autocovariance
600–800	1.0806	0.0473	0.3426
800–1000	1.0851	0.0538	0.3738
1000–1200	1.0879	0.0530	0.4316
1200–1400	1.0889	0.0524	0.4443

):

These computed values show that the flow field is at least weakly stationary.

However, although the flow field has reached stationary conditions it is possible the concentration field has not reached it yet.

Figure 8 presents the pollutant dispersion inside the computational domain for Scenarios A and C at the position with coordinates of X: 7.5 H, Y: 5.5 H, and Z: 0.5 H for the period of 200–1000 (s).

The mean, the variance, and the autocovariance with time delay of 2 (s) values of the pollutant concentration are computed and given in tabulated form below (

Table 2. Mean Value, Variance and Autocovariance for Specific Time Periods.

Averaging Period (s)	Mean Value	Variance	Autocovariance
200–400	0.2449	0.0186	0.9584
400–600	0.2445	0.0198	0.9532
600–800	0.2480	0.0193	0.9548
800–1000	0.2463	0.0195	0.9524

), which demonstrates that a passive scalar attains stationarity earlier in time than the velocity field:

6. Results and Discussion

The air flow field around the cubical geometry is presented in Figure 9 using normalized velocity profiles on the symmetry plane at different streamwise positions upstream, downstream, and at the middle of the roof of the cube (X/H = 5.5) for all three scenarios [44,54,61,62,63,64,65]. Comparisons are also made with Richards and Castro’s experimental data [60] for Scenario A.

Differences in the normalized velocity profile between Scenario C and Scenarios A and B are observed especially in the wake region of the cube, where the mixing layer for Scenario C has moved almost one cube height above the mixing layer of the other two scenarios. This is caused by the strong buoyant forces created in this case.

Figure 10 shows the streamlines of the flow for Scenario A. A detachment of the flow at the top surface of the building and a recirculation of the flow on both the frontal and the rear surfaces of the building are observed.

Different recirculation lengths X_f, X_b, and X_r are defined for the frontal, rear, and roof recirculation of the building, respectively.

Table 3. Recirculation length of the recirculation zones.

	Xf	Xb	Xr
Martinuzzi and Tropea [61]	1.04 H	1.61 H	-
Rodi [62]	0.651 H	2.182 H	0.432 H
Hoxey, Richards [63]	0.75 H	1.4 H	0.57 H
Richards and Norris [64]	0.9 H	1.4 H	0.9 H
Hu, Xuan [65]	-	1.31 H	0.94 H
Scenario A	0.79 H	1.56 H	0.56 H
Scenario B	0.56 H	1.97 H	0.64 H
Scenario C	0.86 H	2.15 H	0.73 H

presents a comparison of the length of different recirculation zones in comparison with available experimental data.

The recirculation zone X_f in front of the upstream surface of the building for Scenario A is 41% higher than the corresponding for Scenario B and 8.14% lower than Scenario C. In addition, the recirculation zone on the roof X_r for Scenario A is 20.81% lower than the recirculation zone of Scenario B and 27.44% lower than for Scenario C. The recirculation zone in the wake region X_b of the building for Scenario A is 13% lower than the corresponding recirculation zone of Scenario B and 23.28% lower than the recirculation zone for Scenario C.

Figure 11 shows the normalized vertical temperature profile for Scenario C at the position with coordinates of X = 7 H and Y = 5.5 H, on the symmetry plane and two heights behind the building.

The present numerical results are in good agreement with the experimental data of Uehara and Murakami [41], who studied the stability of the atmospheric flow inside an urban street canyon placed normal to the wind direction. They found that inside the canyon exists a stable thermal stratification as shown above, which weakens the cavity eddy.

Figure 12 shows the temperature distribution on the vertical and horizontal surfaces of the building for Scenario C, along the intersection lines with the transverse symmetry plane (path No. 3) and two other planes parallel to it at distances 0.25 H (path No. 2) and 0.45 H (path No. 1), respectively.

The surface temperatures of the cube are the result of an energy balance between the incoming thermal radiation that is absorbed and reflected or emitted, and heat convection by the airflow. The amount of radiation that reaches the surfaces of the cube and the visibility of the solar path are determined by Athens’s longitude and latitude. Albedo and emissivity of the surfaces control the shortwave and longwave radiation components reflected and emitted, respectively. The albedo and emissivity coefficients are 0.2 and 0.95, respectively, typical values for concrete material. The thermal heat capacity and thermal conductivity of the cube control its ability to store and conduct heat, respectively.

The temperature on the cube surfaces depends primarily on the convective thermal energy loss. In regions, where turbulence exists high heat exchange takes place. At the lower part of the frontal surface of the cube, where the horseshoe-shaped vortex appears, low surface temperatures appear. On the contrary, on the upper part of the frontal surface of the cube higher surface temperatures are shown. Similarly, on the wake zone of the cube where the arc vortex is present, low surface temperatures are shown because of the high heat exchange rate. On the first part of the roof surface where the flow separates the heat exchange is small, so the surface temperature is high, and on the second part it is reduced where the flow reattaches.

Figure 13 present the vertical mean pollutant concentration profiles computed for the time period from 800–1400 (s), at which all variables are statistically stationary on the symmetry plane for both Scenarios A and C at the Cartesian coordinate positions of X: 7 H, for a, X: 7.2 H for b, and 7.4 H, for c. The height Z extends from 0 to 2 H.

In Figure 13, the highest mean concentration of the pollutant is observed for Scenario A, where a large amount of the pollutant is trapped inside the wake region. In contrast, in Scenario C where because of the high buoyancy forces the maximum pollutant concentration is small in the above region as the pollutant is shifted almost one cube height higher away from the floor, thus creating a recirculation region of larger length and height than that corresponding for Scenario A.

As the pollutant in Scenario A is trapped inside the wake region of the cubical building, it becomes very harmful to the pedestrian’s health. In Scenario C, this problem is exempted to a large extent, as the high levels of pollutant concentrations are at much higher heights than those of the pedestrians. This behavior is the result primarily of the buoyant forces acting on the roof and on the leeward surfaces of the building.

Figure 14 present the distribution of the flux concentration on the symmetry plane because of the x-velocity component for Scenarios A and C at the x-positions: (a) (X: 7 H, Y: 5.5 H), (b) (X: 7.2 H, Y: 5.5 H), and (c) (X: 7.4 H, Y: 5.5 H) for heights Z that extend from 0 to 2 H.

In the wake region, up to the height of the cube, the flux concentration is similar for both cases A and C. At heights higher than the building, the concentration flux for Scenario A is larger of that of C because of the important influence of higher free stream air velocity. Additionally, the wake zone height in Scenario C exceeds the height of the building and the low values of air velocity produce low concentration fluxes.

Figure 15a–c present the distribution of the flux concentration on the symmetry plane because of the z-velocity component for Scenarios A and C at the x-positions: (a) (X: 7 H, Y: 5.5 H), (b) (X: 7.2 H, Y: 5.5 H), and (c) (X: 7.4 H, Y: 5.5 H) for heights Z that extend from 0 to 2 H.

For Scenario A, where the buoyancy forces do not affect the pollutant distribution, the concentration flux is toward the ground floor (Figure 15a,b) within the recirculation zone near the point of reattachment, in accord with the flow field of Figure 10, and away from the floor outside the recirculation zone. On the contrary, for Scenario C, the buoyancy forces lift the pollutants near position X = 7 H (left end of source) where at positions X = 7.2 H and X = 7.2 H (center and right end of source, respectively) the buoyancy forces are counterbalanced by momentum flux forces. It is expected that by placing the pollutant source closer to the leeward face of the building, the lifting forces that act on the concentration would be higher, resulting in higher concentration fluxes away from the ground. Thus, a higher fraction of the pollutants will escape the recirculation zone.

7. Conclusions

In the present work, three different scenarios are studied. The emphasis is on the flow and pollutant dispersion around an isolated cubical building in an open area and with a pollutant source placed on the ground and at a distance of one building height in the symmetry plane behind its rear surface. Scenario A with adiabatic walls of the building and without any buoyancy forces, Scenario B with the same constant temperature of 50 °C on all the surfaces of the building, and Scenario C with convective and radiative conditions imposed by the local microclimate.

Each scenario is simulated using the large-eddy simulation (LES) approach. The numerical results that describe the nature and the behavior of the flow around an isolated building are in good agreement with the corresponding experimental and numerical data that were used for their validation.

It is observed that at higher flow field temperatures, the length and the corresponding height of the recirculation are increased. Thus, the recirculation region for Scenario B has higher length than that for Scenario A and lower than that for Scenario C, for which exist the highest temperatures from the three scenarios studied.

The pollutant in Scenario A is trapped inside the wake region of the cubical building with its highest concentration levels at a height close to what a pedestrian would breathe, thus being very harmful to his health. In Scenario C, the higher temperatures together with the associated buoyancy forces that are developed raise most of the pollutants to a height greater than that of the building, allowing them to escape the recirculation region of the building and be carried away by the wind, without causing any serious harm to the health of a pedestrian.

It is expected that by placing the pollutant source closer to the leeward face of the building, the lifting forces that act on the concentration would be higher, resulting in higher concentration fluxes away from the ground. Thus, a higher fraction of the pollutants will escape the recirculation zone.