A rigorous breakage matrix methodology for characterization of multi-particle interactions in dense-phase particle breakage

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

Contents lists available at SciVerse ScienceDirect

Chemical Engineering Research and Design journal homepage: www.elsevier.com/locate/cherd

A rigorous breakage matrix methodology for characterization of multi-particle interactions in dense-phase particle breakage E. Bilgili ∗ , M. Capece Otto H. York Department of Chemical, Biological, and Pharmaceutical Engineering, New Jersey Institute of Technology, 161 Warren St., 150 Tiernan Hall, Newark, NJ 07102, USA

a b s t r a c t Broadbent and Calcott’s breakage matrix methodology has been used for more than 50 years to model various comminution processes and to determine breakage functions from experimental data. The methodology assumes first-order law of breakage and neglects mechanical multi-particle interactions that are especially prevalent in densephase comminution processes and breakage tests. Although several researchers severely criticized this aspect of the methodology, Baxter et al. (2004, Powder Technol. 143–144:174–178) were the first to modify the methodology toward determining the elements of a feed-dependent breakage matrix. However, no non-linear breakage matrix has ever been constructed from experimental data using the modified approach. In this study, a critical analysis of this modified approach has been performed, and the non-linear breakage matrix was fundamentally derived from a non-linear population balance model. Different approaches were proposed to identify the breakage functions based on the nature of available breakage tests on multiple mono-dispersed feed samples and at least one poly-dispersed sample. Using the derived equations, available experimental data on the breakage of a binary mixture of coarse and fine limestone particles in uniaxial compression test were fitted to quantify the multi-particle interactions. Superior fitting capability of rational approximation to the effectiveness factor was demonstrated. © 2012 The Institution of Chemical Engineers. Published by Elsevier B.V. All rights reserved. Keywords: Particle breakage; Breakage matrix methodology; Multi-particle interactions; Population balance; Parameter estimation; Rational approximation

1.

Introduction

Comminution processes are widely used in the processing of minerals, ceramics, composites, foods, paints and inks, pharmaceuticals, etc. (Prasher, 1987). Industrial comminution processes are energy-intensive, poorly efficient, and costly, thus entailing a quantitative analysis to process design and optimization besides experimentation. Without being too exhaustive, a categorization of the major approaches for modeling particle breakage in comminution processes may be made as follows: • purely data-driven models (e.g., Celep et al., 2011; Pradeep and Pitchumani, 2011) and empirical models such as characteristic particle size-milling time correlations (e.g., Bilgili

∗

et al., 2008; Strazisar and Runovc, 1996; Varinot et al., 1999) and particle size-specific energy correlations (e.g., Austin, 1973; Bilgili et al., 2001; Nomura and Tanaka, 2011) • particle-scale mechanistic models, which explicitly incorporate some material properties to explain particle breakage (e.g., Gahn and Mersmann, 1999a,b; Ghadiri and Zhang, 2002; Vogel and Peukert, 2003) • mechanistic models such as the Discrete Element Method (DEM), Finite Element Method (FEM), or their combination, which account for particle deformation, multiparticle mechanical interactions, and/or collision frequency/energy at the particle ensemble or agglomerate scale (e.g., Ahmadian et al., 2011; Antony and Ghadiri, 2001; Bagherzadeh et al., 2011; Rajamani et al., 2000; Thornton and Liu, 2004; Tsoungui et al., 1999)

Corresponding author. Tel.: +1 973 596 2998. E-mail address: [email protected] (E. Bilgili). Received 9 October 2011; Received in revised form 20 December 2011; Accepted 12 January 2012 0263-8762/$ – see front matter © 2012 The Institution of Chemical Engineers. Published by Elsevier B.V. All rights reserved. doi:10.1016/j.cherd.2012.01.005

1178

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

Nomenclature Notation A first-order (linear) breakage probability constant (–) breakage distribution matrix (–) b cumulative breakage distribution matrix (–) B bij breakage distribution parameter (–) cumulative breakage distribution parameter (–) Bij diag diagonal matrix operator discrete linear DL DNL discrete non-linear F feed (input) mass fraction distribution vector (–) Fc mass fraction of the coarse in the feed with binary mixture (–) mass fraction of the fines in the feed with Ff binary mixture (–) Fi feed mass fraction in size class i (–) functional (effectiveness factor) matrix (–) [ ] effectiveness factor for the coarse during breakc [ ] age of a binary mixture (–) i [ ] diagonal elements of the functional (effectiveness factor) matrix (–) EF experimental mono-dispersed feed matrix with N mono-dispersed feed vectors F(f) (–) EFij elements of experimental mono-dispersed feed matrix (–) EO experimental product matrix containing all of the product vectors P(f) (–) EOij elements of product matrix (–) I identity matrix or unit vector (–) exponent of the first-order breakage probability m (–) M mass fraction distribution vector (–) mass fraction in size class i (–) Mi min minimize command for the optimization n number of breakage events (#) N total number of size classes (#) total number of breakage events (#) ne product (output) mass fraction distribution vecO tor (–) Oc mass fraction of the coarse in the output (–) product mass fraction in size class i (–) Oi P first-order (linear) breakage probability matrix (–) PBM population balance model Pi first-order (linear) breakage probability parameter (–) PSD particle size distribution standard error of the parameter SEP SSR sum of squared residuals rk kth parameter of the polynomial in the numerator of a rational approximation su uth parameter of the polynomial in the denominator of a rational approximation T non-linear breakage (comminution) matrix (–) Tl linear (first-order) breakage matrix (–) W weighting matrix (–) Wiq weighting parameter (–) x particle size (m)

Greek letters ˛ exponent of the weighting parameter in the logistic functional (–) breakage distribution exponent (–) ˇ ␥ breakage distribution exponent (–) degree of non-normalization (–) ı  non-linearity parameter of the logistic functional (m˛ )  breakage distribution constant (–) exponent of the logistic functional (–) ω Subscripts c coarse particles in the binary mixture size class indices i, j f fines (finer particles) in the binary mixture (f) size class index denoting each of the monodispersed feed mass fraction vector F(f) N size class containing the finest particles o belongs to normalizing size of the coarsest particles q size class index for a generic size class intended to represent all classes from 1 to N 1 size class containing the coarsest particles Superscripts (n) nth breakage event 0 initial or feed output or product after a single breakage event 1 −1 matrix inverse

• population balance models (PBMs), which with certain mixing assumptions for the powder/slurry flow (well-mixed, plug-flow, etc.) describe spatial, temporal, or spatiotemporal variation of the particle size distribution during comminution. The readers are referred to Austin (1971) and Bilgili and Scarlett (2005a) for a review of time-continuous PBMs and Bilgili and Capece (2011) for a review of timediscrete PBMs. • models based on computational fluid dynamics (CFD), which describes single and two-phase flows in a comminution equipment (e.g., Blecher et al., 1996; Toneva et al., 2011), and its combination with PBM (e.g., Fan et al., 2004; Rajniak et al., 2008) or DEM (e.g., Teng et al., 2011). Broadbent and Callcott (1956a,b,c, 1957) developed a quantitative understanding of size reduction at the process length scale by proposing a breakage matrix approach, which was associated with the PBM concepts such as selection for breakage probability and breakage distribution (Epstein, 1948). This approach is the major interest of the current study among various modeling approaches. They related product mass fraction distribution O to feed mass fraction distribution F via a linear or first-order breakage (comminution) matrix Tl as follows: O = Tl · F

(1)

O and F are column matrices (vectors) with N × 1 elements denoted by Oi and Fi . Tl has N × N elements denoted by Tlij . In this approach, the particle size domain is discretized into N size classes (usually, but not necessarily, corresponding to sieve cuts), and i and j are the size-class indices that extend

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

1179

from size class 1 containing the coarsest particles to size class N containing the finest particles. Broadbent and Calcott then used the above concepts of breakage probability and breakage distribution to separate Tl into its components as follows: Tl = I − P + b · P

(2)

where I, P, and b are the identity matrix, first-order breakage probability matrix, and breakage distribution matrix, respectively. Note that Tl is lower triangular and b is strictly lower triangular, whereas P = diag[Pi ] and I = diag[1] are diagonal. In index notation, Eq. (1) can be re-written making use of Eq. (2) and the properties of the respective matrices as

Oi =

i  j=1

Tij Fj = (1 − Pi )Fi +

i−1 

bij Pj Fj ,

N≥i≥j≥1

(3)

j=1

With slight modifications in terminology for different applications, the breakage matrix approach described in Eqs. (1)–(3) has been used for over 50 years in the modeling of oncethrough type mills, where residence time of particles in the mill is short, such as roller mills (see, e.g., Austin et al., 1980; Campbell and Webb, 2001; Fistes and Tanovic, 2006), cone mill (Broadbent and Callcott, 1956a), shredder-granulators (Jekel et al., 2007), and single-event fluid energy mills (Teng et al., 2010). In these systems, assumption of a single breakage event seems to be a decent approximation. On the other hand, some researchers (Abou-Chakra et al., 2004; Baxter et al., 2004) further maintain that Eq. (1) even treats breakage in a process with relatively long residence time such as pneumatic conveying, where multiple breakage events are highly likely, in some time-averaged manner. Austin (1971) and Bilgili and Capece (2011) argue against this and prefer a time-continuous population balance model (PBM) for such processes with long residence time (the latter also suggest consideration of multiple breakage events in the context of a time-discrete PBM). Broadbent and Calcott’s breakage matrix approach (1956a) is based on the assumption of first-order law of breakage, which states that the extent of breakage for a given size class is first-order in Fi or F, thus rendering Eqs. (1) and (3) linear. Interestingly, Broadbent and Callcott (1956b) were the first to note “. . .breakage in an assembly is partly due to interaction between particles in the assembly”, while Eqs. (1) and (3) do not reflect this. After a careful review of the wet milling literature, Meloy and Williams (1992) conclude “. . .the assumption of linearity, used in solving the PBM, is certainly wrong in many grinding mills and wrong for most commercial wet grinding circuits”. From a fundamental particle breakage perspective, we assert that Eqs. (1) and (3) inherently assume the independence of particle breakage from the particle size distribution of the surrounding particle population, thus neglecting the impact of mechanical multi-particle interactions on the breakage probability. On the other hand, ample experimental data to this date have suggested that mechanical multi-particle interactions among particles of different sizes at the particle ensemble scale lead to either a population-dependent breakage probability in particle bed compressions tests and high-pressure roller mills (Abu-Nahar et al., 2006; Abu-Nahar and Tuzun, 2007; Gutsche and Fuerstenau, 1999; Hoffmann and Schonert, 1971) or population-dependent specific breakage rates in ball mills (Austin et al., 1990; Meloy and Williams, 1992).

Fig. 1 – A schematic for uniaxial powder bed compression inside a piston–cylinder apparatus: mono-dispersed coarse particles (left) and binary mixture of coarse particles and fines (right). Schematic is not to scale. There is a great need for developing a fundamental understanding of the impact of multi-particle interactions in dense-phase comminution systems where particle concentration is high and mechanical multi-particle interactions (beyond binary interactions) occur through either collisions or enduring soft contacts. Readers are referred to Bilgili et al. (2006) and Bilgili and Capece (2011) for a review of non-first-order effects in comminution processes and their physical origin. As an example of such multi-particle interactions at the particle ensemble scale, let us consider the effect of fines on the breakage probability of coarse particles in particle bed compression tests (see Fig. 1). Gutsche and Fuerstenau (1999) demonstrated the retardation of coarse particle breakage due to presence of fines and attributed it to the redistribution of force flux. DEM (Discrete Element Method) simulations of compacted granular materials (see, e.g., Antony and Ghadiri, 2001; Tsoungui et al., 1999) have predicted the redistribution of force flux on coarse particle surfaces due to increased number of contacts upon addition or presence of fines. The retardation of coarse particle breakage originating from coarse–fines interactions cannot be explained by Eq. (1), which seriously restricts the capability of the breakage matrix approach in characterizing and predicting particle breakage in dense-phase systems. Recognizing the importance and complexity of multiparticle (intra-mixture) interactions and their impact on the population-dependent breakage probability, Baxter et al. (2004) proposed a modified breakage matrix approach. Although they regarded the breakage matrix T as dependent on the feed mass fraction distribution F implicitly for a general poly-dispersed feed, they still used a mathematical form similar to Eq. (1). To be precise, their work correctly implies O = T(F) · F

(4)

unlike the first-order breakage matrix Tl in Eq. (1), which is a constant matrix. Assuming availability of N breakage tests on all mono-dispersed feeds and a single breakage test on a poly-dispersed feed, Baxter et al. (2004) suggested a method for determining Tij for the poly-dispersed feed. Interestingly, no non-linear breakage matrix T has ever been constructed

1180

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

from experimental data using this approach. In this study, a critical analysis of the Baxter et al.’s approach is performed first with a goal of revealing its inadequacies. Then, the nonlinear breakage matrix T is derived from a recently proposed discrete non-linear PBM (DNL-PBM). In view of breakage tests on N mono-dispersed feeds and one poly-dispersed feed at a minimum, various parameter estimation approaches are proposed to identify the breakage parameters explicitly. Using the derived equations, experimental data on the breakage of a binary mixture of coarse and fine limestone in uniaxial compression test (Gutsche and Fuerstenau, 1999) were fitted to quantify the multi-particle interactions through an effectiveness factor. We found a superior fitting capability of rational approximation to the effectiveness factor, which can allow experimenters to quantify multi-particle interactions in an unbiased manner. This paper has demonstrated for the first time how to describe the multi-particle interactions and resulting non-first-order breakage explicitly within the context of the breakage matrix methodology. The rigorous methodology presented in this paper is hoped to revitalize the breakage matrix approach for studying dense-phase comminution systems.

2.

Background

Baxter et al. (2004) proposed a purely empirical methodology with the hope of determining the elements of a breakage matrix that depends on feed particle size distribution (PSD). From N breakage tests on mono-dispersed feeds covering N size classes, one can determine the elements of Tl in the standard breakage matrix approach represented by Eq. (1) (see also Campbell and Webb, 2001). It is well-known that breakage tests on mono-dispersed samples (narrow sieve cuts) allow minimization of non-first-order effects and ensure validity of first-order law of breakage. In addition, using an additional breakage test on a poly-dispersed feed for which non-firstorder effects are not necessarily negligible, one can determine the elements of T. Interestingly, Tuzun and colleagues (AbuNahar et al., 2006; Abu-Nahar and Tuzun, 2007; Baxter et al., 2004) have not constructed any T from experimental data using this approach. Instead, they took the difference between O predicted by Eq. (1) i.e., Tl ·F, and experimentally observed O for a given poly-dispersed (e.g., binary mixture) F as a gauge of multi-particle interactions (intra-mixture interactions in their terminology) (Abu-Nahar et al., 2006; Abu-Nahar and Tuzun, 2007). The interaction was quantified by measuring the coefficient of distance of the predicted O and observed O based on the so-called Canberra Distance approach. Although the use of a single scalar (the coefficient of distance) is certainly useful for quantifying the impact of multi-particle interactions in gross terms, most details of the interactions are left hidden in the elements of T. For example, Abu-Nahar et al. (2006) studied the breakage of mono-dispersed and binary mixtures of various narrow sieve cuts (coarse, medium, and fine) using a compression tester. Upon an increase in either the fines fraction or the size ratio of the interacting sieve fractions, the Canberra coefficient of distance increased, thus suggesting stronger multi-particle interactions. However, the coefficient does not even appear to determine the nature of the non-first-order effect directly (retardation or enhancement or both retardation and enhancement occurring for different size classes).

One major drawback of Baxter et al.’s approach (2004) is that it treats particle breakage inside processing or testing equipment as a “black-box” and that the elements of T are completely empirical. Moreover, all non-first-order effects are coupled with the first-order breakage and breakage distribution, as quantified by the elements of T in a lumped and convoluted manner. The linear breakage matrix Tl is related to P and b through Eq. (2). Austin (1971) reminded this relationship in the context of the standard discrete linear population balance model (DL-PBM), which may be considered as a “gray-box” model. On the other hand, no equivalent explicit relation of T to PBM exists; therefore, non-first-order effects emanating from multi-particle interactions are not explicitly and accurately quantified. Another major drawback of the Baxter et al.’s approach (2004) is that it is purely descriptive, but has no predictive capability. According to this approach, T depends on F implicitly and somewhat arbitrarily. Theoretically, an infinite number of T matrices exist for infinitely possible F vectors. Since the F-dependence of T has not been determined directly and explicitly from Eq. (4), Baxter’s et al.’s approach presents a philosophical modeling dilemma: one needs to know T to predict O; but, T has to be constructed for any given F in the absence of a functional relationship. This deficiency alone casts serious doubt about the usefulness of the breakage matrix approach for dense-phase comminution systems. In fact, a decade earlier, Meloy and Williams (1992) indicated that no general mill (breakage) matrices had appeared or were likely to appear because of the nonlinearity of the selectivity (breakage probability) matrix. In Section 3, an attempt is made to mitigate these drawbacks of the breakage matrix approach by using the recently formulated discrete non-linear PBM (DNL-PBM) (Bilgili and Capece, 2011).

3.

Theory

3.1. Derivation of a non-linear breakage matrix based on DNL-PBM So as to account for multi-particle interactions in retentiontype mills explicitly, Bilgili and colleagues (Bilgili and Scarlett, 2005a; Bilgili et al., 2006) formulated a non-linear (nonfirst-order) PBM framework for rate-based comminution processes. They decomposed the specific breakage rate function into a population-dependent functional [ ] and a first-order, size-dependent breakage rate function. Numerical simulations using this framework predicted most experimentally observed complex and non-intuitive breakage behavior such as cushioning effect, crossing-effect, and acceleration effect in batch (Bilgili and Scarlett, 2005a; Bilgili et al., 2006) and continuous dry milling operations (Bilgili and Scarlett, 2005b, c). Capece et al. (2011a) demonstrated the superiority of the non-linear model over its linear counterpart through a statistical analysis of these models’ fit to the data on batch milling of quartz (Austin et al., 1990). To account for multi-particle interactions in both eventbased and rate-based comminution processes, Bilgili and Capece (2011) have recently derived the following discrete non-linear PBM (DNL-PBM) from the time-continuous

1181

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

non-linear PBM (TCNL-PBM): M(n) = {I − (I − b) · P · [W · M(n−1) ]} · M(n−1) ≡ T(M(n−1) ) · M(n−1)

with M(n=0) = M(0) n = 1, 2, 3, . . . , ne (5)

where M(n) is the mass fraction distribution vector after nth breakage event and [ ] is the population dependent functional matrix. Note that [ ] = diag[i ] is a diagonal matrix. W is a weighting matrix expressing the integral impact of all size class i–q interactions between particles of generic size xq and particles of size xi . In [ ], [ ] stands for the summed (integral) effect of the multi-particle interactions. Motivated by the concept of falsified kinetics in heterogeneous chemical reaction engineering, Capece et al. (2011b) has recently coined the term “effectiveness factor” for the functional [ ], which precisely quantifies the retardation and/or enhancement of particle breakage due to multi-particle interactions. Both terms are interchangeably used throughout the paper. Without [ ], the breakage kinetics in dense-phase systems may be falsified. In Eq. (5), T is the comminution matrix that depends on M(n−1) and the feed PSD as described by M(0) . The comminution process is assumed to take place through ne consecutive breakage events. The term “event” is broadly used here to connote passage, cycle, step, stage, and section that have been alternatively used in a multitude of different comminution applications. Due to the positive-valued nature of the breakage probability and the validity of first-order law of breakage in the limit of mono-dispersed particles, constraints such as 0 ≤ Pi ≤ 1, i [ ] ≥ 0, and i [ ] → 1 as ∀Mq (q = / i, n) → 0 apply to the diagonal elements of P and [ ]. Bilgili and Capece (2011) discussed the accuracy, stability, and applicability of the DNL-PBM. Numerical simulations (Bilgili and Capece, 2011) demonstrated the capability of DNL-PBM to predict experimentally observed non-first-order particle breakage in retention type mills as well as in powder bed compression tests. However, estimation of the DNL-PBM parameters has not been performed; the inverse problem has not been solved yet. In the present paper, the inverse problem is treated within the context of the breakage matrix methodology. Considering a single breakage event ne = 1 and assigning F ≡ M(0) and O ≡ M(1) , we deduce from Eq. (5) O = {I − (I − b) · P · [W · F]} · F ≡ T(F) · F

(6)

Making use of Eqs. (1) and (2), we relate T in Eq. (6) to Tl in Eq. (2) as T = I −  [ ] + Tl ·  [ ]

(7)

Eq. (6) defines the non-linear breakage matrix T in terms of mathematically and physically well-defined breakage functions, whereas Eq. (4) lumps all breakage functions including the multi-particle interactions (non-linearity) in an unknown manner within the elements of T. In Eq. (6), the breakage probability is now given by P·[ ]. Note that Eq. (6) recovers the traditional linear breakage matrix equation (see Eqs. (1) and (2)) in the limit of weak multi-particle interactions, i.e., T → Tl as [ ] → I. Any deviation from the first-order law of breakage emanating from multi-particle interactions is captured by the functional [ ] explicitly. Unlike the Canberra coefficient of distance (Abou-Nahar et al., 2006), [ ] is a diagonal matrix whose

elements i [ ] present a wealth of information: e.g., 0 < i [ ] < 1, i [ ] = 1, and i [ ] > 1 correspond to retardation of breakage, firstorder breakage, and enhancement of breakage, respectively. Through Eq. (6), we have rectified the first major flaw of Baxter et al.’s approach (2004) to breakage matrix methodology (Section 2). Having explicitly and quantitatively defined the multi-particle interactions via the functional or effectiveness factor i [ ], we have related the non-linear breakage matrix T(F) to the breakage functions P, b, and [ ]. Toward fixing the other major flaw, i.e., lack of predictive capability, we need to determine the breakage functions P, b, and [ ] explicitly first. Note the diagonal nature of P and [ ], lower triangular nature of T (with Tij = 0 for i < j), and strictly lower triangular nature of b (with bij = 0 for i < j due to absence of size growth and negligible intra-class breakage). In index notation, Eq. (6) relates P, b, and [ ] to T as

Tij =

⎧ ⎪ ⎨0

if i < j

1 − P  []

i i ⎪ ⎩b P  [] ij j j

if i = j

(8)

if i > j

from which we find the following useful relations for the breakage parameters: Pi i [ ] = 1 − Tii bij =

Tij

(9)

1 − Tjj

3.2. Development of a direct determination method for the breakage parameters Direct identification of all breakage parameters entails a minimum of N + 1 breakage tests: N breakage tests on mono-dispersed feeds and at least an additional test on a polydispersed feed. In a breakage test, a sample (feed) of material with mono-dispersed particles or poly-dispersed particles is subjected to mechanical stresses for a relatively short duration as in a compression tester (Abu-Nahar et al., 2006; Abu-Nahar and Tuzun, 2007; Gutsche and Fuerstenau, 1999; Hoffmann and Schonert, 1971) or in a process equipment such as roller mills (see, e.g., Austin et al., 1980; Campbell and Webb, 2001; Fistes and Tanovic, 2006), cone mill (Broadbent and Callcott, 1956a), shredder-granulators (Jekel et al., 2007), and singleevent fluid energy mills (Teng et al., 2010). Narrow size cuts may be considered mono-dispersed as an approximation since true mono-dispersed feeds are impossible to obtain. First, one performs a breakage test on each and every narrow sieve cut for size classes i = 1, 2, 3, . . ., N (total of N tests). Therefore, there is a feed vector F(f) corresponding to each size class, which is equivalent to a unit vector I(f) due to mono-dispersed nature of the feed. Upon breakage, a resultant product mass fraction vector O(f) is measured by size (sieve) analysis. By placing these column vectors side by side for f = 1, 2, 3, . . ., N, one generates the following experimental mono-dispersed feed matrix EF and corresponding product matrix EO: EF = [ F(1) = [ I(1) EO = [ O(1)

F(2) I(2) O(2)

≡ [ T(I(1) )I(1)

F(3) I(3)

··· ···

O(3) T(I2 )I2

F(N) ] I(N) ] = I

···

(10)

O(N) ]

T(I(3) )I(3)

···

T(I(N) )I(N) ]

Note that both EF and EO have N × N elements because F(f) and O(f) have N × 1 elements. In view of the restriction [ ] → I

1182

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

for mono-dispersed feeds, Eq. (6) results in T(I) = Tl = EO given the information in Eq. (10). Since EO is known experimentally, Eq. (9) with i [ ] = 1 yields the elements of P and b as follows:

If one insists on using Baxter et al.’s approach for determining T directly, which is undesirable, the theoretical analysis performed here allows back-calculation of [ ] from T. Using standard matrix algebra, Eq. (7) is solved for [ ] as

Pi = 1 − EOii bij =

(11)

EOij 1 − EOjj

Having determined the elements of P, b, and Tl from N breakage tests on mono-dispersed samples, we are left with the task of finding the elements of [ ], i.e., i [ ], which contains all information about the multi-particle interactions. If the purpose is to identify all elements of [ ], to a minimum, one needs to perform a single breakage test on a poly-dispersed feed. However, to identify all i–q size class interactions adequately, to make the breakage matrix approach more predictive, and to increase certainty in the estimated parameters, multiple breakage tests on binary–ternary mixtures and poly-dispersed feeds may be needed. Let us assume that at least a breakage test on a poly-dispersed feed with Fi has been performed and the corresponding product mass frac/ I(f) , one tion distribution Oi has been measured. Due to F = cannot use Eq. (10). Instead, starting with Eq. (6), and making use of Eq. (8), we arrive at

Oi =

i−1 

Tij Fj + Tii Fi

and

Tii = 1 − Pi i [ ]

j=1

=

Oi −

i−1

b P j=1 ij j j

[ ] Fj

(12)

Fi

As we have already determined Pi and bij from Eq. (11), the elements of the functional are found from Eq. (12) upon a straightforward algebraic manipulation as

i [ ] =

1 Pi



1−

Oi Fi



1  bij Pj j [ ] Fj Pi Fi i−1

+

(13)

j=1

Equation (13) presents a recursive relation that allows determination of 1 [ ] first for i = 1, then of 2 [ ] for i = 2 with 1 [ ] now known, and so on. Appendix presents the expressions for i [ ] explicitly while we recognize that i [ ] will be most easily generated by performing a nested calculation (a “do-loop”) for given Fi –Oi pair(s) in a computer code. Eq. (13) suggests that i [ ] cannot be determined at all if Fi of the respective size class is zero (no particles in the feed). Hence, the feed mass fraction distribution Fi above should cover all particle size range of interest. For example, the feed size distribution may be close to a uniform distribution; but, this may not be practical experimentally. In that case, one may prefer to use multiple feed size distributions that can collectively cover the entire particle size region of practical interest. Being very flexible, the “least-squares minimization” approach in Section 3.3, along with assumed functional forms of Pi , bij , and i [ ], can handle any type of feed size distributions. Having determined the breakage parameters Pi , bij , and i [ ] from the above N + 1 (or more) experiments, we can construct the non-linear breakage matrix T using Eq. (8). Baxter et al.’s approach (2004) directly aims at finding the elements of T from the same N + 1 experiments while bypassing the fundamental breakage functions Pi , bij , and i [ ] and information about the true non-linear interactions contained within i [ ].

 [ ] = (Tl − I)

−1

· (T − I)

(14)

3.3. Determination of the breakage parameters via least-squares minimization In theory, Eqs. (11)–(13) appear to be sufficient for determining all non-linear breakage parameters from N + 1 breakage tests because the problem is well-posed (equal number of equations and unknowns). There are, however, two major concerns that must be addressed. First, there are too many parameters to be determined: (N2 + 3N)/2. For experiments with 10 sieve size cuts, 65 parameters must be determined! Second, even in the absence of systematic experimental errors, random errors inevitably corrupt experimental mass fraction distributions Fi and Oi . Coupled with the unrealistically large number of parameters to be determined, it is difficult to determine the parameters accurately and precisely. Interestingly, this issue was not considered at all by Baxter et al. (2004) among others (e.g., Teng et al., 2010). Eqs. (13) and (14) can be subject to serious errors as Pi approaches to small values (small particles). These issues are circumvented by introduction of some size-dependent forms of Pi , bij , and i [ ] and identification of the reduced parameter set in the least-squares sense via optimization-based methods.

Possible forms of Pi , bij , and i [ ]

3.3.1.

From a parameter estimation (inverse problem) perspective, one can reduce the model in Eqs. (6) and (12) by introducing some particle size dependent forms of Pi , bij , and i [ ] based on ample experimental evidence from the literature, heuristics, and intuition. In principle, one can enhance the predictive capability of the PBMs by relating the breakage functions to the processing parameters and fundamental material properties, which entails additional experimental and theoretical investigations. Some particle-scale mechanistic models explicitly incorporate fundamental material properties to explain extent of particle breakage (e.g., Gahn and Mersmann, 1999a, 1999b; Ghadiri and Zhang, 2002; Vogel and Peukert, 2003). However, they are restricted under certain theoretical assumptions and have not found widespread use in PBMs. Therefore, a widely used power-law function for Pi (e.g., Jekel et al., 2007; Otwinowski et al., 2007) and a non-normalized cumulative breakage distribution function Bij (e.g., Austin and Luckie, 1972; Austin et al., 1984)

 x m

Pi = A

i

(15)

x0

Bii = 1,

Bij = 

+ 1−

 x −ı x j

x1

x

j

x1

−ı

i−1



xj

xi−1 xj

ˇ i>j

(16)

are considered. Here, A corresponds to the first-order breakage probability of the coarsest particles on taking the normalizing size as x0 = x1 and the exponent m describes the shape or size-dependence of first-order breakage probability Pi . It

1183

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

is reminded that bij is determined from bij = Bij − Bi+1j . If the non-normalization parameter ı is set as zero, the well-known normalized form of Bij will result (Deniz, 2004; Klimpel and Austin, 1977) as a limiting case of the more general function in Eq. (16). The parameters , , and ˇ describe the cumulative mass fraction of the daughter particles in all size classes i > j that formed due to the breakage of particles in size class j only. One should note that the first-order breakage probability parameters Pi (or specific breakage rate parameters in general) are more sensitive to the specific processing/testing conditions than the breakage distribution parameters bij , which are usually assumed invariant to the processing/testing conditions (e.g., Baxter et al., 2004; Berthiaux, 2000; Meloy and Williams, 1992). Therefore, in practice, one can determine the functional dependence of Pi parameters (e.g., A and m) on the processing/testing parameters such as shear rate, compaction stress, etc. depending on the specific processing or testing equipment. Similarly, parameters of i [ ] can be studied empirically under different processing/testing conditions. Choice of the form of i [ ] cannot be completely arbitrary in view of the overwhelming evidence in literature that breakage of mono-sized feeds obeys first-order law of breakage (see, e.g., Prasher, 1987). Therefore, the functional form must / i) → 0. Within this satisfy the restriction: i [ ] → 1 as ∀Fq (q = constraint, theoretically, one can come up with any functional form that can satisfy some experimental data. The role and various forms of i [ ] in non-linear breakage kinetics have been discussed in great length and depth within the context of timecontinuous non-linear PBMs. Readers are referred to earlier work (Bilgili and Scarlett, 2005a; Bilgili et al., 2006) for various forms of i [ ] with varying complexity. Based on experience, intuition, and detailed knowledge of the experimental system, one may suggest various specific, explicit forms of i [ ] within the context of the solution to the inverse problem, where they are discriminated using statistical methods, as illustrated in Capece et al. (2011a). In general, the form of i [ ] is not known a priori; therefore, a general and unbiased approach is desirable for representing i [ ]. Referring to the breakage matrix, Scarlett (2002) asserted “when there is particle interaction, then the terms of the matrix depend upon the ambient size distribution and the equation must include non-linear terms”. In his perspective, non-linearity should be introduced by additional terms to the linear breakage matrix, thus suggesting either a polynomial approximation or a rational approximation to i [ ]. Motivated by his perspective, unlike in all our prior work in the field of non-linear particle breakage (e.g., Bilgili and Scarlett, 2005a; Bilgili et al., 2006; Capece et al., 2011a, 2011b), we seek here a general representation of the functional. Recognizing the various superior features of the rational approximation (rational function) over the polynomial approximation (Hanna and Sandall, 1995), a general form of i [ ] is proposed in this paper as follows:

⎡ i = i ⎣

N 

q=1 or i

⎤ Wiq Fq ⎦ ≡ i [ ] =

2

k

1 + r1 [ ] + r2 [ ] + · · · + rk [ ]

1 + s1 [ ] + s2 [ ]2 + · · · + su [ ]u (17)

where [ ] stands for the summed interaction, rk and su are the unknown parameters (total of k + u) to be determined. It has been found empirically that the best rational approximations usually correspond to k = u or k differs from u by no more than 1. In light of this, one can consider several parameter

estimation cases with a total number of parameters in the range 2–8 with the heuristic |k − u| ≤ 1 and then assess the model’s goodness of fit and certainty of rk and su to decide the optimal k and u (Hanna and Sandall, 1995). The weighting parameter Wiq can be taken as a function of the sizes of the interacting particles, i.e., xq and xi , such as (1 − xq /xi ) (Bilgili et al., 2006) or more preferably surface area differences (1/xq − 1/xi ) (Bilgili and Capece, 2011), both of which automat/ i) → 0 for Eq. ically satisfy the requirement i [ ] → 1 as ∀Fq (q = (17). Having introduced a specific form, especially with a rational approximation, for i [ ], we attempt to rectify the second major flaw to the Baxter et al.’s approach (2004): an explicit mathematical form for multi-particle interactions i [ ] allowing predictive capability. Once determined from experiments, i [ ] will be incorporated into Eqs. (6) and (8), which in turn can be used predictively to determine the product mass distribution for any given feed mass fraction distributions.

3.3.2. Estimation of Pi , Bij , and i [ ] via least-squares minimization Upon introduction of Eqs. (15)–(17), one expects to fit 8–14 parameters depending on the rational approximation chosen for i [ ]. This corresponds to a significant reduction from the 65 parameters that must have been fit originally for 10 size classes. Therefore, the certainty (precision or reliability) of the estimated parameters is expected to increase significantly. Another advantage of the functional forms introduced in Section 3.3.1 is that the number of parameters will most likely remain the same even if the number of size classes increases. Most importantly, this study proposes a rational approximation for i [ ], which should reveal more unbiased information about the multi-particle interactions. Unlike in Section 3.2, we here seek the solution to the inverse problem, i.e., parameter estimation of the PBM, through minimization of sum of squared residuals between model and experimental data, rather than a direct solution. The unknown parameters of Pi , bij , and i [ ] can be determined either sequentially or simultaneously by solving different types of minimization problems. In the first parameter estimation method, the following set of three least-squares minimization problems can be solved sequentially to find the optimal parameters (first Pi and bij , then i [ ]) resulting in minimum sum of squared residuals (SSR):

min A,m

N 

2

(Pi − 1+EOii ) ,

i=1

and min rk ,su

N  i=1

⎡

min

N

N  

,ˇ,,ı

i=1 j=1 i−1 

⎣Fi −Oi −Pi Fi i [ ] +

Bij −Bi+1j −

EOij

2 ,

1 − EOjj

⎤2 (Bij − Bi+1j )Pj j [ ] Fj ⎦

j=1

(18) The residuals were defined using Eqs. (11) and (12) along with bij = Bij − Bi+1j after some rearrangement. Again, EOij and EOii are experimental product mass distribution data that have been gathered from N breakage tests on mono-dispersed feeds (see Eq. (10)). The breakage test data on a single (or multiple if needed) poly-dispersed feed with Fi and resultant product with Oi are also used. The expressions for Pi , Bij , and i [ ] in Eqs. (15)–(17) are inserted into Eq. (18) before the minimization. The sum of squared residuals (SSR) and standard errors of the parameters (SEP) can be calculated with different combinations of rational approximation parameters k and u (see,

1184

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

e.g., Hanna and Sandall, 1995). An optimal k and u set is taken as the one that leads to the lowest SSR and SEP values. Using this method, one may want to perform multiple breakage tests on poly-dispersed feeds if SEP values of parameters are not satisfactorily low. In an alternative method with a single poly-dispersed feed only, Eq. (12) is utilized in

min

A,m,,ˇ,,ı,rk ,su

N  i=1

⎡ ⎣Fi − Oi − Pi Fi i [ ] +

i−1 

⎤2 (Bij − Bi+1j )Pj j [ ] Fj ⎦

j=1

(19) to estimate all parameters simultaneously. This second method avoids laborious and time-consuming N breakage tests on mono-dispersed feeds; however, parameter interactions are highly likely and precisions (certainty) of the parameters can be low. It is certainly prudent to perform multiple breakage tests (less than N tests) on different poly-dispersed feeds generating multiple Fi –Oi data pairs so that reasonably small SEP of the parameters can be obtained. One can use laser diffraction measurement in the determination of Fi –Oi data pairs (>100 size classes) so that errors associated with crude size discretization (due to sieving) can be minimized. A potential issue in parameter estimation using an optimizer is structured parameter non-identifiability (the issue of uniqueness of the optimal solution), which will be briefly addressed here. The readers are referred to Ramachandran and Barton (2010) for a detailed discussion of parameter identifiability. They studied the parameter estimation within the framework of a multi-dimensional population balance model for breakage, consolidation, aggregation, and nucleation and investigated the local identifiability of the associated model parameters. Simulated (computer-generated or artificial) particle size distribution data were obtained with some assumed “true parameter values” (optimal set) in the respective PBM. Following that, various optimization methods were tested by perturbing the parameter values within ±20% of the optimal values, which are known a priori. Using a Fisher Information Matrix (FIM) and its eigenvalues, which were calculated based on the converged parameter set, a breakage–consolidation model was shown to have no parameter identifiability issue (non-singularity of the FIM or no zero eigenvalue). A similar approach was used for pure breakage problem in the context of time-continuous PBMs (Capece et al., 2011a) and time-discrete PBMs (Capece et al., 2011b). Note that the latter type of PBMs was used to derive the PBM in the current work (see Eq. (5)). After the perturbation of the parameters, in the absence and presence of random errors associated with the particle size distribution data, the optimizer converged to the optimal set very closely (case with random errors) or identically (case without random errors). Since the convergence was attained, Fisher information matrix was not directly calculated in Capece et al. (2011a,b). On the other hand, the SEP values were obtained from the parameter covariance matrix, which entailed inversion of the Fisher Information Matrix (FIM). The FIM was not singular for the converged parameter sets in these studies, suggesting that the system was locally identifiable. It appears that none of the above studies indicates parameter identifiability as a serious problem for the PBMs for size reduction. On the other hand, Capece et al. (2011a, 2011b) illustrated the dramatic impact of experimental errors on the SEP, parameter precision, and statistical significance;

hence, they must be determined in a parameter estimation study with genuine experimental data toward discriminating different models.

4.

Results and discussion

Although voluminous experimental data exist on particle breakage in the open literature, non-linear particle breakage phenomenon has not received the attention it deserves, except in few fundamental experimental studies conducted in particle bed compression testers (Abu-Nahar et al., 2006; Abu-Nahar and Tuzun, 2007; Hoffmann and Schonert, 1971; Gutsche and Fuerstenau, 1999) and ball mills (Austin et al., 1990; Meloy and Williams, 1992). One obvious reason for this is that most experimenters purposefully avoid non-linear breakage phenomenon by either performing breakage experiments with multiple mono-sized feeds or performing size reduction for a relatively short period of time, thus reducing the accumulation of finer particles. In fact, some mills are designed to minimize multi-particle interactions. For example, multi-particle interactions were minimized under dilute-phase conditions in fluid energy mills. In general, the lack of a general theoretical framework to treat the nonlinearity mathematically has been an impediment to the study of non-linear breakage. After all, the prevailing 60-year old traditional linear PBM is based on the assumption that the breakage probability of particles is independent of the particle population. Some researchers clearly indicated deviations from the linear particle breakage, while others neglected such deviations (refer to a review of non-linear breakage by Bilgili and Scarlett, 2005a; Bilgili et al., 2006). Mechanical multi-particle interactions cannot be avoided in dense-phase comminution systems and thus non-linear particle breakage is prevalent in such systems (see, e.g., Bilgili and Capece, 2011, and the references cited therein). Considering that few fundamental studies were conducted on non-linear breakage and even a smaller subset of these data were reported with details in the literature, it is hard to fully validate the methodology presented in the paper. However, we take a first attempt here to demonstrate the use of the proposed methodology for characterizing the impact of multi-particle interactions on the breakage probability and hope to render the breakage matrix methodology more rigorous and predictive. Back-calculation of Pi and bij from N breakage tests on mono-dispersed feeds is a standard practice within the context of Broadbent–Calcott methodology (see, e.g., Austin et al., 1980; Broadbent and Callcott, 1956a; Campbell and Webb, 2001; Fistes and Tanovic, 2006; Jekel et al., 2007; Teng et al., 2010) and will not be elaborated any further in the sequel, while the estimation of i [ ] is the major concern to be addressed in the application below. Out of only few experimental data sets available, the uniaxial compression and breakage of a binary mixture of coarse and fine limestone was considered (see Gutsche and Fuerstenau, 1999 for details) because this study presents striking non-linear breakage behavior in a model dense-phase comminution system. Coarse limestone particles (6–8 mesh, 2360–3350 ␮m) were mixed with finer limestone particles (also called fines) in three 3 sieve cuts: 16–28 mesh (600–1000 ␮m), 100–200 mesh (75–150 ␮m), and <400 mesh (<38 ␮m). The binary mixtures containing a fines mass fraction Ff of 0.2, 0.5, and 0.8 as well as purely coarse particles were separately broken under 140 MPa consolidation pressure in a uniaxial

1185

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

c [ ] =

1 Pc



1−

Oc Fc



=

1 Pc

1−

Oc 1 − Ff

 (20)

Here, Fc and Oc denote the mass fraction of the coarse in the feed and in the output after the breakage, and (1 − Oc /Fc ) is the experimentally measured broken fraction of the coarse. In reaching Eq. (20), it was assumed that all broken particles leave the size class of the coarse and the errors associated with the relatively crude size discretization associated with sieving are negligible. It should be noted that these assumptions are due to the relatively small number of data points available from the sieving analysis performed (Gutsche and Fuerstenau, 1999), and thus they are not necessarily inherent to the inverse solution methodology presented in the paper. In general, Eq. (19) can treat any feed size distribution (e.g., data with even more than 100 size classes). With these assumptions, the firstorder breakage probability parameter for the coarse Pc is equal to the broken fraction of the coarse particles when only the coarse was broken (in the absence of fines). Pc was found to be 0.887 for the 6–8 mesh limestone from the broken fraction data. Equation (20) allows us to directly calculate the effectiveness factor without assuming any specific functional form. However, for the PBM to be predictive, the effectiveness factor must be expressed by a functional form. A cursory look at the breakage of the binary mixture suggests a retardation of the coarse particle breakage by the fines. Gutsche and Fuerstenau (1999) attributed this to the redistribution of force flux during uniaxial compression (see Fig. 1). In view of the experimental data, we suggest the following specific functional form (logistic functional): c [ ] =

1 ˛

1 + [((1/xf ) − (1/xc )) Ff ]

ω

(21)

which satisfies the restriction c [ ] → 1 as xf → xc (single sieve cut). Stronger coarse–fine interactions result as xf (size of the fines) differs more and more from xc (size of the coarse). Here,  is a scaling parameter for the non-linearity.  > 0 leads to breakage retardation of coarser particles in the presence of finer feed particles, whereas  = 0 leads to negligible impact

1.0

Effectiveness Factor for the Coarse (6-8 mesh) Particles (-)

powder compression tester similar to that illustrated in Fig. 1. The broken fraction of the coarse (6–8 mesh) alone and in binary mixtures with the fines was measured. Since only one coarse fraction (6–8 mesh) was considered and particle size distribution after breakage was not reported in Gutsche and Fuerstenau (1999), one cannot determine the first-order breakage probability parameters Pi and effectiveness parameters i [ ] for all particle size classes (except those of the coarse Pc and c [ ]) and the breakage distribution parameters bij . For the breakage matrix methodology to be widely applicable, these interactions must be quantified not only for the coarse, but also for particles in all other sieve cuts. This suggests that dense data sets are needed for full characterization of multiparticle interactions. Nonetheless, this is the first time in the open literature that the effectiveness factor c [ ] and multiparticle (coarse–fines) interactions are quantified for particle breakage taking place in a uniaxial compression tester. In a breakage experiment with a binary mixture of the coarse and fines in the feed, the functional or effectiveness factor for the coarse can be calculated approximately by making use of Eq. (13) and Appendix on taking c ≡ i = 1 and f ≡ q as follows:

Experimental Data: Gutsche & Fuerstenau (1999)

0.8

0.6 Effectiveness Factor Calculated from Eq. (20) Fines Added: 16-28 mesh 100-200 mesh <400 mesh DL-PBM Prediction Logistic Functional, Eq. (21) Rational Approximation, Eq. (22)

0.4

0.2 0.0

0.2

0.4

0.6

0.8

Mass Fraction of Fines (-)

Fig. 2 – Effect of various sieve cuts of fines on the effectiveness factor for the coarse limestone particles c [ ] (6–8 mesh) and fit of the effectiveness factor by various models. of multi-particle interactions. ˛ and ω are the exponents of the weighting function and the logistic function, respectively. Another novelty of the present methodology is that we have introduced a rational approximation to the functional in Eq. (17). Considering the nature of the data (see below) and the small number of data points, we assume the following rational approximation: ˛

c [ ] =

1 + r1 ((1/xf ) − (1/xc )) Ff ˛

1 + s1 ((1/xf ) − (1/xc )) Ff

(22)

This rational approximation, with only three parameters similar to Eq. (21), is the simplest possible approximation because an approximation with either r1 = 0 or s1 = 0 cannot explain the monotone-decreasing, approximately concave-down shape of the c [ ] vs. Ff curves (see below). The functional forms in Eqs. (21) and (22) were fitted to the experimental c [ ] values calculated directly from Eq. (20) using Gutsche and Fuerstenau’s experimental breakage data on limestone (Fig. 2). The non-linear optimizer “fmincon”, which is part of the Matlab v7.9 optimization toolbox and based on the interior-point algorithm (Byrd et al., 1999), was used to determine the parameters of the respective models by minimizing the sum of squared residuals (SSR). Fig. 2 depicts the variation of the functional c [ ], the effectiveness factor for the coarse, as a function of the size and mass fraction of the fines added and fits by the logistic functional in Eq. (21) and the rational approximation in Eq. (22), respectively. Although the traditional discrete linear PBM (DL-PBM), which is the basis of the Broadbent–Calcott breakage matrix approach (Eqs. (1) and (2)), cannot account for multi-particle coarse–fines interaction effects, its prediction with c [ ] = 1 and Pc = 0.887 were also shown in the figures as a reference. Table 1 presents the estimated parameters and associated goodness-of-fit as well as parameter certainty (precision) statistics. The effectiveness factor for the coarse c [ ], which was directly determined from the experimental data using Eq. (20) and fitted using the functional forms, was less than 1 suggesting retardation effect of the fines on the breakage of the coarse (see Fig. 2). c [ ] values were much less than 1 suggesting

1186

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

Table 1 – Parameter estimation and associated statistics for two functional forms of the effectiveness factor for the coarse c [ ]. ˛

c [ ] =

Parameter

˛ ω  (␮m˛ ) r1 (␮m˛ ) s1 (␮m˛ ) SSR F-Ratio p-Value

1 ω 1+[((1/xf )−(1/xc ))˛ Ff ]

Estimated value

SEP

0.102 2.058 2.299 – – 8.572×10−3 126.110 <0.0001

0.229 0.103 0.122 – –

0.083 – – –1.342 –0.652 1.336×10−3 761.092 <0.0001

0.101 – – 0.035 0.096

Broken Fraction of the Coarse (6-8 mesh) Particles (-)

0.8

Experimental Data: Gutsche & Fuerstenau (1999)

0.6 Experimental Data Fines Added: 16-28 mesh 100-200 mesh <400 mesh DL-PBM Prediction Logistic Functional, Eq. (21) Rational Approximation, Eq. (22)

0.2 0.2

˛

1 + s1 ((1/xf ) − (1/xc )) Ff

SEP

1.0

0.0

1 + r1 ((1/xf ) − (1/xc )) Ff

Estimated value

that the DL-PBM (basis of the Broadbent–Calcott breakage matrix approach) becomes completely irrelevant when the multi-particles interactions become significant at high fines’ fraction. The effectiveness factor decreased monotonically with an increase in the fines’ mass fraction. A decrease of the fines’ size led to stronger multi-particle interactions, thus leading to a smaller effectiveness factor (stronger retardation). A visual comparison of the fits by the two functional forms (logistic functional and rational approximation) in Fig. 2 and SSR values in Table 1 suggest the rational (function) approximation is superior to the logistical functional although both functional forms have three parameters. The F-ratio values (much greater than 1) and p-value (less than 0.05) suggest that the fitting with the functional forms were statistically significant overall. Except estimated ˛ values, the standard errors of all parameters were small as compared with the estimated parameter values. It should be noted that ˛ is a parameter that scales the effect of multi-particle interactions through size (or surface area) differences. Unfortunately, the experimental data set was obtained with only three different fines’ sizes; thus, with this small data set, it was not possible to identify the “size or surface area effect” through ˛ precisely. Using the estimated parameter values for c [ ], we predicted the broken fraction of the coarse, which is equal to Pc c [ ] (see Fig. 3

0.4

c [ ] =

0.4

0.6

0.8

Mass Fraction of Fines (-) Fig. 3 – Effect of various sieve cuts of fines on the broken fraction of coarse limestone particles (6–8 mesh) and prediction of the broken fraction by various models.

for both functional forms). The profiles are, by definition of the broken fraction, identical to those of the effectiveness factor in Fig. 2 for the respective functional forms. Again, the rational approximation fits the experimental data much better than the logistic functional.

5.

Summary, conclusions, and outlook

Based on the traditional DL-PBM with first-order law of breakage, Broadbent and Calcott’s breakage matrix approach has been used for more than 50 years for various comminution systems. For dense-phase comminution systems, the firstorder law of breakage does not hold due to the presence of significant mechanical multi-particle interactions, which have been experimentally confirmed by researchers conducting powder bed compression tests on binary particle mixtures and various batch milling studies. To improve the approach, Baxter et al. (2004) introduced an implicit feed PSD dependence to the breakage matrix. The modified approach along with Canberra coefficient of distance has been used to quantify the degree of multi-particle interactions in gross terms. However, no non-linear breakage matrix T has ever been constructed from the experimental data. Moreover, this modified approach has accounted for the multi-particle interactions in a lumped/convoluted manner, while offering no predictive capability. In the present theoretical study, after revealing the inadequacies of the Baxter et al.’s approach (2004), we developed a novel breakage matrix approach and proposed various parameter estimation strategies as tools for experimenters. Deriving elements of a non-linear breakage matrix T from a recently formulated discrete non-linear PBM (DNL-PBM) allowed us to rectify the inadequacies. In the history of breakage matrix approach, this study is expected to be the first of a kind to explicitly describe non-linear effects emanating from multiparticle interactions via a non-linear functional i [ ], also referred to as the effectiveness factor. Unlike all other nonlinear PBM work in the literature, a priori form of i [ ] is not required in the proposed approach due to the rational approximation introduced, which will allow experimenters to obtain a relatively unbiased estimate of i [ ] and the non-linear breakage matrix T. The approach was demonstrated on particle breakage during uniaxial compression of the coarse and fine limestone. The retardation effect of the fines was quantified by determining the effectiveness factor for the coarse explicitly. Besides providing an unbiased estimate of the impact of multi-particle interactions, a rational approximation to the effectiveness factor fitted the experimental data much

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

better than a logistic functional. The present theoretical study is hoped to stimulate further interest in the area of non-linear particle breakage. With the DNL-PBM framework (Bilgili and Capece, 2011) and a rigorous methodology for quantifying the multi-particle interactions presented in this study, the complex non-linear particle breakage in dense-phase comminution systems can and should be studied in detail. This paper also suggests that dense sets of particle breakage data, albeit more laborious to generate, are needed for full characterization of all breakage parameters.

Acknowledgments The authors gratefully acknowledge partial financial support from the National Science Foundation Engineering Research Center for Structured Organic Particulate Systems (NSF ERC for SOPS) through the Grant EEC-0540855. This paper is based on Paper No: 392f delivered at the 2010 AIChE Annual Meeting in Salt Lake City, UT.

Appendix. The expressions for 1 [ ], 2 [ ], 3 [ ], and the general recursive relation are given below: 1 P1 1 2 [ ] = P2 1 3 [ ] = P3 .. . 1 [ ] =

i [ ] =

1 Pi





O1 F1  b P F O2 21 1 1 [ ] 1 1− + F2 P F  b P 2 2 F O3 b32 P2 2 [ ] F2 31 1 1 [ ] 1 1− + + F3 P3 F3 P3 F3



1−

1−

Oi Fi



1  bij Pj j [ ] Fj Pi Fi i−1

+

j=1

References Abou-Chakra, H., Chapelle, P., Tuzun, U., Baxter, J., Patel, M., Cross, M., Bridle, I., Bradley, M.S.A., 2004. Sampling issues in assessing particle degradation in pneumatic conveying systems. Part. Part. Syst. Char. 21, 39–46. Abu-Nahar, A., Baxter, J., Tuzun, U., 2006. Understanding intra-mixture interactions in the breakage of dense particulate mixtures. In: 5th World Cong. on Part. Technol., Orlando, FL, (Paper No: 75e). Abu-Nahar, A., Tuzun, U., 2007. The breakage matrix approach in evaluating inadvertent breakage of dense particulate mixtures. In: AIChE Annu. Meet., Salt Lake City, UT, (Paper No: 18e). Ahmadian, H., Hassanpour, A., Ghadiri, M., 2011. Analysis of granule breakage in a rotary mixing drum: experimental study and distinct element analysis. Powder Technol. 210, 175–180. Antony, S.J., Ghadiri, M., 2001. Size effects in slowly sheared granular media. J. Appl. Mech. ASME 68, 772–775. Austin, L.G., 1971. A review: introduction to the mathematical description of grinding as a rate process. Powder Technol. 5, 1–17. Austin, L.G., 1973. A commentary on the Kick, Bond and Rittinger laws of grinding. Powder Technol. 7, 315–317. Austin, L.G., Luckie, P.T., 1972. Methods for determination of breakage distribution parameters. Powder Technol. 5, 215–222. Austin, L.G., van Orden, D.R., Perez., J.W., 1980. A preliminary analysis of smooth roll crushers. Int. J. Miner. Process. 6, 321–336.

1187

Austin, L.G., Klimpel, R.R., Luckie, P.T., 1984. Methods for direct experimental determination of the breakage functions. In: Process Engineering of Size Reduction. SME-AIME, New York (Chapter 9). Austin, L.G., Yekeler, M., Dumm, T., Hogg, R., 1990. Kinetics and shape factors of ultrafine dry grinding in a laboratory tumbling ball mill. Part. Part. Syst. Char. 7, 242–247. Bagherzadeh, A., Mirghasemi, A.A., Mohammadi, S., 2011. Numerical simulation of particle breakage of angular particles using combined DEM and FEM. Powder Technol. 205, 15–29. Baxter, J., Abu-Nahar, A., Tuzun, U., 2004. The breakage matrix approach to inadvertent particulate degradation: dealing with intra-mixture interactions. Powder Technol. 143–144, 174–178. Berthiaux, H., 2000. Analysis of grinding processes by Markov chains. Chem. Eng. Sci. 55, 4117–4127. Bilgili, E., Scarlett, B., 2005a. Population balance modeling of non-linear effects in milling processes. Powder Technol. 153, 59–71. Bilgili, E., Scarlett, B., 2005b. Numerical simulation of open-circuit continuous mills using a non-linear population balance framework: incorporation of non-first-order effects. Chem. Eng. Technol. 28, 153–159. Bilgili, E., Scarlett, B., 2005c. Nonlinear effects in particulate processes. Nonlin. Anal. 63, e1131–e1141. Bilgili, E., Capece, M., 2011. Quantitative analysis of multi-particle interactions during particle breakage: a discrete non-linear population balance framework. Powder Technol. 213, 162–173. Bilgili, E., Arastoopour, H., Bernstein, B., 2001. Pulverization of rubber granulates using the solid state shear extrusion (SSSE) process: part I. Process concepts and characteristics. Powder Technol. 115, 265–276. Bilgili, E., Arastoopour, H., Bernstein, B., Hamey, R., 2008. Some novel applications of grinding and milling technologies: milling of soft materials and nanomilling. In: Yekeler, M. (Ed.), Fine Particle Technology and Characterization. Research Signpost, Kerala, pp. 41–67. Bilgili, E., Yepes, J., Scarlett, B., 2006. Formulation of a non-linear framework for population balance modeling of batch grinding: beyond first-order kinetics. Chem. Eng. Sci. 61, 33–44. Blecher, L., Kwade, A., Schwedes, J., 1996. Motion and stress intensity of grinding beads in a stirred media mill. Part 1: energy density distribution and motion of single grinding beads. Powder Technol. 86, 59–68. Broadbent, S.R., Callcott, T.G., 1956a. A matrix analysis of processes involving particle assemblies. Philos. Trans. R. Soc. (Lond.) A249, 99–123. Broadbent, S.R., Callcott, T.G., 1956b. Coal breakage processes: I. A new analysis of coal breakage processes. J. Inst. Fuel 29, 524–528. Broadbent, S.R., Callcott, T.G., 1956c. Coal breakage processes: II. A matrix representation of breakage. J. Inst. Fuel 29, 528–539. Broadbent, S.R., Callcott, T.G., 1957. Coal breakage processes III. The analysis of a coal transport system. J. Inst. Fuel 30, 13–17. Byrd, R.H., Hribar, M.E., Nocedal, J., 1999. An interior point algorithm for large scale nonlinear programming. SIAM J. Optim. 9, 877–900. Campbell, G.M., Webb, C., 2001. On predicting roller milling performance. Part I: the breakage equation. Powder Technol. 115, 234–242. Capece, M., Bilgili, E., Dave, R., 2011a. Identification of the breakage rate and distribution parameters in a non-linear population balance model for batch milling. Powder Technol. 208, 195–204. Capece, M., Bilgili, E., Dave, R., 2011b. Emergence of falsified kinetics as a consequence of multi-particle interactions in dense-phase comminution processes. Chem. Eng. Sci. 66, 5672–5683. Celep, O., Aslan, N., Alp, I., Tas¸demir, G., 2011. Optimization of some parameters of stirred mill for ultra-fine grinding of refractory Au/Ag ores. Powder Technol. 208, 121–127.

1188

chemical engineering research and design 9 0 ( 2 0 1 2 ) 1177–1188

Deniz, V., 2004. Relationships between Bond’s grindability (Gbg ) and breakage parameters of grinding kinetic on limestone. Powder Technol. 139, 208–213. Epstein, B., 1948. Logarithmico-normal distribution in breakage of solids. Ind. Eng. Chem. 40, 2289–2291. Fan, R., Marchisio, D.L., Fox, R.O., 2004. Application of the direct quadrature method of moments to polydisperse gas–solid fluidized beds. Powder Technol. 139, 7–20. Fistes, A., Tanovic, G., 2006. Predicting the size and compositional distributions of wheat flour stocks following first break roller milling using the breakage matrix approach. J. Food Eng. 75, 527–534. Gahn, C., Mersmann, A., 1999a. Brittle fracture in crystallization processes. Part A. Attrition and abrasion of brittle solids. Chem. Eng. Sci. 54, 1273–1282. Gahn, C., Mersmann, A., 1999b. Brittle fracture in crystallization processes. Part B. Growth of fragments and scale-up of suspension crystallizers. Chem. Eng. Sci. 54, 1283–1292. Ghadiri, M., Zhang, Z., 2002. Impact attrition of particulate solids. Part 1: a theoretical model of chipping. Chem. Eng. Sci. 57, 3659–3669. Gutsche, O., Fuerstenau, D.W., 1999. Fracture kinetics of particle bed comminution—ramifications for fines production and mill optimization. Powder Technol. 105, 113–118. Hanna, O.T., Sandall, O.C., 1995. Computational Methods in Chemical Engineering. Prentice Hall, Upper Saddle River. Hoffmann, N., Schonert, K., 1971. Bruchanteil von glaskugeln in packungen von fraktionen und binaren mischungen. Aufbereitungs Technik 12, 513–518. Jekel, L.J., Asce, M., Tam, E.K.L., 2007. Plastics waste processing: comminution size distribution and prediction. J. Environ. Eng. 133, 245–254. Klimpel, R.R., Austin, L.G., 1977. The back-calculation of specific rates of breakage and non-normalized breakage distribution parameters from batch grinding data. Int. J. Miner. Process. 4, 7–32. Meloy, T.P., Williams, M.C., 1992. Problems in population balance modeling of wet grinding. Powder Technol. 71, 273–279. Nomura, S., Tanaka, T., 2011. Analysis of energy–size reduction relationships in batch tumbling ball mills. Powder Technol. 208, 610–616.

Otwinowski, H., Zbronski, D., Urbaniak, D., 2007. Experimental identification of entropy model of comminution process. Granul. Matter 9, 377–386. Pradeep, P.R., Pitchumani, B., 2011. Effect of operating variables on the production of nanoparticles by stirred media milling. Asia Pac. J. Chem. Eng. 6, 154–162. Prasher, C.L., 1987. Crushing and Grinding Process Handbook. Wiley, Chichester. Rajamani, R.K., Songfack, P., Mishra, B.K., 2000. Impact energy spectra of tumbling mills. Powder Technol. 108, 116–121. Rajniak, R., Dhanasekharan, K., Sinka, C., MacPhail, N., Chern, R., 2008. Modeling and measurement of granule attrition during pneumatic conveying in a laboratory scale system. Powder Technol. 185, 202–210. Ramachandran, R., Barton, P.I., 2010. Effective parameter estimation within a multi-dimensional population balance model framework. Chem. Eng. Sci. 65, 4884–4893. Scarlett, B., 2002. Particle populations—to balance or not to balance, that is the question! Powder Technol. 125, 1–4. Strazisar, J., Runovc, F., 1996. Kinetics of comminution in microand sub-micrometer ranges. Int. J. Miner. Process. 44, 673–682. Teng, S., Wang, P., Zhu, L., Young, M.-W., Gogos, C.G., 2010. Mathematical modeling of fluid energy milling based on a stochastic approach. Chem. Eng. Sci. 65, 4323–4331. Teng, S., Wang, P., Zhang, Q., Gogos, C.G., 2011. Analysis of fluid energy mill by gas-solid two-phase flow simulation. Powder Technol. 208, 684–693. Thornton, C., Liu, L., 2004. How do agglomerates break? Powder Technol. 143–144, 110–116. Toneva, P., Epple, P., Breuer, M., Peukert, W., Wirth, K., 2011. Grinding in an air classifier mill—part I: characterisation of the one-phase flow. Powder Technol. 211, 19–27. Tsoungui, O., Vallet, D., Charmet, J.C., 1999. Numerical model of crushing of grains inside two-dimensional granular materials. Powder Technol. 105, 190–198. Varinot, C., Berthiaux, H., Dodds, J., 1999. Prediction of the product size distribution in associations of stirred bead mills. Powder Technol. 105, 228–236. Vogel, L., Peukert, W., 2003. Breakage behavior of different materials—construction of a mastercurve for the breakage probability. Powder Technol. 129, 101–110.