*2.2. Agglomeration*

The finite volume method developed by Kumar et al. [44] is used to ge<sup>t</sup> the discretized form of the pure agglomeration population balance equation. Taking *<sup>S</sup>*(*x*) = 0 and *b*(*<sup>x</sup>*, *y*) = 0 in Equation (1), the general agglomeration PBE for batch process can be written as

$$\frac{\partial u(\mathbf{x},t)}{\partial t} = \frac{1}{2} \int\_0^\infty \beta(\mathbf{x}', \mathbf{x} - \mathbf{x}\prime) u(\mathbf{x}', t) u(\mathbf{x} - \mathbf{x}\prime, t) d\mathbf{x}\prime - \int\_0^\infty \beta(\mathbf{x}, \mathbf{x}\prime) u(\mathbf{x}, t) u(\mathbf{x}\prime, t) d\mathbf{x}\prime \tag{19}$$

with initial data *<sup>u</sup>*(*<sup>x</sup>*, 0) = *<sup>u</sup>*0(*x*) ≥ 0.

Take the whole domain as [0,*<sup>R</sup>*]. Divide the domain into *L*1 cells with boundaries *xi*−12 and *xi*+ 12 for *i* = 1, 2, ... , *L*1 where *x*1 2 = 0 and *xL*1<sup>+</sup> 12 = *R*. Take the representative of the *i*–th cell *xi*−12 , *xi*+ 12 as *xi* = - *xi*+ 12 <sup>+</sup>*xi*− 12 2 and cell length Δ*xi* = *xi*+ 12 − *xi*−12 . Definetheinitialapproximation *<sup>u</sup>*0(*x*) as

$$u(0, \mathbf{x}) = \frac{1}{\Delta x\_i} \int\_{\mathbf{x}\_{i-\frac{1}{2}}}^{\mathbf{x}\_{i+\frac{1}{2}}} u\_0(\mathbf{x}) d\mathbf{x}.\tag{20}$$

To ge<sup>t</sup> the discretized form of the agglomeration Equation (19), calculate the total birth and death in each cell. For that, collect all those particles from the lower cell which is going to a particular higher cell. This is done by defining some set of indices. Here, the index set is calculated as follows

$$\mathbb{T}^{i} = \left\{ (j,k) \in L\_1 \times L\_1 : \mathbf{x}\_{i-\frac{1}{2}} < \left( \mathbf{x}\_j + \mathbf{x}\_k \right) \le \mathbf{x}\_{i+\frac{1}{2}} \right\}.\tag{21}$$

Integrating Equation (19) over each cell *xi*−12 , *xi*+ 12 , and introducing weights as given by Kumar et al. [44], the discretized form can be obtained

$$\frac{du\_i}{dt} = \frac{1}{2} \sum\_{(j,k)\in\mathbb{I}^i} \beta\_{jk} u\_j u\_k \frac{\Delta x\_j \Delta x\_k}{\Delta x\_i} w\_{jk}^b - \sum\_{j=1}^{L\_1} \beta\_{ij} u\_i u\_j \Delta x\_j w\_{ij}^d. \tag{22}$$

where *wbjk* and *wdij* are the weights responsible for preservation of number and conservation of mass, respectively. These weights are defined as

$$w\_{jk}^b = \begin{cases} \frac{\mathbf{x}\_j + \mathbf{x}\_k}{2\mathbf{x}\_{jk} - \left(\mathbf{x}\_j + \mathbf{x}\_k\right)}, & \mathbf{x}\_j + \mathbf{x}\_k \le R \\\ 0, & \mathbf{x}\_j + \mathbf{x}\_k > R \end{cases} \tag{23}$$

and

$$w\_{ij}^d = \begin{cases} \frac{\mathbf{x}\_{lj}}{2\mathbf{x}\_{ij} - \{\mathbf{x}\_i + \mathbf{x}\_j\}} & \mathbf{x}\_i + \mathbf{x}\_j \le \mathbf{R} \\\ 0, & \mathbf{x}\_i + \mathbf{x}\_j > \mathbf{R} \end{cases} . \tag{24}$$

Here *lij* is the symmetric index of the cell, where the agglomerate particle *xi* + *xj* falls. Rewrite Equation (22) as

$$\frac{d\Delta\mathbf{x}\_i u\_i}{dt} = \frac{1}{2} \sum\_{(j,k)\in\mathbb{I}^i} \beta\_{jk}\boldsymbol{\mu}\_j \boldsymbol{\mu}\_k \Delta\mathbf{x}\_j \Delta\mathbf{x}\_k \boldsymbol{w}\_{jk}^b - \sum\_{j=1}^{L\_1} \beta\_{ij}\boldsymbol{\mu}\_i \boldsymbol{\mu}\_j \Delta\mathbf{x}\_i \Delta\mathbf{x}\_j \boldsymbol{w}\_{ij}^d. \tag{25}$$

Take *gi* = Δ*xiui*, then Equation (25) gives

$$\frac{d\mathcal{g}\_i}{dt} = \frac{1}{2} \sum\_{(j,k)\in\mathbb{I}^i} \beta\_{jk}\mathcal{g}\_j\mathcal{g}\_k w^b\_{jk} - \sum\_{j=1}^{L\_1} \beta\_{ij}\mathcal{g}\_i\mathcal{g}\_j w^d\_{ij}.\tag{26}$$

To ge<sup>t</sup> the particle number at some certain time point τ, the time span [0, τ] is divided into *K* subintervals [*t*, *t* + <sup>Δ</sup>*t*], where *t*1 = 0, *tK* = τ and Δ*t* is the time step.

Now the particle number *g*(*<sup>t</sup>* + Δ*t*) at time point *t* + Δ*t* using transformation matrix *<sup>T</sup>*(*<sup>t</sup>*, Δ*t*) is given by

$$
\mathcal{g}(t + \Delta t) = \mathcal{g}(t) \cdot T(t, \Delta t). \tag{27}
$$

Here *g*(*t*) is a 1 × *L*1 row matrix with *i*-th element *gi*(*t*), *i* = 1, 2, 3, ... , *L*1 at time *t*, and *<sup>T</sup>*(*<sup>t</sup>*, Δ*t*) is a *L*1 × *L*1 transformation matrix, calculated at time *t* to obtain the distribution after time step Δ*t*. The elements of the transformation matrix are given by

$$T\_{ij}(t, \Delta t) = \begin{cases} 1 + \left(\frac{1}{2} \sum\_{\substack{(i,k) \in \mathbb{I}^i \\ (i,k) \in \mathbb{I}^i}} \beta\_{ik} \mathbf{g}\_k(t) w^b\_{ik} - \sum\_{j=1}^{L\_i} \beta\_{ij} \mathbf{g}\_j(t) w^d\_{ij} \right) \Delta t, & i = j \\ \frac{1}{2} \sum\_{\substack{(i,k) \in \mathbb{I}^i \\ (i,k) \in \mathbb{I}^i}} \beta\_{ik} \mathbf{g}\_k(t) w^b\_{ik} \Delta t, & i < j \\ 0, & i > j \end{cases} \tag{28}$$

To improve the accuracy of the scheme (27), the Runge–Kutta method is proposed:

$$
\lg(t + \Delta t) = \lg(t) \cdot T(t, \Delta t) \tag{29}
$$

where

$$\widetilde{T}(t,\Delta t) = \left[T(t,\Delta t) \cdot T\left(t + \Delta t, \frac{\Delta t}{2}\right) + \left(I - T\left(t, \frac{\Delta t}{2}\right)\right)\right].\tag{30}$$

Here *I* is the identity matrix. The method given by Equations (29) and (30) gives the second order accuracy.

Since the scheme is explicit, there would be some restriction on the time step to ge<sup>t</sup> a non-negative solution of the scheme (Equation (29)). Therefore, the Courant-Friedrichs-Lewy (CFL) condition [46] on the time step is given by

$$
\Delta t \le \min\_i \left| \frac{2g\_i(t)}{Q(t)} \right|, \tag{31}
$$

where

$$\begin{aligned} Q(t) &= W\_1(\mathcal{g}(t)) - W\_2(\mathcal{g}(t)) + W\_1(\mathcal{g}^\*(t)) - W\_2(\mathcal{g}^\*(t)),\\ W\_1(\mathcal{g}(t)) &= \frac{1}{2} \sum\_{(j,k) \in \mathbb{I}^i} \beta\_{jk} \mathcal{g}\_j(t) \mathcal{g}\_k(t) w\_{jk}^b, \\ W\_2(\mathcal{g}(t)) &= \sum\_{j=1}^{L\_1} \beta\_{ij} \mathcal{g}\_i(t) \mathcal{g}\_j(t) w\_{ij}^d, \\ g\_i^\*(t) &= g\_i(t) + \Delta t \cdot W\_1(\mathcal{g}(t)) - W\_2(\mathcal{g}(t)) \end{aligned} \tag{32}$$
