Advances in statistical mechanics of rock masses and its engineering applications

Fqun Wu, Ji Wu, Hn Bo, Bo Li, Zhigng Shn, Dhng Kong

a Shaoxing University, Shaoxing, 312000, China

b Rock Innovation (Zhejiang) Co., Ltd., Shaoxing, 312000, China

c School of Highway, Chang’an University, Xi’an, 710064, China

d Huadong Engineering Co., Ltd., PowerChina, Hangzhou, 311122, China

e Department of Geotechnical Engineering, College of Civil Engineering, Tongji University, Shanghai, 200092, China

Keywords: Statistical mechanics of rock masses(SMRM)Jointed rock mass Geometric probability model Failure probability Anisotropic constitutive model Engineering parameters

ABSTRACT To efficiently link the continuum mechanics for rocks with the structural statistics of rock masses, a theoretical and methodological system called the statistical mechanics of rock masses (SMRM) was developed in the past three decades. In SMRM, equivalent continuum models of stress-strain relationship,strength and failure probability for jointed rock masses were established,which were based on the geometric probability models characterising the rock mass structure. This follows the statistical physics, the continuum mechanics, the fracture mechanics and the weakest link hypothesis. A general constitutive model and complete stress-strain models under compressive and shear conditions were also developed as the derivatives of the SMRM theory.An SMRM calculation system was then developed to provide fast and precise solutions for parameter estimations of rock masses,such as full-direction rock quality designation(RQD),elastic modulus,Coulomb compressive strength,rock mass quality rating,and Poisson’s ratio and shear strength. The constitutive equations involved in SMRM were integrated into a FLAC3D based numerical module to apply for engineering rock masses.It is also capable of analysing the complete deformation of rock masses and active reinforcement of engineering rock masses.Examples of engineering applications of SMRM were presented, including a rock mass at QBT hydropower station in northwestern China, a dam slope of Zongo II hydropower station in D.R. Congo, an open-pit mine in Dexing,China, an underground powerhouse of Jinping I hydropower station in southwestern China, and a typical circular tunnel in Lanzhou-Chongqing railway, China. These applications verified the reliability of the SMRM and demonstrated its applicability to broad engineering issues associated with jointed rock masses.

1. Introduction

Engineering properties and mechanical behaviours of rock masses are hot topics that have attracted various scholars in rock mechanics discipline for decades, which have significantly limited the reliability of deformation and stability calculations of engineering rock masses. Earlier rock mechanics theories rely on continuum mechanics and analysis of discontinuities,e.g.the dip rose diagram and the stereographic projection. In fact, these methods cannot comprehensively reflect the mechanical effects of the geological structure of rock mass.To tackle this problem,researches on the geometric features of rock masses have been extensively carried out since the 1970s, leading to the development of geometric probability methods for rock mass structure characterisation. Their work initiated a new stage for the mechanical quantification of jointed rock masses.

Hudson and Priest (1983) and their colleagues developed the scanline method to measure and estimate the probability distribution of the joint spacing and trace length and to obtain the characteristic parameters of joint networks. Subsequently, a sampling window method for estimating the mean trace length of discontinuities was proposed by Kulatilake and Wu (1984). After that, various scholars enriched the study of the sampling window method and developed the theory of statistically homogeneous regions (e.g. Escuder Viruete et al., 2001; Martin and Tannant,2004; Li et al., 2015; Han et al., 2016). To improve sampling workability, non-contact measurement methods such as digital photogrammetry technology(Reid and Harrison,2000;Li et al.,2016)and three-dimensional (3D) laser scanning technology (Slob et al.,2005; Assali et al., 2016) have been applied to structural statistics of rock masses in recent years. Obviously, the quality of measured structural data is a key to successful assessment of jointed rock mass properties,which demands new development of field survey techniques and methods.

Based on the measured field data, structural models of rock masses reflecting the distribution characteristics of joint spacing,trace length,opening,dip direction and dip angle were established.Significant scientific contributions have been made on this subject(e.g. Pan and Xu, 1989; Ye et al., 2012; Wang et al., 2016). These models can be used to obtain the structural parameters, and are conducive to evaluation of the mechanical properties of rock mass.Priest and Hudson (1976) and Sen (1984) proposed a relationship between joint spacing and rock quality designation (RQD) on a broad scale.Zheng et al.(2018)acquired anisotropic values of RQD on the basis of statistical structural parameters and verified the in situ results. The structural parameters such as RQD can also be related to the mechanical parameters of rock mass (Zhang and Einstein, 2004). Therefore, accurate structural information is the basis to assess the mechanical properties of rock masses. For example,the geological strength index(GSI)which is usually used for rock mass quality evaluation is based on the accurate description of the structure of rock mass (Hoek and Brown,1997).

The abovementioned work initiated two sub-directions from the geometric probability, i.e. (i) random simulation of joint networks and its geometric application, and (ii) analytical model of rock mass mechanics based on geometric probability models of the rock mass structure(Fig.1).

For the first sub-direction,since the introduction of Monte Carlo method for two-dimensional(2D)joint network simulation(Priest and Samaniego, 1983), building of joint networks has become a basis for geometric analysis of the rock mass structure. After decades of development, simulation of rock joint networks has been transformed from the original 2D model to the current 3D state.The joint networks are implemented in finite element (Pouya and Ghoreychi, 2001; Wu et al., 2019), discrete element (Ma et al.,2018; Kwok et al., 2020) and some other methods such as integrated finite/discrete element (Elmo and Stead, 2010). The joint network model is capable of calculating the structural parameters of rock masses in an intuitive manner(Zhang et al.,2013).It is also used to calculate the mechanical parameters (Bidgoli et al., 2013;Pan et al., 2019; Vazaios et al., 2019), the joint connectivity (Wang et al., 2016; Li et al., 2018), and the 3D seepage processes of rock masses (Yan and Zheng, 2017). Although random generation of joint networks is fast and convenient, it is still difficult to perform the theoretical analysis, which restricts its application to engineering problems, in particular when rigorous analytical solutions are required.

For the second sub-direction, it is necessary to establish an analytical model for rock masses, in which the geometric characteristics of all elements can be realistically reproduced and the corresponding mechanical models of rock masses can be developed.For simple and regular rock mass structures,their mechanical properties can be quantitatively analysed based on elasticity theories with orthogonal penetrating joints (Gerrard, 1982), but this method does not work for a rock mass with complex structures.The equivalent continuum constitutive models have later been proposed by Oda (1983, 1986) based on the fabric tensor and by Kawamoto et al.(1988)with a damage mechanics method.Using a non-representative volume element (NRVE) approach, the equivalent elastoplastic constitutive model of rock masses with fully persistent joints was established by Pariseau (1999). Min and Jing(2003) developed a method to determine the equivalent elastic properties of jointed rock masses with stochastic joint systems.Agharazi et al. (2012) established a 3D equivalent continuum constitutive model with three random joint sets based on superposition of deformation.Yang et al. (2016) proposed a constitutive relationship considering the displacement and stiffness parameters of the joint planes.Jiang et al.(2017)developed a joint deformation tensor that considered the effects of mechanical and geometric parameters of joints on the deformability of rock mass.Therefore,it is imperative to combine the structural and mechanical parameters during the mechanical analysis of rock mass (Bao et al., 2020a).

In light of the concepts from the geometric probability research on rock mass structures,the parameter models of the joint system were established through 16 theorems(Wu and Wang,2002), and the constitutive and strength models for jointed rock masses were proposed by Wu (1988,1992,1993) and Wu and Wang (2001a, b),resulting in the development of statistical mechanics of rock masses (SMRM). A number of internationally renowned scholars have recognised the achievements of SMRM (e.g. Jing, 2003; Zhao and Zhu, 2003; Smith, 2004; Jimenez-Rodriguez and Sitar, 2006;Jiang et al., 2014; Ferrero et al., 2016; Zuo et al., 2019; Zhou et al.,2019; Zhan et al., 2020), and the SMRM theory has then been widely used in the scientific research and engineering practice.Hu et al.(2011)modified the constitutive model in SMRM by extending the single-parameter negative exponential distributions of joint spacing and trace length into two-parameter functions. Based on the concept of joint roughness coefficient(JRC)and the mechanical properties of joints (Barton,1973; Barton and Choubey,1977), Bao et al. (2019, 2020b, c, d) studied the statistical regularity of JRC and idealised vertical joints,and applied SMRM to determination of the strength and elastic modulus for jointed rock masses.

This article introduces the latest advances in SMRM theories,methods and their engineering applications. These advances mainly include the intelligent acquisition of rock mass structures,the mechanical behaviours of jointed rock masses subjected to shearing, the establishment of complete stress-strain models, the 3D expression of key parameters, and the active reinforcement of engineering rock masses (Fig. 2). It also aims to provide a systematic introduction to SMRM and promote its application to broader engineering issues.

Fig.1. Two sub-directions initiated from geometric probability research.

2. Methodology of SMRM

A theoretical system of rock mechanics should cover essential elements such as structural model of the rock mass, stress-strain relationship, strength criterion, permeability model, and geoenvironmental conditions such as geo-stress and groundwater. A practical methodology for building this system is to transform the non-continuity in structure and mechanical behaviour into continuity.

The SMRM theory is established based on the following idealisations, as shown in Fig. 2.

Fig. 2. Framework of SMRM and its engineering applications.

2.1. Statistical physics

Statistical physics is a classical system to characterise various objects existing in nature. It obtains the macroscopic mechanical behaviour of an object based on random movements of constitutive micro-particles at the molecular level. Analogously, both the distribution and mechanical effect of rock blocks and discontinuities in a rock mass are treated as random elements, and the macroscopic mechanical performance of the rock mass is a reflection of the probabilistic behaviour of the joints and rock blocks. Hence, the methodology of statistical physics can reasonably be used for investigation of rock mechanical behaviours.

2.2. Equivalent continuum method and energy additivity

Rock mass is a discontinuous medium incised by a network of discontinuities. The mechanical effect of discontinuities in a rock mass conforms to the rules of fracture mechanics. However, fracture mechanics is complex and can only address some simple problems theoretically. It is necessary to identify a continuum mechanical approach that can equivalently reflect the mechanical effects of a large number of discontinuities. Fortunately, the continuum mechanics is a rigorous system and the common basis of many applied mechanics, which is particularly suitable for solving this issue.

On the other hand,the external forces acting on a rock mass can cause strain energy both throughout the intact rock and at the crack tip zone due to the stress concentration. The total strain energy of the rock mass can be easily obtained given that the strain energy is a scalar quantity.The constitutive model of the rock mass can then be derived following the continuum mechanics without complex superposition of the stress and strain in the discontinuity system.

2.3. Weak link hypothesis and reliability method

The weak link hypothesis has become the basis of reliability theory. This hypothesis takes a material as a chain made of a number of links, and the failure probability of the chain is determined by the weakest link; the more the number of the links, the lower the overall strength of the chain.

A rock mass is actually a complex system composed of intact rocks and numerous discontinuities. Its strength is a synergistic reflection of the two types of links and is determined by the weakest one.The weak link hypothesis acts as the ideological basis in SMRM in estimation of the strength and the permeability of jointed rock masses.

3. Theoretical models of SMRM

The SMRM system mainly consists of the geometric probability model of the rock mass structure, the stress-strain relationship,the strength and failure probability,and the hydraulics based on the structural model. Step-by-step improvements have been made via a series of integrated theoretical and application studies since the initial frame of SMRM was proposed in the early 1990s.To elucidate the SMRM, we introduced some main theoretical models and application methods in the following sub-sections (Fig. 2).

3.1. Geometric probability model of rock mass structures

Quantitative structural models have been assumed as the geometric basis for all mechanical analysis in SMRM. Herein we only introduced some geometric probability results of rock mass structures for establishment of the mechanical models.

3.1.1. Field data collection and intelligent identification of joints

Traditional surveys at rock faces, such as the scanline method(Priest and Hudson, 1976, 1981; Priest and Samaniego, 1983) and the sampling window method (Kulatilake and Wu, 1984), have been widely used for field data collection. However, these manual methods can be time-consuming, laborious and error-prone,especially in the area with complex geological conditions. These contact measurements are typically restricted to small study areas due to inaccessibility or disaster risks. In recent years, remote sensing techniques(e.g.3D laser scanning,digital photogrammetry and unmanned aerial vehicle (UAV)), which can produce highresolution 3D point clouds, are introduced to overcome these drawbacks (Assali et al., 2016; Kong et al., 2020).

A method based on 3D point clouds was developed for digital field data collection and rapid identification of joint networks(Kong et al., 2020). Fundamental joint parameters, i.e. orientation,number of sets,trace length,set spacing and density,and exact 3D location can be automatically obtained. These original data can be conveniently input into SMRM to carry out geometric probability analysis, establish mechanical models, and apply for rock mass classification system.

The workflow for digital field data collection and intelligent identification of joints is presented in Fig. 3. This method is developed based on machine learning algorithms, and mainly includes the following steps:

(1) Dominant attitude estimation. Using adaptive nearest neighbour searching and iterative reweighted plane fitting algorithms,we calculated normal vectors of every point from geo-referenced point cloud dataset. Then the dominant attitudes of outcrop surface at every point can be estimated based on the exact transformation relationship with normal vectors.

(2) Joint sets clustering.Statistical analysis of the poles(i.e.every point from joint surfaces) in stereographic projection is performed by a density-based algorithm. The principal sets and mean set orientations can be determined by identifying the peak values. The corresponding points to every set are numbered and restored separately.

(3) Joint segmentation analysis. Every single joint is isolated from its family set using density-ratio method, and the best fitting plane of each joint surface is calculated using the random sample consensus algorithm.

(4) Calculation of other parameters. The trace length is extracted as the maximum chord assuming that the area of an exposed joint plane is a convex hull.The set spacing and density can be calculated by setting a virtual scanline on the 3D point cloud model. Based on the centroid of joint cluster, the exact 3D location can be automatically obtained. A comparison of digital extraction results with those from manual measurements was made,and the error(<5%) is acceptable.

3.1.2. Probability functions and indices of geometric elements

The dominant attitude,trace length and radius,normal density,and mean opening of joints are the essential indices for description of the rock mass structure. The treatments for these indices in SMRM are briefly introduced below:

(1) Dominant attitude. The assumption of a normal distribution and the statistical method of the mean value and variance have been the typical approach for processing the orientation data of joint set, which are adopted in SMRM.

(2) Trace length and radius.The negative exponential function is used to fit the full trace length l(m)for a joint set,because it overcomes the limits of the log-normal function, which typically neglects a large number of unmeasured small joints(Wu and Wang, 2002).

Fig. 3. A schematic diagram of workflow for digital field data collection and intelligent identification of joints (modified from Kong et al. (2020)).

(3) Spacing and density.The spacing of joint set,x(m),has been affirmed to have a negative exponential distribution (Priest and Hudson, 1981; Priest and Samaniego, 1983), and the normal density,λ(m-1),is the reciprocal of the mean spacing of joint set.

The area density of joint set, λs(m-2), can be obtained by geometric probability deduction (Wu and Wang, 2002; Wu et al., 2004) and has the same form as the exponential distribution proposed by Priest and Samaniego (1983):

where l′and h（l′）are the semi-trace length of a joint and its distribution function,respectively;and μ is the reciprocal of the mean trace length of a joint set.

The bulk density λv(or bulk number of joints(m-3))of a joint set is derived in a similar fashion (Wu and Wang,2002;Wu et al., 2004):

where a and R are the radii of joint and the integrating range,respectively.

Hu et al.(2011)suggested using the two-parameter negative exponential (TPNE) function for characterising the joint spacing and trace length.

(4) The most probable maximum radius. The most probable maximum radius,am,is usually used for the strength analysis of rock mass under the weakest link hypothesis, i.e. the failure of rock mass will start from the most probable largest joint under a specific stress state.

The most probable maximum radius of joint set, am, is derived as follows (Wu and Wang, 2002; Wu et al., 2004):

where V is the volume of the rock mass studied.

3.2. Elastic stress-strain model

The elastic stress-strain model of jointed rock mass is based on the idea of an equivalent continuum that transforms the discontinuous mechanical behaviour of rock mass into continuous behaviour.This model was continuously improved after it was first introduced (Wu,1992; Wu and Wang, 2001a).

3.2.1. Theoretical basis for establishment of the stress-strain model

The mechanical transformation from a discontinuous rock mass to a continuous medium is realised based on the strain energy additivity,which is much more effective than the complex stress or strain superposition typically used in the damage mechanics methods. The strain energy density of an element of jointed rock mass, u, is calculated by

where u0is the global strain energy density of the intact rock,andis the strain energy density of all the joints.

The stress-strain model is established by calculating the strain energy density caused by the deformation of the intact rock and joints, respectively, and by combining them to obtain the global strain energy density. The elastic strain energy density of intact rocks(Eq.(5)),u0,is directly introduced from the elastic mechanics.The global strain energy density,as a dependent variable,is defined as u0= σijeij/2, where σijis the stress tensor, and eijis the strain tensor.

The strain energy density of the joints, uc, is much more complex but theoretically realisable.The calculation steps are described as follows:

(1) Calculate the normal strain energy of a joint by integrating the products of the normal stress and normal deformation,and then calculate the shear strain energy by integrating the products of the residual shear stress and shear deformation,and finally add the two parts together.This is the key step for the mechanical transformation.

(2) Add up the results of all the joints to obtain the global strain energy density of the joint system.

(3) Substitute u0and ucinto Eq.(5)and reduce the same term of the stress tensor on both sides of the formula to obtain the stress-strain relationship.

3.2.2. Stress-strain model

The stress-strain model for an element of jointed rock mass is finally written in the following form:

where σstis the element of stress tensor acting on an element of the rock mass; C0, Ccand C are the compliance coefficients for intact rock, joint and rock mass, respectively; m is the total number of joint sets;E and ν are the elastic modulus and Poisson’s ratio of the intact rock,respectively;σnand t are the effective normal stress and shear stress acting on a joint,respectively;c and φ are the cohesion and friction angle of the joint,respectively;and niis the directional cosine of the normal of a joint.

Considering that σi= ninsntσstand ti= nsδitσst-ninsntσst,Eq.(6a) can be rewritten as

Although Eq. (6a) seems complex, practical application is not difficult,and its physical meaning is clear.The mechanical effects of the most important geological factors have been rationally reflected in the model:

(1) The strain of jointed rock mass, eij, is larger than that of an intact rock, e0ij, under a specific stress state because of the additional strain of the joint system, ecij= Ccijstσst. This is similar to the equivalent strain in damage mechanics but has a much more reliable mechanical basis.

(2) The strain of the joint system is proportional to the number of joint sets,m,the normal density,λ,and the mean radius,of joint set.Therefore,the more fragmented a rock mass,the greater its compliance, and the less its stiffness.

(3) The additional strain varies with the direction of the normal vector of the joint, which indicates the mechanical anisotropy of the rock mass.

(4) The factors k and h reflect the impacts of the stress state on the joints and the shear strength of the joints in the rock mass, respectively.

The procedure in practical calculation of Eq. (6a) is as follows:(1)calculate the normal and shear stresses,σnand t, acting on the joints and decide the stress state coefficients k and h;(2)calculate the compliance coefficients C0, Ccand C, and then calculate the elements of the strain tensor,eij,corresponding to the stress tensor,σst.

The stress-strain relationship (Eq. (6a)) is an elastic model.Therefore, all the results calculated based on Eq. (6a) are in the elastic state. The failure occurring in the rock mass needs to be calculated with proper strength criterion that will be discussed later.

Obviously, Eq. (6) is an equivalent continuum stress-strain model for an arbitrary element of the jointed rock mass. With the change of the geological feature of the rock mass, the parameters involved in the model will change consequently. By this way, this model can describe the mechanical behaviours of heterogeneous media. Hu et al. (2011) extended the model based on the TPNE distribution of joint spacing and trace length, and intended to improve the accuracy of the data input of the model.

3.3. Strength and constitutive model

The basis of the strength theory for jointed rock masses is the weakest link hypothesis proposed by Weibull in 1939 for brittle materials, which implies that the strength of the rock mass is determined by the weakest element among the intact rocks and all the discontinuities. Additionally, due to the randomness of rock strength, the strength theory has been expressed by not only a failure criterion but also a failure probability (Wu,1992; Wu and Wang, 2001b).

3.3.1. General models for strength and failure probability of rock masses

The strength of the intact rock is defined by the conventional Mohr-Coulomb model, while the strength criterion for the discontinuities is established based on the strain energy release rate of a crack. For a set of joints, the strength is probably determined by the largest crack,which is likely to be the weakest element.Hence,the strength criterion of a set of joints can be derived from the criterion of strain energy release rate in fracture mechanics,written as (Wu,1993):

where KIcis the fracture toughness of mode I crack,and the fracture toughnesses KIIcand KIIIcfor modes II and III cracks have been represented by KIcbased on the criterion of strain energy release rate.

The strength of jointed rock mass is determined by the lowest among the strengths of the intact rock and each joint set. The failure probability of rock mass is written as follows according to the reliability theory (Wu,1993):

where Pband Pcare the failure probabilities of the intact rock and the discontinuity system, respectively; σcmis the unconfined compressive strength (UCS) of rock; k0and m are the fitting coefficients; and acis the critical radius of joint set.

3.3.2. Constitutive model of jointed rock mass

A constitutive model consists of a stress-strain model and a strength model (yield surface). Based on Eq. (6a), the constitutive model of jointed rock mass can be written as (Wu and Wang,2001a):

Apparently, the constitutive model (Eq. (10a)) describes small deformation of rock masses only because it does not consider the crack extension process and large deformation of the rock mass structure.

3.4. Complete stress-strain models

A complete compressive stress-strain model and a complete shear stress-strain model were developed as the extension of Eq.(10a) to estimate the deformation of jointed rock masses.

3.4.1. Complete compressive stress-strain model

Complete rock mass deformation includes two stages, i.e. preand post-peak stages, reflecting a synergistic behaviour of the intact rock and joints following the weakest link hypothesis. The pre-peak stage generally reflects the compaction, elastic deformation and gradual failure of joints;and the post-peak stage is mainly controlled by the sliding/shearing of pre-existing or newlygenerated joints after the failure of intact rocks.

Under the principal stress state, the strain, e11, in the loading(σ11)direction,can be expressed by Eq.(11a).The first formula is an elastic deformation model derived from Eq. (10a) that does not consider the failure of any elements of the rock mass. The second formula is the critical condition for failure of each set of joints and intact rocks, which has been introduced by Eq. (10a).

Several representative complete deformation curves are shown in Fig.4.With the increase of the loading stress σ11under a constant confining stress σ33,more joint sets will fail(be unlocked),and the slope of the stress-strain curve gradually decreases. The peak of the curve is the critical point for the failure of intact rocks.A notable stress drop can be observed at this point as the UCS of rock R0approaches 0,and the strength will then reduce to its residual value controlled by joints only.

3.4.2. Complete shear stress-strain model

As an example of all strain components,the shear strain e13can be written as

where σ and t are the normal stress and shear stress acting on a joint, or acting on the shear failure plane for intact rocks,respectively.

The second formula in Eq.(13)can be written in a specific form as Eq. (14) for any joint sets and intact rocks (Coulomb shear strength criterion) respectively. The first formula of Eq. (14) for a joint set can be derived from Eq. (8).

where f0and c0are the friction coefficient and cohesion of intact rocks, respectively; and k = 1 is for tensional normal stress σ.

The model in Eq. (14) characterises the pre-peak behaviour following elastic mechanics. With the failure of joints (taking KIc= 0 for failed joint sets), in the similar form of the curves in Fig. 4, the slope of the shear stress-strain curve will gradually decrease, and the failure of the intact rock (c0→0) will make the curve bend down from the peak. The post-peak strength was presumed to be identical to the strength of the rock joints.

Fig. 4. Complete deformation curves of jointed rock masses.

4. Key parameter calculation of rock masses based on SMRM

An SMRM calculation platform has been explored for a series of commonly used parameter calculation, rock mass quality rating,and analyses for specialised problems of rock engineering such as slope, underground space, foundation, and rockburst. Several useful tools/methods developed from the SMRM theories are presented below, which are helpful for estimation of commonly used parameters and establishment of related analytical models.

4.1. Structural parameter of rock masses, 3D RQD

The most probable maximum radius, am, bulk density, λv, and RQD are essential parameters defining the structure of rock masses.Particularly, RQD has been a commonly index used for rock mass structure and rock quality evaluations via borehole data since it was proposed by Deere (1964).

A model for RQD calculation along any scanline is deduced by geometric probability for an exponential distribution of intersecting points (Sen,1984):

Actually, RQD can be simply calculated through the geometric probability (Wu and Wang, 2002; Wu et al., 2004):

where the numerator is the sum of joint spacings with x > 0.1 m;and the denominator is the total length of the scanline,L.Then RQD for joint spacing following the negative exponential function is

Clearly, when the length of the scanline L→+∞, Eq. (17)becomes

Eq. (18) is similar to Eq. (15), with a difference of［e-0.1λ - （1 -0.1λ）e-0.1λ］/e-0.1λ= 0.1λ.

Since RQD can be calculated in any appointed scanline direction of joint network,a 3D distribution can be obtained and drawn as a stereographic projection in the following way:

(1) Obtain the global density, λ, by projecting the normal densities of all joint sets on the scanline direction,and the value of RQD can be calculated by Eq. (17) or (18);

(2) The value of RQD can be projected on the equatorial plane of the stereographic projection sphere following the Schmidt regulation. As shown in Fig. 5, any vector from the origin of the coordinate system will intersect the upper surface of the projection sphere. A ray from the intersecting point to the lower pole of the sphere will pass the point at the equatorial plane,i.e.the projection point.The value of RQD is drawn at the projected point in the diagram with a graded colour.

(3) By changing the dip direction and dip angle of the vector at a certain interval (1°or 2°), a 3D stereographic projection of RQD can be obtained (Fig. 6).

To illustrate the directional distribution of RQD through stereographic projection,a set of data shown in Table 1 is collected via scanline method from eight adits in the dam site of QBT hydropower station in north Xinjiang Uygur Autonomous Region,China.The lithology in the region is mainly the thick layered biotite quartz schist with four sets of discontinuities. The stereographic projection diagram is shown in Fig. 6a.

Besides, every single point in the diagram is corresponding to the normal vector of a plane as shown by the arc, and a computer programme has been developed to obtain the values of RQD along the arc and the resulting rose curve is shown in Fig. 6b.

To reflect the structural anisotropy of rock mass,the anisotropy index of RQD is defined as

Fig. 5. Principle of stereographic projection.

where RQDmaxand RQDminare the maximum and minimum values of RQD, respectively.

The representative elementary volume (REV) of RQD is a common topic and can be investigated through the variation tendency of the relative error of RQD. The relative error of RQD with the scanline length, L, can be defined based on Eq. (17) as

Fig. 7 shows the variation in the RQD relative error with the scanline length. The error decreases rapidly with the increase in scanline length.As an example,5 m can be taken as the scale of REV for RQD when λ= 1 m-1and ε = 1%.

4.2. Deformability parameters of rock masses

Deformability parameters of rock mass such as the elastic modulus,Em,and Poisson’s ratio,ν,are commonly used as indices in geoengineering for rock mass quality evaluation.However,they are generally taken as isotropic parameters based on the elastic mechanics, which cannot reflect the directional variation of anisotropic media. Here, these two parameters are re-defined on the basis of the anisotropic stress-strain model, i.e. Eq. (6a).

Fig. 6. 3D distribution of RQD. Stereographic projection (a) and the sectional distribution (b) of RQD correspond to the black point in the stereographic projection.

Table 1 Data of joints collected in adits.

4.2.1. 3D elastic modulus

The elastic modulus of rock mass,Em,is usually obtained from in situ tests, experiential inference and back analysis. The superposition model (Wu, 1988) and the equivalent model of damage mechanics have also been used.The elastic modulus and Poisson’s ratio of rock mass have been proposed as important applications of the SMRM constitutive model (Wu et al., 2004). Scholars have applied them in rock engineering practice and made some improvements in the past years (e.g. Hu et al., 2011; Li et al., 2018).

Fig. 7. Curves of RQD relative error ε vs. the scanline length L.

The elastic modulus in SMRM,Em,is defined from Eq.(6a)by the ratio of the load, σ11, to the strain, e11, in the loading direction,following the theory of elasticity:

In the practical applications,the modulus is usually divided into elastic modulus and deformation modulus due to the change in the slope of the compression curve. The modulus in Eq. (21) presents the former.

As an example,Emis calculated with Eq.(21)based on the joint data in Table 1 and the rock properties and in situ stress data in Tables 2 and 3. These data are provided by China Beifang Investigation, Design & Research Co., Ltd., which is responsible for the geological survey of the QBT project. The stereographic projection and sectional distribution of Emshown in Fig. 8a and b can be obtained in the same way as RQD. The calculated curves of Emvs.confining pressure are shown in Fig. 9.

The discontinuity-induced weakening and anisotropy of the elastic modulus are two important characteristics of jointed rock mass. It is clear from Eq. (21) that E > Emdue to the influence of discontinuities.Here,a reduction index and an anisotropy index of modulus are proposed to quantify these characteristics as follows:

where ζEand ξEare the reduction index and anisotropy index of rock masses,respectively;and Emmaxand Emminare the maximum and minimum values of Emin the full-direction space,respectively.

The calculated results shown in Fig.9 are compared with the in situ testing data conducted by China Beifang Investigation, Design&Research Co.,Ltd.,as shown in Table 4.The in situ tests are done by φ50.5 cm circular plates in horizontal and vertical directions,with the maximum normal load of 10 MPa. Accordingly, the calculated data are selected from the curves with the confining pressure of 0 MPa and normal loads of 3.5 MPa, 7 MPa and 10.5 MPa.The calculated results of Emare in reasonable agreement with the tested results,indicating that the theoretical calculation is reliable.

Table 2 Mechanical parameters of intact rocks.

Table 3 Data of in situ stresses.

4.2.2. Two directional Poisson’s ratio and higher/lower Poisson’s effect

The Poisson’s ratio also differs with directions. The Poisson’s ratios in the directions x2and x3when the loading is acting on the direction x1,i.e.ν21and ν31,are written as follows based on Eq.(6a)and the theory of elasticity (Wu,1993):

where C1111, C2211and C3311are the compliance coefficients derived from Eq. (6b):

Considering a tensile stress state,we have k=1 and h=1,and a compressive stress state, k = 0 and 1 > h > 0, and interesting phenomena called “higher Poisson’s ratio” and “l(fā)ower Poisson’s ratio”of jointed rock mass occur,respectively.Eqs.(23a)and(23b)indicate that under a tensile stress(σ>0),we have ν21<ν,whereas under a compressive stress(σ<0),ν21>ν,and the value of ν21could exceed 0.5. The Poisson’s ratio curves calculated in two directions using the data in Tables 1-3 are drawn in Fig. 10. These curves reveal the significant directional dependence of Poisson’s ratio controlled by the rock mass structure.

4.3. Strength parameters of jointed rock masses

4.3.1. Coulomb compressive strength of jointed rock masses

A preliminary model for the compressive strength of jointed rock mass was proposed by Wu and Wang (2001b). Upon this model, other researchers identified four types of failure models,and established the corresponding strength formulae.

Actually, the Coulomb compressive strength criteria for both intact rock and joints can be written in a unified form as Eq. (24)(Wu and Wang, 2001a). Letting the subscript i = 0 refers to the intact rock and i=1,2,…,m for different joint sets,and considering the weak link hypothesis, the Coulomb compressive strength of rock mass is expressed as

Fig. 8. Stereographic projection (a) and sectional distribution (b) of Em.

Fig. 9. Curves of Em vs. confining pressure.

Table 4 Comparison between the calculated results and the in situ tested data of Em.

where δ is the angle between the normal load σ1acting on the shear plane and the normal direction of the joint.

Following the same treatment for RQD and Em, the Coulomb compressive strength of rock mass can also be drawn as a stereographic projection and a rose diagram, as shown in Fig. 11. Data tabulated in Tables 1-3 were used in the calculation.

The calculated results shown in Fig.11 are compared with the in situ tested data provided by China Beifang Investigation,Design& Research Co., Ltd., as tabulated in Table 5. In situ direct shear tests were conducted in the direction parallel to the downstream of the river, and friction tests were carried out immediately after the direct shear tests. The size of the sample is about 50 cm(length) × 50 cm (width) × 35 cm (height). The frictional coefficient and cohesion of the rock mass and shear plane were obtained.In addition,the UCS of the rock mass calculated according to the direct shear test(σs)and the friction test(σf)were obtained.The actual UCS of the rock mass(σ1)should be between σsand σf,as proved by the theoretical value shown in Table 5, which indicates that the strength of rock mass calculated by SMRM is reasonable.

The notable distinction with the Hoek-Brown criterion is the strength anisotropy influenced by the orientation of joints,which is reflected by Eq. (24). The indices for the strength reduction and anisotropy can also be defined to quantify the strength and its variation as

where σ0is the basic strength of the rock.

The curve of Eq.(24)is calculated along the loading direction,as shown in Fig.12.Notably,for the third group of formulae in Eq.(24),the following cases should be considered: (1) If δ<φ, the discontinuity is “l(fā)ocked up” and loses its mechanical effect; thus, the strength of the rock mass is equivalent to that of intact rock. (2) If δ-φ = 0 or π and δ = π/2 or 3π/2, there are singularities in the equations, and appropriate mathematical treatment is required to satisfy the continuity condition.

It can be seen that the rock mass strength is relatively low under a small confining pressure,because the strength is strongly affected by the mechanical effects of discontinuities.The strength increases fast as the confining pressure increases due to the“l(fā)ocked up”effect of discontinuities, and eventually reaches the strength of intact rock.

In addition, the scale effect is an important aspect of the compressive strength of rock masses,which is characterised by Eq.(24).The strength of rock masses is determined by the weakest link,i.e. the maximum radius of joint sets, am, which is defined within the studied volume as indicated in Eq.(4).It is well known that the UCS of intact rock, i.e. R0(or σc), has a strong scale effect as expressed by the Weibull form of R0as

where σ0is the basic strength of the rock when V = 1 m3and m→+∞.

4.3.2. Coulomb shear strength of jointed rock masses

The Coulomb shear strength of jointed rock masses is typically obtained through in situ shear tests,the results of which are usually scattered with large uncertainty due to the complex geological structure of the tested samples and neglected interrelation between the orientation of the joints and the shear direction.

Actually,the Coulomb shear strength criterion of rock masses in any direction, e.g.σ13, can be expressed as

Fig.10. Curves of two directional Poisson’s ratio.

where σ13is the shear stress towards x1on the shear plane with the normal direction x3; and f and c are the friction coefficient and cohesion for intact rocks and joints consistent with the subscript,respectively.

As a special case, a formula for the shear strength orienting forward and backward to the dip direction of joint set is written as follows:

Fig.11. Stereographic projection (a) and sectional distribution (b) of the compressive strength of rock mass.

where the positive of the sign‘±’indicates forward shearing to the dip direction, and the negative indicates backward shearing.

Clearly, the forward shear strength of rock mass is notably different from that in the opposite direction. The calculated t-σ curves by solving Eq.(28)with one joint set inclined 30°are shown in Fig.13. These curves and Eq. (28) show that: (1) the backward shearing results in a gradually increasing strength along with the increase in normal stress,and the shear strength eventually reaches that of the intact rock due to the locked-up effect of joints where tan（φ+δ）>tanφ0; and (2) for the forward shearing, the shear stress creates a tensional component on the joints where tan（φ-δ）

Following the conventional way for parameter determination,the cohesion of rock mass,cm,can be obtained from the intercept of the curve with t-axis,and the friction coefficient,fm=tanφm,is the slope of the curve.

To verify the theoretically calculated method of the Coulomb shear strength of jointed rock mass, the forward and backward shear tests on jointed rock mass models carried out by Li et al.(2014) under different normal stresses are introduced. The essential parameters of the model for calculation are tabulated in Table 6,where KIcis converted from the tensile strength of the model material according to an empirical relation(Bao et al.,2018).According to the failure phenomena observed from the shear tests, for both forward and backward shearing, the short side of rock blocks(0.02 m in length)in the rock mass model plays a controlling role in the failure process. Therefore, the shear failure along the 0.02 m long joint was mainly considered in the theoretical calculation.The comparison between the calculated and experimental results is shown in Fig.14. The theoretical calculation agrees well with the experimental results.The backward shear strength is ～3 times that of the forward shearing under a prescribed normal stress. More details regarding the shear tests on jointed rock mass models are referred to Li et al. (2014).

Table 5Comparison between the calculated result and the in situ tested data of UCS for a rock mass.

4.4. Rock mass quality rating from the deformability viewpoint

Evaluation of rock mass quality is a key step for rock engineering practices. A great number of rock engineering design and construction practises have been standardised based on the rock mass quality classification.

There have been many systems for rock mass classification,such as RQD, rock mass rating (RMR) (Bieniawski, 1974), Q-system(Barton et al.,1974),GSI(Hoek and Brown,1997),and Chinese basic quality (BQ) system. However, almost all of these systems are designed from the viewpoint of the strength of an isotropic medium,and rarely from the viewpoint of the deformability.They are not suitable for anisotropic rock masses, e.g. layered rock masses,which are frequently encountered in engineering practises.

In SMRM, a rock mass quality rating method based on Emwas developed, which reflects the deformability and anisotropy of rock masses by considering the mechanical effects of geological factors.Empirical formulae proposed by (Serafim and Pereira (1983)) are available to connect the quality grade and elastic modulus of rock mass:

The following two steps are taken for the SMRM rating of rock mass quality:(1)calculate the elastic modulus of rock mass,Em,in the loading direction by Eq.(21),convert it into an RMR score with Eq. (29), and take the result as the SMRM score in the calculated direction,i.e.SMRM=RMR;and(2)determine the grade of the rock mass according to the RMR rating regulations.

Fig.15 shows the stereographic projection of SMRM rating using the data in Tables 1-3. The calculated rating of the rock mass quality varies between 1.9 and 2.47 and the mean value is 2.23.The rock mass quality has been classified as grades II and III in the actual project, which agrees well with the SMRM rating.

Fig.12. Coulomb compressive strength of rock mass under different confining pressures.

Fig.13. Shear strength of rock mass subjected to forward and backward shearing.

Table 6 Parameters of jointed rock mass models (Li et al., 2014).

Due to the anisotropy of Em,the anisotropy index for the SMRM rating can also be defined, which could reflect the directional variation of the deformation behaviour:

The SMRM rating reflects not only the negative influence of geological factors,but also the positive effect via the increase in Em,i.e. a certain magnitude of in situ stress can improve the elastic modulus and in turn the rock mass quality. The contour plots of rock mass quality distributions for slopes and tunnels can also be drawn since the SMRM rating varies with the stress field, which will later be introduced.

Fig.14. Comparison between the calculated and experimental results of the shear strength of jointed rock mass model under different normal stresses:(a)Forward shearing;and(b) Backward shearing.

5. Engineering applications

5.1. SMRM module

Current numerical simulation software for geotechnical purposes is mainly designed for isotropic continuous medium. The Goodman element has been introduced for representing the mechanical effect of discontinuities (Goodman et al.,1968), but it increases significantly the complexity of models and the time of computation. Developers have also introduced some idealised models such as the universal joint model,but it could only consider some certain patterns of rock mass structure.

An SMRM software module for FLAC3D (Itasca, 2006), JointModel, was developed. The anisotropic constitutive model (Eq.(10)) and the failure probability (Eq. (9)) have been incorporated into JointModel, which could easily produce 3D contour maps of elastic modulus and rock mass quality rating that vary with the stress field of engineering rock masses.

On the premise of the acquisition of rock mass structure and the mechanical properties of intact rocks and joints,the JointModel can be directly applied to investigating the engineering properties and mechanical behaviours of rock masses in a slope or underground space. As shown in Fig.16, the JointModel module is designed to help calculate the rock mass parameters, analyse the deformation and failure of rock masses, and determine the required reinforcement for engineering rock masses.This module is presented along with field application examples shown below.

Fig.15. SMRM rock mass rating.

Fig.16. SMRM module and its application.

5.2. Intelligent identification of rock structures

The proposed method of digital field data collection and intelligent identification of joints is applied to Yinshan open-pit copper mine, located in Dexing, northeast of Jiangxi Province in China.

One slope face of the mine with good rock conditions was selected. This slope is sub-vertical (slope angle: 80°-90°), and characterised by several faults and joint sets. Under the effects of vibrations from mine blasting, the slope is vulnerable to rock fall and collapse disasters. UAV photogrammetry integrated with Structure-from-Motion (SfM) technique was employed to produce high-resolution 3D point cloud datasets,and then to interpret and quantify discontinuity parameters(e.g.orientation,number of sets,trace length,spacing, and frequency).

The digital photogrammetric survey was conducted with DJI Phantom 4 Pro UAV and on-board camera. Two hundred and seventy two digital photographs were acquired from a horizontal camera view,and had an overlap and sidelap of about 90%and 80%,respectively.A dense point cloud dataset of ～10 million points was generated.

As shown in Fig.17, a representative rock outcrop of the slope with minimum blast-induced fractures was chosen as the sampling zone. The approximate size of the selected area was 18 m × 6 m. Based on the adaptive iterative algorithm, the normal vector at each point was calculated. A density-based cluster analysis was performed, in which all points of a cluster belong to a set if they have a similar normal vector. Then the mean orientations of discontinuity sets were determined by clustering as shown in the stereographic projection (Fig. 18a).Four sets have been identified, and the corresponding points assigned to discontinuity sets from 3D point cloud are represented by different colours in Fig.17b. The trace map of every set using different colours is also presented. Additionally, the histograms of trace lengths and the fitted log-normal distribution for each set are shown in Fig.18b. The histograms of set spacing and fitted negative exponential distribution for each set are shown in Fig. 18c. Finally, a summary of the mean orientation,number of discontinuities, mean trace length, set spacing and linear frequency is tabulated in Table 7.

Fig.17. (a) A rock slope from Yinshan open-pit copper mine, China. The area selected in the red rectangle is used for analysis.(b)Discontinuity sets represented by different colours. (c) Trace mapping by the proposed method.

5.3. Active reinforcement of rock slopes

To date,reinforcement of rock masses typically follows a passive approach.It commonly aims at resisting the driving force of failure with artificial structures, without considering the self-supporting capacity of the medium. This design philosophy usually leads to a conservative and costly result.

The core idea of the active reinforcement of rock mass is to activate its self-stabilising potential. This method has been developed and verified by many practical projects across the world.The main theory of the method is introduced as follows.

Table 7 Results of grouping and extraction of discontinuity sets using intelligent identification method.

5.3.1. Self-stabilising potential and enhancement

Any medium has a self-stabilising potential,i.e.its strength,and the strength can be enhanced through adjusting the internal stress state. The minimum principal stress σ3can notably increase the strength σ1of the medium according to the Coulomb strength theory:

where σcis the UCS of the medium.

It is widely reported that φ>30°for typical rock masses,hence,tan2θ>3.Therefore,the strength σ1will be at least three times the stress σ3. This is the theoretical basis for the enhancement of the self-stabilising potential of a medium.

5.3.2. Design method of active reinforcement

(1) Determination of the optimal zone to be reinforced

The optimal zone for reinforcement can be determined by calculating Δσ3,which indicates the requirement for reinforcement:

where σ3cis the critical minimum principal stress determined by Eq.(31),and the actual minimum principal stress σ3is calculated by numerical simulation.

It is easy to confirm that the zone needed to be reinforced is located at the lower 1/3-1/2 of the slope surface (Fig.19a).

(2) Stability factor and anchorage force

Fig.18. (a)Stereographic projection(equal area net on lower hemisphere)of the slope part using proposed method;(b)Trace length histogram with fitted log-normal distribution(probability density function, PDF); and (c) Spacing histogram with fitted negative exponential distribution.

Fig.19. Contour map for slope reinforcement: (a) Δσ3; and (b) stability factor K.

Fig. 20. Right dam slope of Zongo II hydropower station, D.R. Congo: (a) Photograph of the cut slope; and (b) Numerical model for calculation.

The stability factor K at any points can be defined via the conventional limit equilibrium method as where τ is the shear stress on the shear plane, and τcis the shear strength of the plane.

The contour map of the stability coefficient K for a rock slope(Fig.19b) is drawn by numerical calculation. Setting a certain K as the required factor of safety for design, the related σ3ccan be obtained as σ3K= σ1-K（σ1-σ3）based on Eq.(31).Considering that Δσ3= （σ1-σ3）（1 - K）, the anchorage force can be designed as

where A is the area covered by each single anchor.

The designed anchor length should reach the area where the stability factor K exceeds the designed factor of safety(Fig.19b).In this case,the anchor dips at 35°.The active reinforcement method has been used in the analysis of deformation and stability for the right dam slope of Zongo II hydropower station in southwestern D.R.Congo.

The slope is 83 m high with an average slope of 40°inclining towards NE60°(Fig.20a).A large number of cracks were produced at the upper part of the slope in early 7-25 June 2014 due to the excavation at the slope toe, and the maximum aperture of cracks reached 28 cm. The active reinforcement has been applied to the slope in early August after a notable displacement was observed.A series of anchoring bolts was installed at the lower part of the slope in the early August to control the deformation.

The slope mainly consists of upper Palaeozoic thick layered sandstone, with an attitude of SW250°∠5°. The mechanical parameters of intact rocks and geometric data of joints are tabulated in Tables 8 and 9, respectively.

Numerical simulation was conducted to investigate the influence of the excavation at the slope toe on the slope stability(Fig. 20b). Fig. 21a-d shows the contour maps of displacement,reinforcement requirement, SMRM rock mass quality rating and displacement after reinforcement, respectively.

The maximum displacement reaches 27 cm at the top of the slope (Fig. 21a) after excavation, which is close to the monitored value.The reinforcement requirement(Fig.21b)is calculated by Eq.(32),which indicates that the bolts should be placed at the toe and the shoulder of the slope surface with a pre-stress value of 0-1 MPa in the horizontal direction. The SMRM quality rating of the rock mass in x-direction,calculated by Eq.(29),varies from 10 to 76,i.e.grade V to III (Fig. 21c). The estimated RMR value is around 40,representing the grade IV rock mass. Apparently, the contour map of SMRM rating provides more information on the spatial variation of rock mass quality. This also indicates that a constant rock mass quality rating is not appropriate for actal applications. Fig. 21d shows the displacement of the slope, which has been effectively restrained by the reinforcement.

The slope deformation was monitored sooner after excavation.Fig.22 shows the measured displacement of point No.1 at the top of the slope. The rate of deformation towards the valley (xdirection) has been significantly retarded after the reinforcement was carried out in the early August, and the increment of displacement in x-direction, around 50 mm, is very close to the prediction in Fig. 21d. This result shows that the active reinforcement has successfully improved the stability of the slope.

Table 8 Mechanical parameters of intact rocks.

Table 9 Geometric and mechanical properties of joints.

5.4. Behaviour of host rock masses of underground spaces

5.4.1. Circular tunnel

Asymmetric deformation in deep tunnel is a common issue along Lanzhou-Chongqing railway in northwestern China due to the inclined layered strata. Significant tunnel supporting systems have to be removed and rebuilt for multiple times, because the standard symmetric design fails to resist the asymmetric deformation.

To understand the influence of the discontinuities on the anisotropic deformation and rock quality of the host rock mass, a representative case study has been conducted for a circular tunnel in the thin layered carbonaceous schist of Trias. The main mechanical parameters of intact rocks are E=10 GPa and ν=0.3.The attitude of the beddings is N0°∠30°that dips towards left in Fig.23,and the parameters of the beddings are a=3 m,λ=60 m-1,cj=0 MPa,and φj= 17°.The tunnel is assumed to be subjected to gravitational stress at a depth of 500 m.

Fig.21. Contour maps of the slope calculated by the SMRM model:(a)Displacement due to excavation(m);(b)Reinforcement requirement(MPa);(c)SMRM quality rating;and(d)Displacement after reinforcement (m).

Fig. 22. Displacement curves of monitoring point No.1.

Fig.23. Contour maps of the anisotropic rock mass for a circular tunnel (three blue lines represent the attitude of the beddings):(a) Displacement (m);(b)Failure probability;(c)Reinforcement requirement (MPa); and (d) SMRM quality rating.

Table 10 Properties of the rock mass model used in simulation.

Fig. 23a shows the asymmetric displacement in the host rock mass due to the inclination of the beddings. Fig. 23b shows the failure probability, which indicates that the failure will most probably take place at the zone near the tangent point between the beddings and the tunnel wall. The map of reinforcement requirement (Fig. 23c) illustrates that the zone with greater failure probability requires stronger reinforcement.

Compared with the estimated RMR rock quality grade IV(20-40 in score), the SMRM quality rating rapidly reduced from 54(grade III) in the rock mass far from the tunnel to 10 (grade V) on the tunnel wall (Fig. 23d). The low quality zone extends from the tangent point towards the rock mass in the normal direction of the beddings.

These results notably reflect the mechanical effect of the anisotropic medium, under a symmetrical stress field and boundary constraints. The results also indicate that an asymmetric supporting system is needed to control the uneven deformation instead of the standard symmetric design.

Table 11 Properties of joints in the host rock mass of an underground space.

Fig.24. Comparison between calculated and simulated results of displacement surrounding a rectangular opening:(a)Photograph of the experimental model;(b)Displacement on the roof; (c) Displacement on the floor; and (d) Displacement on the right lateral wall.

5.4.2. Rectangular opening

Laboratory model tests on mining tunnels situated in layered rock masses were conducted and reported by Liu (2008). Layered rock mass models were established by placing thin plates made of rock-like materials into a two-axial loading vessel at a certain inclination.A rectangular opening was excavated and the loadings on the top (σ1= 1.2-2.9 MPa) and lateral (σ3= 0.4-4 MPa)boundaries were gradually increased until the failure of the opening. The dimensions of the entire model and the opening are 1.6 m×1.6 m×1 m and 25 cm×18 cm×22 cm,respectively.The joints dip at 8°and cut through the model with a spacing of 1 cm.The displacement at multiple locations surrounding the opening was measured. The properties of the rock-like material are tabulated in Table 10.

Comparison between the calculated and measured displacements at the top,bottom and lateral walls is shown in Fig.24.The displacements on all locations increase along with the increasing magnitude of loading. The displacement on the floor is the largest and that on the roof is ～4 times smaller than those at other locations, revealing that the deformation is significantly controlled by the structure of discontinuities. The calculated results by SMRM agree well with the experimental measurements, which validates the applicability of SMRM for deformation analysis.

5.4.3. Rockburst of underground space

Rockburst frequently occurs in underground space within high in situ stress environments. The underground powerhouse and transformer chamber of Jinping I hydropower station in southwestern China are taken as examples.These powerhouse and chamber are of 138-350 m deep in the left bank of Yalong River.During excavation,rockburst was frequently encountered in the site.

The strata consist of thick layered marble, and the mechanical properties of intact rocks are measured as E = 40 GPa, ν = 0.25,c=9.88 MPa,and φ=48°.The in situ stress near the powerhouse is measured as σ1= 35.7 MPa, pointing downwards to the river at 134°∠26°; σ2= 25.6 MPa at 15°∠45°; and σ3= 22.2 MPa at 243°∠34°. The parameters of main joint sets are tabulated in Table 11.

Fig.25a plots the plastic failure state of the elements as shown in the legend,and Fig.25b shows the rockburst in the host rock mass.The value in the legend represents the difference between the stress and the strength, which was calculated through the SMRM module JointModel, indicating the potential of rockburst. These diagrams demonstrate a highly similar characteristic that the failure of rock mass primarily takes place at the two vertical walls and the zone between the two chambers,which has been validated by the actual failure distribution.

Fig. 25c shows the distribution of the maximum stress, which verifies that the rockburst primarily occurs under the tension-shear stress state (Wu et al., 2010). Fig. 25d shows the reinforcement requirement, which indicates that the concentrated zone of the tensile stress requires stronger reinforcement.

Fig.25. Contour maps of failure elements and reinforcement requirement:(a)Plastic failure elements;(b)Rockburst potential(Pa);(c)Maximum stress(Pa);and(d)Reinforcement requirement (Pa).

6. Conclusions

This article presents new developments in SMRM with a detailed introduction to its principal theories, methods and engineering applications.The theoretical system in SMRM was incorporated into a constitutive model composed of stress-strain model,strength criterion and failure probability.All of these models were based on geometric probability models of the rock mass structure. A significant characteristic of theseSMRMmodelsisthatitiscapableofrealistically characterising the mechanical effects of discontinuities and the notable anisotropy of rock masses.

These theoretical models have been extended to a series of methods for parameter calculation of rock masses, including fulldirection RQD, elastic modulus, Poisson’s ratio, Coulomb compressive and shear strengths,and rock mass quality rating.SMRM could also serve as an engineering analytical method for complete deformation, rockburst and active reinforcement calculations of rock masses.A FLAC3D-based numerical simulation module,JointModel,embedded with the SMRM constitutive models has been developed to realise these functions in terms of 3D modelling and analysis.All these efforts aimed at making SMRM easier for engineering application.

Applications of SMRM to the rock masses at the dam site of QBT hydropower station in northwestern China,the dam slope of Zongo II hydropower station in D.R. Congo, and the underground powerhouse of Jinping I hydropower station in southwestern China have verified its reliability. Compared to conventional methods, the SMRM could effectively estimate the anisotropic mechanical behaviours of jointed rock mass and provide spatial distributions of these parameters.These results are helpful for engineering design of large-scale construction projects in strongly anisotropic rock masses,e.g.layered rock masses,subjected to high in situ stresses.

In the present form of SMRM, the medium and minimum principal stresses were presumed to be identical for compressive strength,and the shear strength was only theoretically solved along the joint plane without expanding to arbitrary directions,which are important topics in the future. Meanwhile, the SMRM needs to be applied to more engineering cases of various types to validate and improve the fundamental theories and methods involved.

Declaration of competing interest

The authors wish to confirm that there are no known conflicts of interest associated with this publication and there has been no significant financial support for this work that could have influenced its outcome.

Acknowledgments

The authors are grateful to the financial s upport from the National Natural Science Foundation of China (Grant No. 41831290)and the Key R&D Project from Zhejiang Province, China(Grant No.2020C03092).We are also grateful for the efforts of Dr.Xiuhong Hu and Dr. Xingxing Li in SMRM module programming and the help from Prof. Yumin Chen, Mr. Yunfei Xie, Mr. Ximeng Yue and Mr.Zhaoyuan Wang in the improvement of the software.

List of symbols