{'fraction_non_alphanumeric': 0.044404619673436875, 'fraction_numerical': 0.019115890083632018, 'mean_word_length': 4.466394557823129, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Analytical approach for solving population balances: a homotopy perturbation method

In the present work, a new approach is proposed for finding the analytical solution of population balances for aggregation and fragmentation process. This approach is relying on the idea of the homotopy perturbation method (HPM). The HPM solves both linear and nonlinear initial and boundary value problems without nonphysical restrictive assumptions such as linearization and discretization. It gives the solution in the form of series with easily computable solution components. The outcome of this study reveals that the proposed method can avoid numerical stability problems which often characterize in general numerical techniques related to this area. Several examples including Austin’s kernel, available in literature, are examined to demonstrate the accuracy and applicability of the proposed method. In addition, the analytical solution to two new kernels (the power-law kernel in fragmentation and the Ruckenstein/Pulvermacher kernel in aggregation) are also introduced.


Introduction
Fragmentation, along with aggregation, represents one of the most basic particulate processes. It refers to the process by which a granule, colloidal cluster, polymer chain or other "particle" breaks into a distribution of fragments. This generic process is encountered in a wide range of natural systems and serves as a model of branching processes in general, in which a population is transformed by replacing a member by a distribution of new members. The physical mechanisms that give rise to breakup or fragmentation (the two terms are used interchangeably in this paper), are generally different between systems but a common mathematical description is possible under the rather general assumption that the rate at which a particle breaks up and the distribution of fragments that are produced depend only on the parent particle and are independent by the presence of other particles in the population. This model views breakup as a first order reaction in the concentration of parent particles, with a distribution of products that depend on the parent particle. Moreover, the physical mechanism that give rise to the aggregation process can be described under the general assumption that the rate at which the particles aggregate and the distribution that are depend on the particles which aggregate together to form a new particle.
The breakage as well as aggregation models will be described by two functions. For defining the breakage model, the rate of breakup of the parent particle, a(m), where m is the size of the particle, and the distribution of fragments, k(m|n) that gives the number of fragments of size m produced by a parent of size n are required. If we take size to refer to the mass of the particle, the fragment distribution obeys the following normalization conditions: (1. 2) The first condition gives the average number of fragments,f , and is the same for all parent sizes m; the second condition expresses mass conservation between the parent particle and the fragments. A physically realistic fragmentation model requiresf ≥ 2.
The governing equation for the size distribution of population that undergoes the above fragmentation process is given by the fragmentation equation, ∂c(m, t) ∂t = −a(m)c(m, t) where c(m, t) is the concentration of particles whose mass is in (m, m+dm) (we normalize the mass m by its mean value at time 0).
Moreover, the aggregation model can be defined with the help of aggregation kernel, g(m, m ′ ), where m and m ′ are the sizes of the particles which defines the rate at which particles aggregate. The governing equation to track the change in particle size distribu-tion due to aggregation process is defined by aggregation population balance equation, Here c(t, m) represents the concentration of particles of size m at time t. Here, g(m, n) is the aggregation kernel, which describes the rate at which the particles of sizes m and n coagulate to form a particle of size m + n. It is non-negative and symmetric.
The equations (1.3) and (1.4) have been the subject of several investigations and the analytical solutions have been obtained for a small number of special cases. This work, largely theoretical, (Simha, 1941(Simha, , 1956Tobolsky, 1957;Ziff and McGrady, 1986;Ernst and Szamel, 1993;Singh and Hassan, 1996) seeks analytical solutions to the fragmentation equation and studies the scaling behavior of the size distribution at long times. A review of these developments and a collection of analytical and solutions is given by Ziff (Ziff, 1991).
The analytical solution for aggregation equation can be studied by many researchers for limited number of aggregation kernels. Ranjbar et al. (2010) derived analytical solution for constant kernel using Taylor polynomials and radial basis functions. Moreover, Hammouch and Mekkaoui (2010) uses the Laplace-Variational iteration to obtain approximate series solution for equation (1.4). For additive aggregation kernel, Fernandez-Diaz and Gomez-Gar (2007) proposed analytical solution. However, exact solutions for constant, sum and product kernels are given in Deaconu and Tanré (2000).
The method of homotopy perturbation has been recently developed to obtain analytical solutions from differential and integral equations and has been successfully used to to solve a a wide range of dynamical problems (Ganji, 2006a;He, 2003a;Shahed, 2011). The purpose of this paper is to formulate the homotopy perturbation method (HPM) in a form appropriate for the fragmentation and aggregation equations and to use it for obtaining analytical solutions. The paper is organized as follows. In section 2 we introduce the general methodology of HPM. In section 3, HPM will be adapt it to the fragmentation and aggregation equations. In next sections 4 and 5, we apply the method to a number of known and new fragmentation and aggregation models. Finally in section 6 we examine the transient and scaling solution for a quaternary breakup and aggregation models and relate its behavior.

The Homotopy Perturbation Method
Recently, the homotopy perturbation method (HPM) has been studied by many researchers for solving linear and nonlinear problems He (1999He ( , 2000He ( , 2006; Nazari-Golshan et al. where T is a general differential operator and B is a boundary operator. Usually the operator T can be decomposed into two parts, a linear operator L and a nonlinear operator N, and expressed as NOTE: If the equation is linear then the non-linear part will be zero, i.e., N(u) = 0. According to HPM, we construct a homotopy which satisfies where u 0 is an initial guess to exact solution of (2.1). When p = 0 then L(v(r, 0)) = L(u 0 ) = 0, and when p = 1, then T (v(r, 1)) − g(r) = 0. It is worth noting that as the embedding parameter p increases monotonically from zero to unity, the changing process of p from zero to unity is just that of v(r, p) from u 0 (r) to u(r). This is called deformation, and also L(v) − L(u 0 ) and T (v) − g(r) are called homotopic in topology.
According to the HPM, we can first view the embedding parameter p as a small parameter, and construct the solution as a power series in p, as (2.5) Substituting (2.5) in (2.4) and then letting p = 1, we obtain the solution as The series (2.6) is a convergent for most of the cases and the rate of convergence depends on the nature of the problem He (2000).
The idea of the HPM for different types of models (He, 2003b;El-Shahed, 2005;He, 2006;Ganji, 2006b) will be demonstrated in this section. Firstly, the problem linear integral equation will be solved using HPM. Secondly, non-linear integral equation will be taken into account.
Model 1: Consider the linear integral equation: Here K(u, t, c(u)) should be linear function, i.e., K(u, t, c(u)) = n j=1 K j (u, t)c i,j (u). For solving above equation using HPM method, a homotopy will be constructed (2.8) Suppose above equation has the solution of the following form: where C i (t), i = 1, 2, 3, . . . are the functions which will be determined and C 0 (t) is the initial approximation to the solution which is taken to be a(t). On substituting (2.9) into (2.8) and equating the like powers of p, we have The approximated solution of (2.7) can be obtained by setting p = 1, Model 2: Consider the nonlinear integral equation (2.10) To apply the HPM, we write (2.10) as with solution y(x) and we construct the homotopy H(u, p) as where F (u) is a functional operator with known solution v 0 . We choose a convex homotopy by and continuously trace an implicitly defined curve from a starting point H(v 0 , 0) to the solution function H(f, 1). The embedding parameter p increases monotonically from 0 to 1 as a trivial problem F (u) = 0 is continuously deformed to the original problem L(u) = 0.
According to the HPM, we can first view the embedding parameter p as a small parameter, and construct the solution as a power series in p, as (2.14) Substituting (2.14) in (2.13) and then letting p = 1, we obtain the solution as The series (2.15) is a convergent for most of the cases and the rate of convergence depends on the nature of the problem He (2000).
Substitute (2.14) into (2.11) and equate the terms with identical powers of p, we obtain

Implementation of HPM
In this section, HPM will be implemented to the aggregation and fragmentation equations.
It can be noticed clearly that the breakage equation is linear whereas aggregation equation is non-linear. Therefore, for solving breakage equation, HPM for linear equation will be implemented and for solving aggregation equation, HPM for non-linear equation will be implemented.

Implementation of HPM for Breakage Equation
To apply HPM to fragmentation we first express the population balance equation in the integral form, where c 0 (m) = c(t, m) is the size distribution at time zero. Next we introduce a new function u = u(t, m; p) of time, size, p and define the homotopy of Eq. (3.1) as follows: Thus by continuously varying p from 0 to 1 we obtain a continuous transformation of u(t, m; p) from the initial distribution c 0 to the actual distribution c(m, t) at time t.
Following HPM, we express u(t, m; p) as a power series in p, where the coefficients u k are function of time and size to be determined. According to the homotopy in Eq. (3.2) the solution to the fragmentation equation is obtained by setting p = 1 in the above series: The series is convergent and the rate of convergence depends on the nature of Eq. (3.1).
To obtain the coefficients u k we insert the series (3.3) into (3.2): collect terms in powers of p and set their coefficients to zero. We find: . . . (3.8) The first term, u 0 , is equal to the initial distribution and each subsequent term involves the coefficient of previous order. Thus we may obtain these coefficients recursively. The nth order approximation of the solution is obtained by truncating the series past the nth term.
Here, however, we are interested in cases for which the infinite series can be expressed in closed analytical form. First we reproduce known solutions to the fragmentation equation using HPM, and we also present new results for a case that has not been solved previously.

Implementation of HPM for Aggregation Equation
The homotopy of equation (1.4) is constructed which satisfies Inserting u = ∞ k=0 p k u k in the above equation, Equate the terms with identical powers of p, we obtain the expression of HPM for aggregation equation: . . .
The nth order approximation to the solution is defined as Note that as n → ∞, we have the limit φ = ∞ k=0 u k , which is an analytical solution of (1.4). In the next section 4, the HPM will used to obtain analytical solution of breakage equation for various breakage kernels with different initial conditions.

Case I: Random binary fragmentation with linear selection function
In random binary fragmentation the distribution of fragments is given by (4.1) Here we obtain the solution with linear fragmentation rate, a(m) = m. The general form of the factors u k (m, t) is given by Eq. (3.8), which now becomes Below we obtain the solution for monodisperse and exponential initial distribution.

Monodisperse initial conditions
We start with u 0 = f 0 = δ(m − 1) and obtain the functions u k recursively: where θ is the Heaviside step function. The distribution is given by the infinite series As shown in the Appendix, this is equal to This is the same as the result given by Ziff and McGrady (Ziff and McGrady, 1985).

Exponential initial condition
(4.10) The size distribution in this case is given by This is the same as the solution given by Ziff and McGrady (Ziff and McGrady, 1985).

Case II: Random binary fragmentation with quadratic selection function
In this model the breakup rate is a quadratic function of size, a(m) = m 2 and the fragment size distribution is again given by g(m|n) = 2/n and the initial condition is exponential.
The general recursion for u k (m, t) is (4.12) We obtain solutions for monodisperse and exponential initial conditions:

Monodisperse initial condition
With u 0 (m, t) = δ(m − 1) the functions u k are: The size distribution is This again recovers the solution of Ziff and McGrady (Ziff and McGrady, 1985)

Exponential initial condition
Starting with c 0 (m) = e −m and solving recursively we obtain (4.17) The solution is constructed by computing the infinite series, This is the same as the result obtained by Ziff and McGrady (Ziff and McGrady, 1985) 4.3 Case III: k-nary fragmentation with power-law rate and number of fragmentsf (4.20) By varying α between 1 and 2 the number of fragments ranges from ∞ to 2. With α = 2, in particular, this model defaults to the random distribution of fragments with quadratic breakage rate. Equation (3.8) in this model now is (4.21) Using u 0 (m, t) = e −m the functions u k are (4.22) Inserting this result into the infinite series in Eq. (3.4) we obtain which gives the size distribution at all times. Now, HPM will be apply to solve various models of aggregation equation. The distribution of aggregates for constant aggregation kernel is g(m, n) = 1.
The partial sum of the series solution is defined obtained as As n −→ ∞, we obtain The exact solution is u(t, m) = n 2 (t)e −n(t)m , where n(t) = 2 t+2 as in Ranjbar et al. (2010).
It is difficult to find the generalized term for the solution of sum kernel. So, we approximate the solution by taking the summation of some terms. The comparison of exact solution with the approximate solution is represented in Figure 1. The comparison reveals that by considering the sum of four terms in the series, exact solution can be approached.

Product aggregation kernel g(m, n) = mn
The distribution of aggregates for product aggregation kernel is g(m, n) = mn The general form of the factors u k (m, t) is given by Eq. (3.15), which now becomes The solution will be obtained for two different initial conditions, i.e., for u 0 (m) = e −m and u 0 (m) = e −m m .

Exponential initial condition, u 0 (m) = e −m
Its analytic solution follows Kumar (2006) as We start with u 0 = e −m and obtain the functions u k recursively . . .
The partial sum of the series solution is defined obtained as Now, the solution for another initial condition.

When
The analytical solution can be found in Ranjbar et al. (2010) otherwise.
Here I is the modified Bessel's first kind function: We initialize with u 0 = e −m m and obtain the functions u k recursively The partial sum of the series solution is obtained as As n −→ ∞, we obtain 2, m] 6 Discussion Cases 1 and 2, which describe random binary breakup with linear and quadratic rate functions essentially form the "classical" model in fragmentation. Their solutions have been obtained previously in the literature multiple times by various methods, including probabilistic treatment Ziff and McGrady (1985), and solution of the equivalent discrete problem (Montroll and Simha, 1940). The kernel of Case 3 is a special form of a kernel whose general form is This kernel was introduced by Austin (Austin and Luckie, 1972;Klimpel and Austin, 1977) to model grinding. In cumulative form it is the sum of two power-law distributions, of which one is the dominant distribution for large m (Austin gave this kernel with different coefficients but the two forms can be shown to be the same under appropriate relationships among the coefficients). The same kernel has been the subject of theoretical investigations of the scaling behavior of the fragmentation equation (Ziff, 1991) and of the so-called shattering transition (Ziff and McGrady, 1986;Ernst and Szamel, 1993), a behavior analogous but opposite to gelation that is manifested by the presence of a finite fraction of the mass of the system in particles of zero size. Ziff and McGrady (Ziff and McGrady, 1985;McGrady and Ziff, 1987) gave the solution of the Austin kernel with power-law breakage rate for monodisperse initial conditions. Oukouomi Noutchie and Doungmo Goufo obtained the same solution by Laplace transform (Oukouomi Noutchie and Doungmo Goufo, 2014). In principle, given the linearity of the fragmentation equation the solution for any other initial condition can be reconstructed from the solution to the monodisperse problem. The monodisperse solution, however, employs the confluent hypergeometric function and its manipulation is not trivial. The solution obtained here in Eq. (4.7) for exponential conditions is compact and represents a new result. This solution obeys the scaling derived by Ziff and McGrady (Ziff and McGrady, 1985), as can be easily demonstrated. Dividing Eq. and setting m = (z/t) 1/α we find αc(m, t) m 2 = 1 α t −2/α e t −1/α (−z 1/α )−z αt 1/α z (α+1)/α + αzt 2/α + z 2/α , whose limit at long t is αc(m, t) m 2 → ze −z ≡ Φ(z), and Φ(z) is the scaling function. The scaling limit of the distribution is c(m, t) ∼ αm 2 Φ(tm a ). (6.2) The scaling function given by Ziff and McGrady is Φ(z) = z γ/α e −z Γ(γ/α) (the Filippov model in Table 1 of Ref. (Ziff and McGrady, 1985) ), which for γ = α reduces to the result obtained here.
To illustrate these results we examine the solution for α = 4/3. The fragmentation rate is a(m) = m 4/3 and the distribution of fragments is k(m|n) = 4 3m 2/3 n 3/2 , (6.3) with an average number of fragments f = 4. This distribution is very asymmetric with a large fraction of the population in the small size region (Fig. 2). There is no shattering, however, because the fragmentation rate decreases as the size approaches zero (indeed our solution does include shattering as the condition f ≥ 2 constrains α to be positive).
The size distribution is shown in Fig. 3. After an initial transition the distribution in log-log coordinates reaches an invariant form that is simply translated to lower average sizes. This is a manifestation of the scaling limit.
Finally, to examine the convergence of the series we calculate the solution by retaining only a finite number of terms k. Figure 4 compares truncated approximations of the full solutions for α = 3/2 at t = 10 and k max = 2, 3, 12 and 13 (k max is the maximum order of the term retained in the series). The general feature of the truncation is that the approximation is excellent for 0 < m < m max , where m max is the size beyond which the approximation breaks down and is pushed to larger sizes as the order of the truncation increases. Even-order truncations provide stable approximations approximations of the size distribution and it is interesting to point out that in this particular case, even a very low order approximation such as k max = 2 represents the distribution very well over two decades in c(m).

Conclusions
The homotopy perturbation method has been successfully applied to the fragmentation as well as aggregation population balance equation. The HPM provides a systematic methodology to solve for the size distribution and was shown to reproduce all known solutions and to lead to new results for a special case of the Austin kernel and also reproduce analytical solution for various models of fragmentation and aggregation.