Assessment of excavation intersections’ stability in jointed rock masses using the discontinuum approach

During ore deposit development, vast networks of excavations are designed, and the volume of their intersections reaches 10%. At excavation intersections, the prediction of stress-strain state changes is complicated due to spatial geometry, the cross-coupling effect, construction sequence, etc. Mechanical properties of rocks, joint set parameters and the initial stress field also have a significant impact on stress redistribution. According to studies, up to 40% of the total number of failures in excavations occur during their intersections’ construction or reinforcing. Loads on the intersection support in accordance with existing methods are determined as for single excavations with an equivalent span and defined as the width of the larger excavation. The trend towards the intensification of mining, an increase in depth and the complexity of mining and geological conditions also complicate stress state assessment. Existing approaches need to be revised and updated for a more accurate prediction of the stress-strain state at intersections, and should consider spatial geometry, joint sets and initial stress field parameters. In this research, discrete element numerical modelling in 3DEC is done and the results are compared with existing empirical methods. Numerical models are created in a spatial setting and contain explicit representations of joints in the rock mass. Models are verified based on in-situ data, and the obtained results show a difference of up to 2 times in comparison with empirical results. This indicates that the reliability of the existing empirical methods is low, which may lead to stability loss on an intersection. Therefore, empirical methods should be updated. This can be done on the basis of numerical modelling, which shows sufficient convergence with in-situ data.


Introduction
During an excavation construction in its marginal rock mass, geomechanical processes of various nature occur. These are the processes of deformation, stress redistribution and rock fragmentation, which appear in the form of instrumentally or visually observed displacements, collapses, rock bursts, etc. (Belyakov et al., 2021;Komolov et al., 2020;Trushko et al., 2019). Insitu studies and analysis of numerical modelling results (Reаd et al., 1997;Volokhov, 2003) showed that the following groups of factors have the highest effect on the stability of a rock contour: • the initial stress-strain state field and its change patterns; • mechanical behaviour of the rock mass responding to stress fields changing; • shape and size of excavations; • type and characteristics of the support. The stress-strain state of the rock mass is formed under the influence of the initial stress fields acting in the undisturbed rock mass (Reаd et al., 1997;Singh et al., 2008;Sherpa et al., 2013;Zhang et al., 2007). The primary task before the start of construction is a comprehensive study of stress field parameters and an assessment of their changes under the influence of various factors and conditions, as well as a study of rock properties.
During excavation construction, the equilibrium state of the stress field is disrupted and a stress concentration occurs in the side walls or roof of the excavation. Depending on the initial stress state of the rock mass and geological conditions (Barton et al., 1974;Laghaei et al., 2019;Gospodarikov et al., 2019), the stress concentration coefficient can reach 2-3 (Korsakova et al., 2009), also zones with reduced stresses can be observed. As a result, a number of zones with increased and reduced stresses are formed around the contour of excavation, which may lead to stability loss (Ignatiev et al., 2021). This is specific for rock masses of the Khibiny Massif, which is studied in the current research.
Studies revealed four main joint systems within the Khibiny Massif with different locations, dip angles and directions. They intersect forming a blocky rock mass with average block sizes from 0.1x0.2 to 2x4 m. Deformation of the jointed rock mass consists of two compo-Rudarsko-geološko-naftni zbornik i autori (The Mining-Geology-Petroleum Engineering Bulletin and the authors) ©, 2021, pp. 137-147, DOI: 10.17794/rgn.2022.2.12 nents -deformation of solid rock blocks and displacement of blocks along discontinuities (Bаrton, 1977;Diederichs et al., 1999). The pattern of the displacement's distribution on the excavation contour for the case of insignificant deformations in continuous and discontinuous rock masses is similar. However, within the reduced stress zones, blocks have the possibility to move along the contacts, and values may differ several times with the possible formation of a rock failure zone (Marysyuk et al., 2020).
Statistical analysis of the reasons leading to rock failures in mines of the Khibiny Massif made it possible to establish that the greatest number of failures occur during the construction of excavation intersections and their reinforcing (see Figure 1a).
The dominating development of failures in a rock mass around intersections is not accidental and can be primarily explained by a larger span compared to a single excavation and, thus, by a larger zone of reduced stresses (Rozenbaum et al., 2011;Ogorodnikov et al., 2011). The dimensions of potential failure zones around single excavations and their intersections are inversely related to GSI (geological strength index) (Marinos et al., 2007;Mohammad et al., 2020). For any GSI value, the linear size of the potential failure zone on intersections is almost twice more than the size of a similar zone in a single excavation, as stated in Figure 1b. The unit size of the potential failure zone on Figure 1b, and later, is the ratio of the maximum linear size of the zone to excavation or intersection span. Based on the results of in-situ studies, it can be concluded that the size of the failure zone in the vicinity of an intersection can be almost 1.5-2 times higher than in single excavations.
The present methods for calculating failure zone parameters near intersections in blocky rock masses, accepted within mining enterprises, are based on solving a plane stress problem, in which the size of an intersection span is equivalent and set as the width of the larger excavation. This approach is mainly accepted because of its simplicity and calculation speed, which are important factors for extraction productivity in mines. This approach cannot provide a reliable assessment of displacements occurring in rock masses, as intersections spatial geometry and parameters of joint sets are not considered. It is necessary to solve the problem in a spatial setting and set joints in an explicit form (Gerqek, 1986).
Analytical methods are difficult for solving this type of problem, and empirical approaches are applicable only to specific mining and geological conditions. The most common approach nowadays is numerical modelling (Bobet et al., 2009;Jing et al., 2003), which includes modelling the rock mass as a continuous or discrete medium. In previous research, it has been stated (Barla G. et al., 2000) that continuum or equivalent continuum approaches provide an assessment of the size of the potential failure zone with relatively good convergence with in-situ data. However, to reflect the rock mass in the most realistic way possible, the discontinuum approach should be implemented.
For mining enterprises, it is very hard to implement spatial numerical modelling of every excavation intersection. Thus, it is necessary to establish dependencies to estimate the deformation of the marginal rock mass, which would take into account a structural disturbance and at the same time would be easy to use.
In this research, discontinuum modelling in a 3DEC software package is done to reveal the dependencies of

Methods
In this research, the behaviour of a rock mass in the vicinity of an intersection is studied with the use of the discrete element method. A rock mass is represented as a discrete medium, or as a set of blocks, with the possibility of their movement, rotation, failure, the formation of new contacts, to consider joint contact properties and predict the possibility of failure in the roof. The result of discontinuum modelling is an area in which blocks are unstable at a given stress state. Relatively free displacement of blocks allows for the determination of potential failure zone parameters, where displacements reach critical values. The stress in the blocks can be estimated, however, it cannot be interpreted as critical, since the blocks have the ability to move along contacts.
To obtain precise values, models should be created with respect to in-situ data, but if there is no joint mapping, a structure can be modelled on the basis of average block size.

Model preparation
The input data for modelling are the location of excavation, geometrical parameters (length, width, height, arch height) and the results of in-situ structural mapping (type of joints, dip angle, dip direction, intensity, length of a single joint, joint roughness, parameters of crack filler, crack width, etc.).
Based on structural mapping of excavations located at the Kirovsky Mine on levels in the range +90 -+470 m, the average values of dip angles and dip directions were determined and main joint systems were identified. Figure  2 shows circle diagrams of fracturing for excavations AW 15 on level +470 m and DHW 28 on level +285 m, where AW is airway working, and DHW is drill haulage working.
Numerical discrete models are created as a set of blocks in accordance with the main joint systems, whose parameters are given in Table 1  (1) containing the intersection of two excavations of a certain type (3), and the discrete section of the rock mass (2) presented as a number of blocks (see Figure 3b).
The size of the model -80 m in each direction, -excludes the influence of boundary parameters on the stress distribution in the vicinity of the intersection. Displacements in directions perpendicular to boundary surfaces of the model are prohibited. Dimensions of the intersecting excavations are equal to each other.
Within the mines of the Khibiny Massif, the tectonic distribution of stresses is observed at deep levels: horizontal stresses exceed vertical components (Shabarov et al., 2019;Mazurov et al., 2019;Petrov, 2007). Depending on the depth, the values of horizontal stresses change according to nonlinear laws. In this work, the stress state of the massif is studied with the ratio of the horizontal to vertical components σ h /σ v from 0.5 to 3. Vertical stress is defined as the weight of the overlying rocks.
The model is anisotropic and consider not only intrinsic anisotropy (block structure), but also induced anisot-ropy (new position of blocks as a result of deformation). The behaviour of rock blocks in the discrete element method is described with an elastoplastic model based on Mohr's strength criterion. For the elastoplastic formulation, a minimum number of parameters is required -unit weight of rock, elasticity modulus, Poisson's ratio, cohesion and angle of internal friction. Still, it describes the behaviour of rock mass quite accurately (Labuz et al., 2012). Table 2 shows the values of the physical and mechanical properties of rocks, obtained from the research reports carried out by the Geomechanics Research Center and the Saint-Petersburg Mining University.

Model verification based on in-situ data
The mines of the Khibiny Massif are characterized by complex lithological structures; fractured zones or highly jointed and oxidized rocks are often observed. Most excavations are constructed in a stable ore massif or host rocks, however, areas of highly fractured rocks are quite common. 1 -discontinuous rock mass; 2 -excavation intersection: (a) (b) Visual methods -analysis of rock fracturing, -can be used for analysing and assessing the size of a potential failure zone in the surveyed area. The direction of the main stresses' action is established according to the location of the maximum fractured zone. It is stated ( The multifactor influence of the intersection spatial geometry, roof configuration, initial stress state is presented in the form of a diagram on Figure 4. Data was obtained by measurements of failure zone sizes for different types of intersections within the Kirovsky Mine and statistically processed. Basically, excavations are stable, with the exception of the areas with oxidized zones, or stress concentration (intersection with the chamber or another excavation). However, with a high joint intensity and the predominance of vertical pressure over horizontal pressure in the intersection roof, and larger span, failure can occur. The possibility of rock block failure is provided by insufficient cohesion or lateral pressure.
The sizes of potential failure zones on Figure 4 were identified depending on fracturing intensity and joint contact parameters. The criterion for the zone outline is a crack opening over 0.1 mm. This outlines the zone with the most intense displacements of blocks. The same criterion was used to define the potential failure zone in numerical models.
Prediction accuracy strongly depends on the accepted mechanical model of rock mass and contacts. For a description of the mechanical behavior of continuous rock blocks, Mohr Coulomb criterion was chosen. It demands basic parameters but provides an adequate nonlinear deformation description. This study focuses on the contact model, because it mainly affects the possibility of blocks to form a failure zone and its shape and size.
Discontinuous models with contact description by Mohr Coulomb (Labuz et al., 2012), Barton-Bandis (Barton et al., 1974;Barton, 1982) criteria and with a rough contact surface were created. Both Mohr Coulomb's and Barton-Bandis' criteria describe contact behaviour in a way to identify blocks prone to collapse, but values when using the Mohr Coulomb criterion are overestimated. The Barton-Bandis model shows good convergence with explicit numerical modelling and in-situ data (see Figure 5). In further study, the nonlinear Barton-Bandis strength criterion is used.

Discrete element numerical modelling
Modelling was carried out in a spatial setting with explicitly specified blocks. The behaviour of the model has been studied using sensitivity analysis to consider the maximum variety of possible adjacent conditions. The rock mass is represented as blocky and there are 288 variants with the following limits of parameters change: • excavation depth: from 100 to 600 m; • excavation span: from 2 to 6 m; • roof configuration: arched or flat; • dip angle of the main joint system: from 0 to 80 degrees; • joint intensity, depending on the excavation span: from 0.5 to 10 1/m. Among the main factors influencing the stability of excavation contour are: the stressed state of the rock mass at various stages (single excavation, approaching excavations, construction of the intersection, the approach of stoping), as well as depth; rock strength and joint intensity (tectonic faults, fractured and oxidized zones, minor intrusions); the ratio of average block size to the excavation span; geometric parameters (the shape and size of the excavation, the intersection roof configuration).
Studies showed a difference in shape of failure zones depends on the type of intersection (see Figure 6). Shapes of failure zones illustrated in Figure 6 were ob-tained by an approximation of possible variants of joints systems with a change of dip angle from 0° to 80° and joint intensity from 0.5 to 10 1/m. The joint geometry for each case has been changed with constant average block size l, contours of the potential failure zone were determined by averaging the obtained results. It was found that the position of the point of maximum linear size of the failure zone changes depending on the type of intersection.

Results
Based on numerical simulation, the dependencies of the relative size of a failure zone on pointed main factors were obtained.
The geometrical parameters of an intersection which fully describe its spatial geometry are: span, the angle between the intersection excavations and the roof configuration. The intersection equivalent span should be determined geometrically in accordance with the position of the maximum size of failure zone point, while taking into consideration the curvature of intersecting excavations walls due to spalling and technological factors.
The angle between intersecting excavations affects the change in the size of the failure zone within 30% (Jinkui et al., 2015). The optimal intersection angle is close to 90°. However, if this cannot be ensured, it is necessary to design an arched shape of the intersection roof. This will reduce the size of the zone even with smooth crack walls (see Figure 7a).
The influence of the stress field on the stress-strain state of the rock mass is manifested in the localization of failure zones on the contour of excavation. In a gravitational field, reduced stress zones are formed in the roof and bottom of excavation, and in the case of a flat roof, even tensile stresses can occur. The formation of these zones leads to potential failures due to the relative dis- The shape of failure zones is similar in rock masses with both high and low joint intensity, however, the joint intensity significantly affects its size, especially in rock masses characterized by stress fields close to gravitational (see Figure 7b). Stress distribution is considered by the lateral stress state ratio, which is the horizontal tectonic stress coefficient.
It can be seen that for an intersection situated in a gravitational stress field, the effect of joint intensity on the failure zone size is higher. In stress fields with a predominance of horizontal stresses, failure zones very rarely appear, regardless of the joint intensity.
According to the methodology accepted at the Kirovsky Mine, parameters of intersection support are generally calculated on the basis of determined size of the failure zone h fz . This size is determined considering rock mass fracturing (Q index) (Barton, 1977), the geometry of an intersection with an empirical shape factor s, and the width of the larger excavation as equivalent span (Equation 1): (1) Where: B -span (width) of larger excavation (m); J r -an indicator characterizing the quality of joint contact; Q -Barton index, modified according to mining conditions; s -shape factor of excavation, for single excavation s = 1, for three-sided intersection s = 2.5, for four-sided intersection s = 3.
The methodology was tested for apatite-nepheline deposit massifs, characterized by fracture systems I-II, the geometry in the cross-section is shown in Figure 3b. Systems I and II are of interest for consideration, since the average block size in them is small enough to form a failure zone. Figure 8 shows deformation diagrams for the cross-section of the actual span of a four-sided intersection with an arched and flat roof for joint systems I and II. All excavations have a 4.6 m span, are located at a depth of 300 m in a blocky rock mass characterized by a stress field with horizontal components λ 1 = 1.5; λ 2 = 2.
The contour of the failure zone is determined by blocks, at least two (in a flat setting) or three (in a spatial setting) contact surfaces of which they are separated or the opening of the cracks is more than 0.1 mm, which is highlighted by red and orange colours. Deformation diagrams for joint system I (see Figure 8a and 8b) and for joint system II (see Figure 8c and 8d) show a difference in the nature of deformations depending on the average size of the blocks. It can be stated that the shape of the roof has a significant effect on the size of the potential failure zone: its size increases up to 4 times.
Sizes of failure zones for all intersection types were determined in a similar way and are presented in Table  3. The calculation results presented in Table 3 show that the size of zones according to current methodology can be either overestimated or underestimated, in some cases by up to 1.5 times.

Discussion
The results obtained by numerical modelling quantitatively correspond to in-situ data presented on Figure  4. This is true for the criterion of ultimate crack opening (0.1 mm), which was accepted at the stage of field research to outline the boundaries of the potential failure zone. However, this criterion can be refined in further research. In particular, current research does not consider crack filler heterogeneity at the contact of blocks. This factor may change the nature of the relative displacement and deformation of blocks. The main obstacle is difficulty to assess heterogeneity parameters at the stage of field studies.
Moreover, the research is devoted only to assessing the stability and determining the size of the failure zone in excavations without support. However, an introduction of support, for example roof bolting, would affect the parameters of a marginal rock mass and change the patterns of deformation, especially in overstressed rock masses.
The presented methodology is limited for application in rock and ore massifs, represented by medium and strongly fractured rocks (blocks), in conditions of a stress field that allows the rock units to move relative to each other. If the value of the stress state ratio is more than 3, which is defined as overstressed rock masses, this methodology is not applicable.

Conclusions
In this paper, the stress-strain state of a rock mass around intersections of excavations of different types, located within the Khibiny Massif was studied. The Khibiny Massif is characterized by complex mining and geological conditions, particularly tectonic initial stress fields and a blocky structure of the rock mass.
The conducted field studies showed that, in general, excavations in the mines of the Khibiny Massif are stable. Failure zones are observed mainly in places of stress concentration, for example, on intersections with a chamber or other excavations. The size of the failure zone in the vicinity of an intersection and single excavation can differ by 1.5-2.5 times. Therefore, special attention should be paid to the mechanical and geometric parameters of joint sets in marginal rock masses at intersections.
According to the conducted analysis, the results of this research can be summarized as follows: • It was found that the stress state has the most significant effect on the stability of a marginal jointed rock mass; with the predominance of horizontal stress components, failure in the roof may not occur. In case of gravitational stress distribution, the configuration of the roof plays a significant role, as sizes of failure zones can differ by 2 times. • Continuum mechanics (FEM) numerical simulations predict with reduced stresses, which means failure is possible in such places. Modelling within the framework of continuum mechanics with a setting of equivalent strength and deformation properties makes it possible to estimate the size of the potential failure zone and it can be used for preliminary calculations. • The numerical modelling of a blocky rock mass with specified contact properties made it possible to confirm that in fractured massifs, the most typical form of stability loss is rock block failure, the parameters of which depend on the mechanical properties of the rock and the stress state of the rock mass. A quantitative analysis of factors that have an effect on the stability of the rock mass around intersections showed that the highest impact has an initial stress state and intersection spatial configuration. Thus, with the known joint parameters, if it is necessary to predict the behaviour of a rock mass within the rock mass, it is recommended to carry out numerical modelling with explicit joint representation. However, in spatial modelling of an intersection, the amount of data can be large, which complicates calculations. In this case, it is recommended to specify explicit joints only in the marginal zone of an intersection. The dimensions of this zone can be preliminarily determined based on the results of continuous modelling.
The results can be used in updating existing methodology for assessing the stress-strain state of the rock mass in the vicinity of intersections.