A stochastic model leading to various particle mass distributions including the RRSB distribution

Modern particle size statistics uses many different statistical distributions, but these distributions are empirical approximations for theoretically unknown relationships. This also holds true for the famous RRSB (Rosin-Rammler-Sperling-Bennett) distribution. Based on the compound Poisson process, this paper introduces a simple stochastic model that leads to a general product form of particle mass distributions. The beauty of this product form is that its two factors characterize separately the two main components of samples of particles, namely, individual particle masses and total particle number. The RRSB distribution belongs to the class of distributions following the new model. Its simple product form can be a starting point for developing new particle mass distributions. The model is applied to the statistical analysis of samples of blast-produced fragments measured by hand, which enables a precise investigation of the mass-size relationship. This model-based analysis leads to plausible estimates of the mass and size factors and helps to understand the influence of blasting conditions on fragment-mass distributions.

Although many distributions have been developed, most of them, including the RRSB distribution, are empirical or approximations for theoretically unknown relationships, with little binding with material properties or crushing devices.The problem of finding suitable models that explain these forms is still open.Three notable exceptions are Kolmogorov [23], Brown and Wohletz [24], and Fowler and Scheu [25], who modelled in different extent the fragmentation process.The first introduced the lognormal distribution to particle statistics and the others developed theories that lead to the RRSB and gamma distribution, respectively.
The present paper introduces a general stochastic model that leads to a wide class of particle mass distributions.It does not model the fragmentation process but instead the structure of particle samples, the number and the sizes and masses of particles.The model does not employ ideas of fractal geometry and self-similarity, but it only operates with finite random variables.The model components can 67 Page 2 of 12 be observed and statistically analysed in particle samples.It may help to understand the appearance of the RRSB distribution as well as the Gates-Gaudin-Schuhmann (GGS) power-law distribution.
The model is based on the point process approach introduced by Stoyan and Unland [26].According to this approach, the particles of a sample are assumed to be ordered with respect to size (some quantity of the nature of a length, as for example, a Feret diameter, see Sect. 2) and interpreted as a point process on the real axis.The new idea of the present paper is to mark the points with the corresponding particle masses.These mass marks are random, and their random variability results in a large extent from fluctuations of particle shapes, which lead to random relations between particle size and mass.
This improved approach was inspired by the paper of Zhang et al. [27], which offers data of particle masses for single particles.The data came from blasting experiments in which the fragments larger than 20 mm were individually measured and weighted, while the fragments smaller than 20 mm were sieved.
The present paper first describes the stochastic model, which is closely related to the so-called compound Poisson process, a classical model of applied probability, see Sect. 2. Section 3 shows how known particle mass distributions may be classified as belonging to the model.The data from [27] are analysed statistically in Sect. 4 in the spirit of the model.Finally, the whole approach is discussed thoroughly.

The model including particle masses
This paper develops the theory in Stoyan and Unland [26] by including the particle masses into the theory.It uses the notation of that paper (which is in the same journal and in free access) and starts with its basic idea: to consider the particle sizes as points on the positive real line or x-axis.
Following the classical approach of sieving statistics, "sizes" are here one-dimensional geometrical characteristics of particles.In the example considered in Sect.4, size is the maximum Feret diameter.In [26], where the particles were measured by means of image analysis (CPA), size was the diameter of the circle of equal projected area.And a third choice of size could be the diameter of a volume-equivalent sphere.
If arranged in an increased order, the size-points form a random sequence, which is interpreted as a one-dimensional point process.Following [26] we assume that this point process is an inhomogeneous Poisson process.
The notion of a Poisson process is explained in [26,28,29].An inhomogeneous Poisson process is a random collection of points (in the case here considered: of particle sizes) on the positive x-axis.The numbers of points in disjoint intervals are independent, and the point density is controlled by a so-called intensity function λ 0 (x), which here has the name size intensity.The random number N((a,b)) of points in the interval (a, b) has a Poisson distribution with mean.meaning that the probability p(n) to find n particles in (a,b) is where λ = E(N((a, b))).
The index ''0″ at λ 0 (x) indicates that zero-dimensional numbers are counted.The shape of λ 0 (x) depends on the process by which the particles were produced.In the present paper as well as in [26] the size intensity is decreasing in x since there are more small particles than large ones; of course, other forms are possible.Now, in addition to particle size, also particle mass m (measured in grams) is included in the modelling.First, the cumulated mass of all particles which are smaller than x is considered.The mean of this random variable is denoted in [26] by Λ 3 (x) called cumulative mass intensity function, considered as a function of x.The index "3" indicates that three-dimensional quantities are counted and added.
The mean Λ 3 (x) can be calculated starting with the fol- lowing statistical assumption: the masses of the particles in the point sequence (ordered with respect to size) are independent but, of course, not identically distributed (In point process statistics one speaks of "independent, locationdependent marking").Indeed, for fixed size x the mass can be greatly variable due to variability of shape.It can happen that a particle of large size may have low mass because of a special shape.Fortunately, for the calculation of Λ 3 (x) only the mean mass of particles of size x (and not the corresponding distribution) is needed, which is denoted by μ(x).
This deterministic mean-mass function (x) is an impor- tant model characteristic, with µ(0) = 0 and µ(∞) = ∞.Its form heavily depends on the shape of the particles considered.In the case of spherical particles of identical material, it holds Here is density, e.g., in the unit of g/mm 3 , and x is the diameter.In this special case of identical particle shapes the masses for size x are constant, not random, given by (3).While (x) may monotonously increase with increas- ing x, other forms may also apply, for example when small particles have high density. (1) x 3 6 .

Page 3 of 12 67
Both functions, size intensity 0 (x) and mean-mass func- tion (x) , determine Λ 3 (x) as well as the corresponding derivative 3 (x) with respect to x.Indeed, 3 (x) satisfies the equation.
The product formula (4) can be proved as follows.Let Δx be the length of a small size interval.The probability that in the size interval (x, x + Δx) a particle size is present is approximately 0 (x)Δx .Consequently, the increment of Λ 3 (x) in the size interval is 0 (x)Δx ⋅ (x) and it holds.
Then and division by ∆x and the limit ∆x → 0 yields Eq. ( 4).Note that 3 (x) is fully determined by the mean-value functions 0 (x) and µ(x), which is natural since Λ 3 (x) is also a mean value.

Product-form of particle mass probability density
The cumulative mass intensity function Λ 3 (x) is closely related to the cumulative particle mass distribution function if there is a finite maximum size x max .In the following we mainly consider the infinite case.Note that Q 3 (x) here denotes the theoretical distribution function, whereas in the usual engineering literature it often denotes the empirical distribution function resulting from one (finite) particle sample.Equation ( 4) implies that the corresponding mass density function q 3 (x), the derivative of This has two important statistical consequences: (a) The influences of the two determinants 'frequency' (given by 0 (x) ) and 'size' (given by (x) ) are separated and they appear as factors in ( 6), ( 4) (b) Many forms of 0 (x) and (x) may lead to particle mass distribution functions.Mathematically it is sufficient that the integral over 0 (x) (x) is finite.
Point (a) may help to find relations of Q 3 (x) to material properties and crushing conditions, while (b) may direct the search for forms of Q 3 (x).

Particle mass distributions following the model
Some well-known particle mass distribution functions can be obtained by suitable choices of 0 (x) and (x).This may inspire researchers to find new distributions in their applications.

The RRSB distribution
Assume and with positive parameters b, m, n, u and v in suitable dimensions.Both functions are plausible and nice starting points, with decreasing number of particles and increasing mean mass of particles for increasing size x in a form similar to that in the case of spheres.
Then, if the model is accepted, the relation (6) yields that the particle mass probability density function q 3 (x) is proportional to ax m ⋅ exp (−bx n ), with a = uv, that is, This relation is nearly the same result which Rosin and Rammler had before the interaction by Karl Sperling (the "S" in "RRSB"), see Eq. ( 1) in [2].
The standard form of the RRSB distribution is where r is a positive parameter (called 'uniformity parameter') and x ′ a further positive parameter (called 'character- istic size').The corresponding mass density function The following choices of the parameters a, b, m, and n transform the relation ( 7) into (9): 67 Page 4 of 12

The Gates-Gaudin-Schuhmann distribution
The Gates-Gaudin-Schuhmann distribution has the distribution function where k is a positive parameter.The corresponding probability density function is , where an additional parameter l is introduced.Also, this function has the product form (6), and the term k x k max ⋅ x −l may play the role of 0 (x).
It is similar with the Gaudin-Meloy distribution [30],

The Gamma distribution
The gamma distribution is not so frequently used as RRSB and GGS, perhaps because of the appearance of the gamma function.Its probability density function is with positive parameters c and p.It has the product form (6) with 0 (x) proportional to exp(-cx) and µ(x) proportional to x p .Note that this is not RRSB for r = 1!In order to demonstrate its importance, two quite different applications are mentioned here.It appeared in a statistical analysis of cubic quartz sand particles reported in [31], where the distribution is called Martin-Andreasen distribution.The shape parameter p has the fixed value 4, and the distribution is therefore a one-parameter distribution.
The gamma distribution plays an important role in [25] in which a mass distribution function with variable = −log 2 (x) was empirically found by statistical analysis of fragments from laboratory explosion of volcanic rocks, where x is particle size as in the present paper.This statistical result was then theoretically explained.

A problem of indefiniteness
When the model is accepted, it is natural to know the meanmass function µ(x) and the size intensity function 0 (x) when the density function q 3 (x) is known.Unfortunately, the x p−1 exp (−cx) for x ≥ 0 relation (6) shows that this wish cannot be fulfilled.When µ(x) is multiplied by some positive function z(x) and 0 (x) is divided by the same function, the same q 3 (x) is obtained.This idea can be used to construct suitable functions 0 (x) and µ(x) for standard distributions.Consider, for example, the case of the RRSB distribution with r = 0.5.Then A choice of µ(x) as would be nonsense since μ 0 (x) is decreasing in x, which means 'mass is decreasing in size'.However, a z(x) = sx 3 would yield a µ(x) proportional to x 2.5 and a which makes sense.

Statistics for fragments of rock blasting
All data of particles (called here "fragments") used in this section come from blasting experiments of nine granite cylinders with a diameter of 240 mm and a length or height of 300 mm [27].There was one drill hole charged with explosive PETN at the drill hole bottom.The detailed parameters of the cylinders and the explosive charges are shown in Table 1. Figure 1 shows all fragments of the granite cylinders S1-S9 after blasting.The classical result of particle mass statistics is shown in Fig. 2: the empirical Q 3 distribu- tion functions for the nine cylinders.The accumulated mass passings of S6-S9 are clearly larger than those of S1-S5.Now follows the analysis of the same data in the spirit of the present paper.Only the larger fragments are considered since the data for fragments smaller than 100 mm form unstructured dense clouds of points.Note that these data, size = maximum Feret diameter x (mm) and mass m(g) , were measured by hand [27].Figure 3 presents the empirical relation between m(x) and x. Figure 3a shows the measurement data of mass and size of all fragments for the cylinders in different colours.Very impressive is the great variability of the masses for large fragments.
Before the statistical results are presented, the basic assumptions are considered: (1) the size points form Poisson processes and (2) the mass marks are independent.Both were tested with the best statistical tests available: Poisson with χ 2 goodness-of-fit test as described in [26], Section 4.3.3 and end of 6.3, and independence with the Table 1 The parameters of the rock specimens, explosive, stemming and charge conditions [27] Specimen  phase-frequency test of Wallis and Moore [32].In the latter the up-jumps and down-jumps in the series of masses are considered.While the critical z-value of the test for error probability = 0.05 is 1.96, the test statistics for all cases without S9 were lower than 1.56.Only for S9 the value 2.03 was obtained, which is lower than the critical value for = 0.01 .The Poisson test is described in detail in [26] and carried out exactly in the same way as there, of course with the corresponding intensity functions.While the critical χ 2 values are 2.9 and 19.0, the test statistics are all between 6.2 and 18.5 with exception of S4, where it was 19.4, which is below the upper χ 2 -value for = 0.01 .Thus, the basic assumptions is considered to be satisfied by these data.
Figure 3b shows the measurement data of masses and sizes of the fragments of all rock cylinders S1-S9 and the curve of the corresponding non-linear regression as an estimate of (x) .Note that this is an approximation for the mean function of a stochastic process and not the empirical form of a functional relationship.The corresponding estimation equation is Figure 4 shows all nine empirical mean-mass curves together with the datasets of mass measurement, which are interpreted as estimates of the corresponding mean-mass function µ(x) where the curves end at the largest measured fragment size.The corresponding regression parameters R-squared, RSE and coefficient c are collected in Table 2. Figures 3 and 4 indicate that the fragment masses tend to increase with increasing size, which is consistent with common blast results, i.e., big fragments have larger masses.However, the variability of masses is great, and there is no monotone growth of masses with increasing size.It looks that the variability increases with increasing x.This is reasonable since the size and shape of a large fragment is dependent on several factors such as the material and length of stemming, the geometrical shape of rock sample, the specific charge, etc.For example, as shown in Fig. 4 and Table 2, the curves of four cylinders S6-S9 with full stemming have mostly low value of parameter c (0.00029, 0.00026, 0.00027 and 0.00036 for S6, S7, S8 and S9 respectively), while the curves of five cylinders S1-S5 with partial steel stemming have relatively high value of parameter c (0.00038, 0.00044, 0.00043, 0.00023 and 0.00052 for S1, S2, S3, S4 and S5 respectively).The above description indicates that stemming condition influences the relation between fragment size and fragment weight.In addition, from Fig. 4 and Table 2 it can be found that cylinder S7 has the least parameter c (= 0.00026) and smallest maximum fragment mass (about 2000 g) among all of four cylinders with sand stemming.The reason is that S7 has higher  specific charge than other cylinders, meaning that specific charge effects the relation between fragment size and fragment weight.This result refines the measurement result for the accumulated mass passing vs. fragment size of the nine samples, as shown in Fig. 1. Figure 3 shows that the mean mass of all cylinders increases with increasing fragment size in the form of x 3 .Now the estimates of size intensity 0 (x) are presented.Remember that 0 (x) is the density of size points on the x-axis at size x, note Eq. ( 1).As recommended and explained in [26], a kernel estimator was used for the estimation of 0 (x) .The estimator is.
where the summation goes over all fragment sizes x i .The symbol k(⋅) denotes the kernel function, which is taken as a Gaussian probability density function, The parameter σ controls the smoothness of the estimate, the function 0 (x) .After experimentation we decided to choose σ = 30 mm.
Figure 5a shows the empirical size intensities 0 (x) for x ≥ 100mm ; for smaller x the values are very large because there are many small fragments.As expected, the curves decrease with increasing fragment size x.They look like decreasing exponential functions k ⋅ exp (− x) , with simi- lar in all cases and largest k for S6.Additionally, Fig. 5b shows again the curve for S6 plus the standard deviation of the errors.The standard deviation values were obtained by bootstrapping, following the ideas of [33].This means that Poisson processes with the empirical 0 (x) were simulated 1000 times and their intensity functions were re-estimated.
The above description indicates that for large x the masssize relation of the granite fragments is similar to that of the case of a RRSB distribution.

Discussion
This paper is inspired by a quite unusual form of particle size statistics, where many large particles were measured individually by hand, both in size and in mass.This resulted in a unique data set and led to the new form of modelling the variability of particle samples and to a deeper understanding of particle mass distributions.Of course, it is not recommended to carry out such time-consuming measurements for daily use.For research purposes, they may, if applied to an interesting material, provide valuable information sometimes.Furthermore, when particle samples are sieved, very ( 14) � .
big particles that cannot be sieved should appear in the sieving report.And the total mass of the sample should always be reported.
The model proposed in this paper is closely related to a classical model of applied probability, the so-called compound Poisson process.Figure 6 shows a sample of such a process, which is a standard model in financial mathematics.The points on the x-axis may be instants of insurance claims or of bankrupts, while the jumps are the corresponding values of claim or loan.In the context of the present paper the points are particle sizes and the jump heights are particle masses.The points belong to a Poisson process and the jump heights are independent random variables.The compound Poisson process is then the cumulative sum of jumps as a function of x as shown in Fig. 6

(in classical applications
x is time).The theory of the compound Poisson process is presented in detail in Last and Penrose [29].This book gives formulas that characterize the variability of such processes.However, our Eq.( 3) belonging to the particular model with an inhomogeneous Poisson process and size-dependent mass distributions cannot be found there.
Of a quite different nature is the model in [24], which also leads to the RRSB distribution.It uses physical ideas, in particular branching trees of cracks.There the exponent r is related to the fractal dimension of a self-similar branching tree of cracks.
Note there is an interesting difference between the factors of the RRSB and the gamma distribution.Whereas for gamma the size intensity 0 (x) is proportional to exp( − cx) and the mean-mass function µ(x) is proportional to x p and there is no connection between the two factors, the factors of the RRSB distribution are coupled by the exponent r, which appears in both factors.This may mark the gamma distribution as belonging to a simpler model.
The mean-mass function µ(x) of both distributions, RRSB and gamma, is controlled by a power term, x r and x p , and the empirical µ(x) in Sect. 4 has such a form.This is also the case for spherical particles as in Eq. ( 2).This leads to the recommendation in the case of searching a new distribution for a particle sample to start with a meanmass function proportional to a power of size x with positive exponent.
The size intensity 0 (x) controls in a large extent the shape of q 3 (x) .For example, the existence of a maximum particle size x max is indicated by 0 (x) = 0 for x > x max , see Fig. 5. Furthermore, if 0 (x) has a pole at x = 0 (i.e., 0 (0) = ∞), then q 3 (x) has also such a pole.
The exponent r in the RRSB distribution controls the variability of sizes and masses.Its value is frequently considered to be between 0.7 and 1.5.For the large fragments in Sect. 4 the parameter c plays a role similar to r.The values of c given in Table 2 are between 0.00023 and 0.00052.This large range of c value may depend on multiple factors such as stemming condition, specific charge, rock properties, etc.For example, the value c of S1-S5 is in a large range from 0.00023 to 0.00052, while the value c of S6-S9 is in a smaller range from 0.00026 to 0.00036.Since S4 has shortest stemming length (meaning highest gas ejection) among five cylinders with steel stemming, S4 can be excluded from S1 to S5.If so, it can be found from Table 2 that the value c of S1, S2, S3 and S5 is in a large range from 0.00038 to 0.00052, which is much larger than the value c of S6-S9.This result indicates that the partial steel stemming yields much higher c value than full sand stemming.One of main reasons for this result is that the partial steel stemming resulted in more gas ejection than the full sand stemming [27].Considering that complex stemming conditions were used and only small rock cylinders employed in the blasts of S1-S9, more blasts with measured particle sizes and weights are needed to further develop this stochastic model.Fig. 6 A piece of a sample of a compound Poisson process, in the particle size and mass interpretation.Here x is particle size and M(x) is the cumulative mass of all particles of mass ≤ x .x n and M n are the size and the mass of the n-th particle, respectively.Note that this is not a curve based on sieving results.The jump points are random particle sizes A general stochastic model, which is closely related to an established model of applied probability, the so-called compound Poisson process, is introduced in this paper for particle mass distributions.It explains the appearance of the RRSB distribution, of the gamma distribution and of the GGS and GM distributions.The model has the potential to lead to new particle mass distributions by adapted choices of the factors 0 (x) and µ(x).Future investigations of these functions may shed new light on various fragmentation processes, in particular stemming.

Fig. 5
Fig. 5 The empirical size intensity λ 0 (x) vs. fragment size x for all nine cylinders (a); the standard deviation vs. fragment size x for specimen S6 (b).The errors are not symmetric with respect to the intensity function curve since intensities are always positive

Table 2
The estimated regression parameters a, b, c and residual standard error(RSE)