1. Introduction

The Kyoto Protocol from 1998, was a crucial step towards global conservation of the environment, triggering rules developments for atmospheric pollution mitigation. Currently, pollutant release into the atmosphere follows protocols requiring monitoring and adequacy in the emission limits. Each country has its own regulatory agency, such as CONAMA (Conselho Nacional do Meio Ambiente) in Brazil, which prepares laws and supervises their compliance. These agencies, following treaties and guidelines of international environmental forums, indicate the use of mathematical models as a complementary tool to estimate pollutants concentrations in industry surroundings. However, each model needs to incorporate physical and micro-meteorological characteristics compatible with each region of application.

Understanding the physical phenomenon for later mathematical description is the first step in developing this monitoring tool. Thus, when observing events of pollutants dispersion, the presence of turbulent movements caused by nonlinear flow contributions has a pronounced presence. However, models for pollutant dispersion based on the advection-diffusion equation result from simplifications in order to obtain a deterministic mathematical description. This procedure eliminates the possibility of this model to reproduce turbulent characteristics fundamental for the dispersion phenomenon. Fick’s closure is an example of a formal procedure applied in the advection-diffusion equation, where the nonlinear terms, the turbulent flows, arising from the Reynolds decomposition, are replaced by a mean concentration gradient. Holmes states that nonlinear terms are essential for turbulence, and the elimination by linearization of equations weakens the results.

The lack of a single mathematical definition for turbulence induces the use of some of its characteristics to describe it. The buoyancy parameter, momentum flow, and heat flow are examples of mechanisms that have structured a mathematical characterization of turbulent behavior in the equation. They generate in the dispersion phenomenon movements that form eddies, also identified as coherent structures and are a conceptual tool for reducing turbulence complexity, although there is still no unique theoretical description for them. The redefinition of the closure by a complex turbulent diffusion constant opens pathways to recover at least some of the effects induced by turbulence, i.e. by coherent structures in the present model manifest in the presence of phase differences in the solution.

The vertical diffusion coefficient K_z is the ideal component to receive the phase inclusion in the advection-diffusion model, based on some already consolidated facts in mathematical models. In order to justify the choice of this component we observe the variables of the problem that are responsible for the creation of turbulent effects, such as roughness, pressure and temperature difference and daily cycle of heating and cooling of the Earth’s surface, which drive the vertical component of the phenomenon. There are also theories that attempt to overcome the lack of turbulent effects in the deterministic equations by parametrizing the turbulence for stable and convective boundary layer schemes. These models insert micro-meteorological parameters in the vertical coefficient K_z and create a profile linked to the height of the boundary layer. Thus, a phase term in the vertical component of the equation can be introduced by a complex vertical diffusion coefficient.

This model recovers nonlinear effects in the system through a complex diffusion coefficient K_z introducing phase effects in the sesquilinear concentration distribution, as shown below. In proposing a new approach to advection-diffusion models with inclusion of the phase the authors are aware of the need to explore and understand its relation to the natural phenomenon. We assume that the current turbulence parameterizations applied to the vertical diffusion coefficient need calibration by experimental data in order to be applicable to real scenarios. An initial study focused on the appearance of fluctuations in the concentration distribution referring to the presence of coherent structures. Such a behaviour has never been reported in deterministic models, proving that the present modifications in the advection-diffusion equation and the resulting concentration representation comes closer to the physical description of the phenomenon. We also showed that the variation in the parameter responsible for the phase presence generates different fluctuation patterns in the concentration distributions, besides contributing to the dispersion characteristics of the pollutant. In this work we present a relation of the present model to the data from the Copenhagen experiment, which has low, moderate and high turbulent convective regimes. By simulating the experiment by this model it is possible to show that there are semi-quantitative evidences that the ratio of real to the imaginary part of the turbulence diffusion coefficient relates to the afore mentioned stability regimes.

2. The Advection-Diffusion Equation and New Closure

The pollutants dispersion models is based in the advection-difusion equation [13, p. 131]

(1)

where the variables are represented by averages and fluctuations terms

(2)

Rewriting the equation with Reynolds-averages [14, p. 531]

(3)

the terms and are known as turbulent flows. Solving now the equation analytically requires replacing them here by a modified Fick’s closure [7], i.e. the turbulent flows will be replaced by concentration gradients and a phase is introduced into the vertical diffusion coefficient K_z through a complex paramater in comparison to the traditional Fick closure [15]. Thus the equation of the model, with wind velocity in the x direction is given by

(4)

where C [g/m³] is the pollutant concentration, u[m/s] is the velocity in the direction x and K_x [m²/s], K_y [m²/s] and K_z [m²/s] are turbulent diffusion coefficients in the directions x, y and z, respectively. The initial condition determines null concentration for t = 0 [s], and three point sources are used as boundary conditions aligned to the y coordinate axis.

(5)

Here the source intensity is [g/s], the coordinate of the source along the y axis is y_p [m] and height is given by H_s [m], respectively. The three sources are at the same height and located at x = 0, with coordinates y in 0.1, 0, and −0.1 m. The flux of pollutants is considered null at the boundaries of a domain with dimensions L_x , L_y and L_z. The techniques used to obtain the solution are the separation of variables, using a Sturm-Liouville procedure in the directions y and z and the Laplace transform in x and also t.

The model’s solution is given in sesquilinear form

(6)

with

(7)

(8)

K_x and K_y are constants and K_z is complex given by

(9)

with K_za and K_zb constants. The A_nl represent the coefficients

(10)

(11)

3. Model Validation

The sesquilinear model is used now to simulate the Copenhagen experiment described in detail in reference [12]. The series of experiments provide runs in three stability regimes, classified by a criterion [17] where a ratio between the convective boundary layer height and the Monin-Obukhov length determine the type of test convection. A tracer substance sulphurhexafluoride (SF₆) was released without buoyancy from a 115-meter-high source at a constant flow rate ranging from 2.4 to 4.7 g/s and release time interval of 60 minutes. Measurements were taken at ground level where the terrain roughness is taken into consideration, being an urban region, that is, z₀ = 0.6 meters. Up to three series of samples were collected for each test and positioned between 2 to 6 km from the release location. In this work we analyzed only the maximum concentration values obtained in the measurements of the 9 experiments [18]. The parameter used in the model are presented in table 1.

Table 1. Meteorological data for the Copenhagen experiment [18]

Recalling, the present work focused on the investigation of the contribution of the phase in the turbulent diffusion coefficient and to relate to the micrometeorological characteristics.

To this end the set of parameters K_x, K_y, and K_zb were adjusted by parameter optimisation in order to represent best the experimental values. The optimisation procedure is performed in two steps by the least squares method. First, identical values were considered for K_x, K_y and K_za, with K_zb = 0. This hypothesis is used because for K_zb = 0 one recovers the usual deterministic model with its solution [10]. Ten parameter sets were simulated with numerical values between 10⁻⁷ and 10³ for each experiment. After having determined the optimal values, in a subsequent step K_x and K_y remain fixed whereas K_za and K_zb are varied freely in the afore mentioned range.

The statistical indices (normalised mean square error, fractional bias, fractional variance and correlation coefficient) [19] shown below were calculated to evaluate and compare both cases, with and without inclusion of a phase.

(12)

(13)

(14)

(15)

4. Results

The approach in the advection-diffusion model with modified Fick’s closure allowed to show correlations between a complex turbulent diffusion parameter (the presence of a phase) in the vertical diffusion coefficient. Recalling that the vertical coordinate is highly associated with the existence of turbulent processes in the atmosphere, variation of roughness, momentum, heat exchanges by soil irradiation and temperature gradient justifies the attempt to implement the proposed modification in the vertical component only.

The sesquilinear model concentration distributions do not show fluctuations in cases where the phase is equal to zero K_zb = 0, where the original deterministic model is recovered. Fluctuations in general arise for complex coefficients whose imaginary part is nonzero [10]. The concentration distributions that have fluctuations were evaluated by independent variation of parameters in the real and imaginary part of the vertical coefficient and the result is presented in figure 1.

Figure 1. Scatter plot for predicted (Cp) observed (Co) concentrations for the simulations of the Copenhagen experiment without the phase (Kzb = 0) in black and with the phase (Kzb ≠ 0) in red. The number represents the respective experiment the data are taken from

The variation of K_za does not interfere in the fluctuation patterns and contributes only to the increase of dispersion. However, the variation of parameters K_zb, associated to phase, presented changes in the concentration fluctuation patterns as well as a difference in the pollutant amount present in the distribution [11]. This behaviour pointed to a possible connection between the phase (imaginary to real coefficient ratio) and the planetary boundary layer stability regime. Thus the present study aims to assess the model behaviour through the simulation of the Copenhagen experiments [12].

The least squares technique was applied to calibrate the turbulent diffusion coefficient parameters K_x, K_y and K_z. Note, that the values used so far in the sesquilinear model, for observing the behaviour of the concentration distributions, were empirical. In this study we put focus on two situations, without and with the inclusion of the phase, respectively, K_zb = 0 and K_zb ≠ 0. The best least squares values found for the parameter variation described in section 2 are presented in tables 2 and 3, without and with the presence of the phase, respectively. The direct comparison of both tables shows that the minimum square values are smaller, that is, the model best describes the tests, in 8 of the 9 experiments when the phase is included.

Table 2. Diffusion coefficients for the sesquilinear model, that obtained better values through the least square method (LS) for the Copenhagen experiments (EXP) were classified in high, moderate and low convective regimes (HC, MC and LC)

Table 3. Diffusion coefficients of sesquilinear model, that obtained better value through the least square method (LS) for the Copenhagen experiments (EXP) were classified in high, moderate and low convection regimes (HC, MC and LC) where the standard (P) signifies (Respects - R, Does Not Respect - NR and Indiferent – IND)

By inspection one observes, that the indiferent cases correspond to the moderate convection experiments. The experiment that did not improve the minimum square value has only one measurement (see table 4 which shows the concentration values obtained by both situations). In addition to the minimum square value, performance of the two situations was evaluated through the study of observed (Co) and predicted (Cp) concentrations by the model.

Table 4. Comparison with observed results for the Copenhagen experiment [18] for the sesquilinear model without the phase (Kzb = 0) and with phase (Kzb ≠ 0)

Comparing the statistical indices of the models (see table 5) it is possible to see the fairly good results obtained by including the phase in the model and its potential, since the obtained indices are compatible with consolidated models. They are good results if we consider that little meteorological information was inserted in the model.

Table 5. Comparison of the statistical indices for the simulations of the Copenhagen experiment without the phase (Kzb = 0) and with phase (Kzb ≠ 0)

In table 3 a systematics may be observed that relates the parameters K_za and K_zb to the atmospheric stability in the respective experiment. Under high convection, i.e., when turbulence is more intense, K_zb has values greater than or equal to those of K_za. In the low convection regime there is an inversion, the parameter K_za becomes greater than or equal to K_zb. Seven of nine experiments were evaluated because they were classified with low or high convectivity. Of these 5 met the above described pattern, as can be sen in column P of table 3. Although a semi-quantitative result only, our findings indicate a sensitivity of the phase to stability which was impossible to detected by other models (see also reference [10]).

5. Conclusions

The advection-diffusion model with complex closure and sesqui-linear pollution concentration showed to be able to identify a semi-quantitative pattern for different stability regimes, which added to previously obtained results, confirms again its promising potential. The behavior of the sesqui-linear model shows local fluctuations in the concentration field which is also observed in field experiments. With decreasing stability in the planetary boundary layer the agreement between model predictions and experimental data was achieved by increasing contributions of an imaginary turbulent diffusion coefficient. A model inherent effect is that also the intensity of fluctuations in the density of concentrations increases. These findings were obtained as a result using the minimum square value, confirmed by the scatter plot, and are accompanied by a significant improvement of the statistical indices. Nevertheless, the model needs to be studied and calibrated with micro-meteorological parameters, since we showed sensitivity to changes in the stability regime.

ACKNOWLEDGEMENTS

This work was supported by Linhares Geração SA, Tevisa, CAPES - Brazil and CNPq - Brazil.