Discontinuity frequency and block volume distribution in rock masses

International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

Contents lists available at ScienceDirect

International Journal of Rock Mechanics & Mining Sciences journal homepage: www.elsevier.com/locate/ijrmms

Discontinuity frequency and block volume distribution in rock masses M. Stavropoulou n Department of Dynamic, Tectonic and Applied Geology, Faculty of Geology and Geoenvironment, University of Athens, GR-15771 Athens, Greece

art ic l e i nf o

a b s t r a c t

Article history: Received 26 November 2012 Received in revised form 16 July 2013 Accepted 4 November 2013

The discontinuity spacing density function is theoretically found by applying the first principles of the Maximum Entropy Theory. It is shown that this function is the negative exponential probability density function. Then, the analytical relation between RQD and discontinuity frequency that may be derived, provided that discontinuity sets follow the negative exponential model, is validated against simulation data. It is also found that if the discontinuity spacings follow the negative exponential distribution, then the number of fractures per length measured along scanlines or drilled cores follow quite well the twoparameter Weibull distribution function. Subsequently, following the methodology proposed originally by Hudson and Priest, the closed-form expression of block volume distribution in a rock mass transected by three mutually orthogonal discontinuity sets is found. Also, the left-truncated block volume proportion above a certain block volume size is found analytically. The theoretical results referring to discontinuity frequency and block volume distributions are finally successfully validated against measurements carried out on drill cores and exposed walls, in a dolomitic marble quarry. The methodology presented herein can be applied to rock engineering applications that necessitate the characterization of rock mass discontinuities and discontinuity spacings are reasonably well represented by the negative exponential probability density function. The proposed method for the prediction of marble block volume distribution was applied to data from a quarry from drill cores and scanlines on exposed quarry walls. & 2013 Published by Elsevier Ltd.

Keywords: Discontinuity Scanline Drill core Joint spacing Joint frequency Block volume

1. Introduction There are a number of important practical rock engineering applications in which the knowledge of the geometry of rock discontinuities (or fractures) is of paramount significance for the correct design and construction of surface or underground works in the rock mass. It is noted that ‘fractures’ are ‘discontinuities’ because they form discontinuities in the mechanical continuum. Throughout the paper, the word ‘discontinuity’ is used as a general term encompassing cracks, fissures, joints, shear fractures, slip bands, bedding planes etc. penetrating the rock and characterized by low cohesion, that may adversely affect the strength of the rock mass as well as the quality of an extracted block from a surface or underground quarry. For example, the in situ distribution of rock block sizes formed by the mutually intersecting discontinuities, may be used for evaluating the production capability of a deposit to be mined using the block caving or the sublevel caving mining methods [1] and to assess the requirements and design of material handling systems in the mine. Laubscher [2] has considered natural rock fragmentation as a major factor affecting the design

n

Tel.: þ 30 210 727 4778; fax: þ30 210 727 4096. E-mail address: [email protected]

1365-1609/$ - see front matter & 2013 Published by Elsevier Ltd. http://dx.doi.org/10.1016/j.ijrmms.2013.11.003

process for block caving operations and stated that while all rock masses will cave, the manner of their caving and the resultant fragmentation size distribution need to be predicted if cave mining is to be implemented successfully. In caving operations, fragmentation has a bearing on drawpoint spacing, dilution entry into the draw column, draw control, drawpoint productivity and secondary blasting/breaking costs. In addition, rock mass characterization for the production of large prismatic or irregular blocks used as armourstone, i.e. blocks weighing many tonnes are used for building cover layers to resist wave action [3], constitutes the most important part of the exploration. Also, the success of marble quarrying operations using diamond wire cutting, chain sawing or pre-splitting or combinations of these techniques, depends on the yield of blocks of orthogonal parallelepiped shape with volumes greater than 1 m3 or 2 m3 or more, depending on the market demands for this specific marble. The same remark holds true also for monumental, building or decorative natural stone quarries of other geological origin, like for example limestones, sandstones, granites etc. Further, surface or underground stability analyses employing the rock block distribution as a factor for the quantitative description of the deterioration of intact rock properties [4,5]. For example Barton [6] proposes the ratio RQD/Jn to represent the relative block size, in the Q rock mass classification system, with

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

RQD denoting the Rock Quality Designation and Jn denoting the number of joint sets transecting the rock mass. Furthermore, intact rock properties and the discontinuity structure of a rock mass are among the most important variables influencing blasting results. This influence is considered to be a composite intrinsic property of a rock mass and is referred to as the ‘blastability of a rock mass’. It represents the ease with which the in situ rock mass can be fragmented and displaced by blasting [7]. It could be also mentioned that the in situ block volume distribution influences the pressure transient behavior of wells drilled in naturally fractured reservoirs [8]. The present research was stimulated by the relevant series of milestone papers [9–12]. It aims at improving the approach of prediction of joint spacings and frequencies, RQD and block size distributions based on scanline measurements or borehole loggings. Regarding the prediction of block size distribution it is concerned only with discontinuity sets occurring of parallel persistent planes irrespectively of the size of joints. Initially, in Section 2, discontinuity frequency is discussed in terms of borehole or scanline orientation. Consideration of the possible density function of spacings between discontinuities by virtue of the maximum entropy method [13,14] is given in Appendix A. By virtue of Monte-Carlo simulations it is demonstrated that if the joint spacings are generated via a Poisson process, then the measured joint number per fixed length of a scanline or borehole core follows the Weibull density function. The relation between RQD and joint frequency proposed by Priest and Hudson [9] assuming that the joint spacings follow the negative exponential density function is then validated in the same section against Monte-Carlo simulation data. In Section 3 these ideas are considered along three orthogonal axes, and analytical formulae for the block volume distribution are presented. Such distributions are examined only for discontinuities that occur in sets of parallel planes. On the experimental side and in Section 4, the previous methods are applied for an active open pit dolomitic marble quarry in order to compare the theoretical discontinuity spacing and block size distributions with those that occur in practice. Finally, the main remarks that may be drawn from this work are outlined in Section 5.

following formula x¼

1 N ∑ x Ni¼1 i

As is mentioned in [9], Piteau [15] has used a scanline (measuring tape) survey technique on rock faces and expressed discontinuity intensity as the number of discontinuities per unit distance normal to the strike of a set of sub-parallel discontinuities. It is remarked here that the quantity of discontinuities present in a rock mass is usually expressed by the discontinuity frequency (or otherwise fracture frequency) denoted here by the symbol λ which is defined as the mean number of discontinuities per unit length intersected along a borehole or scanline (measuring tape) set up on a rock exposure. Furthermore, in the succeeding analysis the symbol fλ will denote the number of discontinuities per meter that is measured at small regular intervals (i.e. sample support) along a scanline or borehole core and it is assumed that all discontinuities belonging to the same set or family are mutually parallel to each other. Discontinuities in rock are never uniformly distributed at all orientations, but usually occur in sets. The mean true value of the spacings of a given discontinuity set along a scanline or borehole perpendicular to the discontinuities may be found by the

ð1Þ

where N denotes the number of measured spacings along a scanline and xi is the ith measurement of spacing including the measurement resulting by adding the beginning and end values in a core or interval [11]. For a sufficiently large sample of spacing measurements along a scanline of length L the following approximation holds true: N=L  λ

ð2Þ

The linear frequency value depends on the direction of the line through the rock mass, and there is a maximum value in one direction and a minimum value in another direction. Values of linear frequency in different directions could be crucial in estimating both the true frequency perpendicular to the discontinuities, as well as for the sizes of the formed rock blocks. Let θ be defined as the solid acute angle subtended between the orientation of the borehole or scanline, and the normal to the plane of the fractures. Then, the apparent frequency λ′ along the scanline or borehole will be less than that along the normal because the distance, x, between two successive discontinuities intersected by the normal is increased to x′¼ x/cos θ along the scanline where x′ denotes the apparent or measured spacing. It may be easily inferred that the apparent measured discontinuity frequency is given by

λ′ ¼ λ cos θ

ð3Þ

If the spacing between neighboring discontinuities of the same set is considered to be a continuous Random Variable (RV) and is denoted by X, then its distribution function (known also as cumulative function or cumulative distribution function, cdf) is denoted as F(x)¼ P{Xrx} and designates the probability that a given spacing value X is less than x. For example the density function (or probability density function, pdf) f(x) of spacing values for a negative exponential distribution is given by f ðxÞ ¼ λe  λx

ð4Þ

Then the corresponding cumulative function for a negative exponential distribution is [11] Z x Z x FðxÞ ¼ f ðξÞdξ ¼ λe  λξ dξ ¼ ½  e  λξ x0 ¼ 1  e  λx ð5Þ 0

2. Discontinuity frequency along scanlines or boreholes and its relation with RQD

63

0

The negative exponential pdf of joint spacing values is a frequently encountered function in fractured rock masses. The reason for this, by recourse to the ‘Maximum Entropy’ assumption in heterogeneous rock masses, is demonstrated in Appendix A. The analysis presented there shows that if the density function of spacing values of a given joint set cannot be estimated from joint spacing measurements for any reason, then it could be inferred that the latter obey the negative exponential pdf based on the hypothesis of maximum entropy (lack of information or uncertainty) of a joint system. Hence, the simple measurement of number of joints per length along a scanline or borehole core of sufficient total length gives the mean discontinuity frequency λ and the complete form of the distribution function. This hypothesis is justified from many measurements of joint spacings in the field like in [11] and recently in [25] among others. In the sequel we investigate which type of distribution function is followed by the measured number of fractures per length fλ, that is viewed also as a RV, provided that the joint spacings follow the negative exponential distribution function. For this purpose we use the Monte Carlo simulation method for producing a synthetic sample of large size in the manner proposed by Hudson and Priest [11]. That is to say, spacing values for any number of fracture sets are progressively selected from each of the component distributions

64

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

and as the simulation proceeds, the spacing values of the resultant distribution are generated from the mutual interference of the component distributions. All the component spacing distributions are assumed to follow the negative exponential distribution function. In this simulation procedure the joint frequency values fλ, i.e. the number of fractures per fixed length along the synthetic scanline or drilled core are also counted and stored as well. A typical result of such a simulation in the form of a histogram of number of joints per meter fλ (i.e. in this case we have chosen a fixed interval of length 1 m to count the number of fractures) produced by the simulation of three joint sets with mean frequencies 2, 2, and 0.4 m  1 along a scanline of length 500 m is shown in Fig. 1a. After quite many runs with different combinations of joint sets from one to three and various associated mean joint frequencies λ with this simulation procedure, it was found that the variance (var) of the fλ is approximately equal to the mean frequency value which for the case of negative exponential distribution of spacings is equal to the sum of the mean frequencies of the component spacing distributions. Also, it was found that the Cumulative Distribution Function (CDF) of the fλ values follows quite well the

joint frequencies histogram 120

100 Data histogram

Frequency

b

F w ðf λ Þ ¼ 1  e  ðf λ =aÞ

Eðf λ Þ ¼ αΓ ð1 þb

1

Þ; varðf λ Þ ¼ α2 fΓ ð1 þ2b

1

Þ  Γ ð1 þ b 2

1

Þg

1 b ag g

Γ ðbg Þ

b  1  f λ =ag

f λg

e

;

γ ðbg ; f λ =ag Þ Γ ðbg Þ

ð8Þ

where γ(  ) denotes the incomplete Gamma function. The mean and variance of the gamma distribution function are, respectively,

40

Eg ðf λ Þ ¼ αg bg ; varg ðf λ Þ ¼ a2g bg 20

0

0

2

4

6

8

10

12

14

16

Joint frequency (1/m)

joint frequencies CDF 1 0.9 0.8 0.7 0.6 CDF

ð7Þ

where E(  ) denotes expectation, and Γ(  ) is the Gamma function. Since there is no closed-form solution for the scale and shape parameters of the two Eqs. (6) and (7) with the left-hand-sides of these equations to be E(fλ) ¼ λ and var(fλ) ¼ λ (of course for dimensional homogeneity of both sides of the equation the value of λ has to be multiplied by a constant of 1 m-1), the values of the Weibull parameters are found by nonlinear regression analysis. At a subsequent stage the mean and the variance predicted by the Weibull distribution could be compared with the predicted values, namely E(fλ) ¼ λ and var(fλ) ¼ λ, respectively. Herein along with the Weibull distribution function we have also tested the gamma distribution two-parameter function. The gamma distribution models sum of exponentially distributed RVs and is based on two parameters, namely the scale and shape parameters, denoted here by the symbols ag, bg, respectively. The chi-square and exponential distributions, which are children of the gamma distribution, are one-parameter distributions that fix one of the two gamma parameters. The gamma pdf and cdf are respectively

F g ðf λ Þ ¼

60

ð6Þ

where in a is the scale parameter and b the shape parameter of the distribution. The mean and variance of the Weibull distribution function are, respectively, given by

f g ðf λ Þ ¼

Weibull PDF gamma PDF

80

Weibull distribution function that is given by the equation

0.5 0.4

Weibull CDF gamma CDF

0.3

Data CDF Median Value

0.2

ð9Þ

The hypotheses that the observed frequencies of joint frequencies fλ in the various classes with those which would be given by the Weibull and gamma density functions were tested by the chi-square test that is based on the frequencies in the various bins, the Kolmogorov-Smirnov (K–S) test and the regression coefficient. A typical result of such a fit for the case at hand is shown in Fig. 1a and b, referring to the frequency histogram and the cumulative distribution, respectively. From Table 1 that displays the results of the chi-square test at the significance level of α ¼0.01, it may be observed that while Weibull distribution is passing this criterion, the gamma distribution is not passing it. The p value is the probability, under assumption of the null hypothesis, of observing the given statistic or one more extreme. The K–S test like the chi-square test is one of the agreement between an observed distribution and an assumed theoretical one, but it is based on the cumulative distribution function rather than on the frequencies in the separate classes (bins). As it may be observed from Table 2, the K–S goodness-of-fit value of the Weibull distribution is lower than that of the gamma function but both of them are greater than the critical value for α ¼0.01.

0.1 0

0

2

4

6

8

10

12

14

Table 1 Chi-squared goodness of fit for simulated joint frequencies (critical value at significance level α¼ 0.01 is χ2L¼ 23.2093).

Joint frequency (1/m)

Fig. 1. Simulation results referring to (a) discontinuity frequency histogram and best-fitted Weibull and gamma density functions, and (b) cumulative distribution of observed discontinuity frequencies with Weibull and gamma best-fitted functions.

Distribution

Observed value

p-Value

Performance

Weibull Gamma

11.9492 246.1232

0.71153 1

Accept Reject

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

Table 2 Kolmogorov–Smirnov goodness of fit for simulated joint frequencies (critical value at significance level α ¼0.01 is DL ¼ 0.07243). Distribution

Observed value

Gamma Weibull

0.1556 0.11138

65

Comparison of RQD values 98 Exact Estimated 96

Z LðlÞ ¼ 0

l

λxf ðxÞdx

ð10Þ

If we substitute the expression for f(x) as is given by Eq. (4) into Eq. (10) and integrate we find the result LðℓÞ ¼ 1  e  λℓ ð1 þ λℓÞ

88

86

50

100

150

200

250

200

250

4

3

2

ð11Þ

1

0

-1

-2

-3

-4

0

50

100

150

Distance along scanline [m] maximum absolute error for 95% of measurements

ð12Þ

1 0.9 0.8 0.7 0.6 CDF

where the asterisk as a superscript denotes that this is a theoretical estimation of RQD. Priest and Hudson [9] have derived the same formula and based on a thorough comparison of its predictions with field measurements with scanlines in several sites they have found that the maximum error does not exceed 5%. Herein the validity of this finding was checked against several ‘synthetic’ joint spacing distributions along a scanline produced by recourse to a MonteCarlo simulation program for a combination of three joint sets intersecting a scanline and obeying the negative exponential density distribution. In these simulations the RQD is measured along segments of equal length along the scanline. This length was always at least fifty times the mean discontinuity spacing according to the results presented in [9]. The simulation method is similar to the method proposed in [11]. A typical result of the RQD – with threshold value of 0.1 m – calculated every 10 m of such a simulation, for a scanline length of 250 m that intersects three discontinuity sets which exhibit apparent frequencies λ′1 ¼ λ′2 ¼ 2 m  1 ; λ′3 ¼ 1 m  1 is shown in Fig. 2a. The length of the sampling interval for RQD estimation was selected to be always equal or larger than fifty times the mean spacing of ′ ′ ′ discontinuities that is equal to 1=ðλ1 þ λ2 þ λ3 Þ. The relative error of actual and theoretically estimated RQDn by formula (12), denoted here by the symbol Re, is estimated with the following

0

Distance along scanline [m]

If we are interested in finding the percentage cumulative proportional length that refers to intact core lengths greater than ℓ (commonly this is taken equal to 0.1 m) which essentially refers to the theoretical definition of RQD as has been proposed by Deere [18] with l ¼0.1 m, then from Eq. (11) it is obtained RQDn ¼ 100½1  LðℓÞ ¼ 100e  λℓ ð1 þ λℓÞ

92

90

Relative error [%]

In conclusion, after many simulation runs it may be inferred that in all cases the Weibull distribution function displays a better performance compared to the gamma distribution function. The cumulative length proportion, L(l), is the proportion of scanline or borehole that consists of all spacing values up to a given value, l. This function can be determined once the density function has been identified for a particular set of spacing values. The emphasis in the ensuing presentation is on the negative exponential distribution based on the maximum entropy hypothesis, but the concept can be applied to any distribution of spacing values. Along a scanline containing, on average, λ discontinuity intersection points per unit length, the proportion of spacing values lying in the range x to xþdx is f(x)dx [17,11] where, as before, f(x) is the density function of spacing values. Spacing values within this range make a proportional length contribution of λxf(x)dx to the scanline and hence the cumulative length proportion up to a spacing value ℓ is given by

RQD [%]

94

0.5 0.4 0.3 0.2 0.1 0 0

0.5

1

1.5

2

2.5

3

3.5

4

Absolute error [%]

Fig. 2. (a) Comparison between measured from synthetic experiments and theoretical RQD; (b) distribution of relative error along the scanline and (c) maximum absolute error including the 95% of the measurements.

expression Re ¼

RQDn  RQDs  100 RQDs

ð13Þ

66

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

wherein RQDs denotes the simulated RQD that represents here the measured RQD. For this simulation run, the distribution of the relative error along the scanline is displayed in Fig. 2b. Also, the maximum absolute value of the relative error corresponding to 95% of measurements that is 3.5%, is also shown in Fig. 2c. Based on several such simulations it was found that the formula (12) predicts in a sufficiently accurate manner the RQD of joints following the negative exponential distribution, thus giving additional validity to the observations and conclusions made in [9]. Hence, since RQD measurement on rock cores is a time consuming process prone to measurement errors, the application of relationship (12) can give good results by simply measuring the mean number of discontinuities λ″ fλ encountered in segments along the scanline or borehole of equal length, provided that this length is sufficiently larger than the mean discontinuity spacing.

3. Block volume distribution It is natural to extend the above ideas regarding spacing distribution along lines to block volume frequency and the distribution of block volume values in a 3D domain. Herein we follow the methodology proposed originally by Hudson and Priest [11]. These authors did not eventually find the block volume distribution in a closed form but rather they have presented results based on the Monte Carlo simulation method. In order to achieve the above aim we start here with the simplest model proposed in [11] namely that of the block areas distribution produced by two discontinuity sets having the same strike in a plane perpendicular to their strike which intersect at an acute angle φ in this plane as is shown in Fig. 3. The sketch at the right-hand side in the same figure illustrates how parallelograms with equal areas are generated when the edge lengths, x and y, satisfy xy ¼a for a given value of area, a. In the case of two nonorthogonally intersected discontinuity sets as shown in Fig. 3, y¼y′/sinφ, where y′ is the perpendicular distance between consecutive joints of system 2 corresponding to the true frequency λ2. The probability distribution of areas can be found by integrating the product of the edge length densities for x and y along the appropriate equal area hyperbolae. Alternatively, the probability, P {A ra}, that an area, A, isolated between mutually intersecting cracks will be less than a given area, a, is found by integration of the product of the density function of x and the probability that y is less than a/x for all x values, where x and a/x are spacing values from each distribution along the orthogonal axes [16] Z PfA rag ¼ f ðxÞPfY r ygdx ð14Þ

where y ra/x. Substituting the result that may be found from Eq. (5) PfY r yg ¼ FðyÞ ¼ 1  e  λ2

sin φy

into the integral (14) Z 1 PfA r ag ¼ λ1 e  λ1 x ð1  e  λ2

¼ 1 e  λ2

sin φa=x

sin φa=x

Þdx

ð15Þ

ð16Þ

0

we may find

pffiffiffi pffiffiffi FðaÞ ¼ PfA r ag ¼ 1  ω aK 1 ðω aÞ

ð17Þ

where K1 denotes the modified Bessel function of the 2nd kind and of 1st order, and we have defined the following variable qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ω ¼ 2 λ1 λ2 sin φ ð18Þ The term λ1λ2sin φ appearing in the definition of ω above, indicates the mean number of block areas formed between mutually intersection fractures per unit area of the plane. By differentiation Eq. (16) w.r.t. a we find the following density function for the areas f ðaÞ ¼

pffiffiffi dFðaÞ ω2 ¼ K 0 ðω aÞ da 2

ð19Þ

where K0 denotes the modified Bessel function of 0th order. The same results represented by Eq's. (17)–(19) have been also reached by Hudson and Priest [11] for φ ¼901. Along the same line of reasoning, we may extend the above results, to the case of distribution of volumes produced by three mutually intersecting discontinuity sets in the three-dimensional space. In this case we assume that the strike of the third joint set makes an angle χ with the common strike of the former two sets. This means that the block volumes are parallelepipeds formed by six parallelograms. Hence, the probability P(V rv) that a volume, V, will be less than a given volume, v, is found by integration of the product of the density function of areas given by Eq. (19) and the probability that z is less than v/a for all z values, where v/a are the apparent spacing values along the Oz-axis that is orthogonal to Ox and Oy-axes. Z PfV r vg ¼ f ðaÞPfZ rv=agda ð20Þ in which z rv/a and PfZ r zg ¼ FðzÞ ¼ 1  e  λ3

sin χ z

¼ 1  e  λ3

sin χ v=a

ð21Þ

Combining Eqs. (19)–(21), the following integral representation for the block volume distribution function is derived: Z pffiffiffi 1 1 2 ω K 0 ðω aÞð1  e  λ3 sin χ v=a Þda ð22Þ FðvÞ ¼ P fV r vg ¼ 2 0 The analytical evaluation of the above integral in terms of known functions may be done by recourse to known properties of

Fig. 3. Block area distribution in a plane produced by two persistent discontinuity sets.

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

67

Bessel functions [19] as follows: 1

ð  1Þ3k ðω2 λ3 sin χ vÞk þ 1 ½AðkÞln v þBðkÞ

k¼0

22k þ 3 Γ ðk þ 2ÞΓ ðk þ1Þ2

FðvÞ ¼ ∑

ð23Þ

where in we have defined the following parameters: AðkÞ ¼ 2 2ψ ðk þ 2Þ þ 2 lnðλ3 sin χ Þ  4ψ ðk þ 1Þ þ 4 ln ω  4 ln 2; BðkÞ ¼ 4ψ ðk þ 2Þψ ðk þ1Þ þ2 lnððλ3 sin χ ÞÞ þ π 2

 4 lnðλ3 sin χ Þψ ðk þ 2Þ  4 ln ωψ ðk þ 2Þ

þ ψ ðk þ 1Þ2  ψ ð1Þ ðk þ 2Þ  2ψ ð1Þ ðk þ 1Þ þ 4ψ ðk þ 1Þ2  4 lnðλ3 sin χ Þψ ðk þ 1Þ þ 4 ln 2ψ ðk þ 2Þ þ8 ln ω

þ 8ψ ðk þ 1Þln 2  8ψ ðk þ1Þln ω þ lnðλ3 sin χ Þln ω  4 ln 2 lnðλ3 sin χ Þ  8 ln 2 ln ω þ 8 ln 2

ð24Þ

ψ(x) denoted the digamma function that is defined as d ln Γ ðxÞ 1 dΓ ðxÞ ¼ ð25Þ ψ ðxÞ ¼ dx Γ ðxÞ dx Also, ψ(n)(x) stands for the nth polygamma function, which is

and

the nth derivative of the digamma function [20]:

ψ ðnÞ ðxÞ ¼

n

Fig. 4. Probability density functions of block volumes for three mean block volumes λ1λ2λ3 at hand.

nþ1

d d ψ ðxÞ ¼ n þ 1 ln Γ ðxÞ dxn dx

ð26Þ

By differentiation Eq. (23) w.r.t. the volume, v, we find the following density function for the volumes ð  1Þ3k ðω2 λ3 sin

1

f ðvÞ ¼ ∑

k¼0

χ vÞk þ 1 v  1 ½ðk þ 1ÞfAðkÞln v þ BðkÞg þ AðkÞ 22k þ 3 Γ ðk þ 2ÞΓ ðk þ 1Þ2 ð27Þ

The above infinite series converges rapidly and a finite number of terms n may be retained for accurate results, i.e. it may be approximated with the expression n

ð  1Þ3k ð4vλ1 λ2 λ3 sin φ sin χ Þk þ 1 v  1 ½ðkþ 1ÞfAðkÞln vþ BðkÞg þ AðkÞ

k¼0

22k þ 3 Γ ðk þ 2ÞΓ ðk þ 1Þ2

f ðvÞ  ∑

ð28Þ It is recommended though that several values of the number of terms n in the series expression above should be used to determine the convergent value of frequency for a given λ1λ2λ3sin φ sin χ. The term λ1λ2λ3sin φ sin χ appearing in the expression for the cumulative volume distribution in Eq. (28) above, indicates the mean number of volumes per unit volume of the 3D space occupied by the rock mass and gives the scale of the block volume distribution. To illustrate this, in Fig. 4 the frequencies of block volumes f(v) as they given by Eq. (28) have been plotted for three different values of the product λ1λ2λ3 at hand, assuming mutually orthogonal discontinuity sets (i.e. φ ¼ χ ¼ π/2). In direct analogy to the cumulative length proportion, LðℓÞ, the cumulative volume proportion V(v) is found from the following formula:    1 n Δv VðvÞ ¼ ð29Þ ∑ f ðvi Þ iΔv  V tot i ¼ 1 2 where Δv is the class interval of the volume frequency distribution histogram, f(vi) is the numerical frequency of volume values in the ith class interval of the frequency distribution, and Vtot is the total volume of rock mass. Following the above method, the cumulative volume proportion curves illustrated in Fig. 5 were produced. The vertical axis on the graph gives the percentage volume of rock mass that consists of block volumes less than a volume specified by the value on the horizontal axis. For example, taking values λ1 ¼ λ2 ¼ λ3 ¼1/m, the percentage of a rock mass that consists of blocks with a volume less than 1 m3 is approximately 23% for negative exponential distributions of discontinuity spacings. In marble quarries (a) the volume proportion of blocks larger than a given volume as well as (b) the maximum block

Fig. 5. Cumulative distributions with blocks smaller than given volume for three mean block volumes λ1λ2λ3 at hand.

volume that could be extracted is also of special interest. From the same graph it may be seen that the maximum block size for λ1λ2λ3 ¼10 m  3 is around 7 m3; whereas for the case of λ1λ2λ3 ¼16 m  3 it is approximately equal to 5 m3. It is a usual practice of a decorative stone quarry that only the block volume distribution of extracted marble blocks above a certain volume is assessed (for example vL Z1 m3), since the blocks with lower volume from this threshold are dumped as waste material without further measuring their volumes. In order to compare this experimental block volume distribution with the theoretical one derived previously, we should be able to derive the left-truncated block volume proportion above a certain block volume size. Hence, based on basic statistical principles, given a point of truncation, vL, the left-truncated cumulative volume distribution V(v) denoted as VLTN(v) can be stated as follows ( V LTN ðvÞ ¼

0

vr0

Vðv þ vL Þ  V ðvL Þ ; 1  V ðvL Þ

vZ0

ð30Þ

Fig. 6 displays the truncated distributions of block volumes above v ¼1 m3 for the three mean block volumes at hand.

68

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

4. Case study of a white dolomitic marble quarry 4.1. Basic structural features of the quarry and the method of excavation of blocks Understanding a quarry in terms of its potential for production of decorative or building stones in the form of orthogonal blocks, presents a special challenge for the mining engineer and the engineering geologist. Unlike blasting in aggregates and mining operations, optimization of the extraction process in a marble quarry, for example, has a focus on the potential for production of large orthogonal parallelepiped blocks with volume greater than 1 or 2 m3 that are free of crack-like defects. An actual quarry is considered here that is located in a dolomitic marble formation. The quarry has the form of an open

Fig. 6. Left-truncated cumulative distributions above vL ¼1 m3 with blocks smaller than given volume for the three cases at hand.

surface excavation with vertical benches of 6 m height. There were identified three joint sets transecting the marble, namely the grain, the head-grain that has the same strike and it is almost orthogonal to the former set, and a secondary set with a strike orthogonal to the common strike of the other two families of joints. All the three sets were created during the uplift of the

Fig. 8. Isometric view of the model of a jointed marble bench with excavated panel with dimensions 10  10  6 m3. Diamond wire cutting planes oriented N201W (direction of head-grain planes) and orthogonal to it. Oy-axis points to the North. The trace of the grain and head-grain are indicated in the ZOX plane, whereas the trace of the secondary plane is indicated in the XOY plane.

Fig. 7. Lower hemisphere equal area stereographic projection of joint sets; (a) poles concentration, and (b) great circles of the three principal joint systems.

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

dolomitic marble layer of a thickness of 200–300 m up to the surface. The marble has been initially failed in shear along the weakest bedding planes; at the same time the head-grain planes as well as the secondary joints have been formed in order to accommodate the large shear displacements along the master sliding grain planes. The orientations of the poles of these three sets processed with specialized software [21] are illustrated in the lower-hemisphere stereographic projection diagram of Fig. 7. In the same figure the great circles of these sets are also displayed.

Fig. 9. Method of marking the joints along the drilled core and apparent spacing measurements s1, s2, etc.

69

The marble is excavated by using diamond wire cuts at three mutually perpendicular planes as usual, i.e. one horizontal in a first place 10  10 m2, and then two vertical cuts with dimensions 10 (horizontal)  6 (vertical) m2 oriented along the head-grains and secondary joints, respectively. As was mentioned before the bench height is 6 m according to the usual quarrying practice. The initial dimensions of the panel are 10  10  6 m3. Subsequently four vertical cuts are made along the grain planes at 2–2.5 m apart in order to produce orthogonal parallelepiped sub-panels that are easy to be tilted by the excavator. Such a typical panel of marble with volume 10  10  6 m3 inside the quarry constructed by virtue of a distinct element code [22] is shown in Fig. 8; since this figure is for illustrating the directions of cuts with respect to the joint orientations we have assumed a uniform distribution of spacings of all three joint sets. In order to facilitate the measurements of the spacings between adjacent joints as well as the frequencies of the same set on each photo of a box with cores, the joint traces appearing on the marble cores have been properly identified and then marked carefully with different colors as is shown in Fig. 9, namely: (a) grain planes were marked with green color; (b) head-grain planes were marked with red color, and (c) secondary planes were marked with blue color. 4.2. Experimental results on discontinuity spacing distributions and frequencies

Table 3 Estimated apparent and true frequencies from spacing measurements on the drill cores and on an exposed wall of the quarry. Method of Joint set measurement

Dip angle [deg]

Mean true Number of Apparent frequency, measurements true frequency, λ′ λ [1/m] [1/m]

Drill-core Drill-core Drill-core Exposed quarry wall

40 70 85 85

733 618 41 354

Grain Head-grain Secondary Secondary

1.5 1.3 0.14 –

2 3.7 1.66 1.71

Twenty vertical boreholes from the drilling campaign have been logged for Fracture Frequency (FF) in units of (1/m) by counting the number of joints per meter λf, i.e. the support in the original data is equal to 1 m. ‘Support’ is a geostatistical term indicating the size of a sample. The twenty vertical boreholes penetrating the marble were drilled mainly along EW and NS directions. This approximates a direction perpendicular and parallel to the strike of the bedding or grain planes that is almost coincident with the strike of the head-grains (e.g. Fig. 7). Since the secondary joints are steeply dipping, the joint frequency observed along the drilled cores is mainly due to the grain and head-grain

Fig. 10. Distribution of the four marble qualities expressed as Fracture Frequency (1/m) along the vertical boreholes inside the planned quarry limits (see color bar for the four marble qualities).

70

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

joints (i.e. Table 3). For this reason in order to validate the joint frequency distribution inferred from the drill cores, it was necessary at a later stage to map the secondary planes along a horizontal scanline oriented perpendicularly with the mean strike of these planes. The location of the boreholes and the measured FF's along them at every 1 m apart are shown in Fig. 10. Joint apparent spacing data for the three principal joint sets obtained from the drill core inspection are presented below in the

form of frequency histograms and cumulative distribution plots. Three distribution functions have been best-fitted on each set of data namely, the one-parameter negative exponential distribution function, the Weibull, as well as, the gamma two-parameter distribution functions. Fig. 11a shows the grain spacings histogram deduced from the spacing measurements along the vertical drill cores and the best-fitted density functions at hand. Fig. 11b illustrates the cumulative distributions of measurements and of

joint spacings histogram 800

0.9 0.8

600

0.7

500

0.6

CDF

Frequency

700

joint spacings CDF

1

Data histogram Negative exponential PDF gamma PDF Weibull PDF

400

neg exponential CDF gamma CDF Weibull CDF Data CDF Median Value

0.5 0.4

300 0.3 200

0.2

100 0

0.1 0

1

2

3

4

5

0

6

0

1

2

3

Spacing [m]

4

5

6

joint spacing [m]

joint spacings CDF

joint spacings histogram 1 Data histogram Negative exponential PDF gamma PDF Weibull PDF

600

500

0.9 0.8

neg exponential CDF gamma CDF Weibull CDF Data CDF Median Value

0.6

CDF

Frequency

0.7 400

300

0.5 0.4

200

0.3 0.2

100 0.1 0

0 0

1

2

3

4

5

6

7

8

0

1

2

5

6

7

8

1

Data histogram Negative exponential PDF gamma PDF Weibull PDF

0.9 0.8 0.7

4

neg exponential CDF gamma CDF Weibull CDF Data CDF Median Value

0.6

CDF

Frequency

4

joint spacings CDF

joint spacings histogram

5

3

joint spacing [m]

Spacing [m]

3

0.5 0.4

2 0.3 0.2

1

0.1 0 0

5

10

15

Spacing [m]

20

25

0 0

5

10

15

20

25

joint spacing [m]

Fig. 11. Histograms of measured apparent joint spacings measured along vertical drill cores and best-fitted distribution functions for each joint set, i.e. (a) histogram of grains, (b) cumulative curve of grains, (c) histogram of head-grains, (d) cumulative curve of head-grains, (e) histogram of secondary joints, and (f) cumulative curve of secondary joints.

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

the three theoretical cdf's for the same joint set. In a similar fashion, Fig. 11c and d illustrates the results for the head-grains and Fig. 11e and f for the secondary joints occurring in the quarry. Table 3 presents the main results of the mean apparent (measured) and true (corrected) frequencies of the three main joint sets in marble. Table 4 presents the results of the chi-square goodness-of-fit test of the three distribution functions on each set of joint spacing data at the significance level of α ¼0.01. The following observations could be done from these results, namely: (a) The Weibull and gamma two-parameter distribution functions always exhibit better performance compared to the one-parameter negative exponential distribution function, with the only exception for the case of secondary joint system in which the negative exponential has better performance. This was expected since the negative exponential function has only one parameter while the other two have two parameters. (b) In only one case of the grain system the negative exponential distribution displays a higher observed value than the critical value of the test. (c) The Weibull distribution function always displays a lower observed value and larger correlation coefficient than the other two distribution functions. In order to have a better picture of the performance of the considered three distribution functions against the drill core data, Table 5 illustrates the results pertaining to the K–S test at the same significance level. The main results of these K–S tests may be summarized as follows: (i) The observed values of all theoretical functions are larger than the critical value for the grain and headgrain joint sets, while the Weibull displays the lower observed value. (ii) Regarding the secondary joint set all the distribution functions display a lower observed value compared to the critical one, with the negative exponential exhibiting the better performance compared to the other two.

71

In order to check the validity of the estimated mean frequency of the secondary joints from drill cores, additional measurements of spacings of joints from this set have been carried out along a horizontal scanline of 200 m length on an exposed vertical wall of the quarry oriented in an orthogonal direction with the strike of the almost vertical secondary joints. The results of this additional survey are presented in the form of a frequency histogram and a cumulative distribution in Fig. 12a and b, respectively; whereas the mean joint frequency (in this case the apparent joint frequency is identical with the true) is shown in Table 3. From this table it may be observed that the secondary joints frequency is almost the same for the two sampling methods. Again for this set of data the three theoretical distribution functions have been best-fitted and displayed in the two graphs of Fig. 12a and b. The goodness-of-fit test results referring to the chi-squared and K–S tests are also displayed in Tables 4 and 5, respectively. From Table 4 it may be seen that the negative exponential function exhibits the worst performance, while the Weibull distribution has lower observed value than the gamma distribution function. Also from Table 5 it could be observed that all the three distribution functions display larger observed value compared to the K–S critical value, and the Weibull function exhibits the lowest observed value hence better matches the data. According to the results presented in Section 2, an indirect way to check the validity of the exponential density hypothesis for the joint spacings is by the virtue of the experimental distribution of number of joints measured along fixed length intervals on the drill cores, instead of the time consuming measurement of individual joint spacings. For the case of simply measuring the number of joints of all sets occurring every one meter of drill core extracted from vertical boreholes, labeled as FF and corresponding to the experimental λf values, Fig. 13a to b presents the all 868 data

Table 4 Chi-squared goodness-of-fit for measured apparent joint spacings at significance level α ¼0.01 and correlation coefficient. Joint set

Distribution function

Critical value at α ¼0.01

Observed value

Correlation coefficient

Grain

Negative Weibull Gamma Negative Weibull Gamma Negative Weibull Gamma Negative Weibull Gamma

27.6882 26.217 26.217 30.5779 29.1412 29.1412 44.3141 42.9798 42.9798 30.5779 29.1412 29.1412

43.2453 4.2194 10.5802 15.5321 2.0986 6.018 3.0168 3.3869 3.086 135.0054 3.7125 12.5709

0.98585 0.99409 0.99066 0.98905 0.99403 0.99138 0.99548 0.99552 0.99548 0.95131 0.98871 0.97697

Head-grain

Secondary (drill core)

Secondary (exposed vertical wall)

exponential

exponential

exponential

exponential

Table 5 Kolmogorov–Smirnov goodness-of-fit for measured apparent joint spacings at significance level α¼ 0.01. Joint set

Distribution function

Critical value at α ¼0.01

Observed value

Grain

Negative exponential Weibull Gamma Negative exponential Weibull Gamma

0.059876 0.059876 0.059876 0.065184 0.065184 0.065184

0.13642 0.070794 0.094703 0.10992 0.065866 0.086522

Negative exponential Weibull Gamma Negative exponential Weibull Gamma

0.24904 0.24904 0.24904 0.085993 0.085993 0.085993

0.079552 0.082099 0.08007 0.27336 0.092263 0.13967

Head-grain

Secondary (drill core)

Secondary (exposed vertical wall)

72

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

measurements performed after tilting of a large number of subpanels. It is noted that according to the practice followed in this particular quarry, there are measured volumes of blocks of volume larger than 1 m3, which means a left-truncated distribution of

joint spacings histogram 450

Data histogram Negative exponential PDF

400

gamma PDF Weibull PDF

350

180

250

160

200

140

150

120

Frequency

Frequency

joint frequencies histogram 300

100 50 0

0

2

4

6

8

Data histogram Weibull Distribution PDF gamma PDF

100 80 60

10

Spacing [m]

40 20

joint spacings CDF 1

0

0.9

0

5

10

15

20

25

30

35

Joint frequency [1/m]

0.8 0.7

CDF

0.5 0.4

joint frequencies CDF

1

neg exponential CDF gamma CDF Weibull CDF Data CDF Median Value

0.6

0.9 0.8 0.7

0.3

CDF

0.6 0.2 0.1 0

0.5 Weibull CDF gamma CDF Data CDF Median Value

0.4 0

2

4

6

8

10

joint spacing [m] Fig. 12. Histogram (a) and cumulative distribution (b) of measured secondary joint spacings on an exposed quarry wall aligned perpendicular to the strike of secondary joints and best-fitted distribution functions.

0.3 0.2 0.1 0 0

5

10

15

20

25

30

Joint frequency [1/m]

measured on borehole cores lying inside the planned final excavation boundaries along with the best-fitted Weibull cdf. In order to compare the performance of the Weibull distribution function with another candidate function, we have plotted in Fig. 13a and b the pdf and cdf curves of the gamma distribution function. The goodness-of-fit results performed by virtue of the chisquared and K–S tests are illustrated in Tables 6 and 7, respectively. As it may be observed from Table 6 both distribution functions are passing the criterion with the Weibull function better matching the data since it exhibits a lower observed value. From Table 7 it may be observed that both functions display a higher observed value than the critical one prescribed by the K–S goodness-of-fit test. However, also in this case the observed value of the Weibull distribution is lower than that of the gamma function. 4.3. Block volume distribution Having estimated the mean true joint frequencies of the three main joint sets in the quarry and having measured the volumes of the extracted marble blocks above a certain volume size at the quarry, a comparison could be made between the theoretical model expressed by Eq. (30) and the actual block volume

Fig. 13. Nonlinear regression of measured total number of joints per meter along drill cores with the Weibull and gamma distribution functions; (a) experimental and theoretical histograms and (b) experimental and theoretical cumulative distributions.

Table 6 Chi-squared goodness of fit for measured joint frequencies along the boreholes inside the final excavation boundaries (critical value at significance level α¼ 0.01 is χ2L ¼ 42.9798). Distribution

Observed value

Weibull Gamma

2.5565 2.9936

Table 7 Kolmogorov–Smirnov goodness of fit for counted joint frequencies on the drill cores inside final excavation boundaries (critical value at significance level α¼ 0.01 is DL ¼0.058577). Distribution

Observed Value

Gamma Weibull

0.17712 0.18906

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

Fig. 14. Comparison of the left-truncated actual block volume distribution (circles) with the analytical distribution function given by Eq. (30) (line).

measured volumes. Fig. 14 shows the comparison of the actual left-truncated marble block distribution with the predicted distribution given by the analytical Eqs. (28)–(30) with a mean joint frequency of the three main joint sets (i.e. the sum of the mean frequencies λ1 þ λ2 þ λ3 of the joint sets) of 7.4 1/m displayed in Table 3 and based on exploratory borehole data. The comparison shown in Fig. 14 is very good given the inherent assumptions of the theory and the complexity of the natural rock fragmentation conditions. The maximum measure block volume is 6.5 m3 whereas the theoretical model predicts slightly larger maximum block volume of around 7.5 m3. This can be explained by the fact that at the quarry the original block volume distribution is inevitably affected by the vertical and horizontal diamond wire sawing cuts. It may be also noticed that the theoretical curve is shifted to the left relatively to the experimental one. This can be explained by the inherent assumption of the model that the measured discontinuities along the drill cores are persistent whereas in reality a percentage of them correspond to intermittent joints.

5. Concluding remarks

73

found that if the joint spacings follow the negative exponential distribution, then the measured number of joints per length of drilled core follows the Weibull density function with scale and shape parameters related only to the mean frequency of joints. Following the methodology proposed originally by Hudson and Priest [10] the closed-form expression of block volume distribution has been found. Also, the distribution of block volumes has been found analytically. Furthermore, the left-truncated block volumes distribution has been found in analytical form. The theoretical results are validated against experimental data collected at a dolomitic marble quarry referring to joint spacings and frequencies sampling, as well as marble block volumes. Joint spacings data have been best-fitted by three pdfs namely the oneparameter negative exponential, and the two-parameters gamma and Weibull pdfs. As was expected in most of the cases it was found that the Weibull and gamma pdfs fit better with the data, however the much simpler negative exponential pdf was found to describe adequately the experimental data.

Appendix A. Maximum entropy theory applied in heterogeneous rock masses Methods from Information Theory [23] and entropy theory [13,14] are employed in order to derive the form of density function for the joint spacings in the case of a lack of previous information. This lack of information pertains to: (a) the influence of random factors on stress distribution inside the heterogeneous rock mass that cause high variability of the local stresses from the values which are calculated using the averaged constants of the elastic or plastic solutions, (b) the heterogeneity of rock strength, and (c) the type, succession and intensity of previous tectonic episodes responsible for the current state of fracturing of a rock mass. The density function is derived here from a condition of maximum likelihood of a given state of rock mass fracturing, which corresponds to a maximum entropy of the blocky and fractured rock mass. The condition for the maximum entropy of this process could be written in the following manner Z

1



f ðxÞln ½f ðxÞdx- max

ðA:1Þ

0

The work presented above was stimulated by the relevant series of milestone papers [9–12]. It aims at improving the approach of prediction of joint spacings and number of joints per length, RQD and block size distributions based on scanline measurements or measurements on drill cores at the exploratory phase. Regarding the prediction of block size distribution it is concerned only with discontinuity sets occurring of parallel persistent planes irrespectively of the size of joints. Finally, the results found here are validated against measurements of joints on drill cores taken from a dolomitic marble quarry. The following conclusions may be drawn from this study: The joint density function found in a theoretical manner by applying the maximum entropy theory is the one-parameter negative exponential function. It is noted that a more refined and elaborate model of joint spacings such that presented in Appendix A with appropriate mechanical constraints could be created but this is out of the scope of this paper. The relation between RQD and joint frequency found by Priest and Hudson [9] has been validated against simulation data. Aiming at inferring the mean joint frequency from measurements of number of joints per meter along drill cores instead from the more cumbersome joint spacing measurements, it has been

However, the density function must satisfy two constraints. Rx Since FðxÞ ¼ 0 f ðξÞdξ, and F(1) ¼1, the first constraint is the following well known integral equation Z

1

f ðxÞdx ¼ 1

ðA:2Þ

0

The second constraint may be derived by an energy balance applied to a given volume of the rock V of the internal specific volume energy absorbed by the rock denoted by wV (units of energy divided by the volume of strained rock) and the assumption that all the volume energy is converted into surface energy of cracks wA (units of energy divided by area). Then we assume that the expected value or mean of joint spacings is proportional to the ratio of specific energies as follows Z 0

1

w xf ðxÞdx ¼ k A wV

ðA:3Þ

where k is a proportionality constant. The Lagrangian of the system [24] as usual is the sum of the objective function we want to maximize (i.e. Eq. (A.1)), plus the

74

M. Stavropoulou / International Journal of Rock Mechanics & Mining Sciences 65 (2014) 62–74

constraints i.e (A.2) and (A.3) each multiplied by a Lagrange multiplier, i.e. Z 1 Ζ ðf ; x; λ1 ; λ2 Þ ¼  f ðxÞln ½f ðxÞdx 0 Z 1    Z 1 wA f ðxÞdx þ λ2 k  xf ðxÞdx ðA:4Þ þ λ1 1  wV 0 0 where {λ1, λ2} denote the Lagrange multipliers. From Eq. (A.4) it may be observed that the Lagrangian Z is a function of four variables. Maximizing the Lagrangian, one may obtain the following result ∂Z ¼0 ∂f

3

f ðxÞ ¼ e  ð1 þ λ1 Þ e  λ2 x ;

ðA:5Þ

Finally, inserting the above expression (A.5) for the estimated frequency into Eq. (A.2) one of the Lagrange multipliers may be eliminated as follows: Z 1 f ðxÞdx ¼ 1 3 λ2 ¼ e  ð1 þ λ1 Þ ðA:6Þ 0

Hence, f ðx Þ ¼ λ 2 e  λ 2 x

ðA:7Þ

and it may be noted that the mean frequency of joints is essentially a Lagrange multiplier. Finally, combining the above Eq. (A.7) into (A.3) the frequency parameter λ2 may be expressed as follows Z 1 wA 1 wV λ2 xe  λ2 x dx ¼ k 3 λ2 ¼ ðA:8Þ k wA wV 0 The above formula means that according to the supposed simple model the mean frequency of joints is proportional to the ratio of specific volume to specific fracture surface energies that are responsible for the rock fracturing. References [1] McNearny RL, Abel JF. Large-scale two-dimensional block caving model tests. Int J Rock Mech Min Sci 1993;30(2):93–109. [2] Laubscher DH. Cave mining – state of the art. J S Afr Inst Min Metall 1994:279–93. [3] Latham JP, Van Meulen J, Dupray S. Prediction of in-situ block size distributions with reference to armourstone for breakwaters. Eng Geol 2006;86:18–36.

[4] Palmstrom A. Measurements of and correlations between block size and rock quality designation (RQD). Tunnell Undergr Space Tech 2005;20:362–77. [5] Ellefmo SL, Eidsvik J. Local and spatial joint frequency uncertainty and its Application to rock mass characterisation. Rock Mech Rock Eng 2009;42:667–88. [6] Barton N. Some new Q-value correlations to assist in site characterisation and tunnel design. Int J Rock Mech Min Sci 2002;39:185–216. [7] Latham JP, Ping L. Development of an assessment system for the blastability of rock masses. Int J Rock Mech Min Sci 1999;36:41–55. [8] Cinco H, Samaniego F. Effect of wellbore storage and damage on the transient pressure behavior of vertically fractured well. In: SPE 52nd annual technical conference and exhibition. Denver, 9–12 October 1977, Paper SPE 6752. [9] Priest SD, Hudson JA. Discontinuity spacings in rock. Int J Rock Mech Min Sci 1976;13:135–48. [10] Hudson JA, Priest SD. Discontinuities and rock mass geometry. Int J Rock Mech Min Sci 1979;16:339–62. [11] Priest SD, Hudson JA. Estimation of discontinuity spacing and trace length using scanline surveys. Int J Rock Mech Min Sci 1981;18:183–97. [12] Hudson JA, Priest SD. Discontinuity frequency in rock masses. Int J Rock Mech Min Sci 1983;20(2):73–89. [13] Wilson AG. Entropy in urban and regional modeling. London: Pion Ltd.; 1970. [14] Rietsch E. The maximum entropy approach to inverse problems. J Geophys 1977;43:115–37. [15] Piteau DR. Geological factors significant to the stability of slopes cut in rock. In: South African institute of mining and metallurgy symposium. Planning Open Pit Mines, Johannesburg. 1970; pp. 33–53. [16] Papoulis A. Probability, random variables and stochastic processes. 3rd ed. New York: McGraw-Hill; 1991. [17] Kendall MG, Moran PAP. Geometrical probability. London: Charles Griffin & Co.; 1963. [18] Deere DU. Technical description of rock cores for engineering purposes. Rock Mech Eng Geol 1964;1:17–22. [19] Watson GN. Theory of bessel functions. 2nd ed. Cambridge: Cambridge University Press; 1958. [20] Abramowitz M, Stegun IA, editors. New York: Dover; 1972. [21] Dips6.0. Graphical & Statistical Analysis of Orientation Data, Rocscience, 〈www.rocscience.com〉; 2012. [22] 3DEC4.10. Itasca™, 〈http://www.itascacg.com〉; 2012. [23] Khinchin AI. Mathematical foundations of information theory. New York: Dover Publications; 1957. [24] Lanczos C. The variational principles of mechanics. Toronto: University of Toronto Press; 1949. [25] Annavarapu S, Kemeny J, Dessureault S. Joint spacing distributions from oriented core data. Int J Rock Mech Min Sci 2012;52:40–5.