Crystal Population Balance Formulation and Solution Methods: A

Subscriber access provided by CORNELL UNIVERSITY LIBRARY

Review

Crystal population balance formulation and solution methods: A review Hecham M. Omar, and Sohrab Rohani Cryst. Growth Des., Just Accepted Manuscript • DOI: 10.1021/acs.cgd.7b00645 • Publication Date (Web): 13 Jun 2017 Downloaded from http://pubs.acs.org on June 15, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Crystal Growth & Design is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Crystal Growth & Design

1 Crystal population balance formulation and solution methods: A review

Hecham M. Omar1, Sohrab Rohani1* 1

The University of Western Ontario, Department of Chemical & Biochemical Engineering,

London, Ontario, Canada *

Corresponding author: Sohrab Rohani ([email protected])

Abstract: Crystallization is an important part of many chemical industries. Efforts are being invested to improve the performance of the crystallization process by designing novel crystallizers. An important aspect in the development of new crystallizers is the ability to describe the behavior of such units in terms of rigorous dynamic mathematical models and solving the resulting models efficiently. The current bottleneck in modelling crystallization systems is the complexities associated with the crystal birth, growth and death processes using population balance equations. In this article, various crystal birth, death and growth models are introduced and reviewed. Population balance models as well as solution methods (e.g. analytical, moment methods, discretization (classes/sectional) methods and Monte Carlo methods) are also reviewed and new advances in solution methods are described. Population balance equations are used in other fields and developments from other fields that can be extended to crystal population balance equations are included in this review.

Keywords: Crystallization, Population Balance, Mathematical Modelling

ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 50

2

1. Introduction Crystallization is an important chemical process used to separate and purify solid materials. Crystallization proceeds in two main steps: nucleation and subsequent crystal growth. The driving force for crystallization comes from solute supersaturation. Nucleation can occur via primary or secondary nucleation. Primary nucleation occurs from the supersaturated solution without any preexisting crystals. Secondary nucleation occurs when crystals are already present. The existing crystals lower the energy of nucleation and therefore require less supersaturation when compared to primary nucleation. Secondary nucleation is the common process in industrial crystallizers. Crystallization (nucleation and growth) kinetics are affected by cooling rate, thermal history, mixing performance, impurities, crystallizer volume/geometry and a stochastic aspect of nucleation process.1 These factors are also related to the mixing, mass and heat transfer phenomena, therefore, computational fluid dynamics (CFD) can capture all the complexities involved. Other processes can also create and destroy particles, namely breakage and agglomeration. In particular, models that combine transport phenomena (momentum, mass and energy balances) with the population balance are being developed at an increasing rate.2–6 Using these models, crystallizer designs and operating conditions can be optimized. Another important use of modelling is for the control of crystallization systems.7–10 With better models, more advanced control strategies can be developed and implemented. The transport phenomena involved in the crystallizer models are relatively well-known and developed compared with the population balance. Many review papers have been written on transport phenomena (including CFD).11–14 Theoretically a solution of a population balance equation describes the temporal and spatial variation of the crystal size distribution. In the basic


Page 3 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


3 population balance, two coordinate domains are considered: external and internal coordinate domains. The external coordinates refer to the location in the crystallizer space (i.e. x, y and zdirections). The internal coordinates refer to the properties of the crystal (e.g. crystal length, width, etc.). Apart from the crystallization of inorganic materials and pharmaceuticals, population balance equations have been applied to other fields such as polymer crystallization15– 18

, aerosols19–22 and biological applications.23–26 Many of the insights in other fields can be

adopted to the crystallization. Reviews have been written on the derivation of population balance equations27 and their applications to biological systems.26 Books have also been written on different aspects of population balances and their applications.28–30 This review looks at how population balance equations for crystallization are formulated (different terms in the population balance equation) and the methods to solve the equations. A similar review was written in 198531, however, to the best of our knowledge, a collection of all of the recent developments for solving population balances has not been presented.

2. Population balance models The population balance equation is widely considered to have been derived simultaneously by Hulburt and Katz32 and Randolph.33 The population balance equation can be derived from three different viewpoints: continuum mechanics principles, statistical Boltzmann-like equation and probability principles.27 Population balance equations describe the properties of the particles in space and time. The population balance equation describes not only the movement of the particles but also nucleation, growth, aggregation/agglomeration and breakage of the particles.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 50

4 Population balance equations can be formulated using a Eulerian or a Lagrangian approach (Equations 1 and 2).29 The Lagrangian viewpoint tracks a finite number of particles in a flow field. The Eulerian viewpoint tracks the particles as a bulk continuous phase.34 The Lagrangian approach results in: ∂n + ∇ ⋅ (v e n ) + ∇ ⋅ (v i n ) = B − D ∂t

(1)

And the Eulerian formulation:


Page 5 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


5 ∂n d (log V ) + ∇ ⋅ vin + n =B−D− ∂t dt

∑ k

Qk n k V

(2)

where n represents the population density, v e is the external velocity vector, v i is the internal velocity vector (typically “internal velocity” refers to the crystal growth in units of length per time but can be any state property change in time), B is crystal birth (explained in Section 3.1),

D is crystal death (explained in Section 3.2), V is the crystallizer volume and Qk nk is the volumetric inflow/outflow multiplied by the number density of in the inflow/outflow streams. Population balance equations can be formulated to include more than one particle property. The most common use of more than one dimension is for changes in shape (e.g. changes of aspect ratio of needle-like crystals) where the different particle coordinates (e.g. for a high-aspect ratio crystal: r and z) are the particle properties.35,36 The multiple properties do not necessarily have to be particle size coordinates. The n-dimensions in the population balance equation describe the nindependent properties of the crystals. For example, Gerstlauer et al. tracked both the crystal size and the internal lattice strain of the crystals.37 They used a model to describe the effect of different assumptions of internal lattice strain on the growth of the crystals. They found that steady-state mass density function was strongly dependent on the various assumptions for relaxation of internal lattice strain. The difficulty in solving population balance equations (and the subsequent research into finding effective solution methods) comes from the integro-differential nature of a population balance equation in the presence of agglomeration/aggregation and/or breakage terms (see Equations 5-8). Inputs to the population balance equation (Equations 1 and 2), the crystal birth rate for primary and secondary nucleation, growth rate, particle birth due to agglomeration and breakage,



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 50

6 as well as particle death due to agglomeration and breakage, are given below. Primary nucleation is usually modelled by:

B1 = k N ,1 S α

(3)

where k N ,1 is the primary nucleation constant [1/m3/s], S is the supersaturation [-] and α is an empirical constant. Alternatively, secondary nucleation is usually modelled by:

B2 = k N , 2 N iα M Tβ S γ

(4)

where k N , 2 is the secondary nucleation constant [1·minαm3β/kgβ/m3/s], N i is the impeller speed [rpm], M T is the suspension density [kg/m3], S is the supersaturation [-] and α , β , γ are empirical constants. The most commonly used mathematical models for primary and secondary nucleation take the form of an empirical power law model.38–40 Using a power law (Equations 3 and 4) allows the model to fit to different conditions by estimating empirical parameters using the experimental data. However, this limits their usage to the system for which they were derived. Also, since it is an empirical formulation, no information is given as to the mechanism of nucleation. Equation 4 is sometimes modified by removing the impeller term or by adding extra terms to account for other phenomena (such as attrition). The simplest form of the power law model includes the suspension density and supersaturation effects. Birth and death due to breakage and agglomeration are shown in Equations 5-841. ∞

Bb = ∫ b(t , x, x ')a ( x ')n (t , x ')dx '

(5)

Db = a ( x )n(t , x )

(6)

x

∞

D a = ∫ β (t , x, x ')n(t , x )n (t , x ')dx '

(7)

0


Page 7 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


7 x

Ba =

1 β (t , x − x ' , x ')n(t , x − x ')n(t , x ')dx ' 2 ∫0

(8)

where b(t , x, x ') is the probability density function for the generation of particles (referred to as the breakage kernel), a( x ) and a(x') are the selection functions that describe the rate at which particles with size x and x' are selected to break, respectively, β (t , x − x' , x') is the measure of frequency of collisions between crystals of size x − x' and x' and β (t , x, x') is the measure of frequency of collisions between crystals of size x and x' (referred to as the agglomeration or aggregation kernel). Various forms of the agglomeration/aggregation kernel and breakage kernel are found in Equations 5-8 can be rewritten using length instead of volume using a change of variables.

Table 1 and Table 2. Equations 5-8 can be rewritten using length instead of volume using a change of variables.

Table 1 Agglomeration/aggregation kernel models



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 50

8 Expression β ( x, x ' ) = β 0

β ( x, x ' ) = β 0 ( x + x ' )

β (L, L') = β 0 (L3 + L'3 )

Reference 41 42 43

β (L, L') = β 0 (L + L')(L−1 + L' −1 )

44

β = k β G h B0pτ q

45

Table 2 Breakage kernel models Expression b =1 b = Lα b = exp δL3

Reference 46

( )

b = L3

47

Where β 0 is the constant aggregation kernel, L and L' are the lengths of the particles, k β is the rate coefficient for agglomeration kernel [L/{1·s·(m3/s)h·(1/L/s)p(s)q}], G is the overall crystal volume growth rate [m3/s], h is the exponent of growth rate in the aggregation kernel correlation, p is the exponent of nucleation rate in the aggregation kernel correlation, q is the exponent of mean residence time in the aggregation kernel correlation, B0 is nucleation rate [1/L/s], τ is mean residence time [s] and a is the exponent in the breakage kernel.

3. Crystallization models Crystallization models describe the various ways in which crystals can be generated, destroyed and grown. Many different forms of these equations exist and they are typically derived using a semi-empirical approach. A flow chart showing all the birth/death mechanisms is shown in Figure 1.


Page 9 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


9

Figure 1 Methods of birth and death of crystals (adapted from 48)

3.1 Crystal birth 3.1.1 Primary nucleation Primary nucleation occurs from the solution itself 49 at high supersaturation. The dependence of primary nucleation on supersaturation is highly nonlinear being close to 0 when supersaturation is low and rapidly increases once a specific critical supersaturation is reached.50 Primary nucleation is difficult to control due to poor reproducibility. There are two views to explain the poor reproducibility. One viewpoint assumes that primary nucleation is deterministic and the poor reproducibility is a result of the high sensitivity of the experimental conditions. Therefore, if the experimental conditions were controlled precisely and accurately, nucleation could be controlled. The second viewpoint questions the fundamental assumption that primary nucleation is deterministic and assigns a stochastic process to the nucleation process.51 The stochastic nature of primary nucleation has been reported in literature.52–56 Goh et al. developed a stochastic model of nucleation for both homogeneous and heterogeneous nucleation in droplet-based microfluidic systems.57 The model was capable of being used for any form of nucleation rate expression. The



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 50

10 analytical solution of the model can be used to find nucleation kinetic parameters as was done by Goh et al. Classical nucleation theory provides a hypothesis of how nuclei form from the solution.58,59 However, many experiments differ from predictions made by classical nucleation theory. A theory that is gaining attraction is the two-step nucleation theory.60,61 Atoms or molecules interact in such a way that clusters or aggregates begin to form. Where classical nucleation theory and two-step nucleation theory differ is on what happens after the clusters form. In classical nucleation theory, there is a critical size of the cluster, which depends on the supersaturation, in which the cluster is stable and will grow to a macroscopic size. However, in two-step nucleation theory after the formation of the cluster of solute molecules, the molecules are reorganized into an ordered structure to become a nucleus.62 A comprehensive review of both classical nucleation theory and two-step nucleation theory is presented by Karthika et al.63 Therefore, based on both theories, factors that affect the cluster formation can either aid or inhibit primary nucleation. For example, depending on the intensity of mixing, cluster formation may increase or decrease. Forsyth et al. found that fluid shear had a significant effect on the primary nucleation.64 They found that both primary nucleation rate and growth rate were strongly dependent on the shear rate. They also found that increasing the average shear rate reduced the mean induction time via a power law relationship. At lower shear rates, the nucleation rate increased fairly sharply whilst at higher shear rates, the nucleation rate showed slower increase. Liu and Rasmuson proposed shear-induced molecular alignment, specifically agitation-enhanced cluster aggregation, is the mechanism that should be looked at to explain why shear rate aids nucleation.65 Yang et al. used cluster aggregation to explain the effect of shear on nucleation rate


Page 11 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


11 as viewed through the classical nucleation theory.66 As the shear rate increased a maximum nucleation rate was reached beyond which increasing shear rate decreased the nucleation rate. They suggested that as shear rate increases beyond the maximum, clusters reaching the critical cluster size are disrupted, decreasing the nucleation rate. Fluid shear may cause the clusters to interact at specific energies that favor nuclei formation. Another method to affect cluster formation is the use of acoustic cavitation to stimulate consistent primary nucleation proposed by Virone et al.67 They found that in the presence of acoustic cavitation there was negligible sensitivity to the average supersaturation. The pressure caused by the cavitation induces a high local supersaturation that drives the primary nucleation. With a higher average supersaturation, there is more potential collisions due to a larger number of molecules. Miyasaka et al. performed a similar experiment by using ultrasonic irradiation.68 The study showed that depending on the average supersaturation, ultrasonic irradiation can promote or inhibit nucleation. It was also found that in the presence of ultrasound activated nucleation, the final crystal size decreased and when ultrasound inhibited nucleation, the final crystal size increased. As for the effect of the ultrasound intensity, low levels of ultrasound irradiation inhibited primary nucleation whereas high levels of ultrasound irradiation induced it.69 Primary nucleation can be further divided into two categories: homogeneous and heterogeneous nucleation. Homogenous nucleation represents spontaneous formation of nuclei from the solution. The pathway of homogeneous nucleation is still unclear. This is largely due to the difficulty in direct observation. Classical nucleation theory fails to completely explain the nucleation process. Tan et al. showed that density fluctuations may not be the driving force for nucleation.70 This conclusion was drawn from the observation that the local structures can swing



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 50

12 in and out of the solution and that nucleation did not necessarily start in the densest regions. Heterogeneous nucleation occurs in the presence of impurities (particle, container walls, internals, etc.) that aid the formation of a seed for a nucleus. Heterogeneous nucleation is more common than homogeneous nucleation due to the lower energy required to form the nuclei. Many different factors control the presence of heterogeneous nucleation (e.g. wall roughness, impurity type, etc.). 3.1.2 Secondary nucleation Secondary nucleation occurs in the presence of pre-existing crystals. Primary nucleation precedes secondary nucleation.49 Secondary nucleation begins to dominate when a sufficient number of crystals have formed. However, Kadam et al. found that a single crystal led to a surge of crystals. A single nucleus is formed and grows until eventually attrition occurs leading to the surge of crystals.71 Secondary nucleation occurs via two mechanisms: shear and collision induced secondary nucleation.72 The secondary nucleation rate is a combination of the individual removal mechanisms73 shown in Figure 2 which presents a flowchart of the specific mechanisms that make up the two mechanisms (i.e. shear and collision).


Page 13 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


13

Figure 2 Mechanisms for secondary nucleation (reproduced from 72)

The preexisting (seed) crystals act as the source of nuclei. Attrition, unlike other mechanisms, is independent of the supersaturation. Contact nucleation involves the removal of the surface layer of solute to act as nuclei. Multiple forms of a crystal develop on the surface of the seed crystals. This means that if the layer is removed early in the process, multiple forms of seeds are produced. However, if a long enough time has elapsed so that the solute layer has adopted the configuration of the seed crystal, the crystal form is the same as the seed.72 The number of secondary nuclei produced during contact nucleation depends on the supersaturation and the contact energy.74 Data showed that faster growing substances give higher yields of nuclei.75 Contact nucleation can be used as a way to control the crystal size as shown by Wong et al.76 In this method, the parent crystal is placed on a rotating platform while under applied stress. The size of the secondary nuclei was controlled by the supersaturation and flow rate of the feed solution. Shear nucleation occurs in the same way that contact nucleation occurs. The surface layer of crystals is removed and act as nuclei. Sung et al. found that the nuclei had to be in a region of relatively high supersaturation.77



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 50

14 3.1.3 Breakage/attrition Breakage/attrition is the fracturing of a crystal into fragments induced by a physical force such as crystal-crystal or crystal-impeller collision. These fragments function as nuclei and can grow into larger crystals. It is like contact nucleation in that nuclei are formed from a crystal. However, in contact nucleation, the parent crystal remains intact since only the surface layer is removed. In breakage/attrition, the parent crystal is irreversibly broken into multiple pieces. Time is an important variable on crystal attrition. Mazzarotta et al. looked into the effect of time on attrition in a stirred vessel78. They found that in the first minute, many fragments were produced but immediately afterwards fragment production became negligible. Cracks already present in the crystals were suspected to be the cause. Mazzarotta et al. performed the study with sugar crystals. It would have been advantageous to test different crystals with different mechanical properties to see if the results could be generalized. Also, the sugar was dispersed in a nonsolvent. A nonsolvent ensures no crystallization or dissolution. It would have been worthwhile to try this experiment with various solvents that slightly and significantly dissolve the solute of choice. Bravi et al. looked at the effect of mechanical properties on crystal attrition in a stirred vessel79. Their results showed that there was a clear connection between attrition resistance and the intrinsic mechanical properties of the crystals, although more work is required before predictions can be made about attrition knowing the mechanical properties. Asakuma et al. used micro-hardness parameters (elastic modulus, fracture strength and roundness) to predict attrition behavior.80 The micro-hardness parameters were used to calculate fracture energy and attrition coefficient which could then be used to predict the extent of attrition. Another important consideration is the crystal morphology/shape.81,82 Crystals with a high aspect ratio are easier to break than crystals with a low aspect ratio. Sato et al. looked at


Page 15 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


15 breakage of high aspect ratio crystals.81 Experimental work was performed to obtain the parameters that were then used for a two-dimensional population balance. Crystals were prone to breaking after a certain aspect ratio threshold. The breakage rate was modelled by the product of impact energy, impact frequency and a function that explained the critical aspect ratio dependence. The model that was developed explained the crystal breakage well at different process conditions. Briesen extended the shape dependence further and formulated a twodimensional population balance model that was capable of solving for shape dependent crystal attrition.82 The two-dimensions considered were size and shape of the crystals. In the vast majority of cases, shape dependence is ignored even though it is an important factor, especially in attrition and secondary nucleation. The results showed that the method could model shape dependence in a physically meaningful way but the shape modification functions used had no mechanistic foundation. The next step is to derive the shape modification functions from a physical basis. As stated in Section 3.1, ultrasound under certain circumstances helps the nucleation. Guo et al. looked at the effect of ultrasound on breakage of crystals.83 The effect was found to be more pronounced for larger crystals. Ultrasound affected the crystals in 2 ways: the collision rate between crystals increased and the cavitation contributed to breakage. This could be important in a stirred tank crystallizer such as a mixed suspension mixer product removal (MSMPR) crystallizer where there is a broad size distribution. Ultrasound can be used for nucleation but if prolonged, it could lead to the crystal breakage that is unaccounted for, changing the final crystal size distribution.

3.2 Crystal death 3.2.1 Dissolution



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 50

16 Dissolution occurs when the relative supersaturation of the solution falls below 1 and the solution becomes undersaturated. Dissolution can be affected by crystal growth, crystal shape, crystal size, crystal habit, crystal form (stable versus metastable) as well as the chemical nature of both the crystal and the solvent.84 Dissolution occurs in two steps: (1) the surface reaction and detachment from the surface and (2) the mass transfer of the species from the crystal to the bulk solution though a diffusive layer. Depending on the solubility of the crystal, the rate-limiting step is either the mass transfer in the case of a readily dissolvable species or the detachment of the species from the surface for a sparingly soluble species.38 An expression for dissolution was derived assuming a diffusion controlled process and a spherical geometry.85 Lu et al. extended the work by accounting for shape and size of the particle.86 Shan et al. further extended the dissolution model by including both the mass transfer and the detachment of molecules at the crystal surface.87 The model matched the experimental dissolution characteristics closely. Mangin et al. used a population balance equation to study dissolution in a stirred tank vessel coupled with experimental work.88 They found that near the beginning of the dissolution process, the dissolution is accelerated by disaggregation of the aggregated powder. After the initial accelerated dissolution decreases, the dissolution was found to be controlled by mass transfer. The study should have been repeated with other pharmaceuticals to see whether this could be generalized. Garcia et al. found that at high undersaturation the dissolution process is diffusion-controlled and at low undersaturation, the dissolution process is surface mechanism controlled89. 3.2.2 Aggregation/agglomeration Aggregation and agglomeration are used interchangeably in literature. Nichols et al. proposed a distinct definition of each term90. Aggregation is limited to pre-nucleation situations when


Page 17 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


17 molecules form into supramolecular structures. When particles are involved, the term agglomeration is used. As defined, aggregations may lead to agglomeration. There is further distinction between hard and soft agglomerates where hard agglomerates require considerable force to fracture. Agglomeration occurs in multiple steps. The success/failure of agglomeration depends on both steps. The first step requires 2 particles to come close enough for the growth of an agglomerative bond. The second step is the formation of the agglomerative bond. The bond strengthens over time but can be broken by other collisions. Factors that affect agglomeration include: surface structure (different faces), nature of the agglomerative bond (the bond may not be isotropic/amorphous), hydrodynamic effects, effect of supersaturation, forces involved in aggregation (geometry, speed, orientation of the collision).91 Depending on the size, the agglomerates form differently. With large crystals, agglomerates get larger by the addition of a single crystal. With small crystals, agglomerates get larger with the addition of smaller agglomerates.92 David et al. developed an expression for agglomeration combining fluid mechanics and crystal engineering.93 Using the equation, it was found that agglomeration is not sizeindependent; small crystals can preferentially stick to large crystals compared to crystals of equal size. The model includes the effect of tank size, geometry and stirring power. Yu et al. also found that agglomeration was dependent on crystal size.94 They also found that agglomeration rate showed a maximum as agitation intensity was increased. This was explained by better antisolvent dispersion and the disruption of the agglomeration. Peña et al. developed a new method for formulating agglomeration to be included in population balance equations.95 The method involved splitting the crystals into three groups:



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 50

18 primary crystals, un-agglomerated crystals and agglomerate crystals. With this formulation, the three groups were tracked and compared to the typical formulation where every crystal was treated the same. This also allows optimization of each population separately. Agglomeration can be problematic from an industrial crystallization point of view since agglomeration is difficult to predict on a large scale. Agglomerates can also lead to pharmaceutical stability issues. Larger and more uniform particles (both size and shape) can be produced when the supersaturation is controlled to remain within the metastable zone.96

3.3 Crystal growth Crystal growth is also a function of supersaturation of the solute. The supersaturation can be controlled in various manners: cooling of the supersaturated solution, adding an anti-solvent, via a reaction, via evaporation or by altering the pH. The most common mathematical form of crystal growth is represented using a power law relationship (Equation 9).36,97,98 The crystal growth model can be modified to include additional dependencies/terms.10,38,39 The power law form of the growth rate is given by:

G = kG S α

(9)

where k G is the growth rate kinetic parameter [m/s or µm/s] and α is an empirical constant. By far the most common growth rate used in literature is size-independent growth (McCabe’s ∆L Law 99).43,100–102 A significant factor that can affect the crystal growth rate is the impurities present in the mother liquor. This is an important consideration in industrial crystallization since impurities are inevitable. There are many studies in literature on the effect of impurities on the growth of pharmaceutical crystals.103–105 The growth rate can be increased, decreased or significantly inhibited.106 Additives can also be useful in controlling crystal morphology. Thompson et al.


Page 19 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


19 looked at how additives affected the growth and morphology of paracetamol crystals. Depending on the functional groups of the crystal faces, the structure of the impurities plays an important role.107 Zhang et al. studied the effect of impurities on a two-dimensional population balance model.108 By understanding how impurities affect crystal growth, crystal properties can be controlled. A two-dimensional population balance model can be used to determine the evolution of the shape.

4. Population balance model solution methods Population balance equations can be solved by many different methods. The difficulty in solving the population balance increases exponentially when agglomeration and breakage are included.

4.1 Analytical solution Under certain simplifying assumptions an analytical solution can be found for a population balance equation. Including agglomeration/aggregation or breakage requires numerical solution unless additional assumptions are made.42,109,110 Many analytical solutions for population balance equations come from atmospheric research-related literature.42,111,112 The most common analytical solution is for a mixed suspension mixed product removal (MSMPR) crystallizer. The assumptions include steady-state operation with size-independent growth rate and negligible agglomeration or breakage. Also, it is assumed that the feed has no crystals, the product has the same crystal size distribution as the content of the crystallizer (representative/mixed product removal assumption) and that the MSMPR is well-mixed. Also, the residence time, τ , for the solution and the crystals must be equal. This leads to Equation 10 and its solution Equation 11.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 50

20

dn +n=0 dL

(10)

 L  n = n 0 exp −   Gτ 

(11)

Gτ

where n is the number density [1/m3/µm], n0 is nuclei population density [1/m3/µm], L is the crystal size [µm], G is the crystal growth rate [µm/s] and τ is the residence time [s]. For multistage MSMPR crystallizers, Equation 10 can be extended to yield113:

G1τ 1

dn1 + n1 = 0 dL

(12)

G2τ 2

dn2 + n2 = n1 dL

(13)

Gk τ k

dnk + nk = nk −1 dL

(14)

Equations 12-14 can then be analytically solved to yield Equations 15 and 16. The solution of the kth MSMPR (Equation 14) can be found in Randolph and Larson (1962).114

 L   n1 = n10 exp −  G1τ 1 

(15)

   L   L  0  G1τ 1 L   + n1     n2 = n20 exp − exp − − exp −      G τ   G2τ 2   G1τ 1 − G2τ 2    G1τ 1   2 2 

(16)

Analytical solutions can be combined with moment methods to find the average crystal size to compare with the experimental data.98

4.2 Moments method The moments method was first developed by Randolph and Larson115 and Hulburt and Katz.32 In general, the moments method involves converting the population balance equation into moments equations of the number density. The moments can then be used to find global properties of the


Page 21 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


21 distribution (e.g. total length of all particles, total area of all particles, total volume of all particles) or can be combined to find the average crystal size, and the coefficient of variation of the size distribution. The computational simplicity relative to other methods makes moments method attractive. 4.2.1 Standard method of moments (SMOM) The SMOM is used in literature to solve population balance equations3,32,115–117 and is the foundation for other moments methods. Other moments methods (discussed in later sections) are developed as extensions to handle situations in which the SMOM fails. For example, the SMOM equations cannot mathematically handle size dependent growth rate expressions and breakage and aggregation kernels.47 The definition of moments transformation equation is given in Equation 17 and its general form is given in Equation 18.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 50

22 ∞

m k = ∫ Lk n ( L; t ) dL

(17)

∞ dmk (t ) k = (0 ) B(t ) + ∫ kLk −1 (L )n(L, t )dL + Bk ,a (t ) − Dk ,a (t ) + Bk ,b (t ) − Dk ,b (t ) dt 0

(18)

0

where mk is the kth moment [mk/m3], B (t ) is the birth rate due to nucleation, Bk ,a (t ) is the birth rate due to agglomeration, Dk ,a (t ) is the death rate due to agglomeration, Bk ,b (t ) is the birth rate due to breakage, Dk ,b (t ) is the death rate due to breakage. From calculated moments, the average crystal size and the width of the its distribution can be found and used for data comparison. For example, Zauner and Jones tracked the second moment (total surface area of all crystals) with respect to time.3 They found that in their crystallizer the second moment increased sharply in the beginning and then remained constant for approximately 4 residence times. Falola et al. extended SMOM to size-dependent growth, aggregation/agglomeration and breakage.47 This method was called the extended method of moments (EMOM). EMOM has the flexibility to track any number of moments where as QMOM (discussed in the following section) requires an even number of moments. The EMOM also allows for the calculation of the number density, which may be useful in CFD applications. The EMOM does allow reconstruction of the distribution in limited cases, however in general the unique distribution cannot always be constructed, a priori knowledge of the distribution is required for unique reconstruction in such cases. 4.2.2 Quadrature method of moments (QMOM) The QMOM is an extension of the SMOM by approximating the integral with a quadrature (Equation 19)19:


Page 23 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


23 ∞

N

0

i =1

mk = ∫ Lk n( L)dL ≈ ∑ wi Lki

(19)

The quadrature points are calculated via the product-difference algorithm developed by Gordon.118 The QMOM fixes the closure problem associated with the SMOM. Due to its robustness in handling different types of problems, it is commonly used in literature.119,120 It presents the same advantages as the SMOM (namely computational simplicity) but can handle crystallization processes that the SMOM cannot (e.g. when breakage is included) and is therefore a better choice in which crystallization processes are combined with CFD.119–122 Many of the developments in QMOM were developed for aerosol research but they can be extended to crystallization. The QMOM developed by McGraw is susceptible to errors caused by the moments used for the quadrature points necessitating the use of a small number of moments. Modifications of the QMOM have been developed to improve its performance. The direct quadrature method of moments (DQMOM) was developed to handle multi-coordinate systems and systems with strong coupling between phase velocities (fluid dynamics) and any internal coordinates (e.g. crystal length, volume).22 The QMOM does not always handle these situations well. The DQMOM fixes this issue by tracking the weights and abscissas instead of the moments.123 The DQMOM has the same problem that plagues the QMOM, namely ill-conditioning (unable to handle more than a small number of moments). A general method called the sectional quadrature method of moments (SQMOM) was developed by Attarakih et al.124 The method is characterized by splitting the particle population into two parts: primary particles are used for distribution reconstruction and secondary particles are used for particle interactions such as breakage and agglomeration. What is interesting about this method is that it is general and can be simplified to other sectional methods by assuming one



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 50

24 secondary particle or to the original QMOM or DQMOM by assuming one primary particle. SQMOM was shown to be accurate in reconstructing the distribution. By increasing the number of primary or secondary particles, SQMOM accuracy increases but at the expense of computation. However, since the number of particles can specify the type of QMOM to be used, the general SQMOM can be implemented and then altered to the type of moments method required by changing the number of primary or secondary particles. An extended quadrature method of moments (EQMOM) was developed by Yuan et al.125 The method allows for number density function reconstruction. The framework of the method involves choosing a kernel density function, depending on a single parameter. The single parameter is then varied to fit an additional transported moment to agree with the reconstructed number density function. This method was tested for 13 cases and was shown to be accurate in approximating the moments. The method presented in the paper was used for spatially homogeneous population balance equations but can be extended to spatially inhomogeneous systems using finite volume methods. Another QMOM modification that was developed by Alopaeus et al. is the fixed quadrature method of moments (FQMOM).126 The FQMOM is defined by fixed quadrature points and moment-conserving weights. The relative location of the quadrature points can be optimized for the problem being solved. When optimized, the error is an order of magnitude lower than QMOM. However, there is a computational cost increase associated with the FQMOM. Other quadrature methods such as conditional quadrature method of moments127, automatic differentiation-based quadrature method of moments128 and Jacobian matrix transformation129have been developed, but they have not been used as extensively.


Page 25 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


25 4.2.3 Recovering the distribution



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 50

26 The size distribution can be recovered from moments methods130–133 using deconvolution methods. Ideally, the type of distribution to be recovered would be elucidated from theoretical calculations or experimental evidence.134 This is due to no natural/suitable orthogonal basis set for the distribution fit. Deconvolution of the moments does not define a unique distribution and therefore a distribution form must be assumed a priori. Recovering a unique distribution from moments is still an open area of research in mathematics. The most common method used is to assume a size distribution functional form (e.g. Gaussian, log-normal or a gamma distribution) and to use the moments to fit the function. This method is common due to its relative simplicity and its speed. However, its drawbacks are the need to assume a size distribution in advance and the potential for the form to change as time progresses. See John et al.131 for a thorough testing of this method with different assumed size distribution function forms. A related method was developed by John et al.131 that uses a piecewise polynomial function without requiring an assumed functional form. It was further developed by de Souza et al.132 to include adaptive grid node distribution to optimize the distribution in critical regions. Another common method to find the distribution is the moments inversion method. This is done by equating the moments to the matrix formed from the set of moments equations multiplied by the discrete crystal size distribution µ = An and invert the matrix A to solve for the crystal size distribution n . Using more moments discretizes the distribution into finer intervals; n is discretized into j + 1 intervals requiring j moments. However, as the number of moments and intervals increases, the matrix A becomes ill-conditioned. To improve the condition number, the central moments can be used.135


Page 27 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


27 A third method that is used to reconstruct the size distribution from the moments is the maximum entropy method. This method has received increased attention to reconstruct a distribution from moments. The method maximizes the Shannon Entropy by solving a constrained optimization problem.133 This method has progressed in reconstructing general distributions136,137 and a size distribution in the case of population balance moments.133 It was further extended to multidimensional distributions by Abramov.138,139 Lebaz et al. showed that maximum entropy was superior to beta kernel density function and spline based methods for reconstruction because it required less moments and offered the greatest computational cost reduction and reconstruction accuracy.133 Maximum entropy has only recently been applied to moments for population balance equations and has potential. It is an open area of research in information theory and could have significant application in crystallization modelling. Diemer and Olson developed a reconstruction technique where a set of target moments is matched to a set of predicted moments using nonlinear regression to minimize a weighted rootmean-squared error.130 Hutton et al. came up with a technique where a two-parameter probability density function was assumed and the moments calculated while the two parameters were varied prior to implementation.140 Then by finding where constant moments contour lines intersect the pre-calculated surfaces, the two parameters can be recovered. By knowing the information about the probability density function and the total particle count or solid loading, the distribution can be recovered. This method showed quick and accurate recovery of the distribution. However, as with other a priori methods, this method may show difficulty in handling complex cases and should be tested. This method should be used with crystallization systems where other methods fail (e.g. with multiple peaks, non-symmetric peaks and complicated birth/death terms) to test its robustness.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 50

28

4.3 Discretization methods Discretization methods (also called classes or sectional methods) are, in theory, the preferred method to solve a general population balance because the distribution of the population is preserved. In MOM and QMOM (and other variants), the moments are calculated and the distribution information is lost. However, discretization techniques are more computationally expensive. Discretization methods consist of breaking up the internal coordinates of interest (e.g. crystal size) into finite (discrete) bins. Depending on the range of the distribution, the number of bins can be large and become computationally expensive. For example, if nucleation is considered in the nanometer range and crystal growth up to the micrometer or millimeter range, then to ensure the distribution is properly resolved, the bin size may be in the nanometer range, requiring hundreds (or more) bins. Each one of these bins requires an additional equation and can quickly increase the required computation to solve the population balance. Recently, computational power has increased dramatically and the bottleneck to using discretization methods has gradually decreased. Also, modern discretization algorithms are more efficient than the discretization schemes of the past. 4.3.1 Method of characteristics The method of characteristics involves converting a partial differential equation into an ordinary differential equation. The ordinary differential equation is solved along a hypersurface or characteristic and used to find the solution of the original partial differential equation. In the case of population balance equations, the method of characteristics discovers curves along the L − t plane (in the case of a dynamic formulation) which transforms the partial differential equation into an ordinary differential equation.141 The method of characteristics has some limitations including long calculation time when solving complex and practical systems, careful selection of


Page 29 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


29 time-steps and obligated scalar modelling.142 The method of characteristics is often coupled with other numerical methods10,143,144. In 1997, Kumar and Ramkrishna144 combined a discretization technique they previously developed145 with the method of characteristics. The method can handle nucleation, growth and aggregation. The method has the advantage of their previously developed discretization technique as well as the method of characteristics. In this method, the desired properties of the distribution are obtained to ensure computational efficiency. Mahoney et al. found the characteristics parallel with the size history of the particles under the assumptions of deterministic growth rate and no aggregation.146 The characteristics can then be associated with constant cumulative number densities that can be found from experimental data. From that, the equations are decoupled and can be used to find the growth rate.

Févotte and Févotte used the method of characteristics to solve a population balance

equation where crystallization was performed in the presence of impurities.147 Kubota and Mullin showed that crystallization was suppressed by the presence of impurities.148 The effect was dependant on the time the crystals spent with the impurities. Therefore, the population balance formulation for this situation requires an extra time variable to account for the time spent in contact with the impurities. Févotte and Févotte used the method of characteristics to solve the population balance equation and the results were consistent with industrial observations147. Aamir et al. combined the method of characteristics with QMOM to solve the population balance equation and reconstruct the distribution.10 The method showed promise to be used for offline and online batch crystallization optimization. Fast prediction of the evolution of the distribution with respect to time for a batch crystallization was achieved with this method. This method was further used by Aamir et al. to determine the optimal seed recipe design to control



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 50

30 crystal size distribution in a batch cooling crystallization process.149 This method was also used by Nagy et al. to look at internal fines removal.38 The method of characteristics was further extended to two dimensions by Jiang et al.150 Jiang et al. used the method of characteristics to solve a two-dimensional population balance model to find the growth and dissolution parameters. 4.3.2 Finite difference method The finite difference methods are characterized by approximating differential equations with difference equations. Explicit, implicit and Crank-Nicolson schemes are the most common finite difference methods. John and Suciu compared finite difference schemes of different order using a first order discretization and a third order discretization to discretize a bi-variate population balance equation.151 They found that the higher order method was more accurate. However, there is a trade-off between accuracy and computational cost. Also, finite difference discretization methods may cause non-physical oscillations. John and Suciu took this into consideration by using a monotone first order method that does not produce oscillations. Bennett and Rohani combined the Lax-Wendroff and the Crank-Nicholson methods to solve a population balance equation without agglomeration and breakage.152 The combined method was compared to the individual Lax-Wendroff and Crank-Nicholson methods. When solved until steady-state, the 3 methods produced nearly identical crystal size distributions. However, the Lax-Wendroff and the Crank-Nicholson methods showed some oscillations early in the simulation. When dissolution conditions were present, only the combined method produced 0 population balance; Lax-Wendroff was unstable and Crank-Nicholson was oscillatory. When the population balance was solved dynamically, no oscillations occurred using


Page 31 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


31 the combined method. The individual methods produced oscillations, but did not affect the end distribution. A method that can be used to solve combined CFD-population balance models is to use a traditional discretization algorithm to solve the CFD (e.g. finite element or finite volume), use the finite difference method to solve for the population balance and couple the CFD and population balance together. This approach is effective because CFD simulation has become a mature field and there are many commercial software programs to solve CFD problems. The population balance equation can be simply discretized using finite difference and then added to the CFD code. For example, John et al. used a combined finite element/finite difference approach to solve for the flow field and the precipitation process.153 The finite element method was used for the flow equations whereas the finite difference method was used for the population balance equation. Sheikhzadeh et al. used the population balance equation solved by the finite difference method for model based control of an anti-solvent semi-batch crystallization process.154 The method showed good results at controlling the crystallization process even with internal disturbances and model/process mismatch. 4.3.3 Finite element method In the finite element method, the partial differential equation (or strong form) is converted to the integral form (or weak form). The domain or geometry must then be discretized into a finite number of elements to numerically solve the integral equation and summed together. In the case of a steady-state problem, the partial differential equation is converted into a system of algebraic equations. In the case of a dynamic problem, the partial differential equation is converted into a system of ordinary differential equations. The system of equations (algebraic or differential) is



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 50

32 then solved using any of the many direct or iterative solver methods. The solution of the system of equations is then reconstructed into the solution of the original partial differential equation. The advantage of the finite element method is that unlike finite difference methods, the computational mesh does not have to be uniformly distributed. This provides an important advantage in situations where size ranges need more resolution (e.g. in the case of nucleation producing small crystals) over uniform discretization methods. Nicmanis and Hounslow developed a mixed Galerkin-collocation finite element method to solve steady-state population balance equations.155 This method was able to accommodate adaptive mesh refinement which increases computational efficiency and accuracy by controlling the mesh fineness. The developed algorithm showed ability to predict the solution of the population balance equations and the moments. The method avoided the ill-conditioned matrices that occur during problems with growth. Rigopoulos and Jones developed a first-order finite element scheme to solve dynamic one-dimensional population balance equations.156 The method could handle nucleation, growth, aggregation and breakage phenomena. The method used linear elements, decreasing the number of interpolations when compared to higher-order elements which in theory will decrease the computation time. However, this method may not be as accurate as higher order methods. The time savings can be used to refine the mesh and increase the accuracy. Tsang and Tao used a moving finite element method to solve a multicomponent population balance equation.157 The method showed good results when compared to upwinddifferencing and Smolarkiewicz methods. The method was compared with the analytical solutions and an error function. The advantage of this method is that the region of the mesh that


Page 33 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


33 is most fine can move with the distribution to ensure high resolution without being excessively computationally expensive by discretizing the entire distribution with a fine mesh. 4.3.4 Other discretization types This section includes discretization algorithms/schemes that do not fall in the preceding categories. The main issue with discretization methods is the need for fine discretization to accurately capture the crystallization behavior. Hu et al. developed a method using successive algebraic equations to solve a population balance equation including nucleation, growth, aggregation and breakage which can partially avoid the need for fine discretization by allowing non-uniform discretization101. The method is attractive since algebraic equations are computationally simple to solve and the method showed accuracy. However, this method may show issues when coupled with problems that involve external coordinates. Using external coordinates necessitates the use of partial differentials in the population balance equation which also require discretization by some scheme. This method should be coupled with spatial discretization schemes (e.g. finite difference, element or volume methods) to test its applicability in such situations. Another method to avoid using uniformly fine discretization was developed by Sewerin and Rigopoulos who developed an adaptive technique for discretizing a spatially inhomogeneous and unsteady population balance equation for one particle coordinate.158 The method uses the current solution of the particle number density distribution and current coordinate transformation to determine the future coordinate transformation. The adaptive technique was combined with the Galerkin finite element method, fully upwinded orthogonal collocation finite element method and a high resolution finite volume method and reduced the number of grid points necessary by more than one order of magnitude for the latter two methods.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 34 of 50

34 Another mini-class of methods that has been used in literature is the least squares method. Dorao and Jakobsen developed a general framework using the least squares method to solve a generalized population balance equation.159 The method was tested using cases with analytical solutions to determine its accuracy. It was found that the error in the properties of the distribution (e.g. moments) depended on the order of the expansion; increasing the order of expansion, decreased the error. Zhu et al. further developed the least squares method by including a direct minimization method.160 The addition of the direct minimization allows for solving systems with ill-conditioned matrices (high non-linearity and large scale) that the least squares method cannot solve. However, the direct minimization method does not always produce a symmetric, positive-definite system. Therefore, depending on the situation, these issues must be considered. Solsvik and Jakobsen compared the least squares method with Galerkin, tau and orthogonal collocation methods and found that the least squares method did not perform as well as the other two methods in the solution framework that was presented161. A high-resolution central discretization scheme was developed for nonlinear conservation laws and convection–diffusion equations.162 This discretization scheme was adapted to crystallization and qualitatively validated with existing literature by Woo et al. (2006).5 It was further adapted to an anti-solvent impinging jet crystallizer.6 Kumar et al.163 compared the cell average technique developed by Kumar et al.164 and the finite volume scheme developed by Filbet and Laurençot.165 Kumar et al. found that for fine grids, the finite volume scheme was the superior choice. However, for coarse grids, the cell average scheme gave reasonable results for the number density and the moments. Gunawan et al. extended high resolution finite volume methods that are used in other fields (e.g. aerodynamics, astrophysics, etc.) to multidimensional population balance equations.97


Page 35 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


35 They found that the algorithms gave accurate results for batch and continuous population balance models. They also found that since the developed algorithm had improved numerical accuracy, larger mesh sizes and longer time steps could be used without sacrificing accuracy. This resulted in an order of magnitude reduction in computational cost compared to existing finite difference/volume methods. Another approach that has recently gained attraction to increase computational efficiency for all types of solution methods is using advanced computing strategies (e.g. graphical processing units, parallel and distributed computing).166–170

4.4 Monte Carlo methods Monte Carlo methods have received considerable attention recently due to their adaptability to many different population balance equation types (e.g. nucleation, agglomeration, breakage, multidimensional, etc.). The Monte Carlo method involves inputting randomly generated numbers into the mathematical system and generating a set of solutions. The actual solution (or solutions) will appear more often in the set of solutions. By increasing the number of randomly generated input trials, the solution (or solutions) can be found more accurately. The advantage of the Monte Carlo methods for population balance equations deals with the stochastic nature of the methods, that naturally fit with particle processing (e.g. nucleation, agglomeration).171 Since the Monte Carlo methods use random numbers, they are very robust and closely resemble the type of crystallization processes involved. However, to get an accurate representation of the system, many particles must be tracked and this increases the computational expense. By tracking the particles, this principally becomes like the Lagrangian formulation of the population balance equation. This allows for a greater understanding and accuracy of the crystallization modelling.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 50

36 Ramkrishna discusses the connection between Monte Carlo method simulation and population balance in detail.172 Yu et al. performed analyses on the accuracy and optimal sampling using the Monte Carlo method on population balance equations.173 Implementation requires the choice of initial sample number, number of replicates and number of bins for probability distribution reconstruction. To find the accuracy, the Squared Hellinger Distance (H2) can be used and is related to initial sample number, number of replications and number of bins. It was found that changing the number of bins to increase the resolution, decreases the accuracy, and requires a trade-off by the user. There are different formulations of the Monte Carlo method that have been developed to solve population balance equations. Common ones include the constant-number Monte Carlo method174–176, constant-volume Monte Carlo method177,178, Multi-Monte Carlo179–182 and direct simulation Monte Carlo method.183 The methods can be grouped into two categories: time-driven (Multi and direct simulation) and event driven methods (constant-number and constant-volume). Time-driven methods involve time steps preceding the event, and event-driven methods involve time being calculated from the known event probability. The constant-number Monte Carlo method is performed by expanding/contracting the physical volume to keep the number of particles constant. To relate the sample of particles to the actual volume, two methods were used: constant mass concentration regardless of simulation volume and the constant number concentration regardless of simulation volume. Constantvolume Monte Carlo is characterized by an initial constant volume simulation. The number of particles in the volume changes based on the events. Whenever the number of simulation particles is doubled or decreased by half, the number of simulation particles is halved or doubled,


Page 37 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


37 respectively. Multi-Monte Carlo contains facets of both constant-number and constant-volume Monte Carlo. It involves assigning a statistical weight to each simulation particle that determines how many real particles are represented by each simulation particle. In this method, the statistical weight does not have to be constant for all particles, whereas in the other three methods, each simulation particle represents the same number of real particles. Direct simulation Monte Carlo involves using a period regulation of the number of simulation particles in the simulation box. Every particle is tracked, for each fixed time internal, to determine if the particle participates in an event. Zhao et al. compared the four types of Monte Carlo methods in accuracy in the distribution, accuracy in the number and mass concentration of particles and in the requirement of the CPU.184 It was found that the four methods provided similar accuracy and that the accuracy was primarily controlled by the number of simulation particles. For computation time, event-driven methods were better than time-driven methods since event-driven methods advance time with each event. Reasonably small numbers of particles, O(103), were enough to achieve good accuracy. Zhao and Zheng combined the Monte Carlo method with Eulerian-Lagrangian model for hydrodynamics description.185 In most combined CFD-population balance models, one-way coupling is assumed; the particles do not affect the flow field but in this paper, complete coupling was used. The method is predicated on using a time-step that effectively uncouples the fluid dynamics, particle transport and particle dynamics so that each is simulated successively. This method had a factor of ten decrease in computational time when compared to direct numerical simulation for the fluid dynamics but this method is still computationally expensive in its current form.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 38 of 50

38

5. Conclusions As mathematical modelling becomes more common in crystallization research and crystallizer design/control, efficient and accurate solutions of population balance equations become essential. The different numerical techniques to solve population balance equations were reviewed and compared. The different birth/death terms were also described. Analytical solutions are possible for various situations which often include many assumptions. Moments methods are a popular solution class that involve converting the population balance equations into ordinary differential equations using a moment transformation. Using the moments, the mean particle sizes and the coefficient of distribution can be found. However, information about the size distribution is lost. Some techniques exist to reconstruct the distribution from the moments but this is still a very active area of research. Discretization techniques provide the distribution but come at additional computational cost. As computers become more powerful and numerical algorithms become more efficient, discretization techniques become more valuable. However, in parallel, moments methods have also become more efficient. Moments methods are useful when coupled with computational fluid dynamics and especially in process control applications when quick calculations of the average crystal size are desired. Monte Carlo methods have also been increasingly applied on population balance equations since the methods naturally adapt to solving population balance equations and are robust in dealing with different crystallization sub-processes (e.g. nucleation, growth, aggregation/agglomeration and breakage). However, like discretization techniques, Monte Carlo methods can be computationally expensive to accurately represent a crystallization system.


Page 39 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


39 Research on population balance equations has made tremendous advances. Future developments could include using phenomenological models instead of the semi-empirical (e.g. power law models) that are currently commonly used. Moreover, more realistic aggregation/agglomeration and breakage kernels may prove to provide a better representation of the underlying processes. Even further in the future, quantum mechanics modelling may be included into population balance modelling. Polymorphic distributions may be possible to elucidate from the combined quantum mechanics and population balance models.

Acknowledgements The authors acknowledge the financial support provided by the Natural Sciences and Engineering Council of Canada (NSERC).



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 40 of 50

40

References (1)

(2) (3) (4) (5)

(6)

(7) (8) (9)

(10)

(11)

(12)

(13)

(14)

(15)

(16)

Bogacz, W.; Al-Rashed, M. H.; Lemanowicz, M.; Wójcik, J. A Simple Densimetric Method to Determine Saturation Temperature of Aqueous Potassium Chloride Solution. J. Solut. Chem. 2016, 45, 1071–1076. Piton, D.; Fox, R. O.; Marcant, B. Simulation of fine particle formation by precipitation using computational fluid dynamics. Can. J. Chem. Eng. 2000, 78, 983–993. Zauner, R.; Jones, A. G. On the influence of mixing on crystal precipitation processes— application of the segregated feed model. Chem. Eng. Sci. 2002, 57, 821–831. Kougoulos, E.; Jones, A. G.; Wood-Kaczmar, M. CFD Modelling of Mixing and Heat Transfer in Batch Cooling Crystallizers. Chem. Eng. Res. Des. 2005, 83, 30–39. Woo, X. Y.; Tan, R. B. H.; Chow, P. S.; Braatz, R. D. Simulation of Mixing Effects in Antisolvent Crystallization Using a Coupled CFD-PDF-PBE Approach. Cryst. Growth Des. 2006, 6, 1291–1303. Woo, X. Y.; Tan, R. B. H.; Braatz, R. D. Modeling and Computational Fluid Dynamics−Population Balance Equation−Micromixing Simulation of Impinging Jet Crystallizers. Cryst. Growth Des. 2009, 9, 156–164. Chiu, T. Y.; Christofides, P. D. Robust control of particulate processes using uncertain population balances. AIChE J. 2000, 46, 266–280. Shi, D.; El-Farra, N. H.; Li, M.; Mhaskar, P.; Christofides, P. D., P. D. Predictive control of particle size distribution in particulate processes. Chem. Eng. Sci. 2006, 61, 268–281. Mesbah, A.; Kramer, H. J. M.; Huesman, A. E. M.; Van den Hof, P. M. J. A control oriented study on the numerical solution of the population balance equation for crystallization processes. Chem. Eng. Sci. 2009, 64, 4262–4277. Aamir, E.; Nagy, Z. K.; Rielly, C. D.; Kleinert, T.; Judat, B. Combined Quadrature Method of Moments and Method of Characteristics Approach for Efficient Solution of Population Balance Models for Dynamic Modeling and Crystal Size Distribution Control of Crystallization Processes. Ind. Eng. Chem. Res. 2009, 48, 8575–8584. Norton, T.; Sun, D.-W. Computational fluid dynamics (CFD) – an effective and efficient design and analysis tool for the food industry: A review. Trends Food Sci. Technol. 2006, 17, 600–620. Fimbres-Weihs, G. A.; Wiley, D. E. Review of 3D CFD modeling of flow and mass transfer in narrow spacer-filled channels in membrane modules. Chem. Eng. Process. Process Intensif. 2010, 49, 759–781. Kuriakose, R.; Anandharamakrishnan, C. Computational fluid dynamics (CFD) applications in spray drying of food products. Trends Food Sci. Technol. 2010, 21, 383– 398. Aslam Bhutta, M. M.; Hayat, N.; Bashir, M. H.; Khan, A. R.; Ahmad, K. N.; Khan, S. CFD applications in various heat exchangers design: A review. Appl. Therm. Eng. 2012, 32, 1–12. Crowley, T. J.; Meadows, E. S.; Kostoulas, E.; Doyle III, F. J. Control of particle size distribution described by a population balance model of semibatch emulsion polymerization. J. Process Control 2000, 10, 419–432. Immanuel, C. D.; Doyle III, F. J. Computationally efficient solution of population balance models incorporating nucleation, growth and coagulation: application to emulsion polymerization. Chem. Eng. Sci. 2003, 58, 3681–3698.


Page 41 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


41 (17)

(18) (19) (20) (21)

(22) (23) (24) (25)

(26) (27) (28) (29) (30)

(31) (32) (33) (34) (35) (36)

Kotoulas, C.; Kiparissides, C. A generalized population balance model for the prediction of particle size distribution in suspension polymerization reactors. Chem. Eng. Sci. 2006, 61, 332–346. Kiparissides, C. Challenges in particulate polymerization reactor modeling and optimization: A population balance perspective. J. Process Control 2006, 16, 205–224. McGraw, R. Description of Aerosol Dynamics by the Quadrature Method of Moments. Aerosol Sci. Technol. 1997, 27, 255–265. McGraw, R.; Nemesure, S.; Schwartz, S. E. Properties and evolution of aerosols with size distributions having identical moments. J. Aerosol Sci. 1998, 29, 761–772. Wright, D. L.; McGraw, R.; Rosner, D. E. Bivariate Extension of the Quadrature Method of Moments for Modeling Simultaneous Coagulation and Sintering of Particle Populations. J. Colloid Interface Sci. 2001, 236, 242–251. Marchisio, D. L.; Fox, R. O. Solution of population balance equations using the direct quadrature method of moments. J. Aerosol Sci. 2005, 36, 43–73. Fredrickson, A. G.; Mantzaris, N. V. A new set of population balance equations for microbial and cell cultures. Chem. Eng. Sci. 2002, 57, 2265–2278. Mantzaris, N. V.; Daoutidis, P. Cell population balance modeling and control in continuous bioreactors. J. Process Control 2004, 14, 775–784. Fernandes, R. L.; Carlquist, M.; Lundin, L.; Heins, A.-L.; Dutta, A.; Sørensen, S. J.; Jensen, A. D.; Nopens, I.; Lantz, A. E.; Gernaey, K. V. Cell mass and cell cycle dynamics of an asynchronous budding yeast population: Experimental observations, flow cytometry data analysis, and multi-scale modeling. Biotechnol. Bioeng. 2013, 110, 812– 826. Ramkrishna, D.; Singh, M. R. Population Balance Modeling: Current Status and Future Prospects. Annu. Rev. Chem. Biomol. Eng. 2014, 5, 123–146. Solsvik, J.; Jakobsen, H. A. The Foundation of the Population Balance Equation: A Review. J. Dispers. Sci. Technol. 2015, 36, 510–520. Ramkrishna, D. Population Balances: Theory and Applications to Particulate Systems in Engineering, 1st ed.; Academic Press: San Diego, CA, 2000. Handbook of Industrial Crystallization: Second Edition, 2nd ed.; Myerson, A. S., Ed.; Butterworth-Heinemann, 2002. Yeoh, G. H.; Cheung, C. P.; Tu, J. Multiphase Flow Analysis Using Population Balance Modeling: Bubbles, Drops and Particles, 1st ed.; Butterworth-Heinemann: Amsterdam, 2013. Ramkrishna, D. The Status of Population Balances. Rev. Chem. Eng. 1985, 3, 49–95. Hulburt, H. M.; Katz, S. Some problems in particle technology. Chem. Eng. Sci. 1964, 19, 555–574. Randolph, A. D. A population balance for countable entities. Can. J. Chem. Eng. 1964, 42, 280–281. Ponthot, J.-P.; Belytschko, T. Arbitrary Lagrangian-Eulerian formulation for elementfree Galerkin method. Comput. Methods Appl. Mech. Eng. 1998, 152, 19–46. Ma, D. L.; Tafti, D. K.; Braatz, R. D. High-Resolution Simulation of Multidimensional Crystal Growth. Ind. Eng. Chem. Res. 2002, 41, 6217–6223. Briesen, H. Simulation of crystal size and shape by means of a reduced two-dimensional population balance model. Chem. Eng. Sci. 2006, 61, 104–112.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 42 of 50

42 (37)

(38)

(39) (40)

(41)

(42) (43) (44) (45) (46) (47)

(48)

(49) (50) (51) (52) (53) (54)

(55)

Gerstlauer, A.; Mitrović, A.; Motz, S.; Gilles, E.-D. A population model for crystallization processes using two independent particle properties. Chem. Eng. Sci. 2001, 56, 2553–2565. Nagy, Z. K.; Aamir, E.; Rielly, C. D. Internal Fines Removal Using Population Balance Model Based Control of Crystal Size Distribution under Dissolution, Growth and Nucleation Mechanisms. Cryst. Growth Des. 2011, 11, 2205–2219. Borsos, Á.; Lakatos, B. G. Investigation and simulation of crystallization of high aspect ratio crystals with fragmentation. Chem. Eng. Res. Des. 2014, 92, 1133–1141. Szilágyi, B.; Lakatos, B. G. Model-based analysis of stirred cooling crystallizer of high aspect ratio crystals with linear and nonlinear breakage. Comput. Chem. Eng. 2017, 98, 180–196. Qamar, S.; Warnecke, G.; Elsner, M. P. On the solution of population balances for nucleation, growth, aggregation and breakage processes. Chem. Eng. Sci. 2009, 64, 2088–2095. Scott, W. T. Analytic Studies of Cloud Droplet Coalescence I. J. Atmospheric Sci. 1968, 25, 54–65. Lister, J. D.; Smit, D. J.; Hounslow, M. J. Adjustable discretized population balance for growth and aggregation. AIChE J. 1995, 41, 591–603. Friedlander, S. K. Smoke, dust and haze: Fundamentals of aerosol behavior; Wiley: New York, 1977. Tavare, N. S.; Patwardhan, V. Agglomeration in a continuous MSMPR crystallizer. AIChE J. 1992, 38, 377–384. Marchisio, D. L.; Vigil, R. D.; Fox, R. O. Quadrature method of moments for aggregation–breakage processes. J. Colloid Interface Sci. 2003, 258, 322–334. Falola, A.; Borissova, A.; Wang, X. Z. Extended method of moment for general population balance models including size dependent growth rate, aggregation and breakage kernels. Comput. Chem. Eng. 2013, 56, 1–11. El-Zhry El-Yafi, A. K.; El-Zein, H. Technical crystallization for application in pharmaceutical material engineering: Review article. Asian J. Pharm. Sci. 2015, 10, 283– 291. Viedma, C. Experimental evidence of chiral symmetry breaking in crystallization from primary nucleation. J. Cryst. Growth 2004, 261, 118–121. Garside, J. Industrial crystallization from solution. Chem. Eng. Sci. 1985, 40, 3–26. Maggioni, G. M.; Mazzotti, M. Modelling the stochastic behaviour of primary nucleation. Faraday Discuss. 2015, 179, 359–382. Pino-García, O.; Rasmuson, Å. C. Primary Nucleation of Vanillin Explored by a Novel Multicell Device. Ind. Eng. Chem. Res. 2003, 42, 4899–4909. Jiang, S.; ter Horst, J. H. Crystal Nucleation Rates from Probability Distributions of Induction Times. Cryst. Growth Des. 2011, 11, 256–261. Kadam, S. S.; Kulkarni, S. A.; Coloma Ribera, R.; Stankiewicz, A. I.; ter Horst, J. H.; Kramer, H. J. M. A new view on the metastable zone width during cooling crystallization. Chem. Eng. Sci. 2012, 72, 10–19. Sullivan, R. A.; Davey, R. J.; Sadiq, G.; Dent, G.; Back, K. R.; ter Horst, J. H.; Toroz, D.; Hammond, R. B. Revealing the Roles of Desolvation and Molecular Self-Assembly in Crystal Nucleation from Solution: Benzoic and p-Aminobenzoic Acids. Cryst. Growth Des. 2014, 14, 2689–2696.


Page 43 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


43 (56) (57)

(58) (59) (60) (61) (62) (63) (64)

(65) (66) (67) (68)

(69)

(70) (71)

(72) (73) (74) (75)

Sear, R. P. Quantitative studies of crystal nucleation at constant supersaturation: experimental data and models. CrystEngComm 2014, 16, 6506–6522. Goh, L.; Chen, K.; Bhamidi, V.; He, G.; Kee, N. C. S.; Kenis, P. J. A.; Zukoski, C. F.; Braatz, R. D. A Stochastic Model for Nucleation Kinetics Determination in DropletBased Microfluidic Systems. Cryst. Growth Des. 2010, 10, 2515–2521. Volmer, M.; Weber, Α. Keimbildung in übersättigten Gebilden. Z. Für Phys. Chem. 1926, 119, 277–301. Becker, R.; Döring, W. Kinetische Behandlung der Keimbildung in übersättigten Dämpfen. Ann. Phys. 1935, 416, 719–752. Vekilov, P. G. Dense Liquid Precursor for the Nucleation of Ordered Solid Phases from Solution. Cryst. Growth Des. 2004, 4, 671–685. Vekilov, P. G. The two-step mechanism of nucleation of crystals in solution. Nanoscale 2010, 2, 2346–2357. Erdemir, D.; Lee, A. Y.; Myerson, A. S. Nucleation of Crystals from Solution: Classical and Two-Step Models. Acc. Chem. Res. 2009, 42, 621–629. Karthika, S.; Radhakrishnan, T. K.; Kalaichelvi, P. A Review of Classical and Nonclassical Nucleation Theories. Cryst. Growth Des. 2016, 16, 6663–6681. Forsyth, C.; Mulheran, P. A.; Forsyth, C.; Haw, M. D.; Burns, I. S.; Sefcik, J. Influence of Controlled Fluid Shear on Nucleation Rates in Glycine Aqueous Solutions. Cryst. Growth Des. 2015, 15, 94–102. Liu, J.; Rasmuson, Å. C. Influence of Agitation and Fluid Shear on Primary Nucleation in Solution. Cryst. Growth Des. 2013, 13, 4385–4394. Yang, H.; Yu, X.; Raval, V.; Makkawi, Y.; Florence, A. Effect of Oscillatory Flow on Nucleation Kinetics of Butyl Paraben. Cryst. Growth Des. 2016, 16, 875–886. Virone, C.; Kramer, H. J. M.; van Rosmalen, G. M.; Stoop, A. H.; Bakker, T. W. Primary nucleation induced by ultrasonic cavitation. J. Cryst. Growth 2006, 294, 9–15. Miyasaka, E.; Ebihara, S.; Hirasawa, I. Investigation of primary nucleation phenomena of acetylsalicylic acid crystals induced by ultrasonic irradiation—ultrasonic energy needed to activate primary nucleation. J. Cryst. Growth 2006, 295, 97–101. Kurotani, M.; Miyasaka, E.; Ebihara, S.; Hirasawa, I. Effect of ultrasonic irradiation on the behavior of primary nucleation of amino acids in supersaturated solutions. J. Cryst. Growth 2009, 311, 2714–2721. Tan, P.; Xu, N.; Xu, L. Visualizing kinetic pathways of homogeneous nucleation in colloidal crystallization. Nat. Phys. 2014, 10, 73–79. Kadam, S. S.; Kramer, H. J. M.; ter Horst, J. H. Combination of a Single Primary Nucleation Event and Secondary Nucleation in Crystallization Processes. Cryst. Growth Des. 2011, 11, 1271–1277. Agrawal, S. G.; Paterson, A. H. J. Secondary Nucleation: Mechanisms and Models. Chem. Eng. Commun. 2015, 202, 698–706. Evans, T. W.; Margolis, G.; Sarofim, A. F. Mechanisms of secondary nucleation in agitated crystallizers. AIChE J. 1974, 20, 950–958. Denk, E. G.; Botsaris, G. D. Mechanism of contact nucleation. J. Cryst. Growth 1972, 15, 57–60. Tai, C. Y.; McCabe, W. L.; Rousseau, R. W. Contact nucleation of various crystal types. AIChE J. 1975, 21, 351–358.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 44 of 50

44 (76)

(77) (78) (79)

(80)

(81)

(82) (83)

(84)

(85) (86) (87) (88) (89) (90)

(91) (92) (93)

Wong, S. Y.; Cui, Y.; Myerson, A. S. Contact Secondary Nucleation as a Means of Creating Seeds for Continuous Tubular Crystallizers. Cryst. Growth Des. 2013, 13, 2514–2521. Sung, C. Y.; Estrin, J.; Youngquist, G. R. Secondary nucleation of magnesium sulfate by fluid shear. AIChE J. 1973, 19, 957–962. Mazzarotta, B.; Di Cave, S.; Bonifazi, G. Influence of time on crystal attrition in a stirred vessel. AIChE J. 1996, 42, 3554–3558. Bravi, M.; Di Cave, S.; Mazzarotta, B.; Verdone, N. Relating the attrition behaviour of crystals in a stirred vessel to their mechanical properties. Chem. Eng. J. 2003, 94, 223– 229. Asakuma, Y.; Tetahima, T.; Maeda, K.; Hideo, M.; Fukeu, D. Attrition behavior by micro-hardness parameters in suspension-crystallization processes. Powder Technol. 2007, 171, 75–80. Sato, K.; Nagai, H.; Hasegawa, K.; Tomori, K.; Kramer, H. J. M.; Jansens, P. J. Twodimensional population balance model with breakage of high aspect ratio crystals for batch crystallization. Chem. Eng. Sci. 2008, 63, 3271–3278. Briesen, H. Two-dimensional population balance modeling for shape dependent crystal attrition. Chem. Eng. Sci. 2009, 64, 661–672. Guo, Z.; Jones, A. G.; Li, N.; Germana, S. High-speed observation of the effects of ultrasound on liquid mixing and agglomerated crystal breakage processes. Powder Technol. 2007, 171, 146–153. Blagden, N.; de Matas, M.; Gavan, P. T.; York, P. Crystal engineering of active pharmaceutical ingredients to improve solubility and dissolution rates. Adv. Drug Deliv. Rev. 2007, 59, 617–630. Noyes, A. A.; Whitney, W. R. The Rate Of Solution Of Solid Substances In Their Own Solutions. J. Am. Chem. Soc. 1897, 19, 930–934. Lu, A. T. K.; Frisella, M. E.; Johnson, K. C. Dissolution Modeling: Factors Affecting the Dissolution Rates of Polydisperse Powders. Pharm. Res. 1993, 10, 1308–1314. Shan, G.; Igarashi, K.; Ooshima, H. Dissolution kinetics of crystals in suspension and its application to l-aspartic acid crystals. Chem. Eng. J. 2002, 88, 53–58. Mangin, D.; Garcia, E.; Gerard, S.; Hoff, C.; Klein, J. P.; Veesler, S. Modeling of the dissolution of a pharmaceutical compound. J. Cryst. Growth 2006, 286, 121–125. Garcia, E.; Hoff, C.; Veesler, S. Dissolution and phase transition of pharmaceutical compounds. J. Cryst. Growth 2002, 237–239, Part 3, 2233–2239. Nichols, G.; Byard, S.; Bloxham, M. J.; Botterill, J.; Dawson, N. J.; Dennis, A.; Diart, V.; North, N. C.; Sherwood, J. D. A review of the terms agglomerate and aggregate with a recommendation for nomenclature used in powder and particle characterization. J. Pharm. Sci. 2002, 91, 2103–2109. Brunsteiner, M.; Jones, A. G.; Pratola, F.; Price, S. L.; Simons, S. J. R. Toward a Molecular Understanding of Crystal Agglomeration. Cryst. Growth Des. 2005, 5, 3–16. David, R.; Marchal, P.; Marcant, B. Modelling of agglomeration in industrial crystallization from solution. Chem. Eng. Technol. 1995, 18, 302–309. David, R.; Marchal, P.; Klein, J.-P.; Villermaux, J. Crystallization and precipitation engineering—III. A discrete formulation of the agglomeration rate of crystals in a crystallization process. Chem. Eng. Sci. 1991, 46, 205–213.


Page 45 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


45 (94)

(95)

(96)

(97) (98)

(99) (100) (101) (102)

(103) (104)

(105)

(106) (107)

(108)

(109)

(110)

Yu, Z. Q.; Tan, R. B. H.; Chow, P. S. Effects of operating conditions on agglomeration and habit of paracetamol crystals in anti-solvent crystallization. J. Cryst. Growth 2005, 279, 477–488. Peña, R.; Burcham, C. L.; Jarmer, D. J.; Ramkrishna, D.; Nagy, Z. K. Modeling and optimization of spherical agglomeration in suspension through a coupled population balance model. Chem. Eng. Sci. Fujiwara, M.; Chow, P. S.; Ma, D. L.; Braatz, R. D. Paracetamol Crystallization Using Laser Backscattering and ATR-FTIR Spectroscopy: Metastability, Agglomeration, and Control. Cryst. Growth Des. 2002, 2, 363–370. Gunawan, R.; Fusman, I.; Braatz, R. D. High resolution algorithms for multidimensional population balance equations. AIChE J. 2004, 50, 2738–2749. Morris, G.; Power, G.; Ferguson, S.; Barrett, M.; Hou, G.; Glennon, B. Estimation of Nucleation and Growth Kinetics of Benzoic Acid by Population Balance Modeling of a Continuous Cooling Mixed Suspension, Mixed Product Removal Crystallizer. Org. Process Res. Dev. 2015, 19, 1891–1902. McCabe, W. L. Crystal Growth in Aqueous Solutions. Ind. Eng. Chem. 1929, 21, 112– 119. Hounslow, M. J.; Ryall, R. L.; Marshall, V. R. A discretized population balance for nucleation, growth, and aggregation. AIChE J. 1988, 34, 1821–1832. Hu, Q.; Rohani, S.; Jutan, A. New numerical method for solving the dynamic population balance equations. AIChE J. 2005, 51, 3000–3006. Ridder, B. J.; Majumder, A.; Nagy, Z. K. Population Balance Model-Based Multiobjective Optimization of a Multisegment Multiaddition (MSMA) Continuous Plug-Flow Antisolvent Crystallizer. Ind. Eng. Chem. Res. 2014, 53, 4387–4397. Sangwal, K. Effects of impurities on crystal growth processes. Prog. Cryst. Growth Charact. Mater. 1996, 32, 3–43. Garnier, S.; Petit, S.; Coquerel, G. Influence of supersaturation and structurally related additives on the crystal growth of α-lactose monohydrate. J. Cryst. Growth 2002, 234, 207–219. Gu, C.-H.; Grant, D. J. W. Relationship Between Particle Size and Impurity Incorporation During Crystallization of (+)-Pseudoephedrine Hydrochloride, Acetaminophen, and Adipic Acid from Aqueous Solution. Pharm. Res. 2002, 19, 1068– 1070. Kubota, N. Effect of Impurities on the Growth Kinetics of Crystals. Cryst. Res. Technol. 2001, 36, 749–769. Thompson, C.; Davies, M. C.; Roberts, C. J.; Tendler, S. J. B.; Wilkinson, M. J. The effects of additives on the growth and morphology of paracetamol (acetaminophen) crystals. Int. J. Pharm. 2004, 280, 137–150. Zhang, Y.; Liu, J. J.; Wan, J.; Wang, X. Z. Two dimensional population balance modelling of crystal growth behaviour under the influence of impurities. Adv. Powder Technol. 2015, 26, 672–678. Patil, D. P.; Andrews, J. R. G. An analytical solution to continuous population balance model describing floc coalescence and breakage — A special case. Chem. Eng. Sci. 1998, 53, 599–601. McCoy, B. J.; Madras, G. Analytical Solution for a Population Balance Equation with Aggregation and Fragmentation. Chem. Eng. Sci. 2003, 58, 3049–3051.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 46 of 50

46 (111) Ramabhadran, T. E.; Peterson, T. W.; Seinfeld, J. H. Dynamics of aerosol coagulation and condensation. AIChE J. 1976, 22, 840–851. (112) Brock, J. R. Condensational growth of atmospheric aerosols. J. Colloid Interface Sci. 1972, 39, 32–36. (113) Alvarez, A. J.; Singh, A.; Myerson, A. S. Crystallization of Cyclosporine in a Multistage Continuous MSMPR Crystallizer. Cryst. Growth Des. 2011, 11, 4392–4400. (114) Randolph, A. D.; Larson, M. A. Transient and steady state size distributions in continuous mixed suspension crystallizers. AIChE J. 1962, 8, 639–645. (115) Randolph, A. D.; Larson, M. A. Theory of Particulate Processes: Analysis and Techniques of Continuous Crystallization; Academic Press, 1971. (116) Rohani, S.; Bourne, J. R. A simplified approach to the operation of a batch crystallizer. Can. J. Chem. Eng. 1990, 68, 799–806. (117) Sang-Il Kwon, J.; Nayhouse, M.; Orkoulas, G.; Christofides, P. D. Crystal shape and size control using a plug flow crystallization configuration. Chem. Eng. Sci. 2014, 119, 30– 39. (118) Gordon, R. G. Error Bounds in Equilibrium Statistical Mechanics. J. Math. Phys. 1968, 9, 655–663. (119) Marchisio, D. L.; Pikturna, J. T.; Fox, R. O.; Vigil, R. D.; Barresi, A. A. Quadrature method of moments for population-balance equations. AIChE J. 2003, 49, 1266–1276. (120) Gavi, E.; Rivautella, L.; Marchisio, D. L.; Vanni, M.; Barresi, A. A.; Baldi, G. CFD Modelling of Nano-Particle Precipitation in Confined Impinging Jet Reactors. Chem. Eng. Res. Des. 2007, 85, 735–744. (121) Zheng, D.; Li, J.; Jin, Y.; Zou, D.; Zhu, X. Coupling of CFD with PBM for growth behavior of potassium sulfate in spray fluidized-bed crystallizer. Powder Technol. (122) Rane, C. V.; Ganguli, A. A.; Kalekudithi, E.; Patil, R. N.; Joshi, J. B.; Ramkrishna, D. CFD simulation and comparison of industrial crystallizers. Can. J. Chem. Eng. 2014, 92, 2138–2156. (123) Fan, R.; Marchisio, D. L.; Fox, R. O. Application of the direct quadrature method of moments to polydisperse gas–solid fluidized beds. Powder Technol. 2004, 139, 7–20. (124) Attarakih, M. M.; Drumm, C.; Bart, H.-J. Solution of the population balance equation using the sectional quadrature method of moments (SQMOM). Chem. Eng. Sci. 2009, 64, 742–752. (125) Yuan, C.; Laurent, F.; Fox, R. O. An extended quadrature method of moments for population balance equations. J. Aerosol Sci. 2012, 51, 1–23. (126) Alopaeus, V.; Laakkonen, M.; Aittamaa, J. Numerical solution of moment-transformed population balance equation with fixed quadrature points. Chem. Eng. Sci. 2006, 61, 4919–4929. (127) Yuan, C.; Fox, R. O. Conditional quadrature method of moments for kinetic equations. J. Comput. Phys. 2011, 230, 8216–8246. (128) Kariwala, V.; Cao, Y.; Nagy, Z. K. Automatic differentiation-based quadrature method of moments for solving population balance equations. AIChE J. 2012, 58, 842–854. (129) McGraw, R.; Wright, D. L. Chemically resolved aerosol dynamics for internal mixtures by the quadrature method of moments. J. Aerosol Sci. 2003, 34, 189–209. (130) Diemer, R. B.; Olson, J. H. A moment methodology for coagulation and breakage problems: Part 2—moment models and distribution reconstruction. Chem. Eng. Sci. 2002, 57, 2211–2228.


Page 47 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


47 (131) John, V.; Angelov, I.; Öncül, A. A.; Thévenin, D. Techniques for the reconstruction of a distribution from a finite number of its moments. Chem. Eng. Sci. 2007, 62, 2890–2904. (132) de Souza, L. G. M.; Janiga, G.; John, V.; Thévenin, D. Reconstruction of a distribution from a finite number of moments with an adaptive spline-based algorithm. Chem. Eng. Sci. 2010, 65, 2741–2750. (133) Lebaz, N.; Cockx, A.; Spérandio, M.; Morchain, J. Reconstruction of a distribution from a finite number of its moments: A comparative study in the case of depolymerization process. Comput. Chem. Eng. 2016, 84, 326–337. (134) Wulkow, M. The simulation of molecular weight distributions in polyreaction kinetics by discrete Galerkin methods. Macromol. Theory Simul. 1996, 5, 393–416. (135) Flood, A. E. Thoughts on Recovering Particle Size Distributions from the Moment Form of the Population Balance. Dev. Chem. Eng. Miner. Process. 2002, 10, 501–519. (136) Tagliani, A. Maximum entropy solutions and moment problem in unbounded domains. Appl. Math. Lett. 2003, 16, 519–524. (137) Biswas, P.; Bhattacharya, A. K. Function reconstruction as a classical moment problem: a maximum entropy approach. J. Phys. Math. Theor. 2010, 43, 405003. (138) Abramov, R. A practical computational framework for the multidimensional momentconstrained maximum entropy principle. J. Comput. Phys. 2006, 211, 198–209. (139) Abramov, R. V. An improved algorithm for the multidimensional moment-constrained maximum entropy problem. J. Comput. Phys. 2007, 226, 621–644. (140) Hutton, K.; Mitchell, N.; Frawley, P. J. Particle size distribution reconstruction: The moment surface method. Powder Technol. 2012, 222, 8–14. (141) Peng, Y.; Zhu, Z.; Braatz, R. D.; Myerson, A. S. Gypsum Crystallization during Phosphoric Acid Production: Modeling and Experiments Using the Mixed-SolventElectrolyte Thermodynamic Model. Ind. Eng. Chem. Res. 2015, 54, 7914–7924. (142) Lim, Y. I.; Le Lann, J.-M.; Meyer, X. M.; Joulia, X.; Lee, G.; Yoon, E. S. On the solution of population balance equations (PBE) with accurate front tracking methods in practical crystallization processes. Chem. Eng. Sci. 2002, 57, 3715–3732. (143) Qamar, S.; Mukhtar, S.; Seidel-Morgenstern, A. Efficient solution of a batch crystallization model with fines dissolution. J. Cryst. Growth 2010, 312, 2936–2945. (144) Kumar, S.; Ramkrishna, D. On the solution of population balance equations by discretization—III. Nucleation, growth and aggregation of particles. Chem. Eng. Sci. 1997, 52, 4659–4679. (145) Kumar, S.; Ramkrishna, D. On the solution of population balance equations by discretization—I. A fixed pivot technique. Chem. Eng. Sci. 1996, 51, 1311–1332. (146) Mahoney, A. W.; Doyle, F. J.; Ramkrishna, D. Inverse problems in population balances: Growth and nucleation from dynamic data. AIChE J. 2002, 48, 981–990. (147) Févotte, F.; Févotte, G. A method of characteristics for solving population balance equations (PBE) describing the adsorption of impurities during crystallization processes. Chem. Eng. Sci. 2010, 65, 3191–3198. (148) Kubota, N.; Mullin, J. W. A kinetic model for crystal growth from aqueous solution in the presence of impurity. J. Cryst. Growth 1995, 152, 203–208. (149) Aamir, E.; Nagy, Z. K.; Rielly, C. D. Optimal seed recipe design for crystal size distribution control for batch cooling crystallisation processes. Chem. Eng. Sci. 2010, 65, 3602–3614.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 48 of 50

48 (150) Jiang, M.; Zhu, X.; Molaro, M. C.; Rasche, M. L.; Zhang, H.; Chadwick, K.; Raimondo, D. M.; Kim, K.-K. K.; Zhou, L.; Zhu, Z.; Wong, M. H.; O’Grady, D.; Hebrault, D.; Tedesco, J.; Braatz, R. D. Modification of Crystal Shape through Deep Temperature Cycling. Ind. Eng. Chem. Res. 2014, 53, 5325–5336. (151) John, V.; Suciu, C. Direct discretizations of bi-variate population balance systems with finite difference schemes of different order. Chem. Eng. Sci. 2014, 106, 39–52. (152) Bennett, M. K.; Rohani, S. Solution of population balance equations with a new combined Lax-Wendroff/Crank-Nicholson method. Chem. Eng. Sci. 2001, 56, 6623– 6633. (153) John, V.; Mitkova, T.; Roland, M.; Sundmacher, K.; Tobiska, L.; Voigt, A. Simulations of population balance systems with one internal coordinate using finite element methods. Chem. Eng. Sci. 2009, 64, 733–741. (154) Sheikhzadeh, M.; Trifkovic, M.; Rohani, S. Real-time optimal control of an anti-solvent isothermal semi-batch crystallization process. Chem. Eng. Sci. 2008, 63, 829–839. (155) Nicmanis, M.; Hounslow, M. J. Finite-element methods for steady-state population balance equations. AIChE J. 1998, 44, 2258–2272. (156) Rigopoulos, S.; Jones, A. G. Finite-element scheme for solution of the dynamic population balance equation. AIChE J. 2003, 49, 1127–1139. (157) Tsang, T. H.; Rao, A. A moving finite element method for the population balance equation. Int. J. Numer. Methods Fluids 1990, 10, 753–769. (158) Sewerin, F.; Rigopoulos, S. An explicit adaptive grid approach for the numerical solution of the population balance equation. Chem. Eng. Sci. (159) Dorao, C. A.; Jakobsen, H. A. A least squares method for the solution of population balance problems. Comput. Chem. Eng. 2006, 30, 535–547. (160) Zhu, Z.; Dorao, C. A.; Jakobsen, H. A. A least-squares method with direct minimization for the solution of the breakage–coalescence population balance equation. Math. Comput. Simul. 2008, 79, 716–727. (161) Solsvik, J.; Jakobsen, H. A. Evaluation of the Least-Squares Method for the Solution of the Population Balance Equation. J. Dispers. Sci. Technol. 2015, 36, 68–82. (162) Kurganov, A.; Tadmor, E. New High-Resolution Central Schemes for Nonlinear Conservation Laws and Convection–Diffusion Equations. J. Comput. Phys. 2000, 160, 241–282. (163) Kumar, J.; Warnecke, G.; Peglow, M.; Heinrich, S. Comparison of numerical methods for solving population balance equations incorporating aggregation and breakage. Powder Technol. 2009, 189, 218–229. (164) Kumar, J.; Peglow, M.; Warnecke, G.; Heinrich, S.; Mörl, L. Improved accuracy and convergence of discretized population balance for aggregation: The cell average technique. Chem. Eng. Sci. 2006, 61, 3327–3342. (165) Filbet, F.; Laurençot, P. Numerical Simulation of the Smoluchowski Coagulation Equation. SIAM J. Sci. Comput. 2004, 25, 2004–2028. (166) Gunawan, R.; Fusman, I.; Braatz, R. D. Parallel high-resolution finite volume simulation of particulate processes. AIChE J. 2008, 54, 1449–1458. (167) Prakash, A. V.; Chaudhury, A.; Barrasso, D.; Ramachandran, R. Simulation of population balance model-based particulate processes via parallel and distributed computing. Chem. Eng. Res. Des. 2013, 91, 1259–1271.


Page 49 of 50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


49 (168) Wei, J.; Kruis, F. E. GPU-accelerated Monte Carlo simulation of particle coagulation based on the inverse method. J. Comput. Phys. 2013, 249, 67–79. (169) Xu, Z.; Zhao, H.; Zheng, C. Accelerating population balance-Monte Carlo simulation for coagulation dynamics from the Markov jump model, stochastic algorithm and GPU parallel computing. J. Comput. Phys. 2015, 281, 844–863. (170) Szilágyi, B.; Nagy, Z. K. Graphical processing unit (GPU) acceleration for numerical solution of population balance models using high resolution finite volume algorithm. Comput. Chem. Eng. 2016, 91, 167–181. (171) Hao, X.; Zhao, H.; Xu, Z.; Zheng, C. Population Balance-Monte Carlo Simulation for Gas-to-Particle Synthesis of Nanoparticles. Aerosol Sci. Technol. 2013, 47, 1125–1133. (172) Ramkrishna, D. Analysis of population balance—IV: The precise connection between Monte Carlo simulation and population balances. Chem. Eng. Sci. 1981, 36, 1203–1209. (173) Yu, X.; Hounslow, M. J.; Reynolds, G. K. Accuracy and optimal sampling in Monte Carlo solution of population balance equations. AIChE J. 2015, 61, 2394–2402. (174) Smith, M.; Matsoukas, T. Constant-number Monte Carlo simulation of population balances. Chem. Eng. Sci. 1998, 53, 1777–1786. (175) Lee, K.; Matsoukas, T. Simultaneous coagulation and break-up using constant-N Monte Carlo. Powder Technol. 2000, 110, 82–89. (176) Lin, Y.; Lee, K.; Matsoukas, T. Solution of the population balance equation using constant-number Monte Carlo. Chem. Eng. Sci. 2002, 57, 2241–2252. (177) Kruis, F. E.; Maisels, A.; Fissan, H. Direct simulation Monte Carlo method for particle coagulation and aggregation. AIChE J. 2000, 46, 1735–1742. (178) Maisels, A.; Einar Kruis, F.; Fissan, H. Direct simulation Monte Carlo for simultaneous nucleation, coagulation, and surface growth in dispersed systems. Chem. Eng. Sci. 2004, 59, 2231–2239. (179) Zhao, H.; Zheng, C.; Xu, M. Multi-Monte Carlo method for coagulation and condensation/evaporation in dispersed systems. J. Colloid Interface Sci. 2005, 286, 195– 208. (180) Zhao, H.; Zheng, C.; Xu, M. Multi-Monte Carlo method for particle coagulation: description and validation. Appl. Math. Comput. 2005, 167, 1383–1399. (181) Haibo, Z.; Chuguang, Z.; Minghou, X. Multi-Monte Carlo approach for general dynamic equation considering simultaneous particle coagulation and breakage. Powder Technol. 2005, 154, 164–178. (182) Zhao, H.; Zheng, C. Correcting the multi-Monte Carlo method for particle coagulation. Powder Technol. 2009, 193, 120–123. (183) Liffman, K. A direct simulation Monte-Carlo method for cluster coagulation. J. Comput. Phys. 1992, 100, 116–127. (184) Zhao, H.; Maisels, A.; Matsoukas, T.; Zheng, C. Analysis of four Monte Carlo methods for the solution of population balances in dispersed systems. Powder Technol. 2007, 173, 38–50. (185) Zhao, H.; Zheng, C. A population balance-Monte Carlo method for particle coagulation in spatially inhomogeneous systems. Comput. Fluids 2013, 71, 196–207.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 50 of 50

50

For Table of Contents Use Only Crystal population balance formulation and solution methods: A review Hecham Omar, Sohrab Rohani

TOC graphic:

Synopsis: A review is presented that outlines the different forms the population balance equation can take and reviews the advantages and disadvantages of the different methods (analytical, moment methods, discretization methods and Monte Carlo methods) of solving the population balance equations.


Crystal Population Balance Formulation and Solution Methods: A

Recommend Documents