**1. Introduction**

The management and preservation of air quality require information of the environmental status. Such information involves both cognitive and interpretative features. Networks and measurements for monitoring, together with emission inventory of sources, are of great significance for the building of the cognitive representation, but not the interpretative one. Air quality control demands interpretative tools that are able to extrapolate in space and time the values measured by analytical instrumentation at field sites, while environmental improve‐ ment can only be obtained by means of emissions reduction from a systematic planning by using tools (mathematical models of dispersion) adequate to link the causes (sources) with the effects (pollutants).

The pollutants' transport and diffusion processes are complex and without mathematical models it is impossible to account for them. Such models therefore are an indispensable technical instrument of air quality management.

The theoretical approach to the atmospheric dispersion assumes different points of view. The diffusion in the K approach is proportional to the local gradient concentration from diffused material (at the fixed point in space). As a result, it is fundamentally Eulerian, once it assumes the motion of fluid within a fixed system of reference.

Lagrangian models differ from Eulerian ones in adopting a system of reference that follows atmospheric motions. This category contains all models that consider the pollutant cloud as discrete "elements" (puffs or computer particles). In Lagrangian approach, the dispersion is simulated through the motion of particles whose path allows the concentration field compu‐ tation from the liberated material.

© 2015 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

In this chapter, we consider the Eulerian approach.

#### **2. K models**

The Eulerian models takes into account the resolution of the mass conservation equation of pollutant chemical species [1]:

$$\frac{\partial \mathbf{c}}{\partial t} + \mu \nabla \mathbf{c} = D \nabla^2 \mathbf{c} + \mathbf{S}. \tag{1}$$

In (1), *u* is the wind speed vector, *D∇<sup>2</sup>* c is the molecular diffusion, D is the molecular diffusion coefficient, S is the term source, is an operator, and ∇<sup>2</sup> is the Laplacian.

To obtain the solution of equation (1), it is inevitable to know the wind field, something that is not possible because it is changing in space and time. As a consequence, wind is divided into two parts:

*u*¯, ensemble average

*u* ′ , turbulent fluctuation

Therefore, the wind speed can be expressed as the sum of the two contributions:

$$
\mu = \overline{\mu} + \mathfrak{u}'\tag{2}
$$

with *u*¯ '=0.

Similar hypothesis can be assumed for the concentration:

$$
\mathcal{L} = \overline{\mathcal{c}} + \mathcal{c}' \tag{3}
$$

with *c* ¯′ = 0.

The new *u* and *c* are introduced into (1). Then, the following is obtained:

$$\frac{\partial \overline{\overline{c}}}{\partial t} = -\overline{u} \cdot \nabla \overline{c} - \nabla \cdot \overline{c'u'} + D\nabla^2 \overline{c} + \overline{S}.\tag{4}$$

This equation possesses new unknown variables. This feature leads to a number of unknowns greater than the number of equations. By solving this problem it is possible to parameterize some terms with known quantities.

The widely used approach is the parameterization of second order moments. Such an approach is called as K-theory or flux-gradient theory:

$$
\overline{\mathbf{c'u'}} = -\mathbf{K} \nabla \overline{\mathbf{c}}\_{\prime} \tag{5}
$$

where K is the well-known eddy diffusivity.

A great variety of K-formulations exists [2]. The major part them take into account the similarity theory, and give different results for the same atmospheric stability.

Applying the following approximations:

**a.** Tensor is diagonal

In this chapter, we consider the Eulerian approach.

The Eulerian models takes into account the resolution of the mass conservation equation of

To obtain the solution of equation (1), it is inevitable to know the wind field, something that is not possible because it is changing in space and time. As a consequence, wind is divided

Therefore, the wind speed can be expressed as the sum of the two contributions:

*<sup>t</sup>* (1)

c is the molecular diffusion, D is the molecular diffusion

*uuu* = + ¢ (2)

*ccc* = + ¢ (3)

*<sup>t</sup>* (4)

is the Laplacian.

<sup>2</sup> . ¶ +Ñ= Ñ +

*<sup>c</sup> uc D c S*

¶

coefficient, S is the term source, is an operator, and ∇<sup>2</sup>

Similar hypothesis can be assumed for the concentration:

The new *u* and *c* are introduced into (1). Then, the following is obtained:

<sup>2</sup> . ¶ = - ×Ñ -Ñ× + Ñ + ¢ ¢ ¶ *<sup>c</sup> u c cu D c S*

This equation possesses new unknown variables. This feature leads to a number of unknowns greater than the number of equations. By solving this problem it is possible to parameterize

**2. K models**

330 Current Air Quality Issues

into two parts:

with *u*¯ '=0.

with *c* ¯′ = 0.

*u* ′

*u*¯, ensemble average

, turbulent fluctuation

some terms with known quantities.

pollutant chemical species [1]:

In (1), *u* is the wind speed vector, *D∇<sup>2</sup>*


Then, equation (4) can be written as:

$$\frac{\partial \overline{\overline{c}}}{\partial t} = -\overline{\mu} \cdot \nabla \overline{c} + \nabla \cdot K \nabla \overline{c} + \text{S.} \tag{6}$$

Equation (6) can be solved analytically or numerically, with data for *u*, *K*, and *S* and the initial and boundary conditions for*c*. This equation can be resolved in two ways:

**i.** With analytic methods, obtaining exact solutions

**ii.** With numerical methods, obtaining approximate solutions

In this chapter we take into consideration analytical solutions.

#### **3. Analytical solutions**

Analytical solutions of equations are very important to understand and describe physical phenomena. Moreover, they are fast, simple, and, generally, do not require complex meteoro‐ logical inputs.

There are analytical solutions of two-dimensional advection-diffusion equation [3]:

$$
\mu \frac{\partial \overline{\boldsymbol{c}} \langle \mathbf{x}, \boldsymbol{z} \rangle}{\partial \mathbf{x}} = \frac{\partial}{\partial \overline{\boldsymbol{z}}} \left( K\_z \frac{\partial \overline{\boldsymbol{c}} \langle \mathbf{x}, \boldsymbol{z} \rangle}{\partial \boldsymbol{z}} \right) + S\_\prime \tag{7}
$$

where *u* is longitudinal mean speed, *c*¯is the mean concentration, S is the source term, and *Kz* is the vertical eddy diffusivity. Furthermore, the longitudinal diffusion can be neglected because it is considered less in respect to the advection. Although, very recently a steady-state analytical solution for dispersion of contaminants in low winds by taking into account the longitudinal diffusion in the advection-diffusion equation was formulated [4].

The best-known solution is the so-called Gaussian solution, where both wind and turbulent diffusion coefficients are constant with height. So it is not a realistic solution of the transport and diffusion equation in the atmosphere. In the so-called Gaussian models, the solution is forced to represent real situations by means of empirical parameters, referred to as "sigmas." The different versions of Gaussian models substantially differ in the techniques utilized to calculate the "sigmas" as a function of atmospheric stability and the downwind distance from the emission source. Gaussian models are fast, simple, and do not require complex meteoro‐ logical inputs. For these reasons they are still widely used by all the environmental agencies over the world for regulatory applications. At this point, it is important to mention that there are models based on non-Gaussian analytical solutions.

[5] report a bi-dimensional solution considering cases where the wind speed and vertical eddy diffusivities follow power laws as a function of height:

$$\mathbf{u} \bullet \mathbf{u}\_1 \left(\mathbf{z}/\mathbf{z}\_1\right)^a \tag{8}$$

$$\mathbf{K}\_x = \mathbf{K}\_1 \left(\mathbf{z}/\mathbf{z}\_1\right)^\beta,\tag{9}$$

where z1 is the height when u1 and K1 are computed.

[6] report a bi-dimensional solution for linear profiles of the eddy diffusivity. [7] showed the bi-dimensional equation with u and K power functions of height with the exponents following the law of Schmidt.

[8] shows a solution considering u constant, but K as follows:

$$K\_z = K\_0 z^a \left( h - z \right)^b,\tag{10}$$

where Kzis a constant and *a* and *b* can be expressed as:

$$\begin{aligned} \mathbf{a}^3 \ 0 \ &\mathbf{b} = 0\\ \mathbf{a} \ \mathbf{a} = 0 \ \mathbf{;} \ \mathbf{b} > 0 \ \text{for } 0 \ \mathbf{E} \le \mathbf{z} \le \mathbf{h} \\ \mathbf{a} \ \mathbf{a} = 1 \ \vdots \ \mathbf{b} > 0 \ \text{for } 0 \le \mathbf{z} \le \mathbf{h} \\ \mathbf{a} \ \mathbf{a} = 1 \ \vdots \ \mathbf{b} = 0 \ \text{for } 0 \le \mathbf{z} \le \mathbf{h}/2 \ \vdots \ \mathbf{a} = 0 \ \text{and} \ \mathbf{b} = 1 \ \text{for } \mathbf{h}/2 \le \mathbf{z} \le \mathbf{h} \ \mathbf{h} \end{aligned}$$

where h is the ABL height.

[9] proposed a solution with constant *u* and Kz as:

On the Analytical Formulations for Pollutant Dispersion Simulation in the Atmospheric Boundary Layer http://dx.doi.org/10.5772/60023 333

$$\mathbf{K}\_x \equiv \mathbf{z} \quad \text{for} \; 0 \le \mathbf{z} \le \mathbf{z}\_s \tag{11}$$

$$\mathbf{K}\_x \mathbf{=} \mathbf{K}\_x(\mathbf{z}\_1) \text{ for } \mathbf{z}\_s \preceq \mathbf{z} \le \mathbf{h},\tag{12}$$

where zs is a predetermined height (generally, the height of the surface layer). This solution allows (as boundary conditions) a net flow of material toward the ground:

$$K\_z \frac{\partial \overline{\mathcal{C}}}{\partial \overline{z}} = V\_g \overline{c}\_{\prime} \tag{13}$$

where Vg is the deposition velocity.

because it is considered less in respect to the advection. Although, very recently a steady-state analytical solution for dispersion of contaminants in low winds by taking into account the

The best-known solution is the so-called Gaussian solution, where both wind and turbulent diffusion coefficients are constant with height. So it is not a realistic solution of the transport and diffusion equation in the atmosphere. In the so-called Gaussian models, the solution is forced to represent real situations by means of empirical parameters, referred to as "sigmas." The different versions of Gaussian models substantially differ in the techniques utilized to calculate the "sigmas" as a function of atmospheric stability and the downwind distance from the emission source. Gaussian models are fast, simple, and do not require complex meteoro‐ logical inputs. For these reasons they are still widely used by all the environmental agencies over the world for regulatory applications. At this point, it is important to mention that there

[5] report a bi-dimensional solution considering cases where the wind speed and vertical eddy

( ) α

( ) β

[6] report a bi-dimensional solution for linear profiles of the eddy diffusivity. [7] showed the bi-dimensional equation with u and K power functions of height with the exponents following

= - <sup>0</sup> ( ) ,

a = 1 ; b = 0 for 0 z h/2 ; a = 0 and b = 1 for h/2 z£ h,

£ £ £ £

1 1 u=u z/z (8)

K =K z/z , z1 1 (9)

*<sup>b</sup> <sup>a</sup> K Kz h z <sup>z</sup>* (10)

longitudinal diffusion in the advection-diffusion equation was formulated [4].

are models based on non-Gaussian analytical solutions.

diffusivities follow power laws as a function of height:

where z1 is the height when u1 and K1 are computed.

[8] shows a solution considering u constant, but K as follows:

where Kzis a constant and *a* and *b* can be expressed as:

a = 0 ; b > 0 for 0£ z h a = 1 ; b > 0 for 0 z h

£ £ £ £

a³ 0 ; b = 0

[9] proposed a solution with constant *u* and Kz as:

the law of Schmidt.

332 Current Air Quality Issues

where h is the ABL height.

[10, 11] report bi-dimensional solutions considering elevated sources with wind speed and vertical eddy diffusivity power profiles, but for an unbound atmosphere. That is:

$$K\_z \frac{\partial \overline{\mathcal{C}}}{\partial z} = 0 \text{ at } z \quad = \infty. \tag{14}$$

[12] shows a solution considering identical conditions, but for a vertically limited boundary layer. That is:

$$K\_z \frac{\partial \overline{c}}{\partial z} = 0 \text{ at } \mathbf{z} = \mathbf{h} \text{ .} \tag{15}$$

The solutions [10, 11, 12] are used in the KAPPAG model [13, 14, 15].

[16] using the Monin-Obukhov similarity theory to diffusion, has derived a solution from continuous sources near the ground with u and K follow power profiles. SPM [17] is a model that uses this solution.

[18] report a solution, which was a particular case of [7, 8]. In the sequence, [19] extended the solution to the situation of a growing ABL. [20] extended the solution to the case of nonzero mean vertical wind profiles.

[21] presented a three-dimensional atmospheric dispersion-deposition model for a groundlevel area source.

[22] extended the solution obtained in [12] with boundary conditions considering dry depo‐ sition to the ground.

In the work [23] were reported equations for point sources releases considering the first four moments of the vertical concentration distribution and the magnitude and downwind location of the maximum ground concentration.

[24] obtained a three-dimensional solution (with *u* and Kz constant) for low wind condition where the diffusion coefficients are function of down-wind distance from the source. [25] presented a two-dimensional solution with Kz constant and function of down-wind distance from the source, with a power low wind profile and for a finite boundary layer.

[26] obtained an analytical solution with dry deposition to the ground with any vertical function of wind and eddy coefficients but for a fixed vertical shape of contaminant concen‐ trations.

$$\overline{\sigma}\left(\mathbf{x},\mathbf{z}\right) = F\left(\mathbf{x}\right) \left(1 - \frac{\mathbf{z}}{h}\right)^2,\tag{16}$$

where *F*(*x*) is any function of *x*.

Using ADMM (Advection-Diffusion Multilayer Method) approach, [27] proposed a general solution for any wind and eddy coefficient profiles, but represented by a stepwise function in *z*. Recently, using the ADMM approach, [28] found a solution for two-dimensional steadystate solution considering eddy coefficient profile depending of *x* and *z* variables.

Finally, [29, 30] applying GILTT technique (Generalized Integral Laplace Transform Techni‐ que), found a general two-dimensional steady-state solution for any profiles of wind and eddy coefficient diffusions and a limited ABL.

#### **3.1. GILTT approach solution**

The General Integral Transform Technique (GITT) is a well-known hybrid method that solved a wide class of direct and inverse problems mainly in the area of Heat Transfer and Fluid Mechanics. There is a vast literature about this subject, including papers and books that become impossible to mention all of them. But, for instance, we mention the works of [31, 32, 33, 34], etc. In this chapter, we restrict our attention to the linear problem, because for nonlinear problem the linear result is iteratively used after the linearization of nonlinear transformed equation. Indeed, for the linear problems, the GITT transformed problem is analytically solved by the Laplace Transform technique. Consequently, the GITT solution becomes an analytical solution in the sense that no approximations are made along its derivation. It is important to notice that the GITT solution for nonlinear problems is semi-analytical because, in each iteration, the solution is analytical. This methodology is known as GILTT (Generalized Integral Laplace Transform Technique).

Bearing in mind that the novelty of the GILTT approach relies on the analytical solution of GITT linear transformed problem, in what follows, we restrict our analysis to the ensuing standard GITT time-dependent transformed problem:

$$E\frac{dY(t)}{dt} + FY\left(t\right) = 0 \qquad t > 0,\tag{17}$$

On the Analytical Formulations for Pollutant Dispersion Simulation in the Atmospheric Boundary Layer http://dx.doi.org/10.5772/60023 335

subject to the initial condition,

[24] obtained a three-dimensional solution (with *u* and Kz constant) for low wind condition where the diffusion coefficients are function of down-wind distance from the source. [25] presented a two-dimensional solution with Kz constant and function of down-wind distance

[26] obtained an analytical solution with dry deposition to the ground with any vertical function of wind and eddy coefficients but for a fixed vertical shape of contaminant concen‐

è ø

Using ADMM (Advection-Diffusion Multilayer Method) approach, [27] proposed a general solution for any wind and eddy coefficient profiles, but represented by a stepwise function in *z*. Recently, using the ADMM approach, [28] found a solution for two-dimensional steady-

Finally, [29, 30] applying GILTT technique (Generalized Integral Laplace Transform Techni‐ que), found a general two-dimensional steady-state solution for any profiles of wind and eddy

The General Integral Transform Technique (GITT) is a well-known hybrid method that solved a wide class of direct and inverse problems mainly in the area of Heat Transfer and Fluid Mechanics. There is a vast literature about this subject, including papers and books that become impossible to mention all of them. But, for instance, we mention the works of [31, 32, 33, 34], etc. In this chapter, we restrict our attention to the linear problem, because for nonlinear problem the linear result is iteratively used after the linearization of nonlinear transformed equation. Indeed, for the linear problems, the GITT transformed problem is analytically solved by the Laplace Transform technique. Consequently, the GITT solution becomes an analytical solution in the sense that no approximations are made along its derivation. It is important to notice that the GITT solution for nonlinear problems is semi-analytical because, in each iteration, the solution is analytical. This methodology is known as GILTT (Generalized Integral

Bearing in mind that the novelty of the GILTT approach relies on the analytical solution of GITT linear transformed problem, in what follows, we restrict our analysis to the ensuing

( ) += > ( ) 0 0,

*dt* (17)

*E FY t t*

*<sup>h</sup>* (16)

( ) ( ) <sup>2</sup> , 1, æ ö = - ç ÷

*<sup>z</sup> c xz F x*

state solution considering eddy coefficient profile depending of *x* and *z* variables.

from the source, with a power low wind profile and for a finite boundary layer.

trations.

334 Current Air Quality Issues

where *F*(*x*) is any function of *x*.

coefficient diffusions and a limited ABL.

**3.1. GILTT approach solution**

Laplace Transform Technique).

standard GITT time-dependent transformed problem:

*dY t*

$$Y\begin{pmatrix} 0 \end{pmatrix} = Y\_0. \tag{18}$$

Here, Y(t) is the unknown vector with N components, E and F are constant matrices of order (NxN), and Y0 is a known vector with N components. We begin our analysis recasting equation (17) like,

$$\frac{dY(t)}{dt} + GY(t) = 0 \qquad t > 0,\tag{19}$$

where G = E-1F and E-1 denotes the inverse of the matrix E. Now, applying the Laplace transform technique to equation (19), we obtain,

$$\overline{sY(s)} + G\overline{Y(s)} = Y\_{0'} \tag{20}$$

where *Y* ¯ (*s*) denotes the Laplace transform of the vector *Y* (*t*). Assuming that the matrix G is nondegenerate, we may write,

$$
\mathbf{G} = \mathbf{X} \, \mathbf{D} \mathbf{X}^{-1}.\tag{21}
$$

Here, D is the eigenvalue of the diagonal matrix, X is the matrix of the respective eigenfunc‐ tions, and X-1 its inverse.

Indeed, replacing equation (21) in equation (20) we get,

$$s\overline{Y(s)} + XDX^{-1}\overline{Y(s)} = Y\_{0'} \tag{22}$$

or

$$\overline{Y}\left(sI+XDX^{-1}\right)\overline{Y(s)} = Y\_{0'} \tag{23}$$

where the matrix I is the (NxN) identity matrix. Recalling that X X-1 = X-1X = I, equation (23) reads like,

$$\left(sXX^{-1} + XDX^{-1}\right)\overline{Y(s)} = Y\_{0'} \tag{24}$$

or

$$X\{sI+D\}X^{-1}\overline{Y\{s\}} = Y\_{0'} \tag{25}$$

which has the well-known solution,

$$
\overline{Y(s)} = X\left(sI + D\right)^{-1} X^{-1} Y\_0. \tag{26}
$$

Performing the Laplace transform inversion of equation (26), we have,

$$Y(t) = X \, L^{-1} \left\{ \left( sI + D \right)^{-1} \right\} X^{-1} Y\_0. \tag{27}$$

Here, L-1 denotes the inverse Laplace transform operator. We must notice that the matrix (*sI* + *D*) has the form,

$$\mathbf{A} \begin{pmatrix} sI + D \\ \end{pmatrix} = \begin{pmatrix} s + d\_1 & 0 & \cdots & 0 \\ 0 & s + d\_2 & \cdots & 0 \\ \vdots & \vdots & \ddots & \ddots \\ 0 & 0 & \cdots & s + d\_N \end{pmatrix} \tag{28}$$

where di are the eigenvalues of the matrix G. From the diagonal structure of the matrix G, we can straightly write down its inverse like:

$$\left(sI + D\right)^{-1} = \begin{pmatrix} 1 \\ s + d\_1 \\ \end{pmatrix} = \begin{pmatrix} 1 & 0 & \cdots & 0 \\ 0 & \frac{1}{s + d\_2} & \cdots & 0 \\ \vdots & \vdots & \ddots & \ddots \\ 0 & 0 & \cdots & \frac{1}{s + d\_N} \end{pmatrix} . \tag{29}$$

Performing the Laplace transform inversion of equation (29) by using the standard results of the Laplace transform theory we obtain:

$$L^{-1}\left\{\left(sI+D\right)^{-1}\right\}=E(t)=\begin{pmatrix}e^{-d\_1t} & 0 & \cdots & 0\\0 & e^{-d\_2t} & \cdots & 0\\\vdots & \vdots & \ddots & \ddots\\0 & 0 & \cdots & e^{-d\_Nt}\end{pmatrix}.\tag{30}$$

Finally, substituting the above ansatz in equation (27), we come out with the following solution for problem (17),

On the Analytical Formulations for Pollutant Dispersion Simulation in the Atmospheric Boundary Layer http://dx.doi.org/10.5772/60023 337

$$Y(t) = \mathbf{X} \begin{pmatrix} e^{-d\_1 t} & 0 & \cdots & 0 \\ 0 & e^{-d\_2 t} & \cdots & 0 \\ \vdots & \vdots & \ddots & \ddots \\ 0 & 0 & \cdots & e^{-d\_N t} \end{pmatrix} \mathbf{X}^{-1} \mathbf{Y}\_0. \tag{31}$$

To this point, it is relevant to underline that we are aware that problem (17) has a well-known solution. However, we must point out that the discussed solution is a robust algorithm, under computational point of view, to solve the problem with large N (N of order of 1, 500) and a small computational effort. Furthermore, this methodology can also be readily applied for the solution of equation (17) with boundary condition. This kind of a problem appears in the solution of discrete ordinates equation in a slab by the LTSN approach. For more information, see [35].

Next, we extend this methodology to the solution of more general linear ordinary differential equation having the entries of the matrices E and F varying with time. To reach this goal, let us consider next the problem:

$$E\left(t\right)\frac{dY\left(t\right)}{dt} + F\left(t\right)Y\left(t\right) = 0 \qquad t > 0,\tag{32}$$

subject to the initial condition,

( ) ( ) <sup>1</sup>

( ) ( ) <sup>1</sup> <sup>1</sup>

( ) { ( ) }

1

*s d s d sI D*

0 0

1

*s d*

1 1 1

Here, L-1 denotes the inverse Laplace transform operator. We must notice that the matrix

2

+ è ø

where di are the eigenvalues of the matrix G. From the diagonal structure of the matrix G, we

2

M M OO

<sup>1</sup> 0 0

+ è ø

Performing the Laplace transform inversion of equation (29) by using the standard results of

1

*d t*


*e*

0 0

Finally, substituting the above ansatz in equation (27), we come out with the following solution

<sup>1</sup> 0 0 .

L

L

L

0 0 ( ) .

*d t*

*N*

*s d*

0 0

L L M M OO

æ ö ç ÷

è ø


*d t*

L *<sup>N</sup>*

*e*

<sup>1</sup> 0 0

æ ö ç ÷ <sup>+</sup>

æ ö <sup>+</sup> ç ÷ +

0 0 0 0 ,

L *<sup>N</sup>*

*s d*

L L M M OO

Performing the Laplace transform inversion of equation (26), we have,

( )

( )

{ ( ) }

1

+ = +

<sup>2</sup> <sup>1</sup> <sup>1</sup>


+ ==

*<sup>e</sup> L sI D E t*

*sI D s d*


can straightly write down its inverse like:

the Laplace transform theory we obtain:

for problem (17),

+ =

which has the well-known solution,

(*sI* + *D*) has the form,

336 Current Air Quality Issues

<sup>0</sup> , - *X sI D X Y s Y* + = (25)

<sup>0</sup> . - - *Y s X sI D X Y* = + (26)

<sup>0</sup> . - - - *Y t X L sI D X Y* = + (27)

(28)

(29)

(30)

$$Y\begin{pmatrix} \mathbf{O} \end{pmatrix} = Y\_{\mathbf{o}}.\tag{33}$$

To solve this problem by the Laplace transform technique, we perform a stepwise approxi‐ mation of the matrices E(t) and F(t). So far, bearing in mind the interest of finding a solution for the time ranging from zero to a prescribed time T, we split the interval (0, T) into subinterval (tn-1, tn) for n=1:M, where M denotes the number of sub-intervals. For each sub-interval, we take the averaged values for the entries of matrix E(t) and F(t). It turns out that problem (28) simplifies to the recursive set of problems,

$$E^{m}\frac{dY^{m}\left(t\right)}{dt} + F^{m}Y^{m}\left(t\right) = 0,\tag{34}$$

for t in the sub-interval (tm-1, tm) and m ranging from 1 to M. Here, Em and Fm are now constant matrices. Henceforth, the previously discussed solution can be applied in a straightforward manner. Indeed, the solution for problem (34) generically reads like,

$$Y^{(m)} = X^{(m)} \begin{pmatrix} e^{-d\_1^{(m)}t} & 0 & \cdots & 0 \\ 0 & e^{-d\_2^{(m)}t} & \cdots & 0 \\ \vdots & \vdots & \ddots & \ddots \\ 0 & 0 & \cdots & e^{-d\_N^{(m)}t} \end{pmatrix} X^{(m)-1} Y^{(m-1)} \tag{35}$$

for t in the interval (tm-1, tm) and the initial condition Y(m-1) given by the previously calculated solution at tm-1 for t in the interval (tm-2, tm-1).

Focusing our attention on the task of searching analytical solution for the GITT linear trans‐ formed equation, in the sequel, we report an analytical alternative approach to solve equation (19), skipping the stepwise approximation of matrices E(t) and F(t) entries. To hit this objective, we solve problem (33) by the decomposition method proposed by [36]. In order to apply the decomposition method, let us recast equation (34) as,

$$
\overline{E}\frac{dY(t)}{dt} + \overline{F}Y(t) = A\left(t\right)\frac{dY(t)}{dt} + B\left(t\right)Y(t),\tag{36}
$$

where the matrices A(t) and B(t) are, respectively, expressed as,

$$A\left(t\right) = E\left(t\right) + \overline{E} \tag{37}$$

and

$$B(t) = F(t) + \overline{F}.\tag{38}$$

Here, *E*¯ and *F*¯ are constant matrices properly chosen. Now, expanding the solution of equation (36) in the truncated series,

$$Y(t) = \sum\_{k=0}^{L} \mu\_k(t) \tag{39}$$

and replacing this ansatz in equation (36), we get,

$$\overline{E}\,\frac{d}{dt}\sum\_{k=0}^{L}\boldsymbol{\mu}\_{k}\left(\boldsymbol{t}\right) + \overline{F}\sum\_{k=0}^{L}\boldsymbol{\mu}\_{k}\left(\boldsymbol{t}\right) = A\left(\boldsymbol{t}\right)\frac{d}{dt}\sum\_{k=0}^{L}\boldsymbol{\mu}\_{k}\left(\boldsymbol{t}\right) + B\left(\boldsymbol{t}\right)\sum\_{k=0}^{L}\boldsymbol{\mu}\_{k}\left(\boldsymbol{t}\right).\tag{40}$$

From equation (40), we are in a position to construct the following recursive set of linear ordinary differential matrix equations with constant matrices,

On the Analytical Formulations for Pollutant Dispersion Simulation in the Atmospheric Boundary Layer http://dx.doi.org/10.5772/60023 339

$$
\overline{E}\frac{d\,u\_0(t)}{dt} + \overline{F}u\_0\left(t\right) = 0,\tag{41}
$$

$$\overline{E}\frac{du\_1(t)}{dt} + \overline{F}u\_1(t) = A\left(t\right)\frac{du\_0(t)}{dt} + B\left(t\right)u\_0(t),\tag{42}$$

$$
\overline{E}\frac{d\boldsymbol{u}\_{2}\left(t\right)}{dt} + \overline{F}\boldsymbol{u}\_{2}\left(t\right) = A\left(t\right)\frac{d\boldsymbol{u}\_{1}\left(t\right)}{dt} + B\left(t\right)\boldsymbol{u}\_{1}\left(t\right) \tag{43}
$$

and so forth. Generically, we may write,

( ) 1

*d t*

*m*


*e*

=

338 Current Air Quality Issues

solution at tm-1 for t in the interval (tm-2, tm-1).

and

(36) in the truncated series,

and replacing this ansatz in equation (36), we get,

ordinary differential matrix equations with constant matrices,

decomposition method, let us recast equation (34) as,

where the matrices A(t) and B(t) are, respectively, expressed as,

0 0

( ) 2

*d t m m m m*

*<sup>e</sup> Y X X Y*

*m*

æ ö ç ÷


() () ( ) 1 ( 1)

è ø

0 0

L L M M OO

L

for t in the interval (tm-1, tm) and the initial condition Y(m-1) given by the previously calculated

Focusing our attention on the task of searching analytical solution for the GITT linear trans‐ formed equation, in the sequel, we report an analytical alternative approach to solve equation (19), skipping the stepwise approximation of matrices E(t) and F(t) entries. To hit this objective, we solve problem (33) by the decomposition method proposed by [36]. In order to apply the

> ( ) ( ) ( ) ( ) += + ( ) ( ), *dY t dY t E FY t A t B t Y t*

Here, *E*¯ and *F*¯ are constant matrices properly chosen. Now, expanding the solution of equation

( ) ( ) ( ) ( ) ( ) ( )

From equation (40), we are in a position to construct the following recursive set of linear

*dt dt* (40)

00 0 0

== = = åå å å += + *LL L L kk k k kk k k d d E u t F u t At u t Bt u t*

( ) ( ) =0 <sup>=</sup> å *L k k*

( )

*d t*

*m N*


*dt dt* (36)

*At Et E* ( ) = + ( ) (37)

*Bt Ft F* ( ) = + ( ) . (38)

*Yt u t* (39)

.

(35)

0 0 ,


*e*

$$\overline{E}\frac{d\boldsymbol{u}\_{n}\left(t\right)}{dt} + \overline{F}\boldsymbol{u}\_{n}\left(t\right) = A\left(t\right)\frac{d\boldsymbol{u}\_{n-1}\left(t\right)}{dt} + B\left(t\right)\boldsymbol{u}\_{n-1}\left(t\right),\tag{44}$$

for n = 1:L. Given a closer look at the above equations, we readily realize that the RHS of the recursive set of equations from (42) to (44) are known terms. Therefore, to obtain the solution of these equations, we must solve the following problem with constant matrices and source, namely,

$$
\overline{E}\frac{dY(t)}{dt} + \overline{F}Y(t) = \mathcal{S}\left(t\right) \qquad t > 0,\tag{45}
$$

which has the well-known solution,

$$Y(t) = Y\_h(t)Y\_0 + Y\_h(t) \* S(t),\tag{46}$$

whereas the homogeneous solution Yh(t) is given by equation (35). Here, star denotes convo‐ lution. Therefore, the solution of the generic equation (44) has the form,

$$\begin{aligned} Y(t) &= \mathbf{X} \begin{pmatrix} e^{-d\_1 t} & 0 & \cdots & 0 \\ 0 & e^{-d\_2 t} & \cdots & 0 \\ \vdots & \vdots & \ddots & \ddots \\ 0 & 0 & \cdots & e^{-d\_N t} \end{pmatrix} \mathbf{X}^{-1} \mathbf{Y}\_0 + \\ & \quad \begin{pmatrix} \mathbf{0} & \mathbf{0} & \cdots & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \cdots & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & e^{-d\_2 t} & \cdots & \mathbf{0} \\ \vdots & \vdots & \ddots & \ddots \\ \mathbf{0} & \mathbf{0} & \cdots & e^{-d\_N t} \end{pmatrix} \mathbf{X}(t')^{-1} \mathbf{S}(t - t') \end{pmatrix} \tag{47}$$

and consequently the solution of equation (32) is well determined by equations (39) and (41) to (44).

Finally, regarding the task of solving second-order ordinary differential transformed equation, we must recall that this equation can be transformed into a set of first-order matrix differential equation by changing the variable. Indeed, performing the substitution of,

$$Y(t) = Z\_1(t) \tag{48}$$

and

$$\frac{dY(t)}{dt} = Z\_2(t),\tag{49}$$

in the following transformed problem:

$$\frac{d^2Y(t)}{dt} + A\left(t\right)\frac{dY(t)}{dt} + B\left(t\right)Y(t) = S\left(t\right),\tag{50}$$

we come out with the ensuing result:

$$
\frac{d}{dt} \begin{pmatrix} Z\_1(t) \\ Z\_2(t) \end{pmatrix} + \begin{pmatrix} 0 & -1 \\ B(t) & A(t) \end{pmatrix} \begin{pmatrix} Z\_1(t) \\ Z\_2(t) \end{pmatrix} = \begin{pmatrix} 0 \\ S(t) \end{pmatrix} \tag{51}
$$

which is a set of first-order ordinary matrix differential equation that can be promptly solved by the previous discussed methods for nondegenerated matrix. Finally, to handle transformed problems with degenerated matrix, we proceed by performing the inversion of the symbolic matrix by the Schur decomposition approach. Furthermore, we make the Laplace transform inversion applying the Heaviside expansion formulation for eigenvalues with multiplicity larger than one. Hoping that we have completed the mathematical analysis regarding the task of solving analytically the GITT transformed problems, in the next section, we illustrate the GILTT aptness to simulate pollutant dispersion in the atmosphere by solving the twodimensional steady-state advection-diffusion equation.

#### *3.1.1. An example: Solution of the 2D advection-diffusion equation*

To solve the problem (7) by the GILTT method we recast equation (7) as (*S* = 0):

$$
\mu(z)\frac{\partial \overline{c}}{\partial x} = K\_z(z)\frac{\partial^2 \overline{c}}{\partial z^2} + \left(\frac{\partial K\_z(z)}{\partial z}\right)\frac{\partial \overline{c}}{\partial z} \tag{52}
$$

On the Analytical Formulations for Pollutant Dispersion Simulation in the Atmospheric Boundary Layer http://dx.doi.org/10.5772/60023 341

with the following boundary conditions:

and consequently the solution of equation (32) is well determined by equations (39) and (41)

Finally, regarding the task of solving second-order ordinary differential transformed equation, we must recall that this equation can be transformed into a set of first-order matrix differential

> ( ) <sup>=</sup> <sup>2</sup> ( ), *dY t Z t*

( ) ( ) ( ) ( ) ( ) ( ) <sup>2</sup> + += , *d Y t dY t At Bt Yt St*

( ) ( ) 1 1

which is a set of first-order ordinary matrix differential equation that can be promptly solved by the previous discussed methods for nondegenerated matrix. Finally, to handle transformed problems with degenerated matrix, we proceed by performing the inversion of the symbolic matrix by the Schur decomposition approach. Furthermore, we make the Laplace transform inversion applying the Heaviside expansion formulation for eigenvalues with multiplicity larger than one. Hoping that we have completed the mathematical analysis regarding the task of solving analytically the GITT transformed problems, in the next section, we illustrate the GILTT aptness to simulate pollutant dispersion in the atmosphere by solving the two-

01 0 , æö æö æ ö æö - ç÷ ç÷ + = ç ÷ ç÷ ç÷ ç÷ èø èø è ø èø

*Yt Z t* ( ) = <sup>1</sup> ( ) (48)

*dt* (49)

*dt dt* (50)

*dt Z t Bt At Z t S t* (51)

( )

equation by changing the variable. Indeed, performing the substitution of,

to (44).

340 Current Air Quality Issues

and

in the following transformed problem:

we come out with the ensuing result:

( )

dimensional steady-state advection-diffusion equation.

*3.1.1. An example: Solution of the 2D advection-diffusion equation*

( ) ( ) ( )

*d Z t Z t*

2 2

To solve the problem (7) by the GILTT method we recast equation (7) as (*S* = 0):

*z*

2 2 ( ) () () ¶¶ ¶ æ ö ¶ = + ç ÷ ¶ ¶ è ø ¶ ¶

*cc c K z uz K z x z z z*

*z*

(52)

$$K\_z \frac{\partial \overline{c}}{\partial \mathbf{x}} = \mathbf{0} \text{ at } \mathbf{z} = \mathbf{h} \tag{53}$$

$$
\mu \overline{c} \left( \mathbf{x}, z \right) = \mathbb{Q} \delta \left( z - H\_s \right) \text{at } \mathbf{x} = \mathbf{0} \tag{54}
$$

and we construct the following associated Sturm-Liouville problem:

$$
\Psi\_i''(z) + \lambda\_i^2 \Psi\_i(z) = 0 \text{ at } 0 \le \mathbf{z} \le \mathbf{h} \tag{55}
$$

$$
\Psi'\_{\,\,i}(\mathbf{z}) = 0 \text{ at } \mathbf{z} = \mathbf{0}, \mathbf{h} \text{ .}\tag{56}
$$

which has the well-known solution:

$$\Psi\_{\wedge}(z) = \cos(\mathcal{A}\_{\mathbb{P}} z),\tag{57}$$

where the eigenvalues *λ<sup>i</sup>* are given by *λ<sup>i</sup>* <sup>=</sup> *<sup>i</sup><sup>π</sup> <sup>h</sup>* for *<sup>i</sup>* =0, 1, 2, 3, ….Furthermore, The eigenfunc‐ tions *Ψ<sup>i</sup>* (*z*) satisfy the ensuing orthonormality condition:

$$\frac{1}{\mathbf{N}\_m^{1/2}\mathbf{N}\_n^{1/2}} \Big| \Psi\_m(\mathbf{z}) \Psi\_n(\mathbf{z}) d\upsilon = \begin{cases} 0, \text{ m} \neq \text{ n} \\ 1, \text{ m} = \text{ n}' \end{cases} \tag{58}$$

where N*m* is expressed by:

$$\mathbf{N}\_m = \int\_{\mathbf{v}} \Psi\_m^2 \left( z \right) d\boldsymbol{v} \,. \tag{59}$$

Having solved the Sturm-Liouville problem, we are in a position to construct the GILTT transform formula, which has the form:

$$\overline{c}\left(\mathbf{x},\mathbf{z}\right) = \sum\_{i=0}^{\infty} \frac{\overline{c\_i(\mathbf{x})} \Psi\_i\left(\mathbf{z}\right)}{\mathbf{N}\_i^{1/2}}.\tag{60}$$

Now, replacing the above ansatz in equation (52) we have,

$$
\mu \sum\_{i=0}^{\alpha} \frac{\overline{c\_i'(\mathbf{x})} \Psi\_i(\mathbf{z})}{\mathbf{N}\_i^{1/2}} = \mathbf{K}\_z \sum\_{i=0}^{\alpha} \frac{\overline{c\_i(\mathbf{x})} \Psi\_i''(\mathbf{z})}{\mathbf{N}\_i^{1/2}} + \left(\frac{\partial \mathbf{K}\_z}{\partial \mathbf{z}}\right) \sum\_{i=0}^{\alpha} \frac{\overline{c\_i(\mathbf{x})} \Psi\_i'(\mathbf{z})}{\mathbf{N}\_i^{1/2}}.\tag{61}
$$

Here, we adopt the prime notation for the first derivative and double the prime notation for the second derivative. Multiplying equation (61) by *Ψj* (*z*) N *<sup>j</sup>* 1/2 and integrating from *z* equal zero to *h*, we read:

$$\begin{split} & -\sum\_{i=0}^{\alpha} \frac{\overline{\overline{c\_{i}}^{\prime}(\mathbf{x})} \lambda}{\mathbf{N}\_{i}^{1/2} \mathbf{N}\_{j}^{1/2}} \prod\_{0}^{h} u \boldsymbol{\Psi}\_{i} \boldsymbol{\Psi}\_{j} \, d\mathbf{z} - \sum\_{i=0}^{\alpha} \frac{\overline{\overline{c\_{i}}(\mathbf{x})} \lambda}{\mathbf{N}\_{i}^{1/2} \mathbf{N}\_{j}^{1/2}} \prod\_{0}^{h} \boldsymbol{\Psi}\_{z} \boldsymbol{\Psi}\_{i} \boldsymbol{\Psi}\_{j} \, d\mathbf{z} + \\ & + \sum\_{i=0}^{\alpha} \frac{\overline{c\_{i}(\mathbf{x})} \lambda}{\mathbf{N}\_{i}^{1/2} \mathbf{N}\_{j}^{1/2}} \prod\_{0}^{h} \left( \frac{\partial \mathcal{K}\_{z}}{\partial \mathbf{z}} \right) \boldsymbol{\Psi}\_{i}^{\prime} \boldsymbol{\Psi}\_{j} \, d\mathbf{z} = 0 \end{split} \tag{62}$$

for*j* =0, 1, 2, …, *N* . Rewriting equation (62) in matrix fashion, we obtain,

$$Y'(\mathbf{x}) + FY(\mathbf{x}) = \mathbf{0},\tag{63}$$

where Y(x) is the column vector whose components are cj (x) and the matrix F is defined as *F* = *B* <sup>−</sup><sup>1</sup> *E*. Here, the entries of matrices B and E are, respectively, given by:

$$\boldsymbol{b}\_{i,j} = \bigwedge\_{0}^{h} \boldsymbol{\mu} \boldsymbol{\Psi}\_{i} \boldsymbol{\Psi}\_{j} \, d\boldsymbol{z} \tag{64}$$

and

$$\varepsilon e\_{i,j} = \int\_0^h \frac{\partial K\_z}{\partial z} \Psi\_i^{\prime \prime} \Psi\_j d\boldsymbol{z} - \lambda\_l^2 \int\_0^h \mathbf{K}\_z \Psi\_j^{\prime} \Psi\_j d\boldsymbol{z}.\tag{65}$$

Transforming the boundary condition (54) by a similar procedure, we mean, multiplying the boundary condition (54) by *Ψj* (*z*) N *<sup>j</sup>* 1/2 and integrating from z = 0 to h, we obtain:

$$\int\_{0}^{h} \sum\_{i=0}^{n} \frac{\overline{c\_{i}(0)} \, u \Psi\_{i} \Psi\_{j}}{\mathbf{N}\_{i}^{1/2} \mathbf{N}\_{j}^{1/2}} dz = \int\_{0}^{h} \frac{\mathcal{Q} \, \mathcal{S} \, \left(z \cdot H\_{s}\right) \Psi\_{j}}{\mathbf{N}\_{j}^{1/2}} dz. \tag{66}$$

Using the orthonormality property of the eigenfunctions *Ψ<sup>j</sup>* (*z*), and the integration property of the generalized delta function, we get the following transformed boundary condition:

On the Analytical Formulations for Pollutant Dispersion Simulation in the Atmospheric Boundary Layer http://dx.doi.org/10.5772/60023 343

$$\mathcal{L}\_0(\mathbf{0}) = \frac{\mathbb{Q}\Psi\_0^{\nu}(H\_s)\mathbf{\check{M}}}{\mu(z)\Psi\_0^2(z)dz} \tag{67}$$

for j = 0, and

( ) ( ) ( ) ( ) ( ) ( ) 1/2 1/2 1/2

¢ YY Y ¢¢ æ ö ¶ ¢ = + ç ÷ è ø ¶ åå å *i i i i <sup>z</sup> i i*

. NN N

Here, we adopt the prime notation for the first derivative and double the prime notation for

*z*

*Ψj* (*z*) N *<sup>j</sup>*

2

*Y x FY x* ¢( ) + = ( ) 0, (63)

*ij i j b u dz* (64)

*i j zij*

*u dz K dz*

(61)

(62)

(65)

1/2 and integrating from *z* equal zero

(x) and the matrix F is defined as

(*z*), and the integration property of

00 0

( ) ( )

<sup>λ</sup> N N N N

*c x c x*

*i i i j i j*

å å ò ò

æ ö ¶ + YY = ¢ ç ÷ è ø ¶

*c x <sup>K</sup> dz z*

for*j* =0, 1, 2, …, *N* . Rewriting equation (62) in matrix fashion, we obtain,

*E*. Here, the entries of matrices B and E are, respectively, given by:

0 0

*i j ij i zij*

*e dz K dz*

( ) ( ) 1/2 1/2 1/2

the generalized delta function, we get the following transformed boundary condition:

*i ij s j*

0 . N N <sup>N</sup>

Y Y Y

d

*h h*

2

l. ¶ = YY - YY ¢ ¶ ò ò

Transforming the boundary condition (54) by a similar procedure, we mean, multiplying the

1/2 and integrating from z = 0 to h, we obtain:

*c u Q z-H dz dz* (66)

, 0 = YY ò *h*

*z*

*K*

*z*

0 0 0

*i i j j*

å <sup>=</sup> ò ò *h h*

¥ ¥ = =

> *h i z*

1/2 1/2 1/2 1/2 0 0 0 0

0

*i j*

¢ - YY - YY +

*h h i i i*

*ii i ii i cx z cx z cx z <sup>K</sup> u K*

¥¥ ¥ == =

*z*

the second derivative. Multiplying equation (61) by

( )

1/2 1/2 0 0

where Y(x) is the column vector whose components are cj

,

*Ψj* (*z*) N *<sup>j</sup>*

Using the orthonormality property of the eigenfunctions *Ψ<sup>j</sup>*

¥

=

N N

*i i j*

å ò

¥ =

to *h*, we read:

342 Current Air Quality Issues

*F* = *B* <sup>−</sup><sup>1</sup>

and

boundary condition (54) by

$$c\_i(0) = \frac{\mathbb{Q}\Psi\_{\rangle}(H\_s)\bigvee\_2^h}{\int\_0^h \mu(z)\Psi\_{\rangle}^2(z)dz} \tag{68}$$

for j ranging from 1 to N. To this point, we are in a position to write down the solution of equation (39), subjected to the boundary condition given by equations (67) and (68), using the results discussed in the previous section. Indeed, the solution writes like,

$$\overline{c}\left(\mathbf{x},\mathbf{z}\right) = \sum\_{i=0}^{\circ} \frac{\overline{c\_i(\mathbf{x})} \mathbb{1}\_i(\mathbf{z})}{\mathbf{N}\_i^{1/2}}.\tag{69}$$

Here, ci (x) are the N components of the vector Y(x) given by,

$$Y(\mathbf{x}) = \mathbf{Z} \begin{pmatrix} e^{-f\_1 \mathbf{x}} & 0 & \cdots & 0 \\ 0 & e^{-f\_2 \mathbf{x}} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & e^{-f\_N \mathbf{x}} \end{pmatrix} \mathbf{G}(\mathbf{x}) \mathbf{Z}^{-1} \mathbf{Y} \begin{pmatrix} 0 \end{pmatrix} . \tag{70}$$

where fk, for k ranging from 1 to N, are the eigenvalues of the matrix of eigenfunctions and *Z* <sup>−</sup><sup>1</sup> its inverse. For more details and recent developments, see [29, 30, 37, 38, 39, 40, 41].

#### **4. Remarks**

To deal with more realistic situations we need to change to numerical methods. However, it is helpful to check firstly possible analytical solutions in order to obtain a known framework and test solutions. The analytical solutions are efficient for many applications, for example: supply analysis of contaminant scenarios, leading sensitivity analyses to investigate the effects of parameters included in contaminant transport, extrapolation over large times and distances where numerical solutions may be impractical, as screening models for transport processes that cannot be solved exactly, and for validating numerical solutions.

Today, air pollution problems are not treated in the manner described in the present chapter. There are various air pollution situations that require the use of complex mesoscale models to properly describe the dispersion processes and properly represent the relevant chemistry and emission processes. Complex models, such as the CMAQ model (the Community Multiscale Air Quality model), have been designed to simulate air quality by including state of the art techniques for modeling multiple air quality issues. However, in complex models, increasingly more processes, such as sea breeze circulations, urban heat islands, and waves, are represented. Therefore, these models are often perceived as "closed" that cannot without difficulty or effort report the influence of individual processes on air quality. Finally, for many scientific appli‐ cations, analytical solutions have utility in understanding air dispersion phenomena and some air chemistry phenomena, showing their usefulness in environmental management.
