Computers and Geotechnics 76 (2016) 93–104
Contents lists available at ScienceDirect
Computers and Geotechnics journal homepage: www.elsevier.com/locate/compgeo
Research Paper
A new excavation analysis method for slope design using discrete fracture network based polyhedral modelling Marc Elmouttie a,⇑, Grégoire Krähenbühl a, Ahmed Soliman b a b
CSIRO Energy Flagship, Australia Glencore, Australia
a r t i c l e
i n f o
Article history: Received 18 August 2015 Received in revised form 10 December 2015 Accepted 23 February 2016 Available online 7 March 2016 Keywords: Isoplethogram Discrete fracture network (DFN) Polyhedral modelling Progressive failure Rigid block stability analysis
a b s t r a c t A method to analyse potential excavation geometries, termed excavation isoplethogram, has previously been proposed for use in civil or mining engineering. The results are typically presented in the form of a contour plot identifying orientations of proposed excavation planes that may lead to large volumes (or masses) of unstable blocks and wedges. In this paper, a novel algorithm is proposed to utilise the probabilistic discrete fracture network approach to quantify the risk associated with various slope designs. The method is applied using open cut mine data in order to demonstrate its potential for decision making in mine design and operations. Ó 2016 Elsevier Ltd. All rights reserved.
1. Introduction Accurate understanding of the stability of an excavated or natural slope requires a thorough characterisation of the rock mass. This is normally a complex process involving acquisition, processing and interpretation of data relating to the underlying geology, hydrogeology, and properties of the rock matrix (e.g. strength, porosity) and discontinuities (e.g. permeability, persistence). Depending on these parameters, the complexity of failure mechanisms needing to be understood for proper ground control may range from soft rock failure mechanisms (circular/semi-circular) through to brittle and jointed rock failure mechanisms which are at least partially constrained and controlled by existing structural defects in the rock. For structurally controlled failure mechanisms, characterisation of the rock mass discontinuities associated with the failure modes requires measurement of their surface properties, geometry and the spatial distribution. The 3-dimensional assemblage of these discontinuities constitutes a network of intersecting structures commonly referred to as a fracture network. Modelling methods that explicitly represent the geometry of these structures, typically
⇑ Corresponding author at: CSIRO Energy, PO Box 883, Kenmore, QLD 4069, Australia. Tel.: +61 7 3327 4775; fax: +61 7 3327 4455. E-mail addresses:
[email protected] (M. Elmouttie), greg.krahenbuhl@ csiro.au (G. Krähenbühl),
[email protected] (A. Soliman). http://dx.doi.org/10.1016/j.compgeo.2016.02.014 0266-352X/Ó 2016 Elsevier Ltd. All rights reserved.
as polygons and meshes, are said to be based on the discrete fracture network (DFN) approach. The discontinuities in a rock mass typically consist of faults, joints, bedding, and other planes of weaknesses. The numbers of discontinuities may be large depending on the model size and type of analysis being conducted. The use of statistical techniques to infer the presence of unseen discontinuities based on those that have been sampled through borehole logging or field mapping has been described extensively [1,2] and methods to represent them explicitly as discrete fracture networks for computational analysis are becoming more common [3]. These DFN can consist of vast collections of thousands or more discontinuities with arbitrary properties such as orientation, size and termination with neighbours. The connectivity of the DFN is directly related to global properties of the rock mass such as its permeability and fragmentation. The latter is of particular interest for slope stability analysis of heavily jointed and blocky rock masses. The formation of blocks or polyhedra (the terms are used interchangeably throughout this paper) as a result of the connectivity of the DFN is of interest since accurate block representation in computer modelling of the rock mass improves the reliability of subsequent geotechnical analyses and, of particular interest for this paper, slope stability analysis. A stability isoplethogram is a means of interrogating the stability characteristics of the internal rock mass and was described by Windsor and Thompson [4]. The term isoplethogram was used by these authors but it is more commonly referred to as a contour
94
M. Elmouttie et al. / Computers and Geotechnics 76 (2016) 93–104
diagram. This 2-dimensional contour diagram displays the contour lines (also known as isopleths) associated with the rock mass feature being investigated as a function of excavation plane orientation (Fig. 1). The horizontal and vertical axes normally represent the dip directions and dips, respectively, of planes representative of potential excavation surfaces. Alternatively, other conventions such as strike-dip can be used. A common point is used to define the planes being considered, known as the ‘anchor point’. For the example shown in Fig. 1, we have used the crest of the model (i.e. the top edge of the front face of the model). Windsor and Thompson [4] presented a discussion of the usefulness of this analysis tool based on tetrahedral wedge analysis and only considering daylighting blocks (forming at the free surface). They also used a ubiquitous block model and commented that such an approach is preferred because of its ability to predict more instability on each excavation plane than the block theory approach. In this paper, it will be demonstrated how the DFN based polyhedral modelling approach allows application of isoplethogram analysis based on arbitrarily shaped polyhedra, in the case of both day-lighting and non-daylighting scenarios. It will also be demonstrated that the use of this approach results in more accurate predictions of stability for slope design.
2. Representation of discontinuities for stability modelling Traditionally discontinuities are represented using either equivalent continuum methods or discrete fracture network methods.
Fig. 1. (a) A 3D model of a rock mass with (b) multiple potential excavation planes (shown in blue) being assessed. In this example, planes are anchored to the centre of the crest of the model, and planes corresponding to dips spaced by 30° and dip directions spaced by 45° are shown. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
There is significant literature describing the use of DFNs to represent both the sampled discontinuities using deterministic methods as well as statistical (often termed stochastic) methods to represent the parent population of discontinuities (e.g. [2,5]). Generation of DFN realisations that are statistically equivalent to the real fracture networks requires a multi-stage approach: Acquisition of data through outcrop mapping or borehole methods must address sampling biases. Orientation set clustering and analysis must account for both small population structure sets as well as the isotropic component which is often missed by clustering algorithms. Characterisation of probability distribution functions and parameters for discontinuity sizes, orientations and other parameters. Deterministic and statistical methods to populate the DFN. The primary advantages of the DFN approach include its ability to represent phenomena which can only be understood in the context of connectivity of a subset of discontinuities (e.g. permeability and fragmentation) and that it highlights the uncertainty associated with degree of knowledge of the geometry and properties of discontinuities. Stochastic continuum methods can ‘smear’ and therefore mask this geometric uncertainty and in any case, derivation of the parameters (e.g. mean and standard deviation) required to define the probability density functions (pdf) to represent, for example, anisotropic strength, would require some form of preliminary analysis based on the discrete approach. The primary disadvantage of the DFN approach is that explicit representation of discontinuities in either finite element method (FEM) or discrete element method (DEM) numerical codes for geomechanics analysis imposes significant computational loads in the preprocessing (e.g. mesh generation or particle packing) and processing (e.g. re-meshing, contact detection) stages. If the code is capable of modelling rock mass failure and fracture propagation, these computational loads increase further. In this paper the DFN approach will be used in combination with a rigid block limit-equilibrium stability to form the structural models to be interrogated via the excavation isoplethogram
Fig. 2. 100 m cube containing a 10 m cube. Centroids are coincident and at the origin.
M. Elmouttie et al. / Computers and Geotechnics 76 (2016) 93–104
95
Fig. 3. Top-anchored isoplethogram for the geometry in Fig. 2. The red dashed line shows the theoretical curve associated with all possible planes containing both the anchor point as well as the centroid of the 10 m wide cube. Top figure represents number of failed blocks, bottom represents failed volume (in m3). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
analysis. The polyhedral modeller utilised in this paper represents an evolution of modellers outlined by [6–9] and is described in detail in [10]. In brief, a polyhedral modeller is an algorithm that determines how the polygons and meshes included in the DFN (representing the rock mass discontinuities) divide the model rock mass into distinct blocks or polyhedra. The DFN can be an arbitrary assemblage of these polygons (sometimes referred to as a polygon ‘soup’) and the algorithm will detect all vertices, edges and faces resulting from mutual intersections of the polygons or triangulated surfaces. A face-tracing algorithm is then used to detect polyhedra
within the volume of interest. The polyhedra can be concave and can be composed of an arbitrary number of facets. 3. Algorithm and complexity For the isoplethogram analysis, each potential excavation plane is used to interrogate the polyhedral model. The stability of the polyhedra is assessed and the process is repeated for the next potential excavation plane. The stability assessments used in this study were based on rigid block limit equilibrium analysis where
Fig. 4. Bottom anchored isoplethogram for the geometry in Fig. 2. Top figure represents number of failed blocks, bottom represents failed volume (in m3)
96
M. Elmouttie et al. / Computers and Geotechnics 76 (2016) 93–104
computationally prohibitive and so the following algorithm is defined: For each polyhedron intersected by the plane 1. Check for the existence of internal fractures that intersect the plane. 2. If present, determine the convex hull associated with the vertices defining the internal fractures.
Fig. 5. 100 m cube with unstable blocks at the crest.
each block is assessed in isolation under the influence of gravity and the contact surface forces. An iterative progressive failure analysis was performed. This is defined as: 1. 2. 3. 4.
Generate 3-dimensional DFN realisations. Apply polyhedral modelling to each realisation. Apply limit equilibrium stability analysis to polyhedra. For each realisation. a. Polyhedra deemed to be unstable are removed. b. Newly exposed polyhedra are assessed for stability. c. Repeat the process until no new unstable polyhedra are detected.
Each excavation plane has the potential to form new polyhedra. Regeneration of the polyhedral model for each plane would be
This is a conservative approach in that the volumes of convex hulls will over-estimate the true volumes of the newly identified polyhedra. This modelling does not account for in-situ stress in the rock mass or excavation induced stress changes and therefore rock mass failure (either yielding or fracture propagation). However, it is useful as a first-pass analysis to determine the likelihood of structurally controlled failure in heavily jointed and blocky rock mass. The uncertainty in the discontinuity properties often necessitates the use of the stochastic DFN approach where multiple DFN realisations are generated. Each represents an equi-probable realisation in what is termed a Monte Carlo simulation. By concatenating the geotechnical analyses performed using multiple realisations, confidence intervals for parameters of interest (such as block size distributions) can be derived. 4. Summary of the process The process of generating an isoplethogram for the assessment of current or proposed excavations will now be summarised. Site characterisation considering the structural geology of the setting needs to take place. As discussed previously, definition of the rock mass discontinuities included in the structural model can then occur. Consideration of the uncertainty associated with the sampling methods used is undertaken. Major structures, such as faults and bedding, and minor structures mapped on exposures are typically represented as deterministic structures in the struc-
Fig. 6. Toe anchored isoplethogram for the geometry in Fig. 5. Top figure represents number of failed blocks, bottom represents failed volume (in m3)
M. Elmouttie et al. / Computers and Geotechnics 76 (2016) 93–104
97
Table 1 Data for the open pit mine site used for the analysis.
Bench geometry
Length 70m Face Dip 70º, Strike (grid) 185º Total Wall Height 40m
Poles
Dip º
Dip Dirn º
Length (m)
139
63±8
273±5
11.4±18.9
65
53±53
270±18
2.5±7.9
3
56±20
249±8
0.5
tural model. The population of un-sampled minor structures must be represented stochastically. Careful consideration of the structural domains is required in this process. Multiple stochastic realisations of the structural model can then be generated. For each realisation of the structural model, the multiple excavation planes being considered can be assessed and a stability analysis against each plane is undertaken. In this paper, a rigid block, limit equilibrium stability analysis applied to the polyhedra detected by the polyhedral modeller, is used. However, this same process is applicable to continuum or coupled continuum–discontinuum modelling approaches. Statistics on the stability of each plane (such as volume of failed material) are accumulated. Finally, isoplethograms can be generated displaying the results for each excavation plane and the statistic of interest. 5. Illustrative cases In this section, the use of the excavation isoplethogram will be used to interrogate rock mass structure via some simple illustrative geometries (including the case of blocks confined to certain locations and interpretation of isoplethograms generated at different anchor positions). Consider a 100 m wide cube containing a solitary 10 m wide cube at its centre. The geometry would appear as shown in Fig. 2.
Fig. 7. Site used for the analysis.
The excavation isoplethogram may be generated using any point as the ‘anchor’ for each of the proposed excavation planes. The isoplethogram shown in Fig. 3 has been generated using the centre of one of the top edges of the 100 m cube as the anchor point and generating excavation planes at increments of 5° in dip and dip direction respectively. The figure shows contours which follow a curved with a minimum located at dip direction 270°. The discretisation effects associated with contours are caused by the finite sampling of possible excavation plane orientations (in this case, 5°). The contours indicate which planes intersect this block (i.e. registering a ‘failure’ on the proposed excavation plane). Due to the central location of the internal block and the square geometry of the bounding volume, the dip/dip direction of 45°/270° is detected as a plane that intersects the block as expected. What might be less intuitive is the fact that the contour curves such as to predict failures at all dips greater than this. To illustrate this, the red-dashed curve is shown which represents the theoretical case of a point-sized block located at the centre of the 100 m wide bounding volume. This curve has been generated
Fig. 8. Structural model used in the analysis.
98
M. Elmouttie et al. / Computers and Geotechnics 76 (2016) 93–104
Fig. 9. Volume of interest for the isoplethogram analysis.
by considering all the possible planes containing both the top-crest anchor point as well as the centroid of the internal block. Any excavation planes containing these points will intersect the internal block and register a single failure in the isoplethogram. Note that the more fine the angular resolution used for choosing the excavation planes, the more closely the contour plot will resemble the theoretical red curve. Note that if the centre of the bottom edge of the cube is chosen as the anchor, then due to the symmetrical geometry in this example, the same curved feature is present but with a minimum located at 090°/45° (dip direction/dip) as shown in Fig. 4. Fig. 5 shows the case of a series of unstable blocks at the crest of a slope (rockmass represented as a cube of width 100 m). Fig. 6 illustrates the corresponding isoplethogram analysis with excavation planes anchored at the toe of the bench face. The dominant peak in the isoplethogram corresponds to excavation planes with dips and dip directions intersecting the unstable blocks. The peak values in the isoplethogram shown in Fig. 6 correspond to steeply dipping planes with dip direction around 090°, consistent with the basal sliding planes of the blocks in Fig. 5. Using these exemplars for guidance, more complex geometries, as described in the following two case studies, can be interpreted. 6. Case study: Metalliferous mine This excavation analysis algorithm has been applied to a section of an open cut mine located in Queensland, Australia. The geology of the open pit site consists of four major lithological units of a general dip of approximately 65° to the west (i.e. tubular). A brief geological description of the four units in a westerly-easterly order is as follows: The first unit consists of interbedded metasediments and metavolcanics, metamorphosed volcanic sediments, quartzite and schists. The second consists of carbonaceous, siliceous shale, locally pyritic. The third consists of dolomitic siltstone, shale, laminated sequences with medium/coarse beds of siltstone. Finally, the fourth consists of laminated and graded sequences of siltstones and mudstones, dolomitic, volcanic and graphitic shales. Bed thicknesses range from 0.05 to 300 mm.
Structural data for the open pit was acquired using the photogrammetry technique. A summary of the collected data of the minor structures, in the area of interest, is given in Table 1. The second and third structure sets shown in the table show wide variance in orientation and can be considered an isotropic set in this analysis. Further discussion of potential biases in this data is included in the next section. Fig. 7 shows the portion of the open cut mine that was examined using the excavation isoplethogram analysis algorithm. The figure shows the lower benches of the north-east area of the pit. A structural model of the pit section of Fig. 7 is shown in Fig. 8. The figure illustrates the most likely stable (green1), stable with friction (yellow) and unstable (red) wedges. The polyhedral stability analysis is in good agreement with site observations having identified thin wedge failures at the surface. In order to perform an excavation isoplethogram analysis to analyse the performance of this bench geometry, a volume of interest was defined around this model assuming a theoretical single bench excavation with face angle of 90° and 40 m height. This volume and the associated DFN are shown in Fig. 9. The excavation isoplethogram analysis was conducted using this deterministic DFN. An anchor point for the excavation planes was defined at the top centre of the model (approximate coordinates of [2580E 5000N 3270RL]). The polyhedral model associated with this DFN is shown in Fig. 10a. Polyhedra in contact with the boundary surfaces of the model (namely the base, back and sides) have been excluded from the analysis as shown in Fig. 10b. The isoplethogram plot in Fig. 11 shows a comparison of the volume of removable blocks associated with each excavation plane tested in the analysis. The bright red regions, labelled ‘hotspots A and B’, indicate excavation planes with dips ranging from 50° to 75° as being particularly prone to failure. The actual bench orientation (275°/70°) corresponds to a ‘cool spot’ in the image. The hot spots are approximately located at dips of 60° and dip directions of 220° and 320°. The isoplethogram was recalculated using a finer resolution of 5° for dip and dip direction at these regions of interest to aide interpretation. The results are shown in Fig. 12. The first hot spot is relatively confined in dip direction and the second less so. 1 For interpretation of colour in Fig. 8, the reader is referred to the web version of this article.
M. Elmouttie et al. / Computers and Geotechnics 76 (2016) 93–104
99
6.1. Interpretation of analysis As discussed in the Introduction, analyses using the DFN approach expose uncertainties in the rock mass characterisation process present in all modelling approaches. Interpretation of results should therefore be made in this context. In particular, the influence of sampling biases requires consideration, particularly when the presence of a structure set that is sub-parallel is suspected as being a dominant factor in observed failure mechanisms (as is the case here). The observation in Fig. 11 that the actual bench orientation corresponds to a local minimum in terms of failure volume must therefore be interpreted with this bias in mind. As indicated earlier, the analysis presented in this paper is based on a rigid block stability analysis. The ‘progressive failure analysis’ described in this paper applies this stability analysis iteratively as each unstable polyhedra is removed from the analysis. This approach does not account for fracture propagation or rock mass failure. Therefore, if field observations confirm the likelihood that composite failure mechanisms will constrain the slope design parameters, the excavation isoplethogram analysis must be performed in conjunction with a slope design process capable of capturing these more complex failures. 6.2. Stochastic analysis
Fig. 10. Polyhedral model for the DFN in Fig. 9 showing (a) all polyhedra and (b) excluding polyhedra in contact with boundary surface.
A more detailed analysis was also performed centred on the actual bench orientation (275°/70°) of Fig. 8 and the results are shown in Fig. 13. This confirmed a local minimum in predicted total failed volume.
The prior analysis was based on a deterministic DFN. There are many uncertainties associated with rock mass characterisation and parameter variance associated with properties of the structures (persistence, orientation) can be significant. To investigate this, a stochastic variance (standard deviation of 5°) was applied to the orientations of these structures whilst retaining the structure locations (as determined from their centroids) and assumption of persistence constant. This approach has been termed a quasistochastic analysis [11] because it uses a combination of deterministic and stochastic variables to define the individual structures in the DFN. A Monte Carlo simulation was performed comprising 100 DFN realisations. The excavation isoplethogram describing failed volume averaged over 100 simulations is shown in Fig. 14. Although similar to Fig. 11, there are some notable differences. Hotspot A identified at (220°/60°) is no longer a local maximum although the one at (320°/60°) remains. The peak value in the plot is now 3500 m3, around 1000 m3 less than that predicted using the deterministic analysis. The proposed algorithm using the stochastic approach therefore has significant potential for determining the trustworthiness associated with interpretation of individual features in the analysis.
Fig. 11. Failed volume excavation isoplethogram (in m3) for the data shown in Fig. 10b.
100
M. Elmouttie et al. / Computers and Geotechnics 76 (2016) 93–104
Fig. 12. Failed volume excavation isoplethogram (in m3) using finer resolution (5°) for (a) hotspot A and (b) hotspot B identified in Fig. 11.
Fig. 13. Failed volume excavation isoplethogram (in m3) using finer resolution (5°) centred on the current bench orientation.
Fig. 14. Failed volume excavation isoplethogram (in m3) for the region of interest – quasi stochastic analysis.
7. Case study 2: Coal mine This excavation analysis has also been applied to a coal mine highwall located in Queensland, Australia. This mine has been the subject of another paper describing estimation of in-situ block size distribution and so a detailed description of the data acquisition is described in that paper [12] however a brief description follows. Data was acquired using digital photogrammetry combined with digital mapping methods. Fig. 15 shows the 3D mosaic of
the section of highwall including structural mapping and stereonet of the data. The beds at the site have undergone limited structural deformation with the dominant regional force regime being the E– W extensional force responsible for the development of the depositional basin. The 3D mosaic was constructed and digitally mapped using the Sirovision Photogrammetry System [13]. The data is shown in Table 2.
M. Elmouttie et al. / Computers and Geotechnics 76 (2016) 93–104
101
Fig. 15. Site used for the analysis, including digitally mapped structures. Colours indicate two discontinuity sets have been identified (green and red). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
Table 2 Data for the site used for the analysis.
Highwall geometry
Length 230m Face Dip 68º, Strike (grid) 234º Height ~50m
Bedding layers
Dip 5º, Dip direction 116º,
(persistent)
Elevation 454m and 475m
Poles Dip º
Dip Direction º
K
Persistence (m)
78
71±5 354±12
22 7.2±5.3
100
74±6 295±10
31 6.4±4.4
Fig. 16. An example of a weakness plane (blue) with enough exposure to be mapped. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
102
M. Elmouttie et al. / Computers and Geotechnics 76 (2016) 93–104
Table 3 Properties for the stochastic set representing the weakness planes greater than 1.5 m in extent. Poles
Dip °
Dip direction °
K
Persistence
Poles
100
5±5
115 ± 5
260
2±2
100
The structures were classified into three categories: bedding, major structures and minor structures, and these structural data were imported into the DFN generator. This generator represents the structures as planes (polygons) for use by the polyhedral modeller. Uncertainties associated with the orientations of the structures have been specified. The majority of the structures visible in this section of the wall appear confined to a domain bound by two bedding planes. This confinement was modelled by defining
a domain volume to lie between the two bedding plane structures described in Table 2. There are potentially a multitude of sedimentary layers that can act as planes of weakness should a wedge form. Their presence increases the likelihood of wedge formation but reduces the likelihood of larger block volumes. Stochastic representation of the weakness planes (Fig. 16) corresponding to the lithological boundaries was applied as described in [14]. This method relies on quantifying both the size and frequency of the weakness planes. Sampling windows were used to establish these parameters and the data are shown in Table 3. Planes greater than around 1.5 m in extent were considered. Monte-Carlo simulations were performed using 100 realisations. For each DFN realisation, a polyhedral model was generated with associated stability analysis, and an example shown in Fig. 17.
Fig. 17. (a) Simulated traces for one DFN realisation including stochastic weakness planes. (b) Formation of hazardous blocks in the model.
Fig. 18. Excavation isoplethogram for 100 realisations of the coal mine model. Dashed line indicates dip direction of the current bench face.
M. Elmouttie et al. / Computers and Geotechnics 76 (2016) 93–104
103
Fig. 19. Excavation isoplethogram using a finer resolution (5°) for planes with orientations around the current excavation plane orientation of 324/068. Dashed line indicates dip direction of the current bench face.
Fig. 20. Excavation isoplethogram using a finer resolution (5°) for (a) hotspot and (b) region identified in Fig. 18.
The colour coding indicates the results of a kinematic stability analysis with red indicating kinematically free blocks. The blocks are qualitatively similar to those observed in the field. For each polyhedral model, an isoplethogram analysis was performed and the results are shown in Fig. 18. A finer resolution excavation isoplethogram was regenerated for planes with orientations around the current excavation orientation. The result is shown in Fig. 19, and it indicates that the current bench face dip of 68 degrees is predicted to result around 15,000 m3 of failed material, less than 40% of that predicted for a vertical face (39,000 m3). Two features of interest are also shown in Fig. 20. Both correspond to more theoretical scenarios in which the excavation plane dips in the opposite direction to the current bench face. For completeness, excavation isoplethograms were recalculated using finer resolution for both features and these are shown in Fig. 20. The
hotspot is confirmed as a local maximum however the extended region shows some variability with a local hotspot at approximately 120/045. 8. Conclusions An algorithm has been defined based on a discrete fracture network – polyhedral modelling approach to extend the work conducted by Windsor and Thompson [4] to apply the excavation isoplethogram technique for slope design and operational purposes. A series of illustrative cases has been presented to provide guidance to the reader on how to interpret the analysis. This is of value given our analysis is based not just on daylighting structures but considers non-daylighting blocks and applies a progressive failure analysis for each proposed excavation plane.
104
M. Elmouttie et al. / Computers and Geotechnics 76 (2016) 93–104
Finally, the method has been demonstrated for the case of a metalliferous mine and a coal mine high wall. The presented method not only can be utilised in open pit design, but also in open pit operations in order to improve mine safety, as potential areas of rockfall and/or removable block hazards can be spatially identified. The authors are currently developing a method for application of a similar algorithm to underground excavations and this will be the subject of a future publication. Acknowledgements The authors would like to acknowledge the sponsors of the Large Open Pit Mine Slope Stability Project (LOP), which is managed by the CSIRO, as well as Australian Coal Association Research Program (ACARP) for their support in developing applications to utilise the structural modelling algorithms for the analysis of slope stability phenomena. We also thank Nathan Ferdinands and Dave Edwards for their assistance with the data acquisition for one of the sites. Finally, we thank Binzhong Zhou, Deepak Adhikary and the anonymous reviewers for their improvement of the manuscript. References [1] Baecher GB, Lanney NA, Einstein HH. Statistical description of rock properties and sampling. In: Proceedings of the 18th US symposium on rock mechanics, 5C1-8; 1978
[2] Dershowitz WS, Einstein HH. Characterizing rock joint geometry with joint system models. Rock Mech Rock Eng 1988;21(1):21–51. [3] Jing L. A review of techniques, advances and outstanding issues in numerical modelling for rock mechanics and rock engineering. Int J Rock Mech Min Sci 2003;40:283–353. [4] Windsor CR, Thompson AG. Block theory and excavation engineering. Rock technology white paper; 1996. [5] Meyer T, Einstein HH. Geologic stochastic modeling and connectivity assessment of fracture systems in the Boston area. Rock Mech Rock Eng 2002;35(1):23–44. [6] Lin D, Fairhurst C, Starfield AM. Geometrical identification of three dimensional rock block systems using topological techniques. Int J Rock Mech Min Sci Geomech Abstr 1987;24(6):331–8. [7] Jing L, Stephansson O. Topological identification of block assemblages for jointed rock masses. Int J Rock Mech Min Sci Geomech Abstr 1994;31:163–72. [8] Jing L. Block system construction for three-dimensional discrete element models of fractured rocks. Int J Rock Mech Min Sci 2000;37:645–59. [9] Lu J. Systematic identification of polyhedral rock blocks with arbitrary joints and faults. Comp Geotech 2002;29:49–72. [10] Elmouttie MK, Poropat GV, Krahenbuhl G. Polyhedral modelling of rock mass structure. Int J Rock Mech Min Sci 2010;47:544–52. [11] Elmouttie MK, Poropat GV. Quasi-stochastic analysis of uncertainty for modelling structurally controlled failures. Rock Mech Rock Eng 2012;47 (2):519–34. [12] Elmouttie MK, Poropat GV. A method to estimate in-situ block size distribution. Rock Mech Rock Eng 2012;45(3):401–7. [13] Poropat GV. New methods for mapping the structure of rock masses. In: Proceedings, Explo 2001, Hunter Valley, New South Wales, 28–31 October, 2001; 2001. p. 253–60. [14] Elmouttie MK, Krähenbühl G, Poropat GV, Kelso I. Stochastic representation of sedimentary geology. Rock Mech Rock Eng 2014;47(2):507–18.