Approximate blast movement modelling for improved grade control

Abstract

An optimization algorithm has been developed to model blast-induced displacement of rock in 3-D using gridded pre- and post-blast topography and direct blast movement measurements. The blast movement problem is expressed as an optimization assignment problem. The optimization assignment problem is solved by a heuristic algorithm offering an approximate solution in a reasonable time. The objective function and optimization details are explained. An example demonstrates the steps required to model blast movement. A realistic blast movement model improves grade control by allowing better positioning of dig lines and more correct determination of the correct destination for mined rock. The algorithm is suitable for mapping pre-blast grades, categories, expected profit, or other information onto the post-blast configuration of the muckpile.

Keywords

Open pit mining short-term mine planning geostatistics simulation optimization

Introduction

In many open pit mines, blasting is a necessary step before excavation. It is performed by means of explosives that are placed in dedicated blastholes. The chemical reaction following the initiation of explosives induces a fast release of detonation gases that create an initial impulse on the walls of the blastholes. The gases also perform the subsequent moving and breaking of the rock. Blasting engineers are usually concerned with the harmonic and reasonable distribution of blast energy throughout the bench. This design considers a number of constraints imposed by geology, the configuration of the free face, available explosives and firing systems, environmental and safety regulations, and desired degree of fragmentation among other considerations.

Blasting causes movement of the rock mass depending on the configuration of the blast. A post-blast muckpile is also characterized by a swell factor, that is, an increase in the volume after blasting. Grade control is concerned with decreasing dilution and ore loss due to incorrect classification of the mined rock. Therefore, the blast-induced displacement of rocks and the swell factor should be accounted for in grade control models for more accurate classification (Vasylchuk and Deutsch 2017).

There are two main approaches to model blast movement for grade control: (i) modelling based on the physics of rock breakage, and (ii) modelling based on direct measurements of the blast-induced displacement of rocks.

Early attempts to theoretically model blast movement were limited by computational capability. A classification of early blast movement models is presented by Gilbride (1995, p. 14–24). Cundall (1980) presents a Universal Distinct Element Code (UDEC) that allows modelling the post-blast displacement of rocks represented as discrete blocks. UDEC was later implemented in a 3-D commercial program, 3DEC (Hart et al. 1988). The method accounts for the displacement and rotation of individual blocks. Schamaun (1984) describes two models named BLOCKS and BUMP, where blast movement is represented by discrete blocks and circles, respectively; the dynamic behaviour of the discrete particles for both models is governed by the geological characteristics of mine benches, shapes and sizes of the particles, and the parameters regulating the interactions between the particles. Preece and Taylor (1989) present a Distinct Motion Code that allows incorporating the properties of explosives for modelling the motion of rocks; DMC also uses the principle of moving sphere objects. Some developments in the field of blast modelling are described in Tordoir et al. (2009). A 3-D version of the DMC algorithm is presented in Preece and Silling (2016). It utilizes parallel processing for modelling the movement of millions of discrete particles. The modern version of DMC is able to predict the influence of different initiation schemes on blast movement. Furtney et al. (2016) describe the application of a new numerical model for modelling the blasting process. The model is developed as a part of the Hybrid Stress Blast Model project implemented in software called Blo-Up (Onederra et al. 2013). The program aims to predict the detonation process, fragmentation, blast movement, and a final muckpile. Yang et al. (2017) present a new model simulating the entire process occurring in rocks during blasting including an initial impulse from shock wave propagation, the movement of rock, and fragmentation. It is reported that the simulation results are in good agreement with blasting cylinder and projectile fire tests (Yang et al. 2017).

Modelling the entire blasting process is an interesting approach for predicting blast movement. However, incomplete knowledge of the geological characteristics, fracture locations and orientation, and mechanical properties of the rocks together with uncertainty in blast parameters such as timing, stemming ejection, and the influence of free faces undermines these modelling results. La Rosa and Thornton (2011) discuss the possible economic consequences of having discrepancies between theoretically modelled blast movement and measured blast movement; for two test blasts, the error between a physics-based model and real measurements was from 1 to 7 m, which was estimated to be equivalent to a total loss of from 2.2 to 4.8 million dollars, respectively. La Rosa and Thornton (2011) argue that even a 0.5 m error may make the theoretical modelling less economically effective than actually measuring the blast-induced displacement of rock. Therefore, any theoretical blast models should be calibrated with real pre- and post-blast configurations of a mine bench.

There are two major approaches that have been used to measure displacement vectors: (i) using simple visual markers and (ii) using remote detecting equipment. Zhang (1994) and Taylor (1995) investigate the displacement of rocks during blasting using sand bags and wooden stakes as markers. Gilbride (1995) offers to use remote sensing magnetic target markers, put inside a mine bench, in a combination with magnetic radiometers. La Rosa and Thornton (2011) describe the use of radio frequency tags ID for measuring blast movement. Adam and Thornton (2004) present a new method to measure blast movement remotely using transmitters placed in a protective shell or Blast Movement Monitors. Vasylchuk (2016, p. 42) indicates that there is currently a lack of methods to reconcile pre- and post-blast 3-D models with complex configurations of free faces using sparse displacement vectors.

Taylor and Firth (2003) suggest modelling the directional components of the sparse measurement vectors separately using an interpolation method such as kriging to get displacement vectors at all locations and later model the post-blast dig limit lines. Isaaks et al. (2014) use pre- and post-blast topographic surfaces to create the pre- and post-blast 3-D models of mine benches. The directions of displacement are determined perpendicular to the timing contours of firing patterns. The horizontal displacements are drawn from the distributions of horizontal displacements based on historical data.

This paper presents a new heuristic optimization algorithm for modelling the blast movement of rock in 3-D using topographic surfaces of mine benches before and after the blast. Blast movement measurements can be used as a component of the objective function. The main elements of the algorithm are explained and it is demonstrated with an example.

Methodology

Problem formulation

The movement of rock due to blasting could be considered to be a combinatorial assignment problem. A pre-blast mine bench and a post-blast muckpile can be discretized by a set of blocks where the coordinate location of the centre of each pre- and post-blast block is known. It is assumed that the numbers of pre- and post-blast blocks are the same. A set of post-blast coordinate vectors is defined as follows: where .

A set of pre-blast coordinate vectors is defined as follows: where .

The members of , , should be assigned to the members of , , using a bijective mapping (or permutation) (Burkard and Derigs 2013). Each post-blast block should be assigned to a pre-blast block based on some cost function , where is one permutation out of the set of all possible permutations . Given that is the matrix of cost elements , the cost for every permutation is associated as follows (Burkard and Derigs 2013, p. 1): (1) The objective is to find a mapping such that the cost is minimized: (2) This problem can also be expressed as a bipartite graph in Figure 1.

Due to the nature of complex geological and grade control boundaries, it is desirable to have a high resolution. The resolution of the grade control models is recommended to be around 25% of the sample spacing (Vasylchuk and Deutsch 2017); for example, a 2 × 2 m² grid is recommended for a blasthole pattern with an average spacing of 8 m. The 3-D models of mine benches built according to these recommendations would consist of approximately 10,000–25,000 blocks.

The Kuhn–Munkres or ‘Hungarian’ algorithm (Kuhn 1955) is a popular algorithm for solving linear sum assignment problems in polynomial time if the size of the cost matrix is reasonable. Unfortunately, the current blast movement problem cannot be solved by such an algorithm because some elements of matrix change at each permutation; this is due to one element of the cost function explained below. Therefore, a heuristic algorithm is applied to find an approximate solution in a reasonable time.

Cost function

It is necessary to make assumptions about the blast movement including (i) there is no mixing and collision of particles inside the mine bench during blasting, and (ii) the swell factor is uniform throughout the entire muckpile. Although the rocks mix and collide during blasting it is currently infeasible to model blast movement with correct physics due to computational considerations and uncertainty in initial and boundary conditions. Fortunately, the direct blast movement measurements made by the systems described by Adam and Thornton (2004) allow creating 3-D blast movement vectors that can be used to constrain pre- and post-blast blocks. Topography is a primary source for creating detailed 3-D models. Gridded topographic surfaces of a mine bench before and after blasting may be created using laser scanners (Isaaks et al. 2014) or scanning drone systems.

The proposed cost function consists of four aspects: (i) the distance between the pre-and post-blast blocks, (ii) the distance between each pre- and post-blast block to neighbouring blocks within a corresponding 3-D model, (iii) distance from each pre- and post-blast block to the boundaries of a corresponding model in three orthogonal directions aligned with the free faces, and (iv) the closeness of optimized displacement vectors to real displacement vectors measured with specialized equipment.

Mathematically, the proposed cost function is expressed as follows: (3) where , i = 1,2,3,4 are components (elements) of the cost function dependent on a current permutation ; are from the initial permutation corresponding to an initial guess model; are weights corresponding to the components of the cost function .

Note that the units of all the components are similar, that is, the units of distance. Also, the initial costs are checked to ensure that they are not close to zero.

First element of the cost function

The purpose of this part of the cost function is to assign each post-blast block to the closest pre-blast block based on their corresponding x, y, and z coordinates. Even though the difference in coordinates for the associated blocks of the pre- and post-blast models is expected to be smaller at the back of the mine bench and larger closer to its free face, the general principle is the same. This idea is illustrated in Figure 2 and described mathematically as follows: (4) where is the first element of ; is the coordinate vector of a post-blast block corresponding to the current permutation .

Second element of the cost function

This element of the cost function is calculated as the average quadratic difference between the corresponding off-diagonal elements of pre- and post-blast distance matrices. This element of the cost function represents the assumption that the post-blast blocks do not significantly change their position relative to surrounding blocks if compared to the assigned blocks of the pre-blast model; it accounts for the coherency of the post-blast model. This is the part of the cost function that prohibits using an exact linear sum assignment algorithm; it is expressed mathematically as follows: (5) where is the second element of ; are the elements of the matrix of distances between all pre-blast blocks, ; are the elements of the matrix of distances between all post-blast blocks, .

Third element of the cost function

The next element of the cost function accounts for relative positions of each pre- and post-blast block with respect to the boundaries of a corresponding 3-D model. A set of vectors corresponds to each block of a pre-blast model. Each vector has distances to the boundaries of the pre-blast model as its elements. Similarly, a set of vectors , where corresponds to each block of a post-blast model. The idea behind this element of the cost function is illustrated in Figure 3.

The third element of the cost function is calculated using the following equation: (6) where is the third element of ; is the confinement coordinate vector of a post-blast block corresponding to the current permutation .

Fourth element of the cost function

The specialized equipment for measuring blast movement usually consists of some type of transmitters and detecting equipment (a receiver) (e.g., the equipment described in Adam and Thornton (2004)). The pre- and post-blast positions of the transmitters are used to calculate the directions and magnitudes of displacement for each assigned pair of pre- and post-blast blocks. Unfortunately, it is not economically viable to use many transmitters for a blast. Some researchers report using around 5 transmitters on average for 38 blasts (Yennamani 2010, p. 73–74). Due to possible measurement errors, the post-blast positions of the transmitters should be calibrated using additional information about blast movement such as additional blast movement measurements (using different measuring tools) and/or the elements of the blast design (e.g., firing sequence). The components of the measured pre- and post-blast coordinate vectors of the transmitters should be interpolated and assigned to the closest pre-blast block. After interpolation, the pre- and post-blast positions of the transmitter can be expressed as two sets of vectors and , respectively. The set of measured displacement vectors assigned to each pre-blast block, , is obtained as follows: (7) The set of optimized displacement vectors for the assigned pairs of pre- and post-blast blocks, , is calculated as follows: (8) The fourth element of the cost function accounts for minimizing the average difference between the measured and optimized displacement vectors. This is the element of the cost function that allows calibrating the theoretical blast movement modelling with real measurements without the need for modelling mixing and colliding blocks (Figure 4).

Mathematically, the fourth element of the cost function is expressed as follows: (9) where is the fourth element of ; is an ith optimized displacement vector corresponding to the current permutation .

Optimization algorithm

The primary idea of the optimization algorithm is to find a permutation that approximately minimizes the cost function using random changes to the permutation vector. Each element of the cost function is standardized and has a weight. The algorithm is first implemented on low-resolution (LR) pre- and post-blast models to reduce the number of required permutations; a fine-tuning optimization is then used on high-resolution (HR) models (one quarter of sample spacing). A detailed workflow is provided below:

Build pre- and post-blast 3-D models using topography and polygons.

Determine the grid size for the LR models; there should be at least four blocks in the vertical direction for an adequate representation of the muckpile shape. The grid size is calculated using the following equation: (10) where and are the average heights of the pre- and post-blast topographic surfaces above the bottom elevation, respectively; is the grid size for the LR 3-D models.

The LR grid size defines a new grid definition. The pre-blast gridded topography is re-gridded according to the new definition; the re-gridded topography and the bottom elevation are sufficient for creating the pre- and post-blast LR 3-D models. Due to the swell factor, the post-blast 3-D model has more blocks than the pre-blast one.

Change the grid sizes of the LR post-blast model in the x, y, and z directions such that the numbers of blocks for the LR pre- and post-blast models are the same accounting for the swell factor; an optimization algorithm from Cirio et al. (2002) is used.

After the numbers of blocks for the two models are equal, the pre-blast blocks should be assigned to the post-blast blocks using an initial bijective mapping . The objective for each pair of assigned blocks is formulated separately: (11)

After a post-blast block

, is assigned to a pre-blast block

, this pair is removed from the queue. The procedure above is implemented using a simple search based on sorting.

(5)

Calculate an initial value of the cost function for the LR model corresponding to the initial guess mapping using Equation (3); then, optimize the LR initial guess model by minimizing the cost function using an algorithm based on random restarts and permutations (described below).

(6)

In order to build the second initial guess HR model, the blocks of the optimized LR post-blast model are associated with the blocks of the HR post-blast model. The grid definitions for the HR pre- and post-blast 3-D models are obtained using a grid size provided by the user; similar to the LR models, the pre- and post-blast topography is re-gridded using the new definition and then used to create the HR 3-D models. The grid sizes in the x, y, and z directions of the HR post-blast model are changed to match the number of blocks in the HR models similar to the LR models. Each block of the HR post-blast model is associated with the closest block of the LR post-blast model; naturally, each LR block will have several HR blocks associated with it. Then, each block of the LR pre-blast model is associated with the blocks of the HR pre-blast model in a similar fashion. After Steps 1–5, each block of the LR post-blast model is assigned to a block of the LR pre-blast model; this mapping can be used to assign each block of the HR pre-blast model to a block of the HR post-blast model. It is not known exactly in which order the HR blocks from each of the pre- and post-blast pulls should be assigned to each other; therefore, they are assigned by the closeness to their respective LR blocks in Euclidian space. The main idea behind this step is illustrated in Figure 5.

(7)

Calculate an initial value of the cost function for the HR model corresponding to the second initial guess mapping using Equation (3); then, optimize the HR second initial guess model by minimizing the cost function.

(8)

Start from a random location of the post-blast model; search for a user-specified number of nearby blocks.

(9)

Switch the first starting block with a block in its neighbourhood, calculate a difference in the total cost and save this permutation if the cost is decreased; continue until all the blocks in the neighbourhood are used in a random order and go to the next start location. This allows randomly permuting the blocks of the HR post-blast model assigned to the blocks of the HR pre-blast model.

(10)

Perform the number of random restarts specified by the user.

(11)

Output the optimized model for further processing and obtaining mineable dig limits.

The numbers of restarts for the LR and HR 3-D models could be selected automatically. The algorithm iterates until the improvement to the cost function of a model is less than a tolerance parameter. The improvement to the cost function is calculated and saved as the percentage of the initial cost function value, . The number of neighbours for the LR and HR models should be selected by the user. The authors suggest using no less than 20 and 40 neighbours for the LR and HR models, respectively. The user can also specify the numbers of restarts for the two models manually.

Example

The optimization algorithm described above is implemented in Fortran according to the GSLIB conventions (Deutsch and Journel 1998). An example is developed to demonstrate the steps for blast movement modelling.

A realization of grade is simulated using sequential Gaussian simulation algorithm (Isaaks 1991) implemented in sgsim (Deutsch and Journel 1998) for an area of 70 × 110 m² with a grid size of 1 × 1 m². The reference realization is then sampled at a grid size of 8 × 8 m² with a random error of up to 1 m in coordinates. The variogram model required for the simulation has one spherical nested structure with a nugget effect of 0.2 and the major direction of continuity at 30° from North. The major and minor ranges of continuity are 100 and 35 m, accordingly. Simulated grades in Gaussian units are converted to be lognormal for more realistic results. Then, the simulated grades are converted to expected profit for three destinations using two cutoff grades; the expected profit is generated for the same area with a grid size of 2 × 2 m² to match the resolution of 3-D models.

Pre- and post-blast topographic surfaces are required to build the 3-D models. Also, pre- and post-blast polygons are required for better defining the limits of the models. The artificial pre- and post-blast topographic surfaces illustrating probable blast movement are in Figure 6. Measured displacement vectors are used for calibrating theoretical blast modelling; a set of displacement vectors was created for this example (Figure 7).

After the LR model is optimized (steps 1–5 of the algorithm), it is used to obtain an HR initial guess model (Figure 8(b)) with a corresponding permutation .

Equation (3) is used to calculate the initial value of the cost function corresponding to the permutation.

After the initial guess permutation for the HR model is obtained, the final part the optimization algorithm starts; each new permutation is accepted if it brings an improvement to the cost function. The tolerance parameter was set to 0.005% of improvement. The LR and HR neighbours were set to 20 and 40, respectively. The algorithm performed 58,000 and 100,000 iterations to minimize the cost functions for the LR and HR models, respectively. The post-blast 3-D models after 30% and 100% of the performed number of random restarts for the HR model are illustrated in Figures 9 and 10.

After the optimized post-blast expected profit in 3-D is obtained, it should be used to optimize dig limit boundaries. Figure 11 shows HR pre- and post-blast expected profit maps converted to 2-D HR classification maps. It can be seen that spatial structures of the HR pre- and post-blast classification maps are similar despite apparent differences in the shapes of the underlying 3-D models.

The current example is performed with the third element of the cost function receiving twice less weight than the rest of the elements; this is done to reduce the amount of visual artefacts that might occur due to the artificial nature of the pre- and post-blast topographic surfaces. Giving more or less weight to one of the elements does not necessarily mean it will drastically change the mapping since all the elements are related. However, it might be reasonable to give more weight to the fourth element of the cost function if many transmitters were used to measure blast movement directly. Having very good and detailed pre- and post-blast topography would justify weighting the third elements of the cost function more than the others. Figure 12 shows the graph of improvement to the cost function versus the percentage of random restarts performed.

The improvement to the cost function grows quickly up to around 30% of the specified number of random restarts and then gradually flattens out; this trend is observed for all the elements of the cost function irrespective of the amount of improvement each element brings. A case with more weight to the fourth element of the cost function would be interesting to analyse; unfortunately, only artificial data are available, which may not represent blast movement accurately.

The run-time for this example was about 35 min for 10,149 blocks on an older PC. In the majority of cases, it can be decreased to 10–15 min by reducing the number of random restarts or neighbours; this may introduce some minor visual artefacts in the 3-D models, which, nevertheless, should not significantly impact final results.

Discussion

In order to model post-blast heave, it is necessary to discretize a pre-blast volume of rocks into some units or blocks, which can be transformed or moved to define a post-blast shape. Each pre-blast discrete block can be assigned some information (e.g., grade, category, expected profit, etc.). Blast movement modelling is a procedure that informs on the post-blast positions of each pre-blast block and on the shape and volume of the post-blast heave; this information allows constructing mineable dig limits. The quality of the dig limits has a direct economic impact on a mine's operation.

Pre-blast topography is necessary to estimate a pre-blast volume. Modelling the trajectories of all blocks allows predicting the post-blast positions of all the pre-blast blocks along with the information assigned to them; post-blast topography is not necessary in this case but may be used for defining the post-blast muckpile (Isaaks et al. 2014). This approach also allows modelling the collisions and mixing of blocks. Unfortunately, the uncertainly in the geological composition and fractureness of a mine bench as well as variations in the properties and characteristics of explosives and initiation systems will likely compromise even very good models and lead to significant errors and the loss of profit (La Rosa and Thornton 2011).

The most precise information about blast movement is obtained from real blast movement measurements and detailed post-blast topography. By including measured blast movement vectors, it is possible to approximately estimate the post-blast positions of dig limits. A common approach is the translation of each pre-blast block to its post-blast position and fixing pre-blast polygons based on interpolated blast movement vectors (Thornton et al. 2005). This approach is unlikely to be accurate due to the difference of the pre- and post-blast shapes of a mine bench. Instead, mineable dig limits are obtained based on the post-blast expected profit mapped onto the 3-D model of a post-blast muckpile.

The challenging task of modelling complex processes occurring in a pre-blast bench during blasting is avoided. The trajectory of each pre-blast block is not modelled; the link between each pre-blast block and each post-blast block is predicted based on a series of assumptions and an optimization procedure. The proposed algorithm provides the information that can be used to design mineable dig limits or polygons. The practicality and effectiveness of the proposed algorithm could be based on the amount of lost profit due to misclassification; a proper comparison to existing alternatives would be reasonable.

Conclusions

An approximate optimization-based blast movement modelling workflow is developed. A four part objective function is formulated that considers geometric information from the blasted volume, blast parameters and movement vectors. The stochastic optimization algorithm and all the elements of the objective function are described. A synthetic example is constructed to demonstrate the steps required to model blast movement with the algorithm and use the output information for subsequent selection optimization. The algorithm takes advantage of the most reliable sources of information about the blast movement of rocks in open pits: (i) topography, and (ii) direct measurements.

Comparing to other blast movement software is difficult and knowing an objective ‘truth’ is all but impossible; nevertheless, additional validation and industrial experiment are required. An approximate blast movement modelling would be better than an overly simplistic translation of dig lines.

The proposed blast movement modelling approach based on a heuristic optimization algorithm offers an easy way to solve the blast movement problem. Testing the algorithm with artificial examples demonstrates its ability to map pre-blast profit or grade onto post-blast muckpile in a reasonable time and honouring real information about blast movement. The algorithm could be used independently or, better yet, integrated into a grade control system that optimizes the destination of mined material and accounts for mining constraints.

Disclosure statement

No potential conflict of interest was reported by the authors.

References

Adam

, Thornton

2004. A new technology for measuring blast movement. In: Sabin

, editor. Innovative Mineral developments Symposium; Oct 6; Sydney, Australia. Sydney, NSW: AusIMM Sydney Branch; p. 89–95.

Burkard

, Derigs

2013. Assignment and matching problems: solution methods with FORTRAN-programs (Vol. 184). Berlin: Springer-Verlag Berlin Heidelberg.

Cirio

, Lucidi

, Parasiliti

, Villani

2002. A global optimization approach for the synchronous motors design by finite element analysis. Int J Appl Electrom Mech. 16(1, 2):13–27.

Cundall

PA.

1980. UDEC – a generalised distinct element program for modelling jointed rock. Virginia Water, Surrey, England: Peter Cundall Associates. Report No.: PCAR-1-80.

Deutsch

, Journel

AG.

1998. GSLIB: Geostatistical software library and user's guide, 2nd ed. New York, NY: Oxford University Press.

Furtney

, Andrieux

, Hall

2016. Applications for numerical modeling of blast induced rock fracture. Proceedings of the 50th US Rock Mechanics/Geomechanics Symposium; June 26–29. Houston(TX): American Rock Mechanics Association.

Gilbride

LJ.

1995. Blast induced rock movement modelling for bench blasting in Nevada open-pit mines [master's thesis]. Reno (NV): University of Nevada.

Hart

, Cundall

, Lemos

1988. Formulation of a three-dimensional distinct element model—part II. Mechanical calculations for motion and interaction of a system composed of many polyhedral blocks. Int J Rock Mech Min Sci Geomech Abstr. 25(3):117–125. doi: 10.1016/0148-9062(88)92294-2

Isaaks

EH.

1991. The application of Monte Carlo methods to the analysis of spatially correlated data [dissertation]. Stanford (CA): Stanford University.

10.

Isaaks

, Barr

, Handayani

2014. Modeling blast movement for grade control. Proceedings of 9th International Mining Geology Conference; Aug 18–20: Adelaide, Australia. Melbourne: Australasian Institute of Mining and Metallurgy. p. 433–440.

11.

Kuhn

HW.

1955. The Hungarian method for the assignment problem. Nav Res Logist. 2(1–2):83–97. doi: 10.1002/nav.3800020109

12.

La Rosa

, Thornton

2011. Blast movement modelling and measurement. APCOM 2011. Proceedings of the 35th APCOM Symposium; Sep 24-30; Wollongong, Australia. Melbourne: Australasian Institute of Mining and Metallurgy. p. 297-310.

13.

Onederra

, Furtney

, Sellers

, Iverson

2013. Modelling blast induced damage from a fully coupled explosive charge. Int J Rock Mech Min Sci. 58:73–84. doi: 10.1016/j.ijrmms.2012.10.004

14.

Preece

, Silling

2016. Ore loss and dilution studies of surface mineral blasting with 3D distinct element heave models. Proceedings of the 42nd Annual Conference on Explosive and Blasting Technique, International Society of Explosives Engineers (ISEE); Jan 30–Feb 3; Las Vegas. Cleveland (OH): International Society of Explosives Engineers (ISEE).

15.

Preece

, Taylor

LM.

1989. Complete computer simulation of crater blasting including fragmentation and rock motion. Albuquerque (NM): Sandia National Labs. Report No.: SAND-88-2108C.

16.

Schamaun

JT.

1984. Methods for predicting rubble motion during blasting. Proceedings of the 25th US Symposium on Rock Mechanics (USRMS); June 25–27; Evanston (IL). American Institute of Mining, Metallurgical, and Petroleum Engineers Inc. Permission to Distribute - American Rock Mechanics Association. p. 969–978.

17.

Taylor

, Firth

2003. Utilization of blast movement measurements in grade control. In: Camisani-Calzolari

, editor. APCOM 2003. Proceedings of the 31st international symposium on application of computers and operations research in the minerals industries; May 14–16; Cape Town, South Africa. Johannesburg: The South African Institute of Mining and Metallurgy; p. 244–247.

18.

Taylor

SL.

1995. Blast induced movement and its effect on grade dilution at the Coeur Rochester mine [master's thesis]. Reno (NV): University of Nevada.

19.

Thornton

, Sprott

, Brunton

2005. Measuring blast movement to reduce ore loss and dilution. Proceedings of the 31st Annual Conference on Explosives and Blasting Technique; Feb 6–9; Orlando, FL; Cleveland (OH): International Society of Explosives Engineers. p. 1-11.

20.

Tordoir

, Weatherley

, Onederra

, Bye

2009. A new 3-D simulation framework to model blast induced rock mass displacement using physics engines. In: Sanchidrian

, editor. Fragblast 9. Proceedings of the 9th International Symposium on rock fragmentation by blasting – Fragblast 9; Aug 13–17; Granada, Spain. London: CRC Press; p. 381–388.

21.

Vasylchuk

YV.

2016. Integrated system for improved grade control in open pit mines [master's thesis]. Edmonton, AB: University of Alberta.

22.

Vasylchuk

, Deutsch

CV.

2017. Improved grade control in open pit mines. Trans Inst Min. Metall, Sect A. 127(2):84–91.

23.

Yang

, Xiang

, Chen

, Fang

, Pavlidis

, Latham

J-P

, Pain

CC.

2017. The immersed-body gas-solid interaction model for blast analysis in fractured solid media. Int J Rock Mech Min Sci. 91:119–132. doi: 10.1016/j.ijrmms.2016.10.006

24.

Yennamani

AL.

2010. Blast induced rock movement measurement for grade control at the Phoenix mine [master's thesis]. Reno (NV): University of Nevada.

25.

Zhang

1994. Rock movement due to blasting and its impact on ore grade control in Nevada open pit gold mines [master's thesis]. Reno, NV: University of Nevada.