**Modeling Aerosol Processes**

## Jose Ignacio Huertas

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/62021

The term aerosol refers to an ensemble of liquid or solid particles suspended in a gaseous medium. The study of aerosols is relevant to a variety of fields, including air pollution, combustion, and chemical engineering. Since the 1950s, advances in aerosol science have been driven by investigations into the health effects of radioactive aerosols and industrial aerosols in the workplace and the environment. More recently, much effort has gone into understanding the effect of various natural and anthropomorphic aerosols on global warming. These efforts have largely been aimed at reducing the adverse effects of aerosol particles.

On the other hand, a great deal of knowledge has been recently gained by researchers working with *desirable anthropomorphic aerosols* used to produce nanophase high-technology materials such as ceramic powders, superconducting materials, and optical fibers. In this regard, aerosol processes are commercially employed for large-scale production of particle commodities like carbon black and inorganic powders (titania, zinc oxide, and fumed silica) and high-purity materials (optical fibers, thin films for micro-electronics, and advanced ceramics). The annual U.S. production amounts to hundreds of millions of tons. [23]

## **1. The General Dynamic Equation (GDE) for single-component aerosols**

To understand the detailed behavior of aerosol systems, it is necessary to model the dynamic behavior of the aerosol population. Several mathematical models have been developed to simulate the physical processes that affect aerosols. Such processes include nucleation, coagulation, growth due to condensation of gases, shrinkage due to evaporation, sedimenta‐ tion, and deposition on surfaces. Assuming that a particle in the aerosol can be described by a variable υ to represent its size (e.g., mass, vo1ume, or diameter), the general dynamic equation (GDE) that describes the behavior of a single component aerosol is: [50]

$$\frac{\partial n}{\partial t} + \nabla \cdot \mathbf{V}n = \nabla \cdot \boldsymbol{D} \nabla n + \left(\frac{\partial n}{\partial t}\right)\_{\text{growth}} + \left(\frac{\partial n}{\partial t}\right)\_{\text{coulomb}} + S\_o \mathcal{S} \left(\boldsymbol{\nu} - \boldsymbol{\nu}^\*\right) - \boldsymbol{R} \tag{1}$$

© 2016 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.

where n(υ,χ,t)dυ is the continuous distribution function that specifies the concentration of particles in the size range [υ, υ+dυ] at position *x* and time t. This equation is obtained through a material balance on a control volume fixed in space with respect to laboratory coordinates. The convective term (∇ ⋅*V n*) gives the change in the distribution function due to the bulk flow of the aerosol. The diffusion term (∇ ⋅*D*∇*n*) accounts for changes in *n* due to gradients in particle concentrations. *R* represents the removal processes such as sedimentation or deposi‐ tion on surfaces. *Soδ*(*<sup>υ</sup>* <sup>−</sup>*<sup>υ</sup>* \*) represents the rate of generation of new particles where δ is the delta function. New particles may be generated by homogeneous nucleation, or any other mechanism. For example, a chemical reaction front or a combustion flame can be treated as sources of clusters, particles, and/or gas-phase products.

These four terms represent processes that may be called external since they involve movement across the walls of an elemental volume. However, particles in a given size range may also be modified as a result of internal processes such as growth and coagulation. The growth term in Equation (1) describes the rate of change of the distribution function *n* due to particle growth by condensation of the vapor phase material. The coagu1ation term describes the changes in the distribution function due to coalescing collisions between particles. In Chapters 4 and 5 we will present expressions for the coagulation and the condensation term, respectively.

#### **2. Multicomponent aerosols**

Even though common aerosols, e.g., atmospheric aerosols, are made of several species rather thana single one, little work has been reported on modeling the evolution of multicomponent to model the particle encapsulation process, and the GDE needs to be generalized to account for processes involving several species. To do so, a new variable **Y** is introduced to represent composition. **Y** is a p-dimensional vector where the Y<sup>i</sup> represents the mass fraction of the ith species present in the particle, and p is the total number of species present in the aerosol. Y<sup>i</sup> is given by:

$$Y\_i = \frac{\text{Mass of species i present in the particle}}{\text{Total particle mass}}$$

such that O < Y<sup>i</sup> < 1, and Σ Y<sup>i</sup> = 1. Now, Equation (1) is the appropriate equation for modeling the particle dynamics with n(υ,**Y**,**x**,t)dυd**Y** being the continuous distribution function that specifies the concentration of particles in the size range [υ,υ+dυ] with composition in the range [**Y**, **Y**+d**Y**] at position **x** and time t in the aerosol.

In general, particles of the same size may have different chemical composition and the composition of a particle affects its growth rate. Thus, the computation of size distribution should generally not be performed independent of particle composition. Despite this, it is reasonable to assume that the rate of collision depends only on the sizes of the particles involved and not on their composition. Consequently, the usual expressions for coagulation, diffusion coefficient D, and collision frequency function β can be used. This also implies that interactions such as attractive or repulsive forces developed between particles of different materials are not considered.

Similarly to the case of single component aerosols, the dynamics of particle collisions is reduced to the assumption that the particles stick together with a probability α, which is known as the accommodation or sticking coefficient. For simplicity, α is assumed to be equal to 1. In addition, it is assumed that no chemical reaction occurs among the components in the particles. Consequently, when a collision between a particle of mass υ<sup>i</sup> and composition Yi with another particle of mass υ<sup>j</sup> and composition Yj occurs, a new particle of mass (υ<sup>i</sup> +υ<sup>j</sup> ) and composition Yk = (Yi\*υ<sup>i</sup> + υ<sup>j</sup> \*Yj )/(υ<sup>i</sup> + υ<sup>j</sup> ) is formed.

## **3. Solution methods**

where n(υ,χ,t)dυ is the continuous distribution function that specifies the concentration of particles in the size range [υ, υ+dυ] at position *x* and time t. This equation is obtained through a material balance on a control volume fixed in space with respect to laboratory coordinates. The convective term (∇ ⋅*V n*) gives the change in the distribution function due to the bulk flow of the aerosol. The diffusion term (∇ ⋅*D*∇*n*) accounts for changes in *n* due to gradients in particle concentrations. *R* represents the removal processes such as sedimentation or deposi‐ tion on surfaces. *Soδ*(*<sup>υ</sup>* <sup>−</sup>*<sup>υ</sup>* \*) represents the rate of generation of new particles where δ is the delta function. New particles may be generated by homogeneous nucleation, or any other mechanism. For example, a chemical reaction front or a combustion flame can be treated as

These four terms represent processes that may be called external since they involve movement across the walls of an elemental volume. However, particles in a given size range may also be modified as a result of internal processes such as growth and coagulation. The growth term in Equation (1) describes the rate of change of the distribution function *n* due to particle growth by condensation of the vapor phase material. The coagu1ation term describes the changes in the distribution function due to coalescing collisions between particles. In Chapters 4 and 5 we will present expressions for the coagulation and the condensation term, respectively.

Even though common aerosols, e.g., atmospheric aerosols, are made of several species rather thana single one, little work has been reported on modeling the evolution of multicomponent to model the particle encapsulation process, and the GDE needs to be generalized to account for processes involving several species. To do so, a new variable **Y** is introduced to represent composition. **Y** is a p-dimensional vector where the Y<sup>i</sup> represents the mass fraction of the ith species present in the particle, and p is the total number of species present in the aerosol. Y<sup>i</sup>

> *Mass of species i present in the particle <sup>Y</sup> Total particle mass* <sup>=</sup>

the particle dynamics with n(υ,**Y**,**x**,t)dυd**Y** being the continuous distribution function that specifies the concentration of particles in the size range [υ,υ+dυ] with composition in the range

In general, particles of the same size may have different chemical composition and the composition of a particle affects its growth rate. Thus, the computation of size distribution should generally not be performed independent of particle composition. Despite this, it is reasonable to assume that the rate of collision depends only on the sizes of the particles involved and not on their composition. Consequently, the usual expressions for coagulation,

= 1. Now, Equation (1) is the appropriate equation for modeling

 *<sup>i</sup>* is

sources of clusters, particles, and/or gas-phase products.

38 Montecarlo Simulation of Two Component Aerosol Processes

**2. Multicomponent aerosols**

< 1, and Σ Y<sup>i</sup>

[**Y**, **Y**+d**Y**] at position **x** and time t in the aerosol.

given by:

such that O < Y<sup>i</sup>

From a mathematical point of view, the modeling of multicomponent aerosols and the particle encapsulation process reduces to solving the GDE under the assumptions described in previous sections. The GDE, being a nonlinear integro-differential equation, will be solved numerically. To reach the objectives of this dissertation, the numerical method chosen must:Model coagulation and condensationTreat multicomponent aerosolsInclude the Kelvin effect, particle heating due to condensation, and heat loss from the aerosol

#### **3.1. Solution method for single-component aerosol**

Several works have been conducted to solve the GDE. Much of the effort to solve this equation has involved conditions where only coagulation is important. This information is of particular importance because the nonlinearity of the GDE comes from the coagulation term. Reference 58, which is a comparative review of mathematical models for aerosol dynamics, considers four approaches to the solution of the GDE. They are classified according to their representa‐ tion of the size distribution and differ only in their degree of approximation of the size distribution.


Those techniques that offer any desired degree of approximation are referred to as *continu‐ ous.* A typical example is the method known as the J-Space method. [51] This method is based on the use of the transformation m (**J**) = m(**J**o)\*Exp{ α\*(J-Jo) }, where m is mass, J and Jo are positive integers such that J ≥ Jo and α is an adjustable parameter. The coagulation equation is expressed in terms of this transformation, resulting in a set of integro-differential equations, one for each value of J. The set of equations are solved through the cubic spline interpolation method. The method is useful for numerical solution of continuous distributions over a wide particle size range; however, it has inherent problems associated with mass conservation and numerical dispersion. Numerical dispersion appears as waves ahead of or behind regions of high concentration.

*Parameterized* representations are those methods in which the size distribution is assumed to adhere to a known "self-preserving" form described by a small number of parameters. Solutions found in this way are believed to be asymptotic forms approached after long times. The assumption reduces the problem to solve a set of few ODEs for the parameters representing the distribution function. Within this category is the well known *method of moments.* Here the distribution function is assumed to be lognormal, and the governing equations are solved in terms of the first 3 moments of the size distribution function, representing the total number of particles, mean size, and dispersion of the distribution, respectively. [53] This method is simple but not suitable to model condensation because when condensation occurs, it cannot be assumed that the particle size distribution adheres to a specific distribution function.

In *discrete* representation methods, the size distribution is divided into a discrete number of sections (m), within which all the particles are assumed to be of the same size. Even though this approach is more complex than other methods, it is very popular due to its generality. The *sectional method* [54] falls into this category. In it, the solution is obtained in terms of the first moment for each section ni =) n(υ)dυ. ni represents the total number of particles per unit volume within section *i.* The overall rate of change of ni is reduced to a set of *m* (the number of sections) ordinary differential equations by deriving expressions for the net rates by which particles are (i) added to section *i,* (ii) removed from section *i,* and (iii) remain in section *i.* The degree of approximation of this technique is limited by the number of sections used, i.e., to the maximum number of ODEs that can be solved simultaneously.

The *Monte Carlo* techniques are methods of statistical sampling of events to determine the average behavior of a system. They have been applied to simulate coagulation [55] and agglomeration. [56] The convenience, ease, directness, and expressiveness of the method makes them an attractive alternative; however, they require large amounts of computational time and memory that until recently has been computationally prohibitive.

#### **3.2. Solution methods for multicomponent aerosols**

Due to the diversity of techniques for computing size distribution dynamics of singlecomponent systems, one approach to the multicomponent problem is to use a single variable, such as the mass of one of the components in the particle, to characterize the particle size and composition. Assuming that the characteristic component is conserved during coagulation, the masses of the nonconserved components in the particle may be determined by an auxiliary constraint, typically thermodynamic equilibrium. Such an approach has been reported for H2SO4-H2O aerosols. [57] In this case it is assumed that H2SO4 is the characteristic component that is conserved by coagulation, and the water content of the particle is determined by thermodynamic equilibrium. Therefore, only the distribution of H2SO4 with time needs to be determined. Although an attractive feature of the approach is the reduction of the size of the problem, it is greatly limited to systems for which auxiliary constraints are valid and can be obtained.

In a more general approach, the techniques for single-component aerosols have been extended to include condensation/evaporation of the same or of a second material. Only the sectional method has been considered because the fundamental assumptions of the continuous and parameterized representation approaches are not valid when condensation occurs. Gelbard and Seinfelds8 extended the sectional method to account for composition in multicomponent aerosols. They assumed that all the particles within a section have a constant size and com‐ position and that coagulation and condensation are independent of particle composition. Then the governing equations are reduced to a set of m x z ODEs, where m is the number of sections for the size range and *z* the number of species (components) present in the aerosol. Therefore, this approach neglects that particles of the same size could have different compositions.

particle size range; however, it has inherent problems associated with mass conservation and numerical dispersion. Numerical dispersion appears as waves ahead of or behind regions of

*Parameterized* representations are those methods in which the size distribution is assumed to adhere to a known "self-preserving" form described by a small number of parameters. Solutions found in this way are believed to be asymptotic forms approached after long times. The assumption reduces the problem to solve a set of few ODEs for the parameters representing the distribution function. Within this category is the well known *method of moments.* Here the distribution function is assumed to be lognormal, and the governing equations are solved in terms of the first 3 moments of the size distribution function, representing the total number of particles, mean size, and dispersion of the distribution, respectively. [53] This method is simple but not suitable to model condensation because when condensation occurs, it cannot be

assumed that the particle size distribution adheres to a specific distribution function.

=) n(υ)dυ. ni

within section *i.* The overall rate of change of ni

40 Montecarlo Simulation of Two Component Aerosol Processes

number of ODEs that can be solved simultaneously.

**3.2. Solution methods for multicomponent aerosols**

In *discrete* representation methods, the size distribution is divided into a discrete number of sections (m), within which all the particles are assumed to be of the same size. Even though this approach is more complex than other methods, it is very popular due to its generality. The *sectional method* [54] falls into this category. In it, the solution is obtained in terms of the first

ordinary differential equations by deriving expressions for the net rates by which particles are (i) added to section *i,* (ii) removed from section *i,* and (iii) remain in section *i.* The degree of approximation of this technique is limited by the number of sections used, i.e., to the maximum

The *Monte Carlo* techniques are methods of statistical sampling of events to determine the average behavior of a system. They have been applied to simulate coagulation [55] and agglomeration. [56] The convenience, ease, directness, and expressiveness of the method makes them an attractive alternative; however, they require large amounts of computational

Due to the diversity of techniques for computing size distribution dynamics of singlecomponent systems, one approach to the multicomponent problem is to use a single variable, such as the mass of one of the components in the particle, to characterize the particle size and composition. Assuming that the characteristic component is conserved during coagulation, the masses of the nonconserved components in the particle may be determined by an auxiliary constraint, typically thermodynamic equilibrium. Such an approach has been reported for H2SO4-H2O aerosols. [57] In this case it is assumed that H2SO4 is the characteristic component that is conserved by coagulation, and the water content of the particle is determined by thermodynamic equilibrium. Therefore, only the distribution of H2SO4 with time needs to be determined. Although an attractive feature of the approach is the reduction of the size of the problem, it is greatly limited to systems for which auxiliary constraints are valid and can be

time and memory that until recently has been computationally prohibitive.

represents the total number of particles per unit volume

is reduced to a set of *m* (the number of sections)

high concentration.

moment for each section ni

obtained.

At high rates of condensation, three main problems appear [59]: the first is numerical diffusion. This problem arises because growth by condensation is a hyperbolic differential equation. Numerical diffusion lowers the peak value and broadens the size distribution, as shown in Figure 1. The second problem is conservation of mass of both gas-phase and particles during growth. Often, growth models subtract off the amount of gas removed by aerosol growth. However, subtracting can result in negative gas concentrations, requiring subsequent adjust‐ ments. The third problem is associated with the wide range of timescales for aerosol processes. Since the particle size domain may extend over several orders of magnitude in particle diameter, the characteristic timescales for particle growth may also have large variations. This results in stiff systems of differential equations. Such systems require special solution techni‐ ques that are not robust and are computationally expensive. [61]

monodisperse aerosol whose diameter has increased 100-fold. (b) Deviation of sectional aerosol mass from the analytical solution with time. Q is the total mass in the section and q(x) is the total mass in the section divided by the length of the section. (From Reference 60) **Figure 1.** (a) Condensation growth showing normalized sectional mass distributions for a monodisperse aerosol whose diameter has increased 100-fold. (b) Deviation of sectional aerosol mass from the analytical solution with time. Q is the total mass in the section and q(x) is the total mass in the section divided by the length of the section. (From Reference 60)

The sectional method creates numerical diffusion because when mass moves to larger (condensation) or smaller (evaporation) sections, it distributes itself uniformly throughout the section. As a result, the distributed mass can quickly grow or evaporate to the next larger or smaller section the very next time step. This error can be reduced by drastically refining discretization. However, the computational effort for such approach is prohibitive. Quadrupling the number of sections that span the size range will halve the

condensable vapor concentration, Gelbard61 showed that the fixed-grid approaches cannot overcome numerical diffusion. Applying the *particle ordering principle,* and using a moving grid, the issue was resolved for situations when only condensation or evaporation is important.

<sup>60</sup> By studying the resulting linear differential equations for a constant

Figure 3-1. (a) Condensation growth showing normalized sectional mass distributions for a

*Monte Carlo simulation of two-component aerosols*

numerical diffusion.

8

The sectional method creates numerical diffusion because when mass moves to larger (condensation) or smaller (evaporation) sections, it distributes itself uniformly throughout the section. As a result, the distributed mass can quickly grow or evaporate to the next larger or smaller section the very next time step. This error can be reduced by drastically refining discretization. However, the computational effort for such approach is prohibitive. Quadru‐ pling the number of sections that span the size range will halve the numerical diffusion. [60] By studying the resulting linear differential equations for a constant condensable vapor concentration, Gelbard [61] showed that the fixed-grid approaches cannot overcome numerical diffusion. Applying the *particle ordering principle,* and using a moving grid, the issue was resolved for situations when only condensation or evaporation is important. The particle ordering principle states than when particles grow or shrink by condensation or evaporation, respectively, particles of the same initial size will always preserve their respective ordering in particle size distribution. Consequently, the boundaries of a moving section will never cross. Unfortunately, the particle ordering principle may not be preserved in multicomponent aerosols. For example, in a condensing environment, small hygroscopic particles may surpass in size initially larger but insoluble particles. Furthermore, the moving grid method cannot be used when there is simultaneous condensation and coagulation.

Jacobson and Turco [62] advanced this idea and developed a model that uses a hybrid size grid. They used a moving size grid for condensation/evaporation and a stationary grid for other processes, such as coagulation, and nucleation. With this model, particles are assumed to be composed of an involatile core and a coating of a second condensable phase. The assumption of same size and composition for all particles within a grid still remains. Initially, the size distribution is organized within the sections in order of increasing core mass as opposed to increasing total mass. When coagulation among uncoated cores occurs, the mass of the new formed core is distributed evenly among the sections having cores with masses immediately larger and smaller than the new formed core. For example, if sections A, B, C, and D contain particles with cores of masses 2, 4, 8, and 16 units, respectively, coagulation among particles coming from sections B and C will result in a core of mass 12 units. Then half of the mass of the particle is added to section C and the other half to section D. Consequently, the number of particles in section B is reduced by one, section C is reduced by half, and section D increased by half. On the other hand, the moving grid contains the total mass of the particles (core plus coating). Then, during condensation or evaporation, particles are not transferred between sections; instead, their total masses increase or decrease to their exact sizes. To simulate coagulation among coated cores, it is assumed that when two coated cores combine, they form a third particle of a composition that falls between the characteristic compositions of two adjacent sections, i.e., the methodology for uncoated cores is applied regardless of the composition of the new formed particle. Clearly, this leads to problems in the conservation of mass of the condensable phase.

In summary, it has been shown that the coagulation equation must be solved through numerical methods. Among the many techniques available, the sectional method, although is not the simplest, is well accepted due to its generality. Contrary to the case of coagulation, there is a general analytical solution for the condensation equation, [63] but the numerical techniques available for its solution possess several inherent problems like numerical disper‐ sion, numerical diffusion, and an inability to conserve mass of the components. Only with a moving-grid-type technique can a proper numerical solution be obtained. However, this technique is not compatible or applicable when there is simultaneous coagulation and condensation.

Then there exists a need for a robust numerical technique to solve for coagulation and condensation simultaneously. In addition, the technique should be able to handle multicom‐ ponent aerosols with no restrictions in particle size and composition.

## **4. The solution method proposed**

The sectional method creates numerical diffusion because when mass moves to larger (condensation) or smaller (evaporation) sections, it distributes itself uniformly throughout the section. As a result, the distributed mass can quickly grow or evaporate to the next larger or smaller section the very next time step. This error can be reduced by drastically refining discretization. However, the computational effort for such approach is prohibitive. Quadru‐ pling the number of sections that span the size range will halve the numerical diffusion. [60] By studying the resulting linear differential equations for a constant condensable vapor concentration, Gelbard [61] showed that the fixed-grid approaches cannot overcome numerical diffusion. Applying the *particle ordering principle,* and using a moving grid, the issue was resolved for situations when only condensation or evaporation is important. The particle ordering principle states than when particles grow or shrink by condensation or evaporation, respectively, particles of the same initial size will always preserve their respective ordering in particle size distribution. Consequently, the boundaries of a moving section will never cross. Unfortunately, the particle ordering principle may not be preserved in multicomponent aerosols. For example, in a condensing environment, small hygroscopic particles may surpass in size initially larger but insoluble particles. Furthermore, the moving grid method cannot be

Jacobson and Turco [62] advanced this idea and developed a model that uses a hybrid size grid. They used a moving size grid for condensation/evaporation and a stationary grid for other processes, such as coagulation, and nucleation. With this model, particles are assumed to be composed of an involatile core and a coating of a second condensable phase. The assumption of same size and composition for all particles within a grid still remains. Initially, the size distribution is organized within the sections in order of increasing core mass as opposed to increasing total mass. When coagulation among uncoated cores occurs, the mass of the new formed core is distributed evenly among the sections having cores with masses immediately larger and smaller than the new formed core. For example, if sections A, B, C, and D contain particles with cores of masses 2, 4, 8, and 16 units, respectively, coagulation among particles coming from sections B and C will result in a core of mass 12 units. Then half of the mass of the particle is added to section C and the other half to section D. Consequently, the number of particles in section B is reduced by one, section C is reduced by half, and section D increased by half. On the other hand, the moving grid contains the total mass of the particles (core plus coating). Then, during condensation or evaporation, particles are not transferred between sections; instead, their total masses increase or decrease to their exact sizes. To simulate coagulation among coated cores, it is assumed that when two coated cores combine, they form a third particle of a composition that falls between the characteristic compositions of two adjacent sections, i.e., the methodology for uncoated cores is applied regardless of the composition of the new formed particle. Clearly, this leads to problems in the conservation of

In summary, it has been shown that the coagulation equation must be solved through numerical methods. Among the many techniques available, the sectional method, although is not the simplest, is well accepted due to its generality. Contrary to the case of coagulation, there is a general analytical solution for the condensation equation, [63] but the numerical

used when there is simultaneous condensation and coagulation.

42 Montecarlo Simulation of Two Component Aerosol Processes

mass of the condensable phase.

In the process of selecting the most convenient numerical technique to simultaneously solve the coagulation and condensation equation to study the encapsulation process, it was found that the Monte Carlo method possesses several advantages that make it an attractive route:


The proposed methodology is to simulate coagulation with the Monte Carlo method and couple it through time to a deterministic solution of the condensational growth equation. The coupling is obtained by using the time step for coagulation as the time step for condensation. The following chapter describes the Monte Carlo method adopted. Chapter 5 describes the solution of the equations for particle growth by condensation and the methodology used to couple the two solutions.

## **Author details**

Jose Ignacio Huertas

Address all correspondence to: jhuertas@itesm.mx

Tecnológico de Monterrey, Mexico

#### Montecarlo Simulation of Two Component Aerosol Processes

**Monte Carlo Simulation of Coagulating Aerosols**

44 Montecarlo Simulation of Two Component Aerosol Processes
