Testing the Utilization of a Seismic Network Outside the Main Mining Facility Area for Expanding the Microseismic Monitoring Coverage in a Deep Block Caving

Abstract

In the case of mining in an inclined intrusion using the block caving method, the highest stress is usually concentrated in the seismogenic and abutment zones, especially in the front of the sloping area. In an inclined intrusion of more than 40°, the seismometer network is usually distributed in the facility area where the footwall area is also located. This causes a limitation in microseismic monitoring due to ray coverage. In this study, we conduct a seismometer deployment outside a mining facilities area with borehole seismometers. The study aims to maximize the resolution and minimize the monitoring uncertainty of underground mines. We created two scenarios of seismometer deployment: (i) seismometers are deployed following the intrusion mining level in the mining facility area; and (ii) additional seismometers are deployed in off-facilities areas. Both areas were tested for their raypath responses and sensitivity using the Checkerboard Resolution Test (CRT). The monitoring resolution influenced by the additional borehole seismometers in the off-facilities area can be quantified. The results suggest that the additional seismometers in the off-facilities areas can increase resolution by 30% in the seismogenic and abutment zones.

1. Introduction

Block caving is an underground mining method that allows the profitable extraction of massive and low-grade ore bodies [1,2,3]. The block caving mining method is becoming increasingly popular as it is the underground mining method with the lowest cost and allows large production rates. Of the ten underground mines with the largest production levels, five of them use the block caving method and have a cave depth of more than 1500 m, namely, the copper, gold and silver mines that are part of the Grasberg Operations of PT Freeport, Indonesia; the Cadia Mine of Newcrest Mining in the Cadia Valley, Australia; and the Padcal Mine of Philex Mining in the Philippines; along with the New Afton Mine, New Gold, Canada, and the diamond mine of Udachny, ALROSA, Russia [1]. The block caving method, especially at depths greater than 1500 m, faces challenges that result from the seismic-induced destruction involved [2,4,5,6]. This is due to rock mass removed, as the block caving process causes perturbations to already high stress regions in situ which may induce major earthquakes around cave areas. In such cases, seismic monitoring of block caving mining is required, especially for cave management and safety [1,3].

The main problem in underground mining, especially deep mining at a depth of more than 1500 m, such as the Hongtoushan copper mine (1600 m) in Fushun, China, the Kidd copper–zinc mine in Canada (2800 m) and the Mponeng gold mine in South Africa (4350 m), is a high level of in situ stress caused by the thick overburden [2,7]. The risk of induced seismicity due to the mining process increases as mining depth increases. In general, induced seismicity by mining processes is classified into two categories: (i) seismicity that is directly related to the stress changes caused by mining activities (blasting, hydrofracturing, caving, etc.); and (ii) seismicity that is triggered by geological structural conditions [4,8,9].

Stress accumulation in underground mines occurs in the abutment and seismogenic zones [10]. The highest stress concentrations in these zones are caused by mass extraction when using the caving method [11]. Should ideal conditions be expected to obtain with respect to the mine design, cave propagation will occur in the direction of updip intrusion. However, in the case of incline intrusion, the rock mass in the hanging wall is also induced by the adjacent cave. The copper, gold and silver mines that are part of the Grasberg Operations of PT Freeport in Indonesia and China’s Yanjingou gold mine are examples of mines with incline intrusion models with slopes of 45° and 49°, respectively [12,13].

Geomechanically, the stability challenge lies in the hanging wall area, which experiences compressional stress due to gravity on the inclined plane and the caving process. Should the induced stress be high, the rock mass can collapse, which can affect the direction of cave propagation. Uncontrolled cave propagation can result in an overbreak in the lithological boundary. This is not desirable because it will cause the mixing of high-grade ore with low-grade ore (dilution). To assess the risk of this possibility, it is necessary to monitor the seismogenic zone.

The cave design of inclined intrusive deposits conforms to the intrusion limit to optimize the recovery and grade of ore. In incline intrusion ore bodies, generally, drift access for the haulage shaft path is deployed on the footwall side by considering the stability and confidence level of the rock mass strength value. Generally, seismometer networks in underground mines are deployed in the facility area due to the ease of access regarding deployment and they can be used to monitor the stress distribution in the rock mass around the tunnel development [14]. During development activities and the preparations for mine production, the seismometers installed in the facility area remain relevant. However, once production activities start, cave propagation in block caving mines should also be monitored. Therefore, the seismometer network is not ideal for monitoring microseismic activities throughout the block caving area. The objective of this research is to maximize the resolution of monitoring and minimize monitoring uncertainty in underground mines using the tomography method. Figure 1 shows the design of a seismometer network for an underground mine in an inclined intrusion geological setting.

One of the solutions to overcome areas with poor coverage is to place borehole seismometers outside the main mining facility (off-facility) areas. Borehole seismometers in off-facilities areas have been applied in geothermal and oil and gas reservoirs [15,16,17,18,19]. Takagishi et al. [18] conducted a case study of seismic networks (type and layout of seismometers) for monitoring microseismic events and natural earthquakes around CO₂ injection areas using a (broadband) seismic network on the surface and multi-array borehole seismometers in America and Japan, along with several carbon dioxide (CO₂) injection sites. The results showed that earthquake events with a magnitude of Mw −1 could not be recorded by the seismic network on the surface. The addition of a multi-array borehole seismometer can increase the resolution of microseismic monitoring to levels previously not detectable using seismometers on the surface. Kenedi et al. [16] conducted microseismic monitoring in the Kilauea Lower East Rift Geothermal Zone, Puna, Hawaii, with a borehole seismometer at a depth of 30 to 210 m. The results of the research showed that earthquake events with a magnitude below −0.5 can be well recorded.

To improve the coverage of microseismic monitoring in the abutment and seismogenic zones that are not covered by the typical seismometer network deployed in the facility area, we simulated the deployment of additional seismometers in off-facilities areas to maximize monitoring resolution. We created two scenarios for placing seismometers: (i) a seismometer network following the intrusion mining level with seismometers spread over the mining facility area; and (ii) the seismometer network in (i) plus additional seismometers in off-facilities areas. The results of the seismometer network will be processed using the delay-time tomography method to determine the maximum resolution of the two scenarios, especially in the abutment and seismogenic zones. This study uses the delay-time tomography method for sensitivity testing using the Checkerboard Resolution Test (CRT) [20]. Through this study, the impact of adding a seismometer in off-facilities areas to monitor resolution in block caving areas can be quantified.

2. Materials and Methods

2.1. Research Method

Mining activities change subsurface structures rapidly and continuously. These structural changes include changes in pressure distribution, subsidence and fault shifts. Changes in structure cause earthquakes of small to large magnitude. The hazard of microseismic activity in mining areas is correlated with several factors: pressure, rock mechanical strength and degree of rock heterogeneity [21]. Several physical parameters have been used to analyze underground mines, such as b-value, magnitude, rock stress, P-wave and S-wave velocity, etc. [5,6,22,23]. Changes in velocity are an indication of changes in stress. Therefore, it is necessary to monitor velocity changes as proxies of the stress perturbations due to mining. One method that can be used to monitor the evolution of velocity is time-lapse passive seismic tomography using microseismic data. Microseismic methods have been widely used in underground mining, especially in longwall mining [24,25,26].

Determining the location of a microseismic event is the most important key in local tomographic studies. One of the methods that can be used to accurately locate earthquake events and describe local velocity structures is the double-difference tomography method [27]. This method can perform velocity and location inversion simultaneously [22,23,25,28]. It assumes that should the distance between two adjacent earthquake sources be much smaller than the distance to the station, then the raypaths from the two sources are to be considered almost the same.

We used the double-difference approach to assess the resolution and sensitivity of tomographic inversion in the geological setting of an inclined intrusion. This approach sharpens the results of the velocity structure near the earthquake source and produces a sharper velocity contrast than the standard tomography method [28]. Qian et al. [25] performed a time-lapse double-difference tomography analysis using differential arrival time data from pairs of events in different periods. According to Qian et al. [28], this can reduce some artifacts caused by differences in the distribution and number of earthquakes in a certain time interval. Wang et al. [23] carried out tomography imaging of P-wave velocity in an underground mining area of the Yongshaba Deposit, Guizhou Province, China. They interpreted low-velocity anomalies in their velocity model as empty space, stress release, rock breaking and cracking. Similar results were also obtained from research in an underground nickel mining area in Ontario. An earthquake with a large magnitude causes stress release and reduces the velocity of wave propagation in the area around the earthquake [29].

Sensitivity analysis with synthetic models is widely used in seismic tomography to assess the spatial resolution of the resultant solution [20]. The most commonly used type of synthetic reconstruction test is the CRT. CRT is a standard procedure at the beginning of tomography studies, in which synthetic resolution testing is carried out to determine the accuracy of the parameters used in inversion tomography [20]. The goal is to determine the confidence level of the tomogram obtained from the inversion results based on the distribution of sources and stations. Synthetic data of wave travel time are made from a velocity model, which, like a chessboard, has two patterns of higher and lower velocity. The area where an alternate pattern of high and low velocity can be recovered indicates the well-resolved area.

2.2. Synthetic Model

The location of the seismometer network deployment should consider the concept of stress distribution. The stress concentration of the block caving mining method in the ore body in an incline intrusion geological setting is similar to the stress concentration encountered in the open-stope mining method in smaller veins, where the greatest stress concentration is on the hanging-wall side [30]. This stress imbalance arises due to a combination of compressional forces due to loading and tensional stress due to caving. The stress concentration zone in the base of underground mines is called the abutment. In this zone, if the strength of the rock mass is still able to withstand the stress, the rock experiences elastic deformation (elastic zone). When the elastic limit is exceeded, the rock undergoes crack initiation, which then spreads to fracture. The crack growth process triggers the emergence of microseismic events. Therefore, this zone is called the seismogenic zone [8,10,12,31]. The zone under the seismogenic zone is the yielded zone, where rock is crushed and collapsed. The boundary between these two zones continues to propagate as cave propagation progresses (Figure 2).

The growth and shape of the cave due to mining activities using the block caving method should be monitored and controlled. Growth of the cave that does not conform to the mine design will have an impact on its security, safety and economic aspects. Microseismics is a proven reliable method for monitoring cave propagation [26,32,33]. The conceptual model of stress redistribution for homogeneous rock masses during the undercutting and caving processes has been agreed upon and proven by empirical data based on the results of microseismic monitoring [10,33]. Based on this general model, during the undercutting process, the compressed stress in the abutment and rock in the cave parameter experiences tensional stress towards the air gap at the undercut level (Figure 2).

We built a synthetic model with an incline intrusion geological setting with a slope of 49⁰ consisting of three (3) different mining levels, the first mining level at an elevation of 1820–1880 mean sea level (msl), the second at 1600–1680 msl and the third at 1500–1560 msl (Figure 1). The depth of the block caving was 2000 m. Our abutment zone was set at an elevation of 1620 m and the top cave was at an elevation of 1740 m (Figure 3). The caving mine in this study was designed with a triangular pyramidal shape, with cave heights varying from 10 m on one side and increasing gradually to a height of 120 m on the other side (Figure 3). The abutment of the cave had an area of 50,000 m². Our synthetic model represents a block caving that rises vertically; abundant microseismic events would be observed in the abutment and the top cave. Therefore, imaging the velocity perturbations in the abutment and top cave becomes crucial in monitoring cave stability.

The microseismic sources in this study are illustrated as divided into several clusters concentrated around the abutment area and seismogenic zone. This synthetic model is based on a conceptual model that has been proven by empirical data based on observations and experiences at a number of mining sites using the block caving method [10,33,34,35]. This model has also been validated by physical modeling conducted by previous researchers [32,36].

We set 1800 synthetic microseismic events to represent typical microseismic events distribution in a block caving area. The distribution of microseismic events was generally oriented to northwest–southeast (Figure 4). The data consisted of 400 events at an elevation of 1500–1600 msl (below the abutment zone), 1000 events at an elevation of 1600–1700 msl (in the abutment and seismogenic zones) and 400 microseismic events at an elevation of 1700–1800 msl (in the seismogenic zone). The seismometer network deployed in the mining facility area (blue triangle) consisted of as many as 41 stations and 14 stations distributed in the off-facilities area (green triangle) (Figure 4). Parameterization of the model is the first step to be performed in seismic tomography calculations [37]. The two most important things in model parameterization are the determination of the grid size and the initial velocity model used. The target study area, station distribution and distribution of microseismic events are the main factors to consider in determining the grid size. In areas with high station density and a high number of earthquake events, a relatively small grid size is used.

To test the spatial resolution and sensitivity, a CRT was performed. The synthetic model used in this CRT was a chessboard model with a perturbation of ±10% of the initial velocity model (Figure 5). The initial velocity model used was a 3D model with a P-wave velocity (Vp) of 5.6 km/s and an S-wave velocity (Vs) of 3.175 km/s, representing typical diorite intrusion rock velocities. The initial velocity model of the CRT was made with an alternate high and low model of 10% relative to the initial model of P- and S-wave velocities. The 3D grid node parameters used were divided into two parts: an inner node and an outer node. The inner node is a node on the inside whose size is more detailed in the research area, while the outer node is on the outside and functions as a coverage with a larger size. The parameters of the inner node were 40 m × 40 m × 40 m. The size of the outer node was made uniform, with 1.6 km and 3.5 km (the furthest extent) for the nx and ny directions and 2.5 km for the nz direction. Synthetic travel time data obtained via the forward modeling process using the CRT velocity model were calculated using the location of the hypocenter and station. To obtain absolute velocity data in this study, the forward modeling process and, for inversion, the TomoDD program created by Waldhauser (2000) and developed by Zhang and Thurber (2006) were used.

The results of the forward modeling were synthetic absolute travel times of P- and S-waves. The P and S travel times were then extracted to obtain the input phase format for pairing events using the ph2dt program. As the data used were synthetic, all the absolute travel times of P- and S-waves were recorded by all stations. To represent real seismic observational conditions, random sorting of phase travel time data was carried out with a limit of 1 earthquake event recorded by a minimum of 10 closest stations.

We calculated the direction of propagation of the wave rays and calculated the travel time using the pseudo-bending ray-tracing algorithm. The results of the ray tracing from the forward modeling gave the raypath and calculated travel time (t_cal) of wave propagation along the trajectories of P- and S-waves. The results for the ray tracing of P- and S-waves can be seen in Figure 6.

The synthetic Vp and Vs data were then grouped based on time, event and number of recording stations to be made in phase format as inputs to the ph2dt program. The ph2dt program was used to obtain paired events and the difference in travel time between events and the recording station. The distance from the earthquake to the station and the distance between the hypocenters of neighboring events should be much smaller than the distance from the earthquake to the station to obtain the proper relocation results [28]. The outputs of the ph2dt program were data from earthquakes that were used as pairs and data on differences between travel times for these pairs of earthquakes. These two files were used as inputs for the double-difference relocation method.

3. Results

3.1. Following the Facility Scenario

The distribution of the seismometer network in this scenario was tailored to the mining facilities. Seismometers were deployed at three levels of elevation of the mining facility: 15 with a distribution at elevations of 1840–1880 msl (level 1), 11 at elevation levels of 1640–1680 msl (level 2) and 15 seismometers at 1500–1560 msl (level 3). The results of the forward modeling were synthetic absolute data for P- and S-wave velocities. The data were then prepared to find the event pair using ph2dt. The ph2dt program provided outputs in the form of data from earthquake events that were used as pairs and data on differences in travel times for these earthquake pairs. These files were used as inputs for the double-difference relocation method. The results of the ph2dt program for the scenario following a mining facility with 41 stations were as follows: the total number of P- and S-wave phase pairs was 389,008, the number of P and S phase pairs used was 379,520 (97%), the number of earthquake pairs was 9489 and the average number of records per pair was 79.

The output of the ph2dt program is the catalog difference time (dtct), which is one of the input parameters in the inversion process. As explained in the section on the synthetic model, we also conducted random sorting of phase travel time data with a limit of 1 earthquake event recorded by a minimum of 10 nearest stations. The inversion input parameters were: (1) catalog difference time (dtct), (2) station coordinates, (3) velocity model, (4) sorted phases and (5) earthquake events. Other parameters that need to be determined in performing a TomoDD inversion are damping and smoothing. The smoothing parameter selected for this inversion process was 0. Determination of the damping value used depends on the condition number (CND) value, which can be seen in the output of the TomoDD program. According to Waldhauser, a stable inversion has a CND value between 40 and 80 (Waldhauser, 2001). The damping value used in this study was 350, with the number of iterations limited to 7 iterations. Figure 7 and Figure 8 are Vp and Vs CRT tomograms for the scenario for the facility for each elevation. Areas with good resolution were interpreted as being in areas that have clear alternate positive and negative anomalies.

3.2. Additional Seismometer in the Off-Facilities Area Scenario

We added 14 seismometer stations in the off-facilities area with varying distances of 50 to 300 m from the facility area. In the distribution of the station network in the off-facilities scenario, seven stations were deployed at level 1, which was at an elevation of 1840 msl, and seven stations were deployed at level 3, at an elevation of 1540 msl. As with the scenario following a mining facility, forward modeling was carried out for synthetic data from 55 stations (44 stations in the facility area and 14 stations in the off-facilities area) to obtain absolute Vp and Vs data. The results of the forward modeling of the off-facilities scenario were then processed for the event pairs along with the differences in travel times using the ph2dt program. The ph2dt results consisted of a total of 521,840 P- and S-wave phase pairs; the number of P and S phase pairs used was 379,520 (72%), the number of event pairs was 9489 and the average recording for each pair was 79. The inversion process in the off-facilities scenario was the same as in the following facility scenario.

The results of the inversion tomography analysis, with an additional 14 stations in the off-facilities area of P- and S-wave velocities, are shown in Figure 7b and Figure 8b, respectively.

4. Discussion

Stress accumulation in underground mines generally occurs in the abutment and seismogenic zones [10]. High stress concentrations in the abutment and seismogenic zones are caused by the extraction of rock mass using the caving method [11]. In the case of an inclined intrusion, the stress in the rock mass in the hanging wall is also induced by the adjacent cave. Should the induced stress be extremely high, the rock mass can collapse, which can affect the direction of cave propagation. Uncontrolled cave propagation can result in an overbreak in the lithological boundary. This is not desirable, as it will cause mixing of high-grade ore with low-grade ore (dilution). To determine the risk of this possibility, it is necessary to monitor the seismogenic zone. Monitoring of stress distribution and cave propagation in this zone is necessary to map the zones that are vulnerable to seismic activity. Therefore, it can reduce the potential risk of seismic hazard in underground mines.

In general, the results for the checkerboard inversion show that there are differences in the recovered zones and the zones not covered in each scenario where the ray coverage is relatively consistent in the center area of the model. The results of the Vp and Vs CRT inversions show the dominant resolution coverage to be in the northeast–southwest direction, following the event and station distribution (Figure 7 and Figure 8). According to the results of the Vp and Vs inversion tomograms with the scenario following the facility, the abutment zone is depicted in the area to the southeast of the caving, while the northwest side is not covered (Figure 9a and Figure 10a). The abutment zone on the northwest side is not covered due to the absence of a raypath as a consequence of the absence of a seismometer station. The raypath direction of the earthquake cluster in the northeast abutment zone becomes indistinct towards the station (northeast–southwest), which implies that the area cannot be interpreted (Figure 6a).

To gain an overview of the abutment zone in the northwest area, a simulation of the deployment of a seismometer borehole in the off-facilities area was conducted in this study. The positioning and procedure for placing seismometers in off-facilities areas have been described in the introduction. Tomograms of Vp and Vs in the scenario of adding a borehole seismometer in the off-facilities area are shown in Figure 9b and Figure 10b. Additional borehole seismometers in the northwest area were deployed at elevations of 1560 and 1740 msl to cover the abutment and seismogenic zone. The addition of this borehole seismometer resulted in a better raypath distribution (Figure 6b). Figure 6b shows that the raypath not only leads to one side (northeast–southwest), but also covers the northwest area. Figure 9 shows the difference in the area recovered for the Vp tomogram. The area outside the cave, which is the abutment zone, can be interpreted properly by adding a seismometer to the northwest area. The area that can be interpreted through the addition of the borehole seismometer network is 30% wider than the initial area. The Vs tomogram (Figure 10) shows the same thing as the Vp tomogram, but the recovered area is not as wide as the Vp tomogram; this is caused by a difference in the number of data.

The abutment zones in the northeast and northwest are areas of stress concentration due to caving and incline intrusion geological settings, according to the conceptual model of Duplancic, (2001); Brady and Brown, (2004). According to the caving model described above, in the stress concentration zone there is a seismogenic zone where seismic events occur (Figure 3). The shift of events in the seismogenic zone is an indicator of the cave propagation process. Therefore, the seismometer network coverage has an important role in the monitoring of events, especially in zones with seismic potential. As previously explained, we built a conceptual model of the cave in the form of a triangular pyramid, with the top cave on the east side at an elevation of 1740 m. The height of the cave varied from 10 m to 120 m. To view the resolution above the top cave, we drew Vp and Vs tomograms at elevations of 1780 and 1820 (Figure 11 and Figure 12). The Vp and Vs tomograms at an elevation of 1820 looked similar between the scenarios following the facilities and in the off-facilities area. A difference in resolution was seen at the elevation of 1780 on the northwest side where the coverage area with the off-facilities scenario could provide 5% better coverage area than the seismometer scenario in the facility area. The checkerboard velocity model presented in this study could identify the temporal evolution and geometry of the cave. In general, cave propagation could still be seen and monitored using a network of facilities, but as the caving process moved to the northwest, the addition of a seismometer network in the northwest area was required.

5. Conclusions

We have successfully tested the resolution of deep block caving in the geological setting of incline intrusion. The area that can be interpreted is generally located around the event area, which is surrounded by the seismometer distribution. The resolution test also provides knowledge about the limitations of the tomography data used. Therefore, it is possible to find out the minimum resolution required to detect subsurface structures. The addition of a borehole seismometer around the off-facilities area is proven to cover the area that cannot be covered by the scenario following the facility. An adequate number of station coverage data, the sampling rate, the frequency of the seismometer used, the raypath and the distance of the station to the microseismic source affect the results of the CRT. The borehole seismometer addition in the off-facilities area has been proven to improve resolution and increase the area that can be interpreted. We suggest that borehole seismometers deployed outside the main facility should be broadband and short-period 4.5 Hz seismometers with sampling rates higher than 4 KHz. In the future, borehole seismometers could be combined with Distributed Acoustic Sensing (DAS) cable technology to monitor seismic activity.