SAIMM, SANIRE and ISRM 6th International Symposium on Ground Support in mining and civil engineering construction B Debecker, A Tavallali, A Vervoort
PROBABILITY DISTRIBUTION OF ROCK PROPERTIES: EFFECT ON THE ROCK BEHAVIOUR B. Debecker, A. Tavallali, A. Vervoort K.U.Leuven, Research Unit Mining, Belgium
Abstract When considering risks in a design process, a probability approach is logic and, hence, the introduction of probability distributions for the input-values of the rock properties in numerical simulations is also logic. In this way, the heterogeneous nature of the rock material is quantified. It is shown that in e.g. the simulation of uniaxial compressive tests introducing some weaker elements has a much larger effect than one would expect from a weighted average. For example the introduction of about 20% of weaker elements results in a strength value of about one third of the weighted average strength. Also the deformation and the fracturing process are different in each simulation. A fault or natural discontinuity close to a tunnel is also investigated and the probability of initiation of shear movement is determined, as a function of probability distributions for two relevant properties (roughness and friction angle). Monte Carlo simulations show that for this particular case, the use of mean parameter values instead of parameter distributions results in an underestimation of the risk.
1
Introduction
The most simple numerical models consider a continuous, homogeneous and isotropic material with a linear elastic behaviour. In recent years, a lot of progress has been made in the development of codes implementing non-linear behaviour and/or discontinuous fracturing processes. However, independent of the complexity of the behaviour being integrated into a model, one should never forget that the most typical characteristic of rock material is its heterogeneous nature. This is valid for the intact material, as well as for the natural discontinuities. When considering risks in a design process, a probability approach is logic and, hence, the introduction of probability distributions for the input-values of the rock properties in numerical simulations is also logic. In this way, the heterogeneous nature of the rock material is quantified. This is not always so easy and the main question is if one always understands correctly the influence of taking probability distributions into account. One should neither forget that, apart from the effect on the overall behaviour, the introduction of probability distributions for the rock properties results in the fact that no two simulations are the same, even for the same characteristics of the probability distribution. The exact position within the mesh of for example a relatively strong or weak element can influence in a very significant way the behaviour (e.g. displacement vectors, stress concentration, …) and failure. This paper focuses on two relatively simple cases. First, an uniaxial compressive test is simulated using a continuum code, but whereby weak and strong elements are present in the model. Second, a fault close to a tunnel excavation is investigated, whereby some of the properties of the fault correspond to a specific distribution. The output parameter from the Page 357
SAIMM, SANIRE and ISRM 6th International Symposium on Ground Support in mining and civil engineering construction B Debecker, A Tavallali, A Vervoort
simulations is a type of safety factor corresponding to the initiation of shear movement. This safety factor is of course not a single value, but a distribution, allowing a more detailed interpretation.
2
Simulation of UCS-tests taking heterogeneities into account
The aim of this paragraph is to evaluate the effect of a random distribution of weak elements in a 2D continuum model (FLAC-2D) of UCS-tests. The Mohr-Coulomb failure criterion is incorporated based on a perfect elasto-plastic constitutive relation (Fang and Harrison, 2002). The elastic and plastic properties of strong and weak rock material used in this simulation with the Mohr-Coulomb failure criterion are given in Table 1. The computational model is a simple two dimensional rectangular shape with the dimension of 7 cm × 18 cm. 1 cm thick steel platen similar to laboratory conditions is put on top and bottom of the model. The bottom of the steel plate is fixed in Y direction. The properties of the interface elements that are used between the steel platen and rock specimen are given in Table 2. The mesh that is used consists of elements of 5mm × 5mm. So, the total number of elements is 504. A certain number of the elements are introduced as weak elements. In this study these elements are considered as heterogeneities. Percentages of these weak elements are 5%, 10%, 20%, 30%, 40% and 50% of all elements. The patterns of heterogeneity distribution are random and for each of the above percentages of weak elements, 10 random patterns are analysed. Figure 1 shows one of the random patterns for each of the different percentages of weak elements.
Table 1. Properties of strong and weak elements (Van de Steen, 2001). ____________________________________________________________ Properties strong element weak element Young’s modulus (GPa) 35 35 Poisson’s ratio 0.25 0.25 Tension cut off (MPa) 25 2.5 Cohesion (MPa) 80 8 Friction angle (°) 25 25 _____________________________________________________________
Table 2. Interface properties of rock and steel. ______________________________________________________ Cohesion: 0 Pa Normal stiffness (Kn): 1E15 Pa/m Dilation angle: 0° Shear stiffness (Ks): 2E10 Pa/m Friction angle (ΦI): 25° Tensile strength: 0 Pa ______________________________________________________
Page 358
SAIMM, SANIRE and ISRM 6th International Symposium on Ground Support in mining and civil engineering construction B Debecker, A Tavallali, A Vervoort
5%
10%
20%
30%
40%
50%
Figure 1. Example of random patterns for different percentages of heterogeneity (i.e. weak elements). Black squares are the weak elements. 275 250 225
UCS (MPa)
200 Weighted average
175 150 125 100 75 50 25 0
20
40 60 Heterogeneities (%)
80
100
Figure 2. Uniaxial compressive strength as a function of heterogeneity percentage.
Figure 2 shows the results of the peak stress as a function of heterogeneity percentage for all 10 random patterns. The straight line in Figure 2 shows the equivalent uniaxial compressive strength UCSeq which is simply written as :
Page 359
SAIMM, SANIRE and ISRM 6th International Symposium on Ground Support in mining and civil engineering construction B Debecker, A Tavallali, A Vervoort
H % × UCSW + (100 − H %) × UCS S , 100 where H% is the heterogeneity percentage, UCSW is the uniaxial compressive strength for weak elements (equal to 25.1 MPa) and UCSS is the uniaxial compressive strength for strong elements (equal to 251.1 MPa). Figure 2 shows that the simulation results and the weighted averages differ significantly. The presence of a small amount of weak elements has a relatively large effect on the overall strength. Although the results of 10 random patterns are different from the result of weighted average solution, the difference between the 10 random patterns is small, at least for the final strength (see also further). Therefore one could consider average curves (see Figure 3). These show that the mentioned analytical method for rock as a heterogeneous material could not be acceptable. Many theories have been developed to account for statistical aspects of rock failure (Fang and Harrison, 2002; Whittaker et al., 1992), but rock failure is a complex process and so it is difficult for any simple theory to be sufficiently flexible to cope with all possibilities (Fang and Harrison, 2002; Hudson and Fairhurst, 1969). It may also be that different rock types have different statistical characteristics, which require different functions for their description (Fang and Harrison, 2002; Tavallali et al., 2007). UCS eq =
50
Average for simulations
UCS reduction (%)
UCS (MPa)
275 225 175
Weighted average
125 75
40 30 20 10 0
25 0
20 40 60 80 Heterogeneities (%)
100
0
20 40 60 80 Heterogeneities (%)
100
Figure 3. Effect of heterogeneity percentages on the strength: left, Average of 10 series for UCS; right, Effect on reduction of UCS values in comparison to weighted average.
Different patterns of failure are observed in these series and make it difficult to have a unique interpretation for the results. Heterogeneity is seen to significantly affect the patterns of failure (Yuan and Harrison, 2005). For homogeneous samples, the behaviour is completely symmetric, e.g. when looking at the displacement vectors. However, for the heterogeneous samples the internal deformation is governed by the local presence of the weak elements (see Figure 4). In this Figure a deformed mesh is presented for two different random generations of 5% of weak elements. To be able to compare various results, a displacement index is defined: DI= |Maximum displacement vector| / |Applied displacement|. In the homogeneous models DI=1 and these models always yield in shear. In heterogeneous models different values of DI are obtained and these models sometimes yield in shear and Page 360
SAIMM, SANIRE and ISRM 6th International Symposium on Ground Support in mining and civil engineering construction B Debecker, A Tavallali, A Vervoort
sometimes yield in a combination of shear and tension. When 5% of elements are weak (H=5%), 5 of the models yield in shear and the other 5 yield in a combination of shear and tension (some elements yield in shear and some elements yield in tension). The average value of DI for the five models that yield in shear is 1.5, while the average value of DI for the five models that yield in shear and tension is 6.8. Of course, one should not forget that a continuum model is used and one could question the physical interpretation of displacements after yield, i.e. when in reality distinct fracture(s) would have occurred. Although the applied vertical displacement for both models in Figure 4 is 1.48 mm, the left model in Figure 4.a yields in shear with DI=1.1, but the second model (4.c) yields in shear and tension with DI=8.6. From Figure 4.e it can be concluded that the value and direction of maximum displacement vectors are completely different from the applied displacement vector. Figure 4.e shows that the vectors of displacements divert to the area which contains more weak elements close to the circumference and it is observed that in that area some of the elements yield in tension.
(a)
(b)
(c)
(d)
(e)
Figure 4. Two random patterns for H=5% and their deformed meshes (Black squares are weak elements; applied vertical displacement = 1.48 mm): (a) first model, (b) exaggerated (factor 5) deformed mesh of the first model, yield in shear with DI=1.1 (maximum displacement = 1.63 mm), (c) second model, (d) exaggerated (factor 5) deformed mesh of the second model before failure at the moment that the maximum displacement equals 1.63 mm, (e) deformed mesh of the second model, yield in shear and tension with DI=8.6 (no exaggerated mesh). When 10% of elements are weak, only one of the models yields in shear and the 9 others yield at the combination of shear and tension (see Table 3). The average value of DI for the 9 models that yield in shear and tension is 10.1. When 20% of elements are weak, similar to H=10%, only one of the models yields in shear and the others yield at the combination of shear and tension with the average value of DI=10.5. The condition changes by increasing the number of weak elements. When 30% of elements are weak, only one of the models with DI=6.9 yields in a combination of shear and tension and the 9 others yield in shear with values of DI close to 1. In H=40% and 50% all the models yield in shear.
Page 361
SAIMM, SANIRE and ISRM 6th International Symposium on Ground Support in mining and civil engineering construction B Debecker, A Tavallali, A Vervoort
Table 3. Summery of displacement index (DI) and yielding condition of all the heterogeneous models in this study. Heterogeneity percentage 5 10 20 30 40 50 Total
Yield in shear Amount 5 1 1 9 10 10 36
DI avg 1.5 1 1 1 1 1 1.07
Yield in shear and tension Amount 5 9 9 1 0 0 24
DI avg 6.8 10.1 10.5 6.9 9.43
300 0%
Stress (MPa)
250 200
5% 10%
150
20%
100
30% 40%
50
50% 100%
0 0
2
4
6 mstrain
8
10
12
Figure 5. Stress-strain curve for the models with different percentage of weak elements. For each percentage of heterogeneity, the model which has the UCS-value closest to the average is presented.
Figure 5 shows the stress-strain curve for the models with different percentage of weak elements. The two homogeneous models show the perfect elasto-plastic behaviour. The sharp corner of the curves in the homogeneous models (where the elastic behaviour changes to plastic behaviour) disappears in the heterogeneous models. In the strain range that strong elements are in elastic behaviour but weak elements are in plastic behaviour, the behaviour of heterogeneous models become complicated. Behaviour of heterogeneous models is elastic only before the yield point of the weak elements. Figure 5 also indicates that as the percentage of weak elements increases, a small strain softening behaviour appears after the peak stress. 3 Influence of joint parameter distributions on safety assessment In this paragraph, a fault or natural discontinuity is considered close to a tunnel (see Figure 6). The aim of this exercise is to illustrate the effect of considering probability distributions instead of fixed values. The physical parameter to be studied is the initiation of shear Page 362
SAIMM, SANIRE and ISRM 6th International Symposium on Ground Support in mining and civil engineering construction B Debecker, A Tavallali, A Vervoort
movement along this discontinuity. To examine this risk of initiation of shear movement, the safety factor can be calculated for any given position along this plane of weakness. The factor of safety, SF is calculated as the ratio of the resisting force, Fr and the driving force, Fd : A ⋅ coh + A ⋅ σ n tan φ coh + σ n tan φ F SF = r = = Fd A ⋅τ τ where σn and τ are respectively normal stress and shear stress on the plane, A the area, coh the cohesion and φ the angle of friction. If the SF is smaller than or equal to 1, it is assumed that shearing is initiated. However, since rock is a natural material, the friction angle is mostly not defined by one single value, but rather by a distribution of values over a certain range. Similarly, normal and shear stress are also dependent on the inclination of the plane, pore pressure, the virgin stress state,… which might all vary according to given distributions. The variability of these parameters results in a distribution of the factor of safety. Hoek (2007) defines a ‘considerable’ uncertainty as one where the standard deviation σ is around a quarter of the mean value μ. If this distribution is normal with a mean value of 1, 50% of the safety factors are thus smaller than one, while the other half is larger than 1. However, an uncertainty like this where 50% of the cases would fail is unacceptable. In order to overcome this, the required minimal SF is set to values larger than 1. Thus if one considers a similar normal SF distribution, but this time with a mean value of 1.5 and an identical standard deviation σ (0.25), only 2.3 % of the cases has a SF smaller than 1 (Figure 7). For distributions with smaller uncertainty (thus smaller σ), this percentage is even smaller. When no other information is available, it is commonly assumed that input parameters and safety factors are normally distributed. This assumption is then used to determine the minimal required SF in safety assessment.
Figure 6. Geometry of the case study: circular tunnel with a radius a = 2 m. A horizontal discontinuity is situated at h = 3.2 m from the centre of the circle. Vertical stress p is 5 MPa and horizontal stress K.p = 2 MPa.
Page 363
probability density functon
SAIMM, SANIRE and ISRM 6th International Symposium on Ground Support in mining and civil engineering construction B Debecker, A Tavallali, A Vervoort
NKR N
μ=1.0 σ=0.25
μ=1.5 σ=0.25
50%
MKR M M
MKR
2.3% N
SF
NKR
O
OKR
Figure 7. Two normal distributions of a SF: σ =0.25 and = 1.0 (left) and 1.5 (right). The area of the grey surfaces represents the probability that SF 1)
