Unipolar and bipolar aerosol charging as time continuous Markov processes

The acquisition of charge by aerosol particles is well-known to be stochastic in nature. We review the principles of charging using the conceptually and computationally clear language of continuous time Markov processes. A novel numeric approach is presented that can be used to calculate the time evolution of various particle charging processes. Its modular character makes it easy to implement and allows for quick adaptation to specific problems. We conclude with the application of ergodicity for finite state-space Markov processes in order to determine stationary charge distributions in case of bipolar charging.


Introduction
The stochastic nature of acquisition of electrical charge by aerosol particles is well known [1,2]. In the simplest model, two distinct steps determine the rate of aerosol charging: First, the purely stochastic process of a charged particle meeting an aerosol of a specific initial charge. And second the charging process itself. The latter is governed by diffusion of electrically charged particles through a spherical shell around the aerosol particles in the background of the electrical field produced by the initial charge of the aerosol. The rate of charging can be modeled by acquisition coefficients depending on k B T as well as on the size and charge of the aerosol [3,4,5,6].
In this work, we focus on the stochastic process determining the first step mentioned above. Whereas in [1], integration of a Boltzmann-type equation resulted in ordinary differential equations for the change of concentration of aerosols of different charge, here we first identify the stochastic process of charging as a specific time continuous Markov process. The idea is to regard aerosol charging as a Poisson process whose time constant changes whenever charge is captured. We also allow for an arbitrary initial distribution of charges among the aerosols as well as acquiring different amounts of charge in one step. Using this approach, the equations given in [1] are a direct consequence of the Markov process: The Kolmogorov forward equations. The purpose of our approach is logical transparency which allows a conceptually clear way to generalize to more complicated situations such as timevarying acquisition coefficients. Furthermore we develop a numerical MatLab routine to implement different situations in a convenient way. Finally, the Markov process formalism allows for universal statements in aerosol charging, such as the existence and uniqueness of steady-state distributions. These follow directly by applying existing ergodic theorems for Markov processes.
The paper is organized as follows. In section 2 we give a brief review of basic notions of Markov process theory. This is textbook knowledge and we provide our personal choice of literature. In section 3 we identify the type of Markov process relevant for aerosol charging. In order to keep the notation simple, we discuss the bipolar case in detail and mention the unipolar and general multipolar cases as useful applications. Finally the concept of ergodicity is applied to the bipolar and multipolar cases to show the existence of steadystate charge distributions. In section 4 we apply a numerical study of the derived concept to show the advantage of the formalism in order to establish realistic modeling of aerosol charging. An appendix containing a MatLab code is provided.

A brief review of time continuous Markov processes
Time continuous Markov processes are among the most studied stochastic processes and are widely used in many fields as diverse as probability theory, statistical mechanics and economy. They describe parameterized families of events (random variables) whose future only depends on their current state. This is called Markov property. In order to set up the basic notation, we briefly review elementary aspects of Markov processes needed to formalize aerosol charging. For a systematic introduction to the theory we refer to the literature on the subject, we used in particular the textbooks [7,8] and references therein. Definition 1. Let X t be a family of random variables, parameterized by t ∈ [0, ∞) and taking values in a countable set E. X t is called time continuous Markov process, if for all t > 0, s > s n > · · · > s 1 ≥ 0, with i, j, k 1 , . . . , k n ∈ E and a stochastic matrix Π t .
The countable set E is usually termed state space. Only the conditional probability P (X s+t = j|X s = i) is relevant, and we assume that it only depends on the time difference t (time homogeneity). The quadratic matrix Π t is stochastic, i.e. its row elements are either zero or positive and add up to 1. Its dimension is determined by the number of elements of the state space E. It is straightforward to see that one has the semi-group property of time evolution Π s+t = Π t Π s and initially Π 0 = 1, where 1 is the identity matrix. Important information about the process is obtained by looking at its evolution in infinitesimal time: Definition 2. Let X t be a continuous Markov process. The infinitesimal generator of the process is defined by By writing 1(i, j) we assume the rows and columns of the identity to be indexed by the elements of E and 1(i, j) is the matrix entry at position (i, j) ∈ E × E. The following properties of the infinitesimal generator are of importance for our work. Their proofs are straightforward and can be found in textbooks on the subject.
• j G(i, j) = 0, i.e. the sums over row elements of G vanish.
• The Kolmogorov equations hold, d dt Π t = GΠ t and d dt Π t = Π t G.
• The Kolmogorov equations are solved by Π t (i, j) = e tG (i, j), where as usual the matrix exponential is defined by the power series e tG (i, j) = ∞ n=0 t n n! G n (i, j) .
Once the state space and generator are identified, various probability distributions of the process are calculated using the solution to the Kolmogorov equations. Let X t be the random variable with values in E, the probability of X t = j when starting at X 0 = i is given by It is easily seen that, given an initial distribution x 0 = (x k 0 , x k 1 , . . . ), for k l ∈ E and k l ∈E x k l = 1, the distribution at time t > 0 is calculated by applying the exponential of the generator: Using this, the Kolmogorov equations translate in a set of first order linear ordinary differential equations: d dt Typical examples of Markov processes appear in queuing theory, birth-death processes or population dynamics. In case of discrete time, random walks on lattices are characteristic examples. A very basic example of a continuous time Markov process is the Poisson process, which we now define and review in the terminology of Markov processes. This is the starting point for generalizations and applications to aerosol charging in the next section. Let (L i ) i≥0 be a series of independent random variables, exponentially distributed with common parameter µ ∈ R, i.e. with distribution function ρ L i (t) = µ e −µt . Furthermore, let T k := k i=1 L i be the time of the k'th event happening. One can then show that the number of events X t that happened until some fixed time 2 t is Poisson distributed with parameter µt, i.e.
As for this process the infinitesimal rate of change from the state N t = k to N t = k + 1 is µ, it is clear how to construct the generator of the corresponding Markov process. As we count events starting from 0, the state space is E = Z + 0 and Indeed, using equation (3) for the Markov process the probability of k events happening, starting from k = 0 is coinciding with (6). The probability of j events happening, starting at event i is computed analogously for i, j ∈ E. As the rate µ is not changing in the whole process, this model clearly is too easy to describe aerosol charging, but conceptually it is now immediate how to generalize. Finally we state a basic result of ergodic theory for continuous Markov chains. Intuitively, if for sufficiently large times every state can be reached from any other state in a finite state space, there is hope to reach a steady state distribution. The mathematical condition for the latter is irreducibility of the generator. G is called irreducible, if for all i, j ∈ E one can find an N ∈ N and k 0 , k 1 , . . . k N ∈ E with k 0 = i and k N = j such that N l=1 G(k l−1 , k l ) = 0. For irreducible generators one has the following basic result: Theorem 2.1. Let E be a finite state space and G an irreducible generator of a continuous Markov process and Π t (i, j) = e tG (i, j). Then lim t→∞ Π t (i, j) = α(j) independent of i ∈ E. Furthermore, α is the only probablility distribution on E which satisfies α G = 0, or equivalently α Π s = α for all s > 0.
For the proof we refer to the textbooks on the subject, e.g. [7,8]. Realistically, any charging process for aerosol particles has a finite state space, as it is not possible to infinitely charge an aerosol particle. In case of an irreducible generator theorem 2.1 then guarantees the existence of a steady state distribution and gives a simple recipe to calculate it.
3 Unipolar and bipolar aerosol charging as Markov processes Using the notation introduced in the previous section, we identify aerosol charging as continuous Markov processes. We start with the general case of bipolar charging with different concentrations of positive and negative ions and treat unipolar charging as a special case. After commenting on various possibilities for generalizations we apply ergodicity to prove the existence of a steady state charge distribution for the bipolar case.

Bipolar charging with different ion concentrations
We begin with the idealized situation that arbitrary positive and negative charging of aerosol particles by the acquisition of positively as well as negatively charged ions is possible. To simplify the identification of the Markov process, we assume that there are only ions with one positive or one negative elementary charge. The concentration as well as mobility of the two types of ions can be different. Finally we assume that the aerosol particle can not loose ions, therefore its charge increases by acquiring a positive ion and decreases by the acquisition of a negative ion. Hence, in case of arbitrary possible charging of an aerosol, the state space is given by integer numbers E = Z. In contrast to the Poisson process, the rate of acquisition of an ion depends on the current charge state. Furthermore there are different rate constants (µ i ) i∈Z and (η i ) i∈Z for positive and negative ions, respectively. By adjusting the generator (7) for the Poisson process accordingly we get for the generator of bipolar charging: The rows and columns of the generator extend over the integer numbers. The rate constants for acquiring a positively charged ions appear in the upper off-diagonal and the rate constants for acquiring a negatively charged ion in the lower off-diagonal, respectively. To implement the fact that positive and negative ions come in different concentrations, we take advantage of the fact that G bi can be written as the sum of two matrices G bi = η + µ, for Implementing different initial concentrations n − and n + for positive and negative ions, respectively, is done by multiplying the corresponding matrix, i.e. G bi = n − η + n + µ.
In practice, the existence of arbitrary high positively or negatively charged aerosol particles is highly unlikely. It is more realistic to assume a maximally positive charge and a minimally negative charge a particle can acquire. It turns out that typically charge numbers from −20 to 20 appear and the probabilities to charge more positive or more negative can be assumed to be zero. This means that the state spaces will be finite (e.g. E = {−20, −19, . . . , 19, 20}, the matrix µ have a zero column for maximally positive charge and the matrix η will have a zero row for maximally negative charge. As an example, in case fo E = {−2, −1, 0, 1, 2} and n + = n 1 = 1, without loss of any generality the generator for bipolar charging is It is clear how the structure looks for general finite state spaces E = (−n, −n + 1, . . . , n − 1, n). It is not necessary to have the states symmetrically distributed around 0 ∈ E, but for convenience we will focus on this case. For completeness, the Kolmogorov equations in the form (5) for a given initial probability distribution x 0 = (x −n , x −n+1 , . . . , x n−1 , x n ) are given below, using Newton's notation for time derivatives: . . .
Systems of ODE's of this type were studied earlier in the literature for various scenarios of particle charging [1]. The Kolmogorov equations thus establish the link between the Markov process formalism to earlier work on the subject. We will use matrices of the form (11) combined with concentrations n − , n + for detailed numerical studies in section 4. Before doing so we will look at variations of the formalism to treat other situations and investigate the steady state behaviour of bipolar charging.

Other scenarios: Unipolar charging and charging with multiple valued ions
An advantage of the formalism used so far is its conceptual simplicity. In particular it allows for an easy implementation of different situations frequently used in aerosol physics. A case of interest is unipolar charging, where there are only ions of one type. It is easy to see that setting either n − = 0 or n + = 0 gives the desired results, as this sets either the contribution of η or µ to the generator to zero. In case of a state space E = {−2, −1, 0, 1, 2} this means We note that the zero row in G uni will result in a row (0, 0, 0, 0, 1) in the corresponding stochastic matrix Π t = e tG . This is what is expected: In case of charging without loosing ions and a maximal positive charge value, the probability of staying at the element −2 ∈ E is 1. An element of E of this kind is sometimes also called absorbing state.
The formalism in addition allows for the extension to charging by multiply valued ions. To illustrate this without overloading the notation, we assume again a state space of E = {−2, −1, 0, 1, 2} and the existence of ions which carry one or two positive charges and one or two negative charges. The respective charging rates are dentoted by µ + i , µ ++ i as well as η − i , η −− i , respectively. In case of all concentrations equal, we therefore get for the generator where the Σ i are chosen in such a way to satisfy the properties of a generator, i.e.
We observe that e.g. η −− −1 does not exist due to our assumptions on the state space. Similarly for µ ++ 1 and so forth. Furthermore, the implementation of different concentrations of the differntly charged ions is done analogously to the bipolar case: The matrix (14) can be decomposed in matrices containing the rates for the different types of ions: G mult = η −− + η − + µ + + µ ++ . Multiplying each matrix with the corresponding concentration gives the generator for the entire process: As a result, our formalism easily allows to treat aerosol charging with ions carrying different charge values and existing in differnt concentrations and with or without a maximal or minimal total charge of the aerosol particle.

Ergodicity for the bipolar and multipolar cases
The ergodic theorem for finite state space continuous Markov processes reviewed at the end of section 2 allows to investigate whether there are steady state charge distributions in the various scenarios studied so far. Clearly, for unipolar charging, the generator is not irreducible, as one can not move backwards (i.e. have charge loss) in the charging process. More precisely, starting with a charge k and ending in a charge l, for k > l, in the product of matrix elements needed for irreducibility, at least one matrix element has to be taken from below the diagonal, but all these elements vanish and hence the product vanishes. However, the generator in the bipolar and multipolar cases has the irreducibility property. As an example, in the bipolar case, using the notation of theorem 2.1 and section 2, first take k 0 = k and k N = l for k, l ∈ Z and k < l. Then for N ∈ N such that k 0 = k, k 1 = k + 1, . . . , k N = l one has N n=1 G(k n−1 , k n ) = G(k, k + 1)G(k + 1, k + 2) · · · G(l − 1, l) Second, in case k > l, taking N such that k 0 = k, As we assumed all µ k and η k nonvanishing (i.e. nonvanishing charging rates), the products in the last two expressions are nonzero, for all k, l ∈ Z. Thus, G is irreducible and ergodicity is guaranteed. Note that this is even more evident in the multipolar case. In this case even if some of the η's or µ's may vanish, in case another path using other (multiple charge) rates exists, there still exist stationary charge distributions. All of them can be computed rigorously by solving the corresponding eigenvalue problem mentioned in theorem 2.1. We will present a numerical study of the results in the next section.

Numerical results
To investigate the applicability of the formalism outlined above, we modeled the charging of particles with a diameter d of 300nm at moderate temperature (T = 300K) in a bipolar ion cloud with constant properties using the Markov chain approach developed above. To apply the equations listed in chapter 3, the rate coefficients µ i and η i for the uptake of an additional positive or negative elementary charge at any relevant charge number q i have to be known. Without any loss of generality, the simple analytic formulae derived by Gunn [3] are used to calculate µ i (19) and η i (20) although more precise derivations have been published in the past ( [4,9]).
with λ = e 2 4π 0 d k B T . Here e represents the elementary charge, m ± the ion mobility of the positively or negatively charged ion, respectively, 0 the dielectric permittivity, d the particle diameter (300nm), k B the Boltzmann constant and T the absolute temperature (300K).
For neutrally charged aerosols, the rate coefficient equations listed above are not valid. For this special case, the rate coefficients are given by In analogy to equations (9)-(11), we constructed G bi (i, j) using µ i and η i for any integer starting charge numbers ranging from −20 to 20. Assuming a bipolar environment of 10 7 singly charged ions per cm 3 of each polarity as well as ion mobilities of m + = 1.35·10 4 V m/s and m − = 1.6 · 10 −4 V m/s for positive and negative ions (taken from [10]), respectively, we calculated the probability distribution matrix Π t = e tG using equation (3)   The i-th row of each probability distribution matrix represents the charge probability distribution that would be expected for the final state after charging aerosols with start charge number q si for t seconds. The probability distribution matrix at time point 0 (Π t=0 ) is a diagonal matrix indicating that the charge distribution is not altered after zero time. With increasing time, the probability distribution matrix gradually approaches a state, where the probability distribution is independent from the start charge q si of the aerosol.
The charge distribution after charging an aerosol with given initial aerosol charge distribution x 0 for t seconds is calculated by applying equation (4). Technically, this procedure corresponds to generating a charge distribution matrix by weighting the rows of the probability distribution matrix with the respective components of the initial charge distribution of the aerosol and summing along its columns. Applying this approach for several time points, we calculated the time evolution of the charging process in a bipolar environment for an arbitrary initial charge distribution (Fig. 2). The arbitrary initial charge distribution shown in Fig. 2A contains a distinct aerosol component loaded with nine positive elementary excess charges. The charging process in the bipolar environment leads to a slow change in charge distribution until an equilibrium steady state charge distribution is reached after approximately 10 −1 s (Fig. 2B). The final charge distribution after a charging time of 1s is shown in Fig. 2C (black line). This distribution is equivalent to the ergodic limit of the steady state charge distribution that is given by the kernel of the matrix G as outlined in theorem 2.1 (see Matlab code example in the Annex for details).
Both, the final charge distribution after 1s and the ergodic limit of the steady state charge distribution are consistent with the analytic solution for the steady state charge distribution suggested by [11] (Fig. 2C, green X).
The charging process occurs on several timescales. From the Kolmogorov equations, i.e the system of linear ODEs (5), (12) follows that the time evolution of the charge number population can be expressed by a weighted sum of exponentials whose rate constants are certain combinations of the entries in the generator matrix. More precisely, these rate constants are given by the eigenvalues of the transposed generator matrix (see Fig.  2D) and are therefore independent from the initial charge distribution. The eigenvalue "0" (corresponding to 0/s) describes a time independent component, namely the final state (Fig. 2C). The next higher rate constant describes the slowest process involved in charging and is therefore a measure of the velocity of the overall charging process. The corresponding lifetime of this slowest process (given by the negative inverse rate constant) is approximately 40 ms for the example discussed here. This means that for the conditions of the assumed bipolar environment (i.e. a bipolar neutralisator), a charging time of less than a second would be sufficient to reach the final charge distribution.
In analogy to results published earlier using Fuchs theory [10], size dependent final bipolar charge distributions can be calculated by varying the parameter d in the Markov chain approach (Fig. 3). By using more precise attachment coefficients than done in our study, comparing theoretical values with experiment becomes feasible as well. Figure 3: Bipolar charge distribution for different particle size. The overall shape of the curves resemble the ones published in [10], more realistic charge distributions can be derived by replacing the simple attachment coefficients by more precise ones that take additional physical properties into account.

Conclusion and outlook
Using time continuous Markov processes, we studied unipolar, bipolar and multipolar aerosol charging processes. This conceptually simple and numerically powerful viewpoint potentially allows for calculating charging probability distributions known in the literature by using experimentally determined charging rates. One of the advantages of the formalism is its structural clarity and the use of powerful general results such as the ergodicity theorem for finite state space Markov chains. The latter e.g. allows for a numerically fast computation of limiting steady state charge distributions.
Another advantage of the formalism is its extendability to more complicated (realistic) charging processes. Aerosol clouds with particles of multiple diameters and varying initial charge distributions can be treated as well as arbitrary mixture of multiply charged ion types, each coming with its own concentration. Further generalizations not treated in this work are time-varying concentrations usually found in situations where the number of ions continuously decreases when there is no infinite supply [9] or ions are created by an experimental setup. In our treatment it is clear how to implement such cases. The mathematical theory behind these processes is more complicated, as the stochastic matrix used in the Markov property will become time-dependent in the sense that not only differences play a role, but also the time of the event one is looking at. The theory of time-inhomogenious continuous Markov processes provides the correct setup to treat these cases, including more intricate ergodicity statements.