Numerical simulation of different pollution sources in an urban environment

. A model for representing different pollution sources in an Urban Environment is presented in this paper. An in-house incompressible flow solver called caffa3d.MBRi was used and modified in order to include the sources from chimneys, cars and houses to compute the dispersion of different pollutants simultaneously. Moreover, different urban domains with a wide variety of sources were simulated reaching promising results.


Introduction
Air pollution in urban environments has become a public concern due to the great amount and wide variety of pollution sources that can be found inside a city. Moreover, the presence of buildings generates a characteristic microclimate that may decrease the pollution dispersion, causing high concentration areas that may present a health risk to inhabitants. With this in mind, different approaches have been considered to study and solve or at least diminish the problem [1], being one of them the use of numerical simulations to diagnose and predict the concentration field in a certain domain. Many different studies use Computational Fluid Dynamics (CFD) to predict the effects of certain sources under different wind conditions and city geometries [2]. Taking this into account, this paper aims to present a CFD tool developed in order to study the emission of different types of sources at the same time in Urban Environments.
In order to simulate the fluid dynamic and the pollution dispersion a CFD in-house open source code called caffa3d.MBRi [3,4] was used, adding some capabilities in order to take into account the different sources included in the study domain. To do so, the emission sources were classified into three categories depending on the emitting entity. The first of them are the point sources that correspond to industrial and commercial chimneys, the second type simulates linear sources that represent the vehicular emissions on a given street or road, and the last one refers to superficial sources that include the residential emissions due to domestic heating. Once the sources were categorized the corresponding emission rates were obtained from bibliography and measurement campaigns [7,8,9,10]. Last but not least, the different sources were included in the aforementioned solver so that they can be computed at the same time considering different pollutants with different emission rates.
Chapter number 2 presents the solver used and the models proposed to represent each different kind of source with multiple pollutants. Later on, on section 3 it is shown how the different input data was acquired, emphasizing the methods to determine the emission rates for each type of source. Then on chapter 4 different results and applications are shown, to then reach the conclusions on section 5.

Numerical model
The CFD caffa3d.MBRi solver is an open source inhouse simulation engine that uses the finite volume method, with second order accuracy in space and time. It allows the user to simulate three-dimensional incompressible flows in block structured curvilinear grids. The program counts with MPI parallelization which allows improving computation times reaching bigger domains and higher resolution. Moreover, it counts with a number of features such as different turbulence models, being LES (Smagorinski) the one used in this case; an Immersed Boundary Conditions module [5] used to represent complex geometries with simple meshes; and last but not least a Scalar Transport model that allows the computation of the pollutants dispersion treated as non-reactive gases. Apart from this features, there are many others that were not used for this particular study case. As far as the iteration scheme is concerned, it works by computing the Momentum Balance, then the Mass Balance and lastly the Scalar Transport for each time step as it can be seen in the scheme shown in Figure 1. Furthermore, the SIMPLE method for pressure-velocity coupling is used [6]. In this scheme, the scalar transport is done to one field at a time, this was used to allow the computation of multiple fields in order to represent different pollutants at the same time. Before a pollutant transport is to be computed, the solver calls a subroutine to handle the tracer's information and determine the boundary conditions and source terms to be taken into account when computing the scalar transport for each pollutant. In order to understand the computation of each type of source their models are presented below.

Point sources
This kind of sources correspond to particular chimneys that are either industrial or commercial ones to be found in the study domain. Given that the emission is located around a single point, the computation of this sources is the simplest one, the emitting cells ought to be the ones whose centre is inside a sphere centred in the emission point, above the chimney, with a given diameter, dependent on the size of it. So, every cell is checked to see if it should emit, so that once they are determined, the total source emission is divided among them proportionally to each cell volume, resulting in an equal emission concentration for all those cells.
Nonetheless, as the total volume is needed in order to distribute the source among the respective cells, the same procedure to determine those cells is done at the initialization and the total volume of each source is computed and stored. In consequence, once the source term is to be computed to a particular cell, corresponding to a particular chimney, one already has all the needed information: the emission rate of it, the total volume of the cells corresponding to that source and the volume of the particular cell, resulting in an emission rate for that cell equal to the one computed in equation 1. (1) In the next Figure the emitting cells above a chimney are shown, while the structure in black is the representation of the buildings of a Power Plant located within the domain obtained through the Immersed Boundary Conditions Module, (results for those emissions are shown later in this paper).

Fig. 2. Emitting cells (green) for a point source.
Black structure is the representation of the buildings.

Linear sources
These sources are meant to represent the vehicular emissions along a street or road, by computing the emission along straight lines. This is done by defining a rectangular prism around that line that goes from a given point to another, both at the same height above ground. Once the prism is defined, the method is similar to the one used in point sources, the cells whose centres are inside the prism will be the emitting ones, and the corresponding emission rate for them will be proportional to the cell volume. Similarly, the total volume of the cells included in the prism is computed at the initialization process, so as to have it at the time of computing the source term for each cell, which is done in the same way as point sources. In Figure 3, the emitting cells can be seen in light blue or yellow. Even though there may seem only cells of a certain height are emitting, the prism includes several layers that are between 0.33 m and 0.93 m above ground. Moreover, given the prisms overlap, some cells (yellow) will have an emission rate corresponding to the sum of both sources, so one has to be careful when working with sharp angles. The difference between these sources and the point ones, is that the latter has a known amount of pollutant emitted from a single point, whereas in the linear ones there is an emission along a straight line, with a known length. So in order to determine the amount of pollution emitted by it, one has to know the emission per length, then, the total emission is determined with the length of the line that defines a particular source of this kind.

Superficial sources
Unlike the other sources, this type represents a sum of sources along a surface, specifically the sum of the residential heating emissions. The reason each house chimney is not taken into account as a point source is that they are usually too many of them and are not N S E W necessarily emitting at the same time, so statistical values are used. Therefore, for these sources an area of emission is defined with a certain height so that the cells whose centres are inside that area and their height is immediately above a fixed one will emit.
In order to compute this, each superficial source is defined by a convex polygon, in order to check whether a column of cells is included in it, and an emission height so as to choose which cell of that column will be the emitting one. Apart from the nature of the emission, these sources differ from the previous ones because the emission rate of each cell will be proportional to the top area of the cells and not to its volume. This is done so that the total emission of a superficial source will be proportional to its surface area. Therefore, the emission rate needed for each superficial source will need to be expressed in mass per time per square metre. In Figure 4, the polygons used to define a superficial source above a residential domain are shown, in black one can see the built blocks to be found in the domain, represented by the Immersed Boundary method. It can be appreciated that, as the polygons have to be convex, two different sources had to be used, although they represent the same kind of emission.

Tracer transport procedures
In order to work with the different pollutants and sources a module was generated that includes the different procedures used. Two different type of functions can be identified, the data collection and initialization one, and the sources and boundary conditions implementation one.
The first kind of functions are called when starting the solver and consist basically on the field initialization, the sources input data reading, and the computation of the sources volume for the point and linear ones. This last one is a key element in order to divide the emissions correctly among the cells in a source. An important aspect of this procedure is that MPI communication is needed as one source may include cells from two different processors, so a reduction for the volume of each source has to be done.
The second kind of functions are in charge of the boundary conditions for the scalars and for the computations of the different sources described above. They are called for each pollutant when their transport is to be computed, and the source term for each cell will be the sum of the source terms corresponding to every source of every kind assigned to that cell. By doing this, a cell can be assigned multiple emission from different type of sources and different sources of the same kind.

Statistics computations
Apart from the modifications done in order to work with the pollutant dispersion, a module in order to compute the mean fields was added. This represents a key feature in pollutant dispersion as most of the pollution is measured in time averages, and the legal limitations are usually expressed in those terms. This module works by updating the mean field for each tracer on every time step, allowing the solver to print that data for post processing. Later on, results showing the mean concentration fields will be shown.

Numerical model input
As presented in the previous chapter, each type of emission source to be computed needs a special set of input data. With this in mind, two distinct set of data can be identified for each source, firstly location data in order to compute which cells should be emitting representing the particular source. Secondly, there has to be emission rates that indicates how much is emitted by each source, and that rate can have different units depending on the emission kind, however, the program will end up with a source term for each cell with units of mass per time. Given that each type of emission has different data to acquire, the methods used to determine such data will be explained for each source separately.

Point sources
For this sources the geometrical data needed to find the emitting cells is very simple, as all that's needed is a point and a diameter. Once this information is processed and the cells to emit are determined, what is needed is how much mass of each pollutant is emitted per second. With this in mind, those emission rates were determined for each chimney in the domain. First of all, for the industrial ones the emission rates were asked to each of the industries given that the information cannot be determined in any other way [13]. S E W N On the other hand, one can find commercial chimneys in the domain that correspond mainly to Restaurants and Bakeries. Given that this emission are due to the combustion of firewood in stoves for cooking, the amount of firewood used would determine the emission rates. So using a recent census on the use of fuels in the city, the amount used by Restaurants and Bakeries (separately) was obtained, thus, getting the wood consumption for each business in a year [12]. Furthermore, using the Atmospheric Emission Inventory [7], the total emission for this chimneys was determined.

Linear sources
These sources need more information than the previous ones, especially for the determination of the emitting cells. First of all, two points in the horizontal plane are needed to define the line of emission, then a height, a width and a distance from the ground are needed to create the prism, and evaluate whether the cells fall inside it or not. As far as the emission rates are concerned, what is needed for each source is an emission rate per metre, so in order to do that a vehicular count was done in key points in the city enabling the determination of the vehicular flow in routes and avenues inside the domain. This campaign was done in different points and at different moments during the day allowing the user to choose one to do the simulation for a specific scenario. Once the vehicles passing were counted, they were divided in three categories, trucks and buses, cars, and motorcycles, and a specific emission rate from bibliography was used for each of them. Table 2 shows the values used. With the collected data, the emission rates for the considered streets were achieved by multiplying the emission factors in Table 2 by the amount of vehicles passing per hour.

Superficial sources
In order to determine the emission cells for this sources what is needed is a set of straight lines (in the horizontal plane) that define a convex polygon to the lines left hand (when going from the first defined point to the second one), this data plus the domain limits will define the area of emission while an input height will determine the emitting cell for each column inside the previously determined area.
When it comes to determining the emission rates per square metre needed for this type of source, the procedure is fairly similar to the one used for the commercial chimneys. Statistical data from a census was used to determine the amount of wood used by houses for heating and cooking purposes and the way it is burnt, whether in stoves or boilers [12]. With this information and emission data for stoves and boilers using wood [7], the emissions per house were obtained and are shown in Table 3. Table 3. Emission factor for different vehicles in grams per km travelled [8,9,10]. With this information, what is left is to determine the amount of houses in the study area and the emission rate will be the emission factors in Table 3 times the houses in the area, divided by the considered area, therefore, reaching an emission rate in kg per year per square metre.

Results and discussion
Once the models were implemented and showed consistent results, different cases were simulated in order to validate the model and study different scenarios. Two different domains were studied one being a Residential Neighbourhood (Case 1), and the other one being the surroundings of a Power Plant (Case 2). In both cases the meteorological conditions are the same with a Northern wind and a Neutral Atmospheric Boundary Layer.
Firstly, Case 1 was simulated in order to validate the model [11]. Particularly the study case for a residential neighbourhood of Montevideo was studied, showing coherent results with measurements taken in place. In this scenario, a superficial source of PM10 was considered above the houses (see Figure 4) getting interesting results, such as the one shown in Figure 5, where the concentration field for PM10 at ground level is shown. Secondly, simulations were done for Case 2, consisting of a power plant that is located next to an avenue with heavy traffic, so both linear and point sources were represented in the domain. Different results were obtained, showing the effects caused by the presence of building in the domain and again the concentration of pollutants at ground level that allow air quality evaluation and possible human health risks presence.
In Figure 6, the concentration field of SO2 at a height of 60 m above ground is shown, different concentration peaks can be seen associated with vortices caused by the presence of the structure at that height. Moreover, given that the concentration field is available in the whole domain, concentration isosurfaces (plumes) can be computed showing the distribution of the pollutant along the domain, as it is shown in Figure  7. In this case a plume for a concentration of 15000 µg/m³ is depicted, and shows a clear dispersion of the pollutant downstream. On the other hand, the impact of the linear sources was studied by getting the concentration field at ground level for CO. In Figure 8, that field can be appreciated showing that the instant concentration peaks can be found downstream of the sources (black rectangles). Nonetheless, when the mean concentration field over 10000 seconds is studied (Figure 9), the maximum concentration decreases and the area of impact is clearly depicted.  Another interesting result to study is the interaction of the sources present in the domain. For instance, NOX is emitted by both linear sources and point sources. So a plume for a concentration of 1500 µg/m³ was plotted in Figure 10. This Figure clearly shows the presence of both sources, as a plume coming from the chimney can be seen and small plumes at ground level are also present. Even though both kind of sources are appreciated it can be seen that the point one has a higher impact, but as the pollutant does not reach ground level at such concentration the linear emissions are still important.

Conclusions and future work
As far as the implemented models are concerned, they represent the sources coherently and as expected, while the dispersion of the pollutants is consistent with what would be expected with the velocity fields to be found in the domain. All in all, the concentration fields obtained seem to be a close representation of reality, allowing both quantitative and qualitative studies of different problems. For instance, one could evaluate the possible health risks and compliance of the air quality standards. Moreover, the location of concentration peaks and low dispersion areas can be identified, allowing one to use this tool to design measurements campaigns or redesign structures to reduce possible pollutant accumulation spots.
Furthermore, as there is a predefined set of input data, the use of the model became fairly easy, and allowed emission changes without the need of recompilation. All in all, the generated model allows the user to easily represent different emission sources to be found in an Urban Environment, in order to compute the dispersion of many pollutants at the same time. So with this model and some of the procedures to collect the emission rates shown above, it is easy to evaluate the air quality conditions of the domain.
Despite the satisfactory results obtained with this model, many improvements and additions are to be made. One of the main ones that is currently being studied is the passage of this capabilities to a Heterogeneous-Computing (CPU-GPU) program, so as to allow the representation of bigger domains and better resolution, a common problem when working in Urban Environments. Another issue to be addressed, is the possibility to use moving sources, to represent vehicles explicitly.