Statistical Mechanics of Discrete Multicomponent Fragmentation

We formulate the statistics of the discrete multicomponent fragmentation event using a methodology borrowed from statistical mechanics. We generate the ensemble of all feasible distributions that can be formed when a single integer multicomponent mass is broken into fixed number of fragments and calculate the combinatorial multiplicity of all distributions in the set. We define random fragmentation by the condition that the probability of distribution be proportional to its multiplicity and obtain the partition function and the mean distribution in closed form. We then introduce a functional that biases the probability of distribution to produce in a systematic manner fragment distributions that deviate to any arbitrary degree from the random case. We corroborate the results of the theory by Monte Carlo simulation and demonstrate examples in which components in sieve cuts of the fragment distribution undergo preferential mixing or segregation relative to the parent particle.


I. INTRODUCTION
Objects disintegrate into fragments via impact, detonation, degradation of structure or cleavage of the bonds that hold the structure together. Objects may range from molecules and cells to macroscopic masses to celestial bodies and the physical mechanisms by which fragmentation occurs are as diverse. This diversity of scale and physics is united by a common mathematical description based on the size distribution of the fragments that are produced. Fragmentation has been studied extensively in particular with respect to modeling the evolution of the size distribution [1], a problem that commonly arises in polymeric, colloidal and granular materials as well as in biological populations where fragmentation refers to the splitting of groups of animals, for example [2]. It is a rich mathematical subject that in addition to its practical relevance to scientific and engineering problems exhibits under certain conditions a remarkable transition, shattering, a process akin to a phase transition [3][4][5][6]. Several analytic results have been given in the literature in the context of population balances under various models in the discrete and continuous domain [7][8][9][10][11][12][13]. In these studies fragmentation is nearly always binary and the particle consists of a single component. In one-component fragmentation the distribution of fragments is a function of a single variable, "size," which usually refers to the mass of the fragment. A multidimensional case was considered by Krapivsky and Ben-Naim [14] who considered the fragmentation of objects in d-dimensional space. A treatment of multi-nary fragmentation was given by Simha [8,9] who studied fragmentation in the context of polymer degradation. With the exception of Simha who used combinatorial arguments to obtain the distribution of fragments in linear [7] and branched [9] chains, most studies adopt either binary random aggregation, which produces a uniform distribution of fragments, or * txm11@psu.edu some empirical distribution of fragments that allows arbitrary number of fragments and non-random distribution of sizes (see for example Ziff [13], Austin et al. [15]).
Multicomponent fragmentation with no limitation on the number of fragments is the rule rather than the exception in most problems of practical interest. This cannot be obtained by simple extension of the one-component problem. In addition to the size distribution of the fragments we must also tackle the distribution of components and provide rules (a "model") for apportioning components to the fragments. We must offer a definition of what is meant by "random fragmentation" when both size and composition are distributed and provide the means for constructing models that deviate from the random case to any desired extent.
The purpose of this paper is to formulate the statistics of a single fragmentation event in the discrete domain for arbitrary number of fragments and components and to provide the means of connecting this mathematical formalism to physics. The main idea is this: We start with a multicomponent particle that contains discrete units of multiple components, subject it to one fragmentation event into fixed number of fragments and construct the set of all fragment distributions that can be obtained. We calculate the partition function of this ensemble of random fragments, assign probabilities in proportion to the multiplicity of each distribution and obtain the mean distribution in terms of the partition function. We then introduce a bias functional that biases the distribution away from that of random fragmentation. We present results from Monte Carlo simulations to corroborate the theory and show that components may preferentially mix or unmix in the fragments depending on the choice of the bias functional.

A. One-component random fragmentation
In discrete fragmentation, a particle composed of m integer units breaks up into N fragments, {f 1 , f 2 · · · f N } that satisfy the mass balance conditions We construct the distribution of fragments n = (n 1 , n 2 · · · ), such that n i is the number of fragments that contain i units of mass. We suppose that N is fixed but n is not; that is, if the fragmentation event is repeated with an identical parent particle the distribution of fragments may be different but the total number of fragments is always N . We refer to this process as Nnary fragmentation. All fragment distributions produced by this mechanism satisfy the following two conditions: The first condition states that the number of fragments is N ; the second that their mass is equal to the mass of the parent particle. Conversely, any distribution that satisfies the above two equations is a feasible distribution of fragments by N -nary fragmentation of mass m.
Thus the set E m,N of all distributions that satisfy Eqs.
(2) and (2) forms the ensemble of fragment distributions produced from m. We will call the process random fragmentation if all ordered lists of N fragments produced by the same mass are equally probable. This views the ordered list of fragments, which we call configuration, as the primitive stochastic variable in this problem.

Probability of random fragment distribution
Theorem 1. The probability of distribution n produced by random N -nary fragmentation of mass m is where n! is the multinomial coefficient of n = (n 1 , n 2 · · · ) Proof. First we note that the number of ordered lists (in all permutations) that can be formed by breaking integer m into N fragments is This is the number of ways to partition integer m into N parts and can be shown easily as follows [16]: thread m balls into a string and partition them into N pieces by cutting the string at N − 1 points (Fig. 1). There are m − 1 points where we can cut and must choose N − 1 of them. The number of ways to do this is the binomial factor on the RHS of Eq. (6). If all ordered lists of fragments are equally probable, the probability of ordered list m = (m 1 , m 2 · · · m N ) is There are n! ordered lists with the same distribution of fragments n. Accordingly, the probability of n is where n is the fragment distribution in m. This proves the theorem.
The multinomial factor n! represents the multiplicity of distribution n, namely, the number of configurations (ordered lists of fragments) represented by n. Using ω(n) = n! to notate the multiplicity of distribution, the probability of distribution is expressed as and Ω (1) satisfies n ω(n) = Ω (1) m;N .
The summation over all distributions n ∈ E m,N , namely, all distributions that satisfy Eqs. (2) and (3). Accordingly, Ω (1) m;N is the total multiplicity in the ensemble, equal to the number of ordered configurations of fragments that can be produced from integer mass m breaking into N fragments. We refer to Ω (1) m;N as the partition function of the one-component ensemble of fragments.

Mean fragment distribution
Each distribution n appears in the ensemble of fragment distributions with probability P (n); the mean distribution of fragments is their ensemble average: n = n n P (n) (11) with P (n) from Eq. (4) and with the summation going over all distributions that are produced by N -nary fragmentation of integer mass m. Theorem 2. The mean distribution in N -nary random fragmentation is with k = 1, · · · m − N + 1.
Proof. First we write the probability of distribution in the form with α i > 0 and note that this reverts to Eq. (4) when α i = 1. We will retain the factors α i and will set them equal to 1 at the end. The normalization condition on the probability P (n) reads The derivative of log Ω (1) m;N with respect to n k is where n k is the mean value of n k in the ensemble of fragments. We also have The summand in the expression in the middle amounts to removing one fragment of mass k from all distributions of the ensemble; accordingly, the quantity in braces is the partition function Ω Equation (12) was previously obtained by Montroll and Simha [7] via a combinatorial derivation. Notably it is the same distribution as in discrete binary aggregation (the reverse process of binary fragmentation) with constant kernel, derived by Hendriks et al. [17] who also credit older unpublished work by White. For large m the fragment distribution becomes This is the continuous limit of random fragmentation of a straight line into N segments and is an elementary result that has been derived multiple times in the literature. The earliest report known to us is by Feller [18] who corrected an earlier approximation by Ruark [19].
B. Two-component random fragmentation

Representations of bicomponent populations
a. The bicomponent distribution We now consider a particle that is made of two components.
The set E mA,mB ;N of all distributions that satisfy the above conditions constitutes the set of feasible distributions in bicomponent fragmentation. b. Color-blind distribution The color-blind size distribution or simply "size distribution" n A+B = {n k } is the distribution of the mass of the fragments k = a b regardless of composition: and satisfies the conditions These are the same as Eqs. (2) and (3) in the onecomponent case for a particle with mass m A + m B . Accordingly, the feasible set of the color-blind distribution is E mA+mB ;N .
c. Sieve-cut distribution The sieve-cut distribution n A|k = {n a|k } is the number of fragments with size k that contains a units of component A: n a|k = n a,k−a , (a = 1 · · · k, k = 1 · · · ∞). (24) and satisfies the normalizations Normalizing the sieve-cut distribution by the number of fragments of size k we obtain the compositional distribution of component A within fragments of fixed size k, The compositional distribution is normalized to unity and may be interpreted as a conditional probability: it is the probability to obtain a fragment with a units of A given that the fragment has mass k. The bicomponent distribution may now be expressed in terms of the color-blind distribution n A+B and the compositional distribution c a|k in the form If we divide both sides by the total number of fragments the result reads as a joint probability: The probability n a,k−a /N to obtain a fragment with mass k that contains a units of component A is equal to the probability n k /N to obtain a fragment of mass k times the probability c a|k to obtain a fragment with a units of component A given that the mass of the fragment is k.

The ensemble of random fragment distributions
Random fragmentation in implemented by analogy to the one-component case: we line up the unit masses in the particle into a string and cut at N − 1 places. Every cut is equally probable and so is every permutation in the order of the beads. Theorem 3. The probability of fragment distribution n is where n! is the multinomial coefficient of the bicomponent distribution, and Ω (2) mA,mB ;N is the two-component partition function, given by Proof. First we count the number of ordered sequences of fragments (configurations). Configurations are distinguished by the order the fragments and by the order of components within fragments (Fig. 1). We color the components and place them in a line in some order. There are m A A's and m B B's; the number of permutations is mA+mB mA . Each permutation produces Ω (1) mA+mB ;N configurations with Ω (1) given in Eq. (6). The total number of configurations therefore is their product and proves Eq. (29).
Since all configurations are equally probable the probability of fragment distribution n is proportional to the number of configurations with that distribution. This is equal to the number of permutations in the order of the fragments and in the order of components within fragments. The number of permutations in the order of fragments is given by the multinomial factor of bicomponent distribution in Eq. (28). The number of permutations of components within a fragment that contains a units of A and b units of B is a+b a and since there are n a,b such fragments, the total number of internal permutations in distribution n is The probability of distribution n is equal to the product of Eqs (28) and (30) divided by the total number of configurations, given by Eq. (29): This proves the theorem.
As a corollary we obtain the multiplicity of the bicomponent distribution, Thus we write with Ω (2) = n ω(n). An alternative equation for P (n) is obtained by expressing the bicomponent distribution n in terms of the color-blind distribution n A+B and all sieve-cut distributions n A|k . The result is and is based on the identity and n A|k ! is the multinomial coefficient of the sieve-cut distribution,

Mean fragment distribution
Theorem 4. The mean distribution of fragments in random bicomponent fragmentation is Proof. The proof follows in the steps of the onecomponent problem. We express the multiplicity and the partition function in the form With α a,b = a+b a we recover the result for random fragmentation but for the derivation we treat α a,b as a variable. Following the same procedure that led to Eq. (12) we now obtain To arrive at this result we note that differentiation of the partition function with respect to α a,b by analogy to Eq. (16)   Alternative proof The mean distribution can be obtained by a mean-field argument beginning with the observation that the mean color-blind distribution is the same as in the one-component case. This follows from the fact that the choice of the points at which the string o beads is cut is independent of the compositional makeup of the particle (Fig. 1). Thus n k is given by Eq. (12) with m = m A + m B : We obtain the compositional distribution by the following construction. Imagine that all possible distributions are stacked vertically to form a table so that column 1 contains the first fragment in all distributions, column 2 contains all second fragments and so on. All columns are permutations of each other (this follows from the construction of the fragments illustrated in Fig. 1) and since all permutations are equally likely (this follows from the condition of random fragmentation), all columns have the same fragment and compositional distribution, therefore we only need to consider one of them. The equivalent problem now is this: count the number of ways to select a beads from a pool of m A A's and k − a beads from a pool of m B B's and take its ratio over the total number of ways to pick k beads: The mean distribution then is the product of the size and compositional distributions: or It is straightforward algebra to show that this is equivalent to Eq. (38). For m A , m B ≫ 1, a + b ≪ m A + m B the compositional distribution goes over to the binomial: with As a means of a demonstration we show the results of a Monte Carlo simulation, which are seen to be in excellent agreement with theory. The binomial distribution, also shown for comparison, is only in qualitative agreement because the fragment masses are small and the conditions for asymptotic behavior are not met in this case.

C. Any number of components
Extension to any number of components follows in a straightforward manner from the bicomponent case. Suppose the parent particle consists of K components A, B · · · and contains m A units of A, m B units of B and so on. The distribution of fragments is now expressed by the K-dimensional vector n = {n a,b··· } that gives the number of fragments that contain a units of component where m Z is the mass of component z = a, b · · · in the parent particle. The set all distributions that satisfy the above conditions constitutes the ensemble of all distributions that are produced by the fragmentation of the parent particle into N integer fragments. Random fragmentation is once again implemented as shown in Fig. 1: Given a string of colored beads we cut it at N − 1 random points to produce N fragments. All permutations of the beads are equally probable. Accordingly, all ordered lists of fragments (configurations) are equally probable. The number of configurations is where m = m A +m B +· · · is the total mass of the particle and m! is the multinomial coefficient of the compositional vector m = (m A , m B · · · ) of the parent. The multinomial coefficient of the compositional vector is the number of permutations of the string of beads and the binomial factor is the number of ways to cut it into N pieces. The multiplicity ω(n) of distribution n is the number of configurations with that distribution and is given by the number of permutations in the order of fragments and in the order of components within each fragment: (50) where c! is the multinomial factor of the compositional vector c = (a, b, c · · · ) of the fragment. The probability of distribution n is its multiplicity normalized by the total number of configurations with The normalizing factor Ω (K) mA,mB ··· ;N is the partition function of the ensemble of distributions that satisfy the constraints in Eqs. (47) and (48).
The mean distribution of fragments is and is the generalization of (38). Alternatively, the mean distribution can be expressed by analogy to Eq. (44) as the product of the color blind distribution with a mean compositional distribution: The mean color-blind size distribution n k /N is the same as in one-component fragmentation, with m = m A + m B + · · · and c a,b···|k is the conditional probability that the compositional vector of fragment size k is c = (a, b · · · ): This is the generalization of Eq. (43).

III. NON RANDOM FRAGMENTATION
In random fragmentation we produce permutations of the configuration of fragments at random and accept them with uniform probability 1. We now bias the acceptance of the permutation by a functional W (n) of the fragment distribution such that the probability of fragment distribution n in is with ω(n) from Eq. (50). These are the general forms of Eqs. (51) and (52) in the bicomponent case and are easily extended to any number of components. Here ω(n) is the intrinsic multiplicity of n in the ensemble of fragments, while the product ω(n)W (n) . =ω(n) is its apparent (biased) multiplicity as weighted by the bias functional and distinguished by the tilde. Similarly, the partition func-tionΩ is the summation of the apparent (biased) multiplicities of all distributions in E mA,mB;N . With W = 1 we recover the random case. In this sense "random" and "unbiased" both refer to the case of uniform bias.

A. Linear ensemble
The bias functional W will remain unspecified. This allows us to choose the bias so as to produce any desired distribution of fragments. A special but important case is when W is of the factorizable form where w a,b depend on a and b but not on the fragment distribution itself. The log of the bias is then a linear function of n: The result states that the log of the bias is homogeneous functional of n with degree 1, i.e. log(λn) = λ log W (n) for any λ > 0. We refer to this case as linear bias with the understanding that linearity actually refers to the log of W .

Theorem 5. The mean distribution of fragments under the bias in Eq. (59) is
Proof. Using Eq. (59) the apparent multiplicityω(n) of distribution n is with The result then follows directly from Theorem 4.

B. Composition-independent bias
If the bias factors are of the form w a,b = g(a + b), where g(x) is a function of a single variable, the acceptance probability of a configuration of fragments depends on the mass k = a + b of the fragment but not on its composition. This leads to a simple expression for the mean distribution by the following argument. With reference to Fig. 1, fix the points where the string is cut; this amounts to fixing the size distribution of the fragments. All permutations of the colored beads are equally probable because they have the same value of W . Accordingly, the compositional distribution is the same as in the random case and is given by Eq. (43). The size distribution on the other hand is biased and is the same as when the same bias is applied to one-component distribution. The final result is where n k is the one-component size distribution under bias w a,b = g(a + b), Except for special forms of g(k) the partition function will not be generally available in closed form. Table I summarizes three cases for which exact results are possible. All three cases are associated with distributions encountered in binary aggregation [20]. The partition functions in cases 1 and 2 refer to the constant and sum kernels, respectively, and are exact; Case 3 is associated with the product kernel and gives the asymptotic limit of the partition function for m, N ≫ 1, m/N < 2, conditions that refer to the pre-gel state [21]. In the general case w a,b depends on both a and b explicitly and affects both the size and compositional distributions. This case will be demonstrated by simulation in the next section.

IV. SIMULATION OF BIASED FRAGMENTATION
Except for certain special forms of the bias the mean fragment distribution cannot be calculated analytically and the only recourse is stochastic simulation. Here we describe a Monte Carlo (MC) algorithm for sampling the ensemble of distributions. We will then use this method to demonstrate result for two cases of biased bicomponent fragmentation.
A. Monte Carlo sampling by exchange reaction Suppose C = ((a 1 , b 1 ), · · · (a N , b N )) is a configuration of N fragments such that fragment i contains a i units of component A and b i units of component B. Since all configurations with the same distribution have the same probability, given in Eq. (57), the probability of configuration is where W (C) = W (n) and n is the distribution of the configuration. If the bias functional is of the linear form in Eq. (59), its value on C is Suppose that two fragments i and j exchange mass according to the reaction under the conservation conditions a i + a j = a ′ i + a ′ j and This amounts to a transition between configurations, C → C ′ , with equilibrium constant This has the form of a chemical equilibrium constant for the exchange reaction in (70) with w playing the role of the activity of the species. The ensemble of fragment distributions may then be sampled via Monte Carlo by direct analogy to chemical reactions [22]. Given a configuration of fragments, pick two at random and implement an exchange of mass by a random amount. The resulting configuration is accepted by the Metropolis criterion: accept if rnd ≤ K C→C ′ , where rnd is a random number uniformly distributed in (0, 1); otherwise reject. With W = 1 every exchange reaction is accepted, which amounts to random fragmentation. We implement the random exchange between fragments as follows. We represent fragments as a list of 1's (representing component A) and 0's (component B). We pick two clusters i and j and random. We merge them into a single list, randomize the order of components, and break into two new fragments by picking a break point at random. The randomization of the order of components in the merged list ensures that all permutations are considered with equal probability.

B. Two examples
In random fragmentation (w a,b = 1) the compositional distribution is given by Eq. (43). We may choose the bias functional so as to produce deviations in either direction relative to the random case. It is possible to produce positive deviations (preferential segregation of components in the fragments relative to random mixing) or negative deviations (more intimate mixing than in random mixing). We demonstrate both behaviors using the two examples below:

Case I (positive deviations)
2. Case II (negative deviations) In Case I the fragment bias w a,b is an additive function of the amounts of the two components. Considering that a + b is constrained by mass balance, the fragment bias is large for fragments that are rich in either component but small for fragments that are relatively mixed. This ought to favor the formation of fragments in which the components are relative segregated. The fragment bias in Case II is a multiplicative function of the amounts of the two components. It is large in fragments that contain both components but quite small if one component is present in excess of the other. This form ought to produce fragments that are better mixed than fragments produced by random fragmentation.
We test these behaviors in Fig. 3 which shows results for α = 4. In this example the particle contains an equal number of units of each component, m A = m B = 18, and breaks into N = 4 pieces. In both cases the size distribution deviates from that in random fragmentation. Indeed the bias influences both the size and the compositions distribution unless its is of the special form w a,b = g(a + b). Compositional distributions are shown for sieve-cut masses k = 2, 4 and 8. The additive bias (Case I) produces distributions that are more spread out relative to the random case. For k = 2, in particular, the compositional distribution is inverted relative to the random case, indicating strong segregation as the majority of fragments contains pure component A or B and only few fragments in this size contain both components. As the fragment size increases the separation of components is less strong but always present, as indicated by the fact that the random distribution is always narrower. The opposite behavior is observed in Case II: distributions are narrower than those in random fragmentation, especially at the smaller fragment sizes.
As a general trend in both cases, small fragments are less mixed while large fragments progressively approach the distribution of random fragmentation. This is because there are not enough units of each component to produce large fragments that consist predominantly of one component. This limitation is not present when the fragment size is small.

V. DISCUSSION AND CONCLUSIONS
We have presented a treatment of multicomponent fragmentation on the basis of random fragmentation in combination with an appropriate functional that biases the ensemble of feasible distributions. The two key notions in this treatment are the set of feasible distributions and the multiplicity of distribution within this set as established by the rules that define "random" fragmentation. In the random-fragmentation ensemble distributions are proportional to their multiplicity. This problem is analytically tractable and we have presented its solution for any number of component and number of fragments. A third key notion is that of the bias functional that modulates the probability of distributions of feasible distributions and allows us to obtain fragment distributions other than that of random fragmentation.
It should be pointed out that the random case is not endowed with universal physical significance but applies in certain cases such as the linear chain in Fig. 1. In this particular case selecting the bonds to break at random might be a reasonable physical model, as Montroll and Simha [7] explain. The primary significance of random fragmentation is mathematical. Similar to the "fair coin" or the "ideal solution," it provides an analytically solvable baseline ("reference state") from which to calculate deviations. The mathematical tool that quantifies these deviations is the bias functional. This functional permits the systematic construction of distributions that exhibit any degree of deviation from the random case. This is main result of this formulation. The fragment distribution per fragmentation event is one of two elements required in order to build a population balance model of a fragmentation process. The other element is the rate at which particles break up. This question is not addressed here beyond the generic observation that this rate must be a function of the compositional vector m = (m A , m B · · · ) of the particle.
In single-component fragmentation the quantity of interest is the mean size distribution of the fragments. In multicomponent systems we are additionally concerned with the compositional distribution. This introduces a new dimension to the problem and raises questions of mixing and unmixing of components. Do fragments inherit the compositional characteristics of the parent particle? Do they become progressively more well mixed or less? Both behaviors are possible and are quantified via the bias functional W . This functional is where the mathematical theory of fragmentation presented here makes contact with the physical mechanisms that lead to the disintegration of material particles. To make this connection quantitatively, one must begin with the a physical model of fragmentation that assigns probabilities to all possible distributions of fragments that can be generated. This is a major undertaking and is specific to the particular problem that is being considered. The point we wish to make is that the formulation presented here offers an entry point to physics via the bias functional.
Lastly, the connection to statistical mechanics should not be lost. We have constructed an ensemble whose fundamental element ("microstate") is a the ordered configuration of fragments; its total number in the ensemble is the partition function. The higher-level stochastic variable (the observable) is the distribution of fragments and its probability is determined by its multiplicity in the ensemble. The form of the probability in Eq. (13), also known as Gibbs distribution [23], is encountered in time reversible processes as well as in population balances of aggregation and breakup [17,[23][24][25][26]. The derivation of the mean distribution in the random case follows in the steps of the Darwin-Fowler method [27]. Additionally, the compositional distribution in random breakup is given asymptotically by the binomial distribution in Eq. (46). This establishes a reference for compositional interactions analogous to that of the ideal solution in thermodynamics. In fact, the Shannon entropy of the binomial distribution is the ideal entropy of mixing when two pure components coalesce into a single particle that contains mass fraction φ A of component A. These connections are not coincidental. Biased sampling from a distribution generates a probability space of distributions and when the base distribution is exponential, this ensemble obeys thermodynamics [22]. In fragmentation the base distribution is a multicomponent exponential: the size distribution in Eq. (12) goes over to the exponential distribution when m, N ≫ 1. In this limit the ensemble of fragments becomes mathematically equivalent to a thermodynamic ensemble of two components with interactions that lead to positive or negative deviations relative to ideal solution.