Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Investigations into the Shear Strength Reduction method using distinct element models 2008

You don't seem to have a PDF reader installed, try download the pdf

Item Metadata

Download

Media
ubc_2008_fall_fournier_mathew.pdf
ubc_2008_fall_fournier_mathew.pdf [ 36.29MB ]
ubc_2008_fall_fournier_mathew.pdf
Metadata
JSON: 1.0052558.json
JSON-LD: 1.0052558+ld.json
RDF/XML (Pretty): 1.0052558.xml
RDF/JSON: 1.0052558+rdf.json
Turtle: 1.0052558+rdf-turtle.txt
N-Triples: 1.0052558+rdf-ntriples.txt
Citation
1.0052558.ris

Full Text

Investigations Into the Shear Strength Reduction Method Using Distinct Element Models by Mathew Fournier B.Sc., Queen’s University, 2006 A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIREMENTS FOR THE DEGREE OF Master of Applied Science in The Faculty of Graduate Studies (Geological Engineering) The University Of British Columbia (Vancouver) October, 2008 c©Mathew Fournier 2008 Abstract This thesis reports a detailed investigation into the use of the Shear Strength Re- duction (SSR) method to determine factor of safety values in discontinuum mod- els using the Universal Distinct Element Code. The SSR method depends on the definition of failure within the model and two different criteria were compared: the numerical unbalanced force definition and a more qualitative displacement- monitoring based method. A parametric study was first undertaken, using a simple homogeneous rock slope, with three different joint networks representing com- mon kinematic states. Lessons learned from this study were then applied to a more complex case history used for validation of the SSR method. The discontinuum models allow for the failure surface to propagate based on con- stitutive models that better idealize the rockmass than simpler methods such as limit equilibrium (e.g. either method of slices or wedge solutions) and even nu- merical continuum models (e.g. finite difference, finite element). Joints are ex- plicitly modelled and can exert a range of influences on the SSR result. Simple elasto-plastic models are used for both the intact rock and joint properties. Strain- softening models are also discussed with respect to the SSR method. The results presented highlight several important relationships to consider related to both nu- merical procedures and numerical input parameters. The case history was modelled similar to how a typical forward analysis would be undertaken: i.e. simple models with complexities added incrementally. The results for this case generally depict a rotational failure mode with a reduced fac- tor of safety due to the presence of joints within the rockmass when compared to a traditional limit equilibrium analysis. Some models with large persistence of steeply dipping joints were able to capture the actual failure surface. Soften- ing models were employed in order to mimic the generation and propagation of joints through the rockmass in a continuum; however, only discontinuum models using explicitly defined joints in the model were able to capture the correct failure surface. ii Table of Contents Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ii Table of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iii List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vi List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vii Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xii 1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.1 Problem Statement . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Research Objectives . . . . . . . . . . . . . . . . . . . . . . . . 2 1.2.1 Sensitivity Analysis . . . . . . . . . . . . . . . . . . . . 2 1.2.2 Case History . . . . . . . . . . . . . . . . . . . . . . . . 3 1.3 Thesis Structure . . . . . . . . . . . . . . . . . . . . . . . . . . 3 2 Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2.1 Definition of Failure . . . . . . . . . . . . . . . . . . . . . . . . 5 2.2 Use of Factor of Safety . . . . . . . . . . . . . . . . . . . . . . . 6 2.2.1 Method of Slices . . . . . . . . . . . . . . . . . . . . . . 9 2.2.2 Numerical Models . . . . . . . . . . . . . . . . . . . . . 11 2.2.2.1 Continuum Modeling . . . . . . . . . . . . . . 13 2.2.2.2 Distinct-Element Modeling . . . . . . . . . . . 13 2.2.2.3 Shear Strength Reduction Method . . . . . . . 13 3 Methodology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 3.1 Geomechanical Behavior of Rockmasses and Intact Rock . . . . 15 3.1.1 Elasticity . . . . . . . . . . . . . . . . . . . . . . . . . . 17 3.1.2 Plasticity . . . . . . . . . . . . . . . . . . . . . . . . . . 19 iii Table of Contents 3.1.3 Yield Criterion . . . . . . . . . . . . . . . . . . . . . . . 20 3.1.3.1 Mohr-Coulomb . . . . . . . . . . . . . . . . . 20 3.1.3.2 Hoek-Brown . . . . . . . . . . . . . . . . . . 22 3.1.4 Rockmass Classification Schemes . . . . . . . . . . . . . 23 3.1.4.1 Rock Mass Rating . . . . . . . . . . . . . . . . 26 3.1.4.2 Geological Strength Index . . . . . . . . . . . 26 3.2 Behaviour of Discontinuities . . . . . . . . . . . . . . . . . . . . 27 3.3 Discontinuity Network . . . . . . . . . . . . . . . . . . . . . . . 28 3.4 Limit Equilibrium . . . . . . . . . . . . . . . . . . . . . . . . . 30 3.5 Universal Distinct Element Code . . . . . . . . . . . . . . . . . 30 3.5.1 Constitutive Models . . . . . . . . . . . . . . . . . . . . 35 3.5.1.1 Intact Rocks . . . . . . . . . . . . . . . . . . . 36 3.5.1.2 Joints . . . . . . . . . . . . . . . . . . . . . . 37 3.6 Shear Strength Reduction Using UDEC . . . . . . . . . . . . . . 38 3.6.1 Displacement-Based SSR . . . . . . . . . . . . . . . . . 41 4 Parametric Study . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 4.1 Cycling Mode . . . . . . . . . . . . . . . . . . . . . . . . . . . 50 4.2 Element Size . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58 4.3 Element Types . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 4.4 Insitu Stress . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 4.5 Elastic Properties . . . . . . . . . . . . . . . . . . . . . . . . . . 67 4.6 Joint Shear Stiffness . . . . . . . . . . . . . . . . . . . . . . . . 72 4.7 Strain-Softening Models . . . . . . . . . . . . . . . . . . . . . . 79 4.8 Joint Network Cases . . . . . . . . . . . . . . . . . . . . . . . . 81 4.8.1 Daylighting Joints . . . . . . . . . . . . . . . . . . . . . 81 4.8.2 Discontinuous Daylighting Joints . . . . . . . . . . . . . 85 4.8.3 Slope Parallel Joints . . . . . . . . . . . . . . . . . . . . 90 4.9 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94 5 Case History . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 5.1.1 Geology . . . . . . . . . . . . . . . . . . . . . . . . . . 97 5.1.2 Slope Failure . . . . . . . . . . . . . . . . . . . . . . . . 98 5.1.3 Structure . . . . . . . . . . . . . . . . . . . . . . . . . . 100 5.1.4 Geotechnical Properties . . . . . . . . . . . . . . . . . . 104 5.1.5 Modelling Approach . . . . . . . . . . . . . . . . . . . . 106 5.1.6 Selection of Critical Sliding Surface . . . . . . . . . . . . 109 iv Table of Contents 5.2 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117 5.2.1 Limit Equilibrium . . . . . . . . . . . . . . . . . . . . . 117 5.2.2 UDEC: Unjointed Slope Models . . . . . . . . . . . . . 120 5.2.3 UDEC: Continuous Bedding Models . . . . . . . . . . . 128 5.2.4 UDEC: Discontinuous Joint Models . . . . . . . . . . . . 129 5.2.4.1 Strain-Softening Behaviour . . . . . . . . . . . 138 5.2.4.2 Material Property Variations . . . . . . . . . . 144 5.2.4.3 Water Influence . . . . . . . . . . . . . . . . . 148 5.3 Continuous Cycling Solution Mode . . . . . . . . . . . . . . . . 152 5.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 153 6 Discussion and Conclusions . . . . . . . . . . . . . . . . . . . . . . 157 6.1 Framework . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 157 6.2 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 159 6.3 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 163 6.3.1 Guidelines for Use . . . . . . . . . . . . . . . . . . . . . 164 6.3.2 Recommendations for Further Work . . . . . . . . . . . 166 Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 169 Appendices A Displacement-based SSR Code . . . . . . . . . . . . . . . . . . . . . 173 B Example Input Files for Parametric Study . . . . . . . . . . . . . . 185 C Case History Joint Networks . . . . . . . . . . . . . . . . . . . . . . 192 v List of Tables 2.1 Statics satisfied and interslice forces in various methods modified from (Krahn, 2003) . . . . . . . . . . . . . . . . . . . . . . . . . 10 4.1 Slope geometries used in parametric study . . . . . . . . . . . . . 45 4.2 Discontinuous and slope parallel joint model parameters . . . . . 49 4.3 Daylighting joint model parameters . . . . . . . . . . . . . . . . 49 4.4 Number of cycles for 10−5 force ratio limit condition in 10 m discontinuous daylighting joint case . . . . . . . . . . . . . . . . 51 4.5 UDEC-SSR results for two different ratio limits . . . . . . . . . . 56 4.6 Summary of D-SSR FOS values for element size variation . . . . 59 4.7 Input parameters for Poisson’s Ratio sensitivity models . . . . . . 71 4.8 Peak-residual parameters used for strain-softening sensitivity . . . 80 4.9 SSR results for different daylighting joint spacings . . . . . . . . 85 4.10 Discontinuous joint generation settings . . . . . . . . . . . . . . . 86 4.11 Selected results from UDEC-SSR for slope parallel joint case . . . 91 4.12 Distance of critical sliding surface from crest for Slope parallel cases . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92 4.13 Summary of parametric study results . . . . . . . . . . . . . . . . 95 4.13 Summary of parametric study results (cont.) . . . . . . . . . . . . 96 5.1 Summary of geotechnical parameters used in the case history . . . 105 5.2 Model and cycle stages of yielded elements . . . . . . . . . . . . 109 5.3 Summary of LEM results, surfaces can be seen in figure 5.17 . . . 118 5.4 Summary of unjointed SSR results . . . . . . . . . . . . . . . . . 121 5.5 Summary of continuous bedding SSR results . . . . . . . . . . . 129 5.6 Joint generation parameters for jointed case history model . . . . 131 5.7 Summary of discontinuous SSR results with a wide range of joint networks, revealing small range in factor of safety values . . . . . 132 5.8 Summary of continuous displacement-based SSR for jointed model cases . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 152 vi List of Figures 2.1 Continuum and discontinuum model examples . . . . . . . . . . . 12 3.1 Scale effects of rockmass failures and properties, modified from (Hoek, 2007) . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 3.2 Stress strain relationships of varying rockmass quality, from Hoek (2007) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 3.3 Graphical representation of the Mohr-Coulomb failure envelope, modified from (Hudson and Harrison, 2005) . . . . . . . . . . . . 21 3.4 Approximate M-C failure envelope from H-B parameters, modi- fied from (Hoek et al., 2002) . . . . . . . . . . . . . . . . . . . . 24 3.5 Disturbance factor guidelines for slopes, modified from (Hoek et al., 2002) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 3.6 Computational cycle used in UDEC, from (Itasca, 2004c) . . . . . 31 3.7 Difference between actual jointing and modelled joints . . . . . . 33 3.8 Nomenclature for JSET joint generator in UDEC, from (Itasca, 2004a) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 3.9 Flowchart for UDEC-SSR showing stability criteria and bracket- ing search, modified from (Itasca, 2004c, 2008) . . . . . . . . . . 40 3.10 Flowchart for general displacement-based SSR method . . . . . . 42 3.11 Conceptual figure showing displacement response and selection of strength reduction factor indicating instability . . . . . . . . . . 44 4.1 Example of scoping of joints (a) and meshing (b) within slope for computational efficiency . . . . . . . . . . . . . . . . . . . . . . 47 4.2 An example of the three joint networks used in the parametric study. Note that the view shown is zoomed to the area of interest compared to the entire model. . . . . . . . . . . . . . . . . . . . 48 4.3 Results of displacement-based SSR for varying number of cycles . 52 4.4 Results of displacement-based SSR for varying force ratio limits per strength reduction step . . . . . . . . . . . . . . . . . . . . . 53 vii List of Figures 4.5 Velocity plots showing the same activated portions of slope (seen by the presence of large velocity vectors) for 10 m and 15 m dis- continuous joint cases, at the onset of instability . . . . . . . . . . 55 4.6 Example of sawtooth horizontal displacement results during D- SSR analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57 4.7 Results of displacement based SSR for element size sensitivity . . 59 4.8 Plastic state of zones in a single block showing detail of shear band localization defined by elements yielding in shear (x yield in shear, o yield in tension) . . . . . . . . . . . . . . . . . . . . . . 61 4.9 Comparison between cycled and force-equilibrium results for 7.5 m grid size . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 4.10 D-SSR results for element type sensitivity for discontinuous day- lighting joints . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64 4.11 Failed volume indicators for 7.5 m triangular elements at SRF stage = 1.76 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 4.12 Failed volume indicators for 7.5 m quadrilateral elements at SRF stage = 1.74 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 4.13 Horizontal stress contours for varying modelled insitu stress ratios 68 4.14 Vertical stress contours for varying modelled insitu stress ratios . . 69 4.15 D-SSR results reduction results for varying horizontal to vertical stress ratios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70 4.16 Results of displacement based SSR for Poisson’s Ratio sensitivity 72 4.17 Results of D-SSR analysis for joint shear stiffness sensitivity . . . 74 4.18 Number of cycles for 10−5 ratio limit for ks = 8 x 108 and 8 x 1010 Pa/m cases . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75 4.19 Results of displacement-based SSR for varying ratio limits when ks = 8 x 108 Pa/m . . . . . . . . . . . . . . . . . . . . . . . . . . 76 4.20 Results of displacement-based SSR for varying ratio limits when ks = 8 x 1010 Pa/m . . . . . . . . . . . . . . . . . . . . . . . . . 77 4.21 Comparison between active portions of the slope for different ratio limits when ks = 8 x 108 Pa/m . . . . . . . . . . . . . . . . . . . 78 4.22 D-SSR for strain limits for strain-softening cases . . . . . . . . . 80 4.23 D-SSR results for different daylighting joint spacings . . . . . . . 82 4.24 Horizontal displacement and shear along joints for varying day- lighting joint spacings (m) . . . . . . . . . . . . . . . . . . . . . 83 4.24 Horizontal displacement and shear along joints for varying day- lighting joint spacings (m) (cont.) . . . . . . . . . . . . . . . . . 84 viii List of Figures 4.25 D-SSR results for discontinuous daylighting joint cases with vary- ing joint network properties . . . . . . . . . . . . . . . . . . . . . 88 4.26 Joint network differences between 7 m and 8 m spacing cases . . . 89 4.27 Results of displacement based SSR for slope parallel joint spacing cases . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90 4.28 Plasticity plot showing crushing at the toe, tension cracks at the crest, and location of critical discontinuity for 15 m case (x yield in shear, o yield in tension) . . . . . . . . . . . . . . . . . . . . . 93 5.1 Cross-section of pre-failure pit, showing major geological units . . 99 5.2 View of the south wall showing the failure, taken July 13th 2005. . 100 5.3 Stereonet showing faults in the failure area . . . . . . . . . . . . . 101 5.4 Stereonet showing bedding in the failure area . . . . . . . . . . . 102 5.5 Stereonet showing joints in the failure area . . . . . . . . . . . . . 103 5.6 Example of telescoping of joints (a) and meshing (b) within slope for the case history . . . . . . . . . . . . . . . . . . . . . . . . . 107 5.7 Yielded elements at 10−6 ratio limit for UDEC-SSR giving FS = 1.12 (x yield in shear, o yield in tension) . . . . . . . . . . . . . . 110 5.8 Yielded elements at 5,000 cycles for D-SSR at SRF = 1.11 (x yield in shear, o yield in tension) . . . . . . . . . . . . . . . . . . . . . 110 5.9 Yielded elements at 10−6 ratio limit for D-SSR at SRF = 1.11 (x yield in shear, o yield in tension) . . . . . . . . . . . . . . . . . . 111 5.10 Yielded elements at 10−6 ratio limit plus 5,000 cycles for D-SSR at SRF = 1.11 (x yield in shear, o yield in tension) . . . . . . . . . 111 5.11 Yielded elements at 10−6 ratio limit plus 15,000 cycles for D-SSR at SRF = 1.11 (x yield in shear, o yield in tension) . . . . . . . . . 112 5.12 Yielded elements at 10−6 ratio limit plus 30,000 cycles for D- SSR at SRF = 1.11, showing the ultimate failure surface (x yield in shear, o yield in tension) . . . . . . . . . . . . . . . . . . . . . 112 5.13 Yielded elements at 10−6 ratio limit plus 50,000 cycles for D- SSR at SRF = 1.11, showing the ultimate failure surface (x yield in shear, o yield in tension) . . . . . . . . . . . . . . . . . . . . . 113 5.14 Example of selecting sliding surface from UDEC velocity vector output . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114 5.15 Velocity plot derived surfaces from models states listed in table 5.2 115 5.16 Velocity plot derived surfaces from models states listed in table 5.2 with yielded elements overlain . . . . . . . . . . . . . . . . . 116 5.17 Critical LEM failure surfaces predicted by SLIDE . . . . . . . . . 119 ix List of Figures 5.18 Unjointed model geometry after bench removal . . . . . . . . . . 121 5.19 Results of displacement-based SSR for models U1-U5 showing a critical SRF range of 1.12 – 1.17 . . . . . . . . . . . . . . . . . . 122 5.20 Comparison of LEM surfaces and D-SSR SRF = 1.13 yielded el- ements for model U1 . . . . . . . . . . . . . . . . . . . . . . . . 124 5.21 Comparison of LEM surfaces and UDEC-SRF velocity-derived surfaces for model U1 . . . . . . . . . . . . . . . . . . . . . . . . 125 5.22 Comparison of LEM surfaces and D-SSR SRF = 1.14 yielded el- ements for model U2 . . . . . . . . . . . . . . . . . . . . . . . . 126 5.23 Comparison of LEM surfaces and UDEC-SRF velocity-derived surfaces for model U2 . . . . . . . . . . . . . . . . . . . . . . . . 127 5.24 Plot showing shear displacements, in red, on bedding surfaces at the point of failure in UDEC-SSR for model C3 . . . . . . . . . . 129 5.25 Results of D-SSR for jointed case history models showing range of critical SRF values from 1.06 – 1.15 . . . . . . . . . . . . . . . 133 5.26 Velocity surfaces representing critical sliding surfaces, from UDEC- SSR, that provide a good match to the actual failure surface . . . . 134 5.27 Velocity surfaces representing critical sliding surfaces, from UDEC- SSR, that provide a good match to LEM Surface Two . . . . . . . 135 5.28 Velocity surfaces representing critical sliding surfaces, from UDEC- SSR, that do not provide a good match and predict a deeper seated failure . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 136 5.29 Yielded elements showing internal shearing of the failing mass in the J4 joint network for a SRF stage just prior to catastrophic failure137 5.30 Results of displacement-based SSR for jointed strain-softening case history models . . . . . . . . . . . . . . . . . . . . . . . . . 139 5.31 Velocity surfaces representing critical sliding surfaces for the J12 joint network for elastic perfectly-plastic behaviour and strain soft- ening behaviour, at the point of instability as defined by D-SSR . . 140 5.32 Results of displacement-based SSR for models jointed strain soft- ening models R1-R3 with low residual properties to capture frac- ture generation within a continuum . . . . . . . . . . . . . . . . . 142 5.33 Velocity surfaces representing critical sliding surfaces for R1–R3 strain-softening models . . . . . . . . . . . . . . . . . . . . . . . 143 5.34 Failure surface and factor of safety sensitivities with respect to initial parameters before strength reduction . . . . . . . . . . . . 145 5.35 Results of displacement-based SSR for models ks sensitivity for joint network model J8 . . . . . . . . . . . . . . . . . . . . . . . 146 x List of Figures 5.36 Variations in failure surface location with initial joint strength pa- rameters for joint network model J8 . . . . . . . . . . . . . . . . 147 5.37 Estimated static water table location, in blue, for case history . . . 149 5.38 Results of displacement-based SSR for models W8 and W12 with static water tables . . . . . . . . . . . . . . . . . . . . . . . . . . 150 5.39 Plasticity indicators showing clear shear band after removal of last bench in W12 model . . . . . . . . . . . . . . . . . . . . . . . . 150 5.40 Velocity plot showing slope failure along critical discontinuity be- hind the crest after removal of last bench in W8 model . . . . . . 151 5.41 Results of displacement-based SSR for models CD1–CD8 show- ing identical behaviour between continuously-cycled and initial state displacement-based SSR . . . . . . . . . . . . . . . . . . . 153 5.42 Photograph roughly perpendicular to failure slope showing pres- ence of a steep, continuous structure . . . . . . . . . . . . . . . . 155 6.1 Typical risk assessment framework for rock slopes, modified from Duzgun and Lacasse (2005) . . . . . . . . . . . . . . . . . . . . . 158 xi Acknowledgments I would like to express my gratitude to my supervisor Dr. Erik Eberhardt. His door was always open to my questions, no matter how trivial. I would like to also thank Dr. Oldrich Hungr for his advice during my research. This document would not have been produced without their patient help. I would like to thank the Natural Sciences and Engineering Research Council of Canada (NSERC) and Golder Associates (particularly the Burnaby office) for the generous funding through the Industrial Partnership Scholarship. I would also like to thank Darren Kennard and Al Chance of Golder Associates for the use of their time and facilities at the Burnaby office. Rocscience was kind enough to donate a graduate license for some of their software and it was appreciated. Piteau Associates and Barrick Gold Corporation provided me with case history information and guidance through the latter part of the thesis. Nick Rose, in par- ticular, needs to be thanked for supplying the case history data and model guid- ance. Many thanks to Doug Stead for time and effort on my behalf as my external ex- aminer. Thank you also to my helpful second, third, and fourth pairs of eyes: Becca Thomas, Katie Grubbe and Katie Jones, for donating their time to slogging through some less-than-ideal early drafts. xii Chapter 1 Introduction 1.1 Problem Statement Problems arising from slope stability are increasing as population densities in- crease throughout the world in mountainous terrain, as open pit mining is delving deeper, and as engineering works are being created in poorer quality rockmasses in order to keep up with the growing demand for natural resources. Increasingly, large slopes are being instrumented and monitored in order to assess and estab- lish trigger limits for early warning of unsafe slope conditions or to indicate when rehabilitation work needs to be done. Coupled with this are advances in numer- ical modelling techniques, and increased computational power, including shifts towards more complex analyses that incorporate both the behaviour of intact rock and the joints within the overall rockmass. The monitoring data provides an excel- lent reference against which these advanced numerical models can be tested. Slope stability is of prime concern to geological engineering from both an eco- nomic standpoint and with respect to the risk large slopes can pose to mine person- nel and the general public. Traditional methods such as limit equilibrium (LEM) are commonly used for soil and planar joint or wedge stability assessments to obtain a factor of safety. However these methods suffer from a number of limi- tations and assumptions, particularly when applied to rockmasses instead of soils or persistent daylighting structure. These include: a large disconnect between modelled and real behaviour, assumptions of a priori failure surfaces, incorrect stresses in the slope, need for assumptions related to the inclinations and loca- tions of interslice forces, and the inability to include deformation in the definition of failure. Complex rock deformation, including both intact rock failure and structurally con- trolled movement, is difficult to model and relies on advanced continuum and discontinuum codes. Results of such modeling typically do not include a value 1 Chapter 1. Introduction of safety factor, which relegates these models to determining possible failure mechanisms or volumes of potential failures rather than deterministic analysis. A technique called the Shear Strength Reduction (SSR) method has been pro- posed and used in simple, continuum finite-element and finite-difference models and compares favorably to results from the limit equilibrium analysis for simple geometries. Most advanced numerical modelling packages, including 2D and 3D distinct-element modelling codes, include SSR routines; however there is little experience in using this method in these models which include additional consti- tutive models for explicit definition of joints and modeling complexities related to these joints, compared to the simpler continuum models it was originally applied to. Research was undertaken through this thesis to investigate the use of the SSR in distinct-element models to determine sensitivities and differences between the possible implementations of the SSR method through the use of a set of parametric studies based on different slope configurations. Insights gained through the para- metric study were then used in the back analysis of a complicated case history. Ultimately, this thesis aims to determine the validity of the SSR in discontinuum models, and to propose guidelines for the engineering practitioner as to how to properly apply it. 1.2 Research Objectives The purpose of this research is to investigate the sensitivities of the SSR method in distinct-element models, the differences between the stability criterion of SSR, and to apply the SSR method to a complicated case history in order to compare results to a traditional limit equilibrium analysis. Much of the work focuses on the use of basic constitutive models with parameters that are derived from a typical site investigation without the use of specialized testing strictly for model input. The specific objects of the research are described below. 1.2.1 Sensitivity Analysis 1. Develop numerical models based on three basic configurations (a) Continuous daylighting joints 2 Chapter 1. Introduction (b) Slope parallel joints (c) Non-continuous daylighting joints 2. Develop a user-defined displacement-based SSR (D-SSR) stability criterion to compare to a numerical force ratio based stability criterion (the latter is that used by the commercial code UDEC; referred herein as UDEC-SSR). 3. Investigate the sensitivity of factor of safety to various discontinuum con- figurations and possible input parameters on the different results returned from the two SSR stability criteria. 4. Develop general guidelines for application of the SSR method to jointed rockmasses using basic constitutive models 1.2.2 Case History 1. Determine the limit equilibrium safety factor and least stable slide surface. 2. Determine the range of possible joint networks involved. 3. Develop numerical models to apply both the displacement-based SSR and UDEC-SSR criterions, which increase in complexity in terms of detail of joint network, water, and constitutive models. 4. Compare SSR results with that of a traditional LEM solution. 5. Refine ambiguous sensitivity cases from the parametric study by way of comparison to the LEM solution. 6. Refine guidelines for application of the SSR method to jointed rockmasses. 1.3 Thesis Structure The problem statement and research objectives are presented above constitute Chapter 1. Chapter 2 introduces various failure and modelling terminology used throughout the analysis. An overview of slope stability and slope stability modelling methods is presented. 3 Chapter 1. Introduction The methodology behind the SSR method is discussed in Chapter 3. It includes specifics to the discontinuum modelling code utilized for the thesis, as well as a discussion on the different implementations of SSR used. Joint networks, be- haviour of rockmasses, generation of model inputs, and constitutive models used in the modelling are also discussed. The parametric study is discussed in detail in Chapter 4. This includes three com- mon slope configurations and the sensitivities to the SSR methods due to model setup and inputs. Chapter 5 focuses on the back analysis and comparison to LEM results of a case history. This includes an introduction to the problem with geology, a description of the discontinuity network, events leading up to the failure, and results as model complexity increases. A general discussion on the SSR method with reference to its applicability as a slope analysis tool is presented in Chapter 6. Limitations and uncertainties about the modelling work and the SSR method in general are reviewed. Conclusions and practical guidelines for use are also given. 4 Chapter 2 Background 2.1 Definition of Failure Failure of geological rock and soil materials can be complex due to heterogeneity, the range of materials involved, structural interactions, and site-specific condi- tions. Because of this, in order for failure to be understood in an engineering manner, the failure needs to be defined in terms of the mechanism and properties before any analysis can take place. The mechanism may be a simple one, or may be a complex multi-staged mecha- nism with time dependence. However, the sequence of events is particularly es- sential. Material properties are important, not just geotechnical properties derived from lab testing, but qualitative site investigation insights as well (Scott, 1987). Large scale failure requires the use of experience and qualitative investigations in order to expand the narrowly focused, small-scale results from lab testing to something meaningful for the larger portions of the slope involved. Only once these two stages have been completed can a proper analysis, which incorporates both the mechanism and the properties, be completed. This analysis first begins with the application of known theory, principles, and methods in order to see if the failure can be classified in a traditional sense. If these methods fail, then some new procedure may need to be developed in order to classify a poten- tial new mechanism, consequence, or material behaviour. Thus, failures can be classified by the combination of results, difficulties, and uncertainties as outlined by (Scott, 1987): Type A Standard measurement of material properties with an associated analysis that results in an accurate representation of the failure mechanism involved. Material properties can be varied to study their effects on the failure as the mechanism, properties, and analysis are clear. 5 Chapter 2. Background Type B The mechanism is easily identified but the analysis does not match the failure using the material properties from accepted methods. Type C Natural variability in materials or key defects, or an ill-understood mech- anism are involved. Key information has been destroyed by the failure, re- sulting in difficult or impossible analysis. Type D Material properties are difficult to obtain due to the material involved, scale effects, or natural variability. Despite the mechanism being clear, anal- ysis is impossible due to problems obtaining proper material inputs. Type E The failure mechanism is difficult to define requiring the use of numerical models for potential insight. This results only in speculation of what may have actually happened instead of a solid analysis. In order to fit into one of the above categories, a minimal threshold for failure has to be established. In a strictly engineering sense, a rock sliding out of an open pit bench constitutes failure; however, in an operational sense, it is not a failure if there is no consequence. Scale and consequence are important and divide the concept of failure into an engineering sense and an operational sense. Failures which involve little consequence may be noted and recorded but largely ignored. They may be commonly expected and designs may incorporate protective mea- sures (e.g. a wide bench to retain bench-scale sized failures) as preventing the failure may be subject to economic constraints. Conversely, large deformations may develop but not lead to a catastrophic failure. Yet if they are large enough to affect or impede operations then it may be classified as failure. Thus, the definition of failure is fluid and is subject to perspective and site specific conditions. 2.2 Use of Factor of Safety Slope stability analysis to determine failure is generally limited to the assessment of the factor of safety (FOS). The concept of the factor of safety has existed for a long time in engineering nomenclature, and consequently has many champions, critics, and specific definitions. In it’s simplest form, the factor of safety can be defined as the ratio of the material’s ability to resist load (strength) to the actual or calculated load. It was originally designed as a safeguard against errors in predict- ing the load, variations in material properties, uncertainties in material properties, and differences between modelled idealized behaviour and true behaviour (El- 6 Chapter 2. Background ishakoff, 2004). Factor of safety can either be deterministic or probabilistic in nature. It is important to note that in designs giving a single, deterministic value of factor of safety that the important factor of variability if often overlooked. Es- sentially, this means that factor of safety can be an arbitrary number derived from engineering judgement and analytical methods and that it is not necessarily the probability of failure (Elishakoff, 2004). In slope stability problems with geoma- terials, it is possible to apply probability theory for some methods. Other tech- niques can only be used with deterministic methods due to current computational time constraints or the model mechanics involved. Also, the use of single values may only be possible as the data required to properly define parameter variability in the form of a probability distribution function may not be available, may be too highly variable for a meaningful statistical distribution, or subject to economic constraints in order to fully investigate. There is familiarity in industry with deterministic factors of safety due to limit equilibrium methods being used, in practice, for a long time. There is less confi- dence with newer and more complex numerical modeling schemes, many of which do not return a simple factor of safety or find it via indirect methods with com- plicated assumptions and methodologies whose implications are yet to be fully understood. Modeling of any system requires good boundary conditions, reason- able knowledge of material input parameters, and selection of proper constitu- tive models. This is difficult in geological engineering practice due to rockmass complexity and the relatively small amount of material sampled compared to the system being modelled (i.e, scale effects). Compounding this problem is that in- situ measurements are highly variable, include some sort of testing bias, and are difficult to constrain. Some models include parameters that are difficult to mea- sure in the field, let alone generate a statistically meaningful distribution for, and must be selected carefully using sound engineering judgment and past experience. There are a range of techniques available , each with inherent strengths and weak- nesses: 1. Analytical techniques using simple closed-form solutions and/or the heavy use of design charts for design slope angles. 2. Limit equilibrium methods, which encompass the following: (a) Structurally controlled (e.g. closed-form wedge solution). (b) Stress-controlled (e.g. method of slices based on force and/or moment balance). 7 Chapter 2. Background 3. Numerical modelling, both 2D and 3D for: (a) Continuum i. Finite-element method ii. Finite-difference method (b) Discontinuum i. Discontinuous deformation analysis ii. Distinct-element method The selection of modeling code depends on the site conditions and expected range of failure mechanisms. For example, limit equilibrium is appropriate for use in a structurally controlled failure mechanism such as a sliding wedge or for strictly homogenous materials in a slope failing through shear. For all methods, predic- tion of instability requires high quality input data. As the quality and quantity of data decreases, modeling results shift from prediction of instability to providing a limited insight into the underlying mechanisms that affect the stability of the slope. Increased modelling capabilities and computer power has allowed for more com- plicated rockmass problems to be modeled. This has introduced the temptation to model very detailed rock behavior that may not be useful for general rock engi- neering practice. Model verification is an important aspect of sound engineering that becomes increasingly difficult as model complexity increases. At the same time, a model may not capture the fundamental failure mechanism if it is oversim- plified. As factor of safety approaches are applied to more complex codes, there is an increased variability of factor of safety values returned due to additional model complexities that are not present in simpler models. Because of this, it is im- portant to understand how the fundamental behaviour of the modelling code and constitutive models of materials affects the single value of factor of safety re- turned. This is particularly true in continuum and discontinuum modeling where a priori failure surfaces are not assumed as it can be difficult to constrain the fac- tor of safety to a particular failure mode without some other failure becoming the dominate mechanism as analysis continues. These difficulties are the result of the following three problems, related to the use of numerical modelling: 8 Chapter 2. Background • Generally long solution times compared to LEM methods, as is the case with discontinuum models. • Difference in fundamental behavior and constitutive models that makes it difficult or impossible to compare the results from one analysis to the other, particularly when switching between a homogenous, continuum-like model to a discontinuum model. • Lack of experience, published case histories, back analysis, and established practices in using deterministic analysis beyond LEM methods. 2.2.1 Method of Slices The method of slices using limit equilibrium, referred throughout this thesis as LEM, is the most commonly used method for assessing the stability of a slope. It requires the least amount of material parameter inputs and offers no selection of governing constitutive models other than the static stability method presented. The limit equilibrium approach has been in engineering practice for a long time, and several advances in the method have occurred regarding assumptions of force equilibrium, internal shear, and moment equilibrium. There are several different possible answers obtained when the different LEM solution methods are applied. Each method divides the slope into a number of vertical slices, except for the Sarma method which allows for inclined slices. A test failure surface, circular or non-circular depending on the method, is assumed. Forces such as the shear and normal force along the base of the slice and between slices are calculated. The factor of safety is determined by the ratio of the resisting forces in the slope over the driving forces. Table 2.1 below reproduced from Krahn (2003) summarizes the main assumptions of some of the different method of slices solutions. The different solutions require assumptions about force and moment equilibrium, interslice forces, and the shape of the failure surface. Each method treats these assumptions differently and it is common to have a different factor of safety for the same trial surface and/or different critical surfaces for the same limiting value of factor of safety between the methods. The general limit equilibrium equations are presented below (Fredlund and Krahn, 1977): 9 Chapter 2. Background FSm = Σ(c′βR + (N − uβ)R tanφ′) Σ(Wx−Nf)±Dd) for moment equilibrium (2.1) FSF = Σ(c′β cosα + (N − uβ) tanφ′ cosα) Σ(N sinα)−D cosw for force equilibrium (2.2) where the slice base normal force is: N = W + (Xr −Xl)− c′β sinα+uβ sinα tanφ′FS cosα + sinα tanφ ′ FS (2.3) Table 2.1: Statics satisfied and interslice forces in various methods modified from (Krahn, 2003) Method Moment Horizontal force Interslice Interslice Inclination of X/E equilibrium equilibrium normal (E) shear (X) Ordinary Yes No No No No force Bishop’s simplified Yes No Yes No Horizontal Janbu’s simplified No Yes Yes No Horizontal Spencer Yes Yes Yes Yes Constant Morgenstern-Price Yes Yes Yes Yes Variable Corps of Engineers -1 No Yes Yes Yes Inclination of a line from crest to toe Corps of Engineers -2 No Yes Yes Yes Slice top ground surface inclination Lowe-Karafiath No Yes Yes Yes Average of ground surface slope and slice base inclination It is important to note that the entire analysis is done using static forces on as- sumed failure surfaces. LEM methods do not take into account any deformation, strain development, or strain softening within the rockmass. This results in the LEM model not having interslice and slip surface stresses that correctly represent the actual field conditions. The forces are instead false and by their formulation only ensure that force or moment equilibrium is maintained (Krahn, 2007). This results in the model assigning erroneous forces such that each slice has the same global factor of safety (Krahn, 2003). Despite these assumptions and disconnect from actual field and failure conditions, the method works well and gives reason- able factors of safety. 10 Chapter 2. Background 2.2.2 Numerical Models Numerical models are computational methods that attempt to capture the mechan- ical response (e.g. deformation) of a rockmass using a set of material parame- ters, initial conditions, and boundary conditions. Typical output from a numerical model include displacements, strains, and stresses which are then interpreted to determine whether or not the result is in a collapsed state or an equilibrium state (Wyllie and Mah, 2004). Different formulations may allow for numerical defini- tions of equilibrium or failure to be defined based on specific numerical factors like maximum gridpoint velocity or unbalanced force ratio. In order to provide a solution, Pande et al. (1990) states that the numerical model must satisfy the following: 1. Equilibrium 2. Strain compatibility 3. Stress-strain relations of the rockmass (via constitutive models) 4. Boundary considerations Two broad categories exist for numerical models in rock mechanics: continuum models where zones are connected, or discontinuum models where explicitly lo- cated surfaces can separate the zones. An example of both of these can be seen in figure 2.1 showing continuum zoning (a) compared to a discontinuity network (b). Constitutive models are applied to a material (or discontinuity) in order to describe the idealized stress/strain response of the material. These constitutive models can vary in complexity and proper selection is important. Stressed rock materials can exhibit a range of possible behaviours: elastic, plastic, and viscous deformation of intact material. Also, joint surfaces may open, close, and slip allowing individual blocks to undergo rigid body displacement (Bawden, 1993). 11 Chapter 2. Background    UDEC (Version 4.00) LEGEND    28-Aug-08  14:25  cycle        0  time    0.000E+00 sec zones in fdef blocks  4.250  4.350  4.450  4.550  4.650  4.750  4.850  4.950  5.050  5.150 (*10^2)  4.000  4.100  4.200  4.300  4.400  4.500  4.600  4.700  4.800  4.900 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada          (a) Continuum zoning example with connected zones    UDEC (Version 4.00) LEGEND    28-Aug-08  14:24  cycle        0  time    0.000E+00 sec block plot  4.250  4.350  4.450  4.550  4.650  4.750  4.850  4.950  5.050  5.150 (*10^2)  4.000  4.100  4.200  4.300  4.400  4.500  4.600  4.700  4.800  4.900 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada          (b) Discontinuum example showing discrete blocks Figure 2.1: Continuum and discontinuum model examples 12 Chapter 2. Background 2.2.2.1 Continuum Modeling Continuum models assume that the material is continuous throughout the entire model. There can be different materials ascribed to different zones but each is interconnected and fully influences the neighbouring zones. This method is ar- guably more applicable for soils than for rocks; however, it is also used in cases involving weak rockmasses with no dominant structure. The basis of the mod- els also rely on using ‘equivalent continuum’ properties for the material, such as those derived from rockmass characterization procedures, to account for the un- modelled discontinuities. There are no explicitly modeled discontinuities; how- ever, strength anisotropy can be introduced to model preferential weaknesses such as those due to bedding through the use of ubiquitous joint models. Otherwise, equivalent continuum assumptions only account for the role the discontinuities play in weakening the rockmass ignoring discontinuity control on deformation and failure kinematics. If the slope is unstable in a continuum model and the proper failure mechanism is captured, there may be no reason in going to a more complicated discontinuum model unless the failure mode is not what was expected (Wyllie and Mah, 2004; Bawden, 1993). 2.2.2.2 Distinct-Element Modeling In distinct-element modeling, the rockmass is represented as a collection of dis- crete blocks. Interfaces between the blocks (e.g. joints) are considered to be boundary conditions. As such, unlike with continuum methods, discontinuities are explicitly modelled. Movement of the blocks occurs in response to a distur- bance applied at the boundary. A force-displacement law is applied at the contacts to find contact forces from displacements. Newton’s second law resolves the mo- tion of the discrete blocks from the forces acting on them (Hart, 1993; Nichol et al., 2002). Because of this, both movement and energy conservation is achieved within the model. This technique will be described in full detail in Chapter 3 as it is the technique employed in this thesis study. 2.2.2.3 Shear Strength Reduction Method The Shear Strength Reduction (SSR) method is a modeling approach used to re- duce the strength parameters of a material in a numerical model until some stabil- 13 Chapter 2. Background ity threshold is overcome (Dawson et al., 1999; Griffiths and Lane, 1999). This stability threshold can be rather arbitrary, as it has to be defined numerically. In continuum analysis, it is often taken as the point at which the model no longer converges (Griffiths and Lane, 1999; Hammah and Yacoub, 2006; Diederichs et al., 2007). In discontinuum methods, there are several different approaches that have been put forth for stability criteria: maximum unbalanced force, max- imum displacements, maximum number of fallen blocks, etc (Diederichs et al., 2007). A factor called the Strength Reduction Factor (SRF) is applied to Mohr-Coulomb strength parameters as follows (Dawson et al., 1999): c′trial = 1 SRF c′measured (2.4) φ′trial = tan −1 ( 1 SRF tanφ′ ) (2.5) τmax = c ′i trial + σn tanφ ′ trial for the i th material (2.6) where the factor of safety is the SRF when the stability threshold has been over- come. This method can be applied to both continuum and discontinuum models. It has also been applied to more advanced constitutive models such as the Hoek- Brown failure criterion (Hammah et al., 2005b). This allows for the factor of safety to be determined using more realistic modeling methods than the assumptions in a typical LEM analysis allows. The main benefit of this method is the elimination of an a priori failure surface in the model and the inclusion of strain and deformation behavior, via constitutive models that better represent the rockmass, allowing for improved treatment of support interactions and progressive failure (Diederichs et al., 2007). There has been a lot of recent study in the area of SSR applied to finite-element and finite-difference continuum approaches. However, there has not been much study into SSR applied to discontinuum models particularly in cases involving complex failure modes or weak, highly structured rockmasses exhibiting complex behavior requiring the use of selective discontinuities and equivalent continuum blocks. Diederichs et al. (2007) cautions that it is not clear that the SSR technique is applicable to problems that do not have a single sliding surface. 14 Chapter 3 Methodology 3.1 Geomechanical Behavior of Rockmasses and Intact Rock Prediction of slope deformation and instability using numerical technique requires knowledge of the rockmass response to stress-induced changes. This response is heavily site-specific and is a function of the stresses in the rockmass, deformation characteristics, the intact rock strength, mechanical behavior of the intact rock, and any discontinuities. Different modeling methods or stability criterion will ex- plicitly include the discontinuities or downgrade the intact rock strength by some method to account for the weaknesses in order to create an ‘equivalent continuum’ material to represent the fractured rockmass. Rockmasses have three basic modes of failure that can occur separately or be com- bined into a more complex failure. First, the purely structurally-controlled failure, involves planar or wedge surfaces sliding on discontinuities. On the other end of the spectrum are very weak rock masses that fail in a stress controlled, rotational, manner. Shear stresses that develop within the weak slope are the dominant con- trols. The third category of failure mode exists between these two end members and includes stress controlled failure of intact rock bridges between joints to cre- ate a planar, wedge, or step-path sliding failure, and the crushing of intact rock at the base of a high slope. These are also scale dependent (figure 3.1) depending on which portion of the rockmass is being included. Large scale failures may in- clude the rockmass as a whole while bench scale failures may only involve smaller volumes sliding on one or two discontinuities. Stress redistribution occurs rapidly in the creation of a new slope (open-pit mine, rock cut, etc.) whereas a natural slope may be in stress equilibrium until some new external force such as construction or an earthquake introduces new stresses. 15 Chapter 3. Methodology Figure 3.1: Scale effects of rockmass failures and properties, modified from (Hoek, 2007) 16 Chapter 3. Methodology Some amount of stress will result in recoverable elastic strain and the rest will be taken up in the rockmass as irrecoverable plastic strain or movement along discon- tinuities. Structurally controlled failure is a bench-scale problem in the open pit environment with larger, deeper-seated failure involving intact rock only occur- ring when the pit height reaches a height that is capable of imparting significant stresses to the intact rock. In highway slopes, this is often not the case as the slope heights are not as high as their open pit counterparts; however, if the rock is weak, there is potential for even the moderate slope heights encountered in mountainous highway construction to enable deeper seated failure. Idealized constitutive relationships have been developed to describe the range of stress-strain behavior observed in rockmasses under varying stress conditions. Different relationships exist for the elastic strain, yielding, and plastic strain occur- ring as the rock deforms. The proper selection of controlling constitutive model is the first step of modeling a system. This is done with the aid of laboratory testing and fitting linear and non-linear models to the actual laboratory data. The im- portant parameters in laboratory testing of intact rock is the behavior in the elastic range, yielding stress levels, shear-banding and fracturing of the material, and post peak-strength behavior of the failed rock (Brady and Brown, 1993). Figure 3.2 be- low show the different response between brittle behavior (hard rock) and ductile behavior (weak rock). The two extremes are often represented by the difference in behavior between crystalline rocks exhibiting brittle behavior and clay-rich rocks exhibiting ductile behavior, with the rest being some combination of both (Brady and Brown, 1993). 3.1.1 Elasticity A material is said to behave in an elastic manner if the strain that occurs during deformation is completely recoverable during unloading. Equations of elasticity have been developed for the stress-strain behavior of intact rock. Stress and strain are related to each other, for isotropic solids, in a linear elastic manner by an inversion of Hookes Law (Brady and Brown, 1993): [σ] = [D][] (3.1) where [σ] is the stress matrix, [D] is the elastic stiffness matrix, and [] is the resulting strain. This equation can further be expanded for isotropic materials 17 Chapter 3. Methodology Figure 3.2: Stress strain relationships of varying rockmass quality, from Hoek (2007) 18 Chapter 3. Methodology to:  xx yy zz γxy γyz γzx  = 1 E  1 −v −v 0 0 0 −v 1 −v 0 0 0 −v −v 1 0 0 0 0 0 0 2(1 + v) 0 0 0 0 0 0 2(1 + v) 0 0 0 0 0 0 2(1 + v)   σxx σyy σzz σxy σyz σzx  (3.2) Both Youngs Modulus (E) and Poissons ratio (ν) are needed to solve the above matrix which uses Hookes Law. This is the basis for models for stress and strain behaviour in the elastic region. 3.1.2 Plasticity Plastic strain occurs once the material has reached its yield point (or elastic limit) and can no longer sustain any additional stress without incurring irrecoverable deformation. Once this point is reached, there are many different behaviours that can be followed including: perfectly plastic, strain hardening, strain softening, visco-plastic, etc. The theory of plasticity has been around since Coulomb studied plastic behaviour in 1773. Centuries of further development has occurred, with a strong emphasis on the plasticity of soils in the fifties and sixties. Simple rheological models have been developed which involve masses, elastic springs and visco-plastic dashports (for time-dependent problems) to develop mathematical models to mimic real ma- terial behaviour. Pande et al. (1990) lays out the four basic tenets of plastic the- ory: 1. A stress-strain relationship before yielding: typically elastic, but not always 2. Yield criterion: The yield point in a uniaxial case defined by: F (σ) = 0 (3.3) Where σ is the vector of stresses such that any stress situation within the surface F (σ1, σ2, σ3) < 0 (3.4) 19 Chapter 3. Methodology is an elastic state, whereas stresses on the surface: F (σ1, σ2, σ3) = 0 (3.5) are in a plastic state. 3. Flow rule, which determines the direction of plastic strain increment. This flow rule can either be associated where the plastic potential function is the same as the yield function, or non-associated where it is not. Typically, associated flow rules predict greater dilation angles than are found in the real world; however, associated flow rules cause stress and velocity to co- incide allowing for a better comparison to traditional LEM solutions than a non-associated flow rule (Pande et al., 1990; Griffiths and Lane, 1999; Zienkiewicz et al., 1975). UDEC uses a non-associated plastic flow rule, for which dilation angles are set to zero as default. 3.1.3 Yield Criterion Yield criterion differs from failure criterion by being more than just a method to directly compare strengths and load. Instead of defining an envelope that separates admissible (unfailed rock) and inadmissible (failed rock) stress states, a yield cri- terion is used to define in the model when plastic yielding occurs. The assumption is made that rock fails in shear, thus most yield criterions are plotted in 2D and that the highest shear stress occurs between σ1 and σ3, with some form of tensile cutoff. However, each site-specific rock has it’s own specific mathematical rep- resentation of it’s strength surface. The yield criterions presented below are an attempt to generalize these strength parameters to be used on a wide range of ma- terials from accepted material strength tests done on samples from site (Franklin and Dusseault, 1989). 3.1.3.1 Mohr-Coulomb Mohr-Coulomb is a simple linear failure criterion used for soils and rock. How- ever, it is one of the least accurate representations of rock behaviour. It is typically used in soil mechanics (Franklin and Dusseault, 1989). A material is a Coulomb material where the shear strength of the sliding surface is expressed by two terms: 20 Chapter 3. Methodology cohesion (c) and internal friction angle (φ) (Coulomb, 1773). The peak shear stress strength is defined by the following equation: τ = c+ σntanφ (3.6) where σn is the normal stress acting on the failure plane The criterion describes a yield envelope, assumed to be linear, that is tangential to a series of Mohr circles, as seen in figure 3.3. These Mohr circles are created from different configurations of principal stresses at failure from triaxial testing. Tensile strength is not sufficiently covered by the linear envelope, and is overestimated when fit to the data extending from the compressive region. This prompts the use of some minimal tension cutoff value (Franklin and Dusseault, 1989). Failure will take place on any plane when the shear stress acting across that plane reaches some value of critical shear strength. Figure 3.3: Graphical representation of the Mohr-Coulomb failure envelope, mod- ified from (Hudson and Harrison, 2005) 21 Chapter 3. Methodology 3.1.3.2 Hoek-Brown The assumption that the failure envelope is linear is not always true. Hoek and Brown developed an empirical criterion in which shear strength is represented as a curved Mohr envelope. Hoek and Brown derived the criterion from Griffith crack theory and many observations of rock behaviour in laboratory settings. Fur- ther refinements seeking to link the empirical criterion to geological observations resulted in the Generalized Hoek-Brown Strength Criterion (Hoek et al., 2002), which when expressed in terms of major and minor principal stresses is: σ′1 = σ ′ 3 + σci ( mb σ′3 σci + s )a (3.7) where σci is the uniaxial compressive strength (UCS),mb is a reduced value of the empirical material constant mi for intact rock given by: mb = miexp ( GSI− 100 28− 14D ) (3.8) where mi is found from charts which give typical back analysis values, D is a disturbance factor, and s and a are constants for the rock mass given by: s = exp ( GSI− 100 9− 3D ) (3.9) a = 1 2 + 1 6 (e− GSI 15 − e−20−3 ) (3.10) GSI is the Geological Strength Index discussed later in section 3.1.4.2. The UCS for the rockmass, σc above, can be found by setting the principal stress (σ3) to zero in the above equation giving: σc = σci ( s )a (3.11) Similarly, the tensile strength (σt) can be found by setting σ′1 = σ ′ 3 = σt in the same equation giving: σt = −sσci mb (3.12) 22 Chapter 3. Methodology Geotechnical software is generally written in terms of Mohr-Coulomb failure cri- terion, necessitating the need to be able to determine equivalent Mohr-Coulomb parameters from Hoek-Brown parameters. The process involves fitting an aver- age linear relationship to the Hoek-Brown curve over a specified range of minor principal stresses. The following equations are used to find φ′ and c′, to generate an approximate envelope as seen in figure 3.4: φ′ = sin−1 [ 6amb(s+mbσ ′ 3n) a−1 2(1 + a)(2 + a) + 6amb(s+mbσ′3n)a−1 ] (3.13) c′ = σci[(1 + 2a)s+ (1− a)mbσ′3n](s+mbσ′3n)a−1 (1 + a)(2 + a) √ 1+(6amb(s+mbσ ′ 3n) a−1) (1+a)(2+a) (3.14) where σ3n = σ′3max/σci The disturbance factor (D) ranges between 0 for undisturbed rockmasses and 1 for very disturbed rockmasses. The disturbance is a means to correct for blast damage and stress relaxation. In the construction of large slopes the rockmass may be exposed to heavy blast damage and stress relief due to the removal of overburden to create the benches and ultimately, the pit. Hoek cautions against using D = 0 in these situations as it results in non-conservative strength properties (Hoek et al., 2002). However, selection of this parameter can be difficult as it is heavily subject to user experience and site conditions. A guide, figure 3.5, provides a starting point for values of disturbance factors based on a variety of slope conditions. Numerous authors (namely Carvalo et al. (2007); Mostyn and Douglas (2000)) have pointed out limitations in the Hoek-Brown failure criterion regarding: poor matches to the behaviour of weak rocks, incorrect estimation of the tensile limit, and the mi and σc parameters not being material constants particularly when the a parameter is commonly fixed at 0.5. 3.1.4 Rockmass Classification Schemes Rockmass classification schemes provide a way to characterize, both qualitatively and quantitatively, the rockmass in order to estimate their strength and deformabil- ity. Different schemes can vary widely for assumptions and inclusions of certain properties or considerations. Two main schemes are used in slope stability: the 23 Chapter 3. Methodology Figure 3.4: Approximate M-C failure envelope from H-B parameters, modified from (Hoek et al., 2002) 24 Chapter 3. Methodology Figure 3.5: Disturbance factor guidelines for slopes, modified from (Hoek et al., 2002) 25 Chapter 3. Methodology Rock Mass Rating system (Bieniawski, 1976, 1989) and the Geological Strength Index (Hoek et al., 1995, 2002). 3.1.4.1 Rock Mass Rating The Rock Mass Rating (RMR) system was developed empirically by Bieniawski (1976). Like most classification schemes, it has evolved over time to be quite different than the original system. Two different versions of RMR are commonly used: RMR76 and RMR89. The differences between the two are mainly due to joint descriptions and ratings (Bieniawski 1989). The original scheme was developed using data mainly from civil engineering excavations in sedimentary rock (Brady and Brown, 1993; Bieniawski, 1976). The rating is based off six parameters: 1. Intact rock strength 2. Rock Quality Designation (RQD) 3. Discontinuity spacing 4. Discontinuity conditions 5. Groundwater conditions 6. Orientation factors For various ranges in input parameters, different rating values are assigned. All the ratings are summed to give the RMR value. Not all parameters have the same influence on the total rating. The rating system charts can be found in Bieniawski (1976). 3.1.4.2 Geological Strength Index The Geological Strength Index (GSI) is a rating system based on a range of struc- tural conditions and joint surface conditions of the rockmass. It was developed to be used with the Hoek-Brown failure criterion discussed earlier in section 3.1.3.2. The two simple parameters provide field engineers with a method to describe a wide range of rock mass types and structural conditions (Hoek et al., 1998) with a focus on geological controls. However, the classification does not work for all conditions likely to be encountered such as extremely weak and laminated/sheared rockmasses. GSI charts can be found in most papers published by Hoek on the 26 Chapter 3. Methodology Hoek-Brown failure criterion, such as Hoek et al. (1998). GSI is based on the same first four parameters as RMR and thus RMR can be used to estimate the GSI using one of the following: RMR′76 = GSI (3.15) RMR′89 = GSI− 5 (3.16) This is possible only if the other parameters in RMR are set to specific values, as outlined in (Hoek et al., 1998). Namely, the groundwater conditions are neglected and set to the maximum rating value, and no adjustment is made for joint orien- tation. It is expected that these influences will be treated separately in a stability analysis. There is criticism that the GSI method requires extensive engineering experience in order to properly apply it in the field. More quantitative methods have been applied to the GSI method (such as those presented in Cai et al. (2004)); however Hoek cautions against this as these quantifications do not work well in situations where the structural fabric is destroyed or difficult to interpret (Marinos et al., 2005). It is important to note that GSI is based on the assumption that the rock- mass behaves as an isotropic mass. The GSI system should not be applied to rockmasses with a clear dominant structural control. Also, the GSI system re- quires rockmass exposure and the use of GSI with borehole data is discouraged. Numerous studies have been done to extend GSI to be used for deep, strong rocks (Carter et al., 2007) and low strength rocks (Marinos et al., 2005; Carvalo et al., 2007). 3.2 Behaviour of Discontinuities In slope stability problems involving rockmasses, most failures will have some form of movement along discontinuities. The failure may even be entirely struc- turally controlled. Conditions for movement along these discontinuities are con- trolled by their shear strengths, which can be tested much like intact rocks, to provide a failure envelope. Shear and normal stiffnesses can also play a role in distribution of stresses and displacements in a rock mass (Brady and Brown, 1993). 27 Chapter 3. Methodology Joint roughness is an important factor with rougher joints typically being stronger than smooth joints due to the dilative response needed to move up and over asper- ities along the joint surfaces or the force required to shear through any asperities on the joint surface. After some amount of strain, a joint may show peak residual behaviour as asperities are sheared off and mechanical grinding has made the joint smoother. However, determining roughness is difficult as roughness is a scale de- pendent attribute. Also, testing of field scale joints (those controlling stability within the rockmass) is exceptionally difficult. True shear performance of the dis- continuity is controlled by the surface roughness, rock strength at joint surface, applied normal force, and amount of shear displacement (Wyllie and Mah, 2004). Because of this, various empirical strengths have been created to account for dila- tency due to surface roughness, such as the one by Barton (1973), who proposes a peak shear strength defined by the following equation: τjoints = σ ′tan ( φr + JRClog10 ( JCS σ′ )) (3.17) Where φr is the residual friction angle, JRC is the joint roughness coefficient ranging from 1 (smoothest) to 20 (roughest), JCS is the compressive strength of the rock at the fracture surface, and σ′ is the effective normal strength. However, capturing this dilative response in a numerical model is difficult and requires many model inputs and thus is rarely used. Infilling of the discontinuity is also important as it can modify the shear properties of the joint. The peak strength envelope, for many infilled discontinuities, is found between that for the filling and that for the clean discontinuity (Brady and Brown, 1993). Also, if the thickness of the infilling is more than 25-50% of the amplitude of asperities in the discontinuity, then the properties of the discontinuity will be the properties of the infilling (Wyllie and Mah, 2004). Joint behaviour can be quite complicated, with strengths changing with confinement, displacement, varying stiffnesses, dilative effects, etc. 3.3 Discontinuity Network The stability of large rock slopes is heavily influenced by the structural geology of the rockmass. Thus, it is essential that the geometrical properties of the discon- 28 Chapter 3. Methodology tinuity network are carefully described, in addition to their geotechnical strength properties. The goal of the geological mapping program is to define the set (or sets) of controlling discontinuities for a particular area. One issue for a geological mapping program is to determine how many disconti- nuities need to be mapped, and at what scale, in order to properly define the design input sets. A variety of factors influence this decision which requires a wide range of properties to be examined. For example, a typical approach may involve: • Block size, RQD, or some other measure of blockiness of the rockmass • Number and type (joint, shear, fault) of joint sets, with the following for each set: – Orientation - strike and dip or dip/dip–direction – Spacing - between discontinuities of the same set – Persistence - a measure of the continuous length of the discontinuity – Roughness - JRC or some other metric – Infilling/aperture Properties – Seepage conditions of the joint sets as well as the overall slope Some of these parameters are easy to note on surface outcrops. Persistence, how- ever, can be difficult to measure as only a small part of the discontinuity may be present on a particular face (Wyllie and Mah, 2004). Persistence is important as it is coupled with orientation and spacing to define the fracture network indicating the size of the blocks that can slide from the slope face and the amount of intact rock bridges which need to break in order for failure to occur. Intact rock bridges defined by the persistence can have a significant strength influence on the stability of the slope (Wyllie and Mah, 2004). Some, but not all, of these parameters can be found from orientated drilling; however, scale effects do apply. Structural geology is also complex with much natural variation. Dip and dip di- rections can vary significantly within a single set, as can spacing and persistence. This makes the definition of the true, insitu, fracture network impossible due to the lack of information, even with detailed geotechnical drilling. Probabilistic meth- ods have been applied to these parameters and complicated statistical methods have been developed to account for the likely range of values to be expected. 29 Chapter 3. Methodology This data may be sourced from a field mapping program such as line or window mapping programs coupled with survey data. This allows for the spatial distribu- tion of joints to be plotted on standard geological maps and contoured pit maps. Stereographic projections can be used to determine the general variability and number of sets for a particular area. Coupled with observations of persistence and spacing, larger features can be found by trace mapping. 3.4 Limit Equilibrium Limit equilibrium work in this thesis has been carried out using SLIDE version 5.0 by Rocscience. It is a 2D method of slices limit equilibrium slope stability package with probabilistic and deterministic factor of safety modes, as well as the ability to include steady-state groundwater analysis using a finite-element solver. It provides factor of safety using the method of slices by the typical methods: ordinary, Bishop simplified, Janbu simplified, Janbu corrected, Spencer, Corps of Engineers #1 and #2, Lowe-Karafiath, and Morgensten-Price. The method and it’s formulation were described in detail in section 2.2.1. 3.5 Universal Distinct Element Code UDEC (Universal Distinct Element Code) by Itasca Inc. (2007) is a 2D modeling code which combines two different modeling approaches into one commercial package (Havenith et al., 2003). The code has been developed to include the distinct element method for modeling the behavior of distinct interacting blocks, and a finite-difference method to determine the deformation of the intact blocks (Hart, 1993). The version used for modelling in this thesis is version 4.0. A problem is divided into a number of blocks that are separated by joints in the model. These joints (interfaces) are treated as a boundary condition in addition to the boundary conditions set at the extents of the model. Figure 3.6 shows the computational loop and equations used during solving in UDEC. The contacts between each block are considered to be soft-contacts to treat the relative normal displacements at the block contacts (Hart, 1993; Nichol et al., 2002). 30 Chapter 3. Methodology Figure 3.6: Computational cycle used in UDEC, from (Itasca, 2004c) 31 Chapter 3. Methodology Both force and moment equilibrium are achieved through the use of a force- displacement law at joint contacts and Newtons second law of motion on the blocks. The equations of motion are solved using an explicit time-marching code with the caveat that the time step is sufficiently small such that disturbances do not propagate between elements. Derivations of the equations of motion and energy can be found in Itasca (2004c). A joint is represented in the model as a contact between two blocks. This contact surface is comprised of data elements that represent point contacts. There is a numerical problem with this representation of joints that can cause blocks to be- come locked within the model. This is because the blocks cannot be broken by the creation of new fractures. All fractures must intersect to form discrete blocks, as seen in figure 3.7. This makes it difficult to capture the true fracture network as persistence and continuity of joints is rarely this fully developed. In addition, cor- ners, which concentrate stress, are not crushed but remain as points. Rounding the corners of the blocks allows for the blocks to slide past one another (Hart, 1993). The contacts between blocks are updated when block motion occurs. When the displacements are large, the code ensures that forces are the same when contacts are added or deleted as the block moves to interact with a new block. Joint stiffnesses (ks and kn) should be selected carefully as they control an impor- tant mechanical response. The time it takes to complete a cycle in the model is also controlled by ks and kn. Generally, ks and kn are selected to be less than ten times the equivalent stiffness of the stiffest zone, or: ks and kn ≤ 10.0 ( max [ K + 4 3 G ∆zmin ]) (3.18) where K and G are bulk and shear moduli, and ∆zmin is the smallest zone in the model (Itasca, 2004c). Joints in the model can be set by hand by specifying end points for a line or using the UDEC JSET joint generator with the following syntax where m is the mean and d is the standard deviation: JSET am, ad tm, td gm, gd sm, sd x0, y0 32 Chapter 3. Methodology    UDEC (Version 4.00) LEGEND     3-Sep-08  10:48  cycle        0  time    0.000E+00 sec block plot  4.530  4.550  4.570  4.590  4.610  4.630  4.650  4.670  4.690 (*10^2)  4.250  4.270  4.290  4.310  4.330  4.350  4.370  4.390  4.410 (*10^2) JOB TITLE :  UDEC parametric studies - 10 m dip parallel University of Brish Columbia Vancouver, BC Canada (a) Real joint network, with joints terminating within a block    UDEC (Version 4.00) LEGEND     3-Sep-08  10:49  cycle        0  time    0.000E+00 sec block plot  4.530  4.550  4.570  4.590  4.610  4.630  4.650  4.670  4.690 (*10^2)  4.250  4.270  4.290  4.310  4.330  4.350  4.370  4.390  4.410 (*10^2) JOB TITLE :  UDEC parametric studies - 10 m dip parallel University of Brish Columbia Vancouver, BC Canada (b) Same jointing as modelled in UDEC with in- complete joints joints removed Figure 3.7: Difference between actual jointing and modelled joints 33 Chapter 3. Methodology Where a is the angle of the joint, t is the trace length of the joint, g is the gap length between segments, s is the spacing normal to the joint and x0, y0 is the start of one joint trace (figure 3.8). Figure 3.8: Nomenclature for JSET joint generator in UDEC, from (Itasca, 2004a) The blocks in the model can either be rigid or deformable. If deformable, each block is discretized with triangular or quadrilateral constant strain elements. Quadri- lateral elements provide a better solution to plastic deformations (Itasca Inc., 2007). The block constitutive model, such as elastic or elasto-plastic, then con- trols the block behavior. Vibrations and force imbalances are dissipated by interblock sliding or internal block deformation. However, this does not occur in the elastic portions of the model which introduces the need for artificial damping in the model. There are two types of damping that may be used separately or in combination (Hart, 1993): • Mass-proportional damping (viscous damping) applies a force which is pro- portional to the mass velocity but in the opposite direction. This, typically acts on lower frequency vibrations and is described by the following equa- 34 Chapter 3. Methodology tion: {d} = −α[M ]{u} (3.19) where {d} is the mass damping force, [M ] is the mass matrix, {u} is the velocity matrix, and α is the mass proportional damping factor • Stiffness-proportional damping applies to higher-frequency vibrations and is described by the following equation: {s} = β[k]{u} (3.20) where {s} is the stiffness damping force, [k] is the stiffness matrix, {u} is the velocity matrix, and β is the stiffness proportional damping factor Selecting the damping factors (α, β) is not trivial for a complicated multi-degree of freedom system. The factors are chosen to be a fraction (usually 2–5%) of the critical damping: λ = 1 2 ( α wn + βwn ) (3.21) Alternatively, the λ value is adjusted as a fraction of the rate of energy displace- ment and rate of change of kinetic energy in the model. There is also the modeling dilemma of how many discontinuities to include. More joints in the model adds more calculations per timestep, which can significantly increase solution times. It is impractical and often impossible to model every single joint in the rockmass. Instead, some tradeoff between a practical amount of explicitly defined joints and the size of the ‘equivalent continuum’ in between must be made. However, the relationship between the modeled joint density and equivalent continuum strength degradation of intact rock is not well defined. 3.5.1 Constitutive Models This section outlines the constitutive models used in the analysis throughout the thesis. It is unclear how the SSR method would be applied to complicated consti- tutive models that include joint closure, dilation, continuously yielding materials, etc., only simple constitutive models were used that require simple strength pa- rameters such as joint and intact friction angle, cohesion, and tension. 35 Chapter 3. Methodology 3.5.1.1 Intact Rocks As mentioned earlier, the Mohr-Coulomb constitutive model is used for the yield criterion. In UDEC, this is implemented as a Mohr-Coulomb shear criterion with a tension cutoff. The plastic flow rule is non-associated and the tensile flow rule is associated (Itasca, 2004b). The failure envelope is defined as: f s = σ1 − σ3Nφ + 2c √ Nφ (3.22) with a tensile failure envelope: f t = σt − σ3 (3.23) where: Nφ = 1 + sinφ 1− sinφ (3.24) and: σtmax = c tanφ (3.25) The non-associated flow rule uses a shear potential function defined as: gs = σ1 − σ3Nψ (3.26) where: Nψ = 1 + sinψ 1− sinψ (3.27) An associated flow rule for the tensile potential function such that: gt = −σ3 (3.28) A full implementation of the Mohr-Coulomb constitutive model can be found in Itasca (2004b), along with plastic corrections made to the model and calculation of strain increments for each timestep. 36 Chapter 3. Methodology 3.5.1.2 Joints In nature, a large part of the rockmass behaviour is controlled by discontinuities. There are many different possible behaviours for these discontinuities, and thus many possible different joint constitutive models to be used. In this study, only elasto-plastic models of joints have been considered, as the SSR formulation for more complicated joint models with many different non-linear parameters has not been developed yet. In these types of models, the elastic regime is captured by the normal and shear stiffnesses (kn and ks respectively). Peak strength and dilative behaviour is con- trolled by the use of a failure criterion and flow rule (Pande et al., 1990). The joint constitutive model used in this thesis is a Coulomb-slip yield criteria which includes a joint tensile strength (T ). If the joint tensile strength is exceeded, then the model sets Σn to zero. The shear response of the joint is controlled by two different responses depending on the amount of shear stress (τs) on the joint such that if: |τs| ≤ c+ tanφ′ = τmax (3.29) then ∆τs = Ks∆u e s (3.30) else if |τs| ≥ τmax (3.31) then τs = sign(∆us)τmax (3.32) where ∆ues is the elastic component of incremental shear displacement, and ∆us is the total incremental shear displacement (Hart, 1993). Stress-displacement is in the normal direction and assumed to be linear with a stiffness kn such that: ∆σn = Kn∆un (3.33) 37 Chapter 3. Methodology where ∆σn is the effective normal stress increment and ∆un is the normal dis- placement increment (Hart, 1993). Additionally, joint dilation can occur. A basic displacement-weakening response can be achieved by setting the tensile and cohesive strength to zero when the ten- sile or shear strength is overcome (Hart, 1993). There are drawbacks to the simple coulomb-slip model presented. The problem is that the value of cohesion and fric- tion angle selected depend on joint confinement which can vary widely within the rockmass, making the selection of a single value for these two parameters diffi- cult. More advanced models such as peak-residual behavior, Barton-Bandis joints, and user-defined joint constitutive models are also available in UDEC. However these are not discussed as the application of the SSR method using non-linear constitutive models is outside the scope of this thesis. 3.6 Shear Strength Reduction Using UDEC The choice of criterion to describe a failed or unstable slope in SSR in discon- tinuum models is important. There are purely numerical schemes such as using the ratio of maximum unbalanced force, to more interpretative methods such as checking for abrupt changes in displacement as strengths are reduced. In this the- sis, both the cohesion and friction have been factored by the same amount, at the same time. This assumes that the failure mechanism is shear only and does not incorporate more complex brittle failure of intact rock. Diederichs et al. (2007) mentions in passing that the validity of equal factoring of these two parameters needs to be debated, particularly in light of work by Martin and Chandler (1994) indicating that frictional strengths are only mobilized after a significant loss of cohesive strength. However, determining the validity of equal factoring of these two parameters is beyond the scope of this thesis. In this thesis, all six frictional strength properties (friction angle, cohesion, tension for both intact rock and dis- continuities) are all included in the SSR routine, and all factored equally. The built-in UDEC factor of safety calculator (UDEC-SSR) only allows for solv- ing by use of a ratio limit which is defined as the average unbalanced mechanical force magnitude divided by the average applied mechanical force magnitude for all gridpoints in the model (Itasca, 2004a). This built in method for strength reduc- tion is only applied to materials using a Mohr-Coulomb based constitutive model. It allows for the user to select which strength parameters (friction, joint friction, 38 Chapter 3. Methodology etc.) and/or which materials will undergo strength reduction. The method brack- ets the factor of safety by using whether or not a given SRF stage is stable or unstable to decide which SRF stage to investigate next. It uses the following test for stability within UDEC (Itasca, 2004c): • A representative number of steps (Nr) is found which is the number of steps it takes the model to return to equilibrium after increasing the cohesion value and changing the internal stresses • For a given SRF, the model properties are set and Nr steps are executed • The force ratio is checked – if it is less than 10−5, then the system is considered to be in equilibrium and the loop exits as stable (to the next SRF stage if applicable) • if it is more than 10−5, then another Nr steps are executed and the average mean force ratio over the current Nr steps is compared with the mean force ratio over the previous Nr steps – If the difference between these two averages is less than 10%, the sys- tem is in non-equilibrium and the loop exits as unstable – If the difference is greater than 10%, blocks of Nr steps are run until one of the following: 1. The force ratio is less than 10−5 (stable) 2. The difference is less than 10% (unstable) 3. Six blocks of Nr have been run (unstable) Note that the above force ratio limit of 10−5 can be changed to any user selected value. The automated process uses a bracketed search to find the minimum SRF value that leads to non-equilibrium, in a scheme similar to divide and conquer (quicksort) algorithm in computer science which can be seen in the flowchart in figure 3.9 39 Chapter 3. Methodology Stage 1 Stage 2 Stage 3 Initial model geometry Elastic gravity loading Plastic properties Removal of benches SSR initial state Determine representative number of steps (N) for the characteristic response time of the system Set SRF = 1.0 and keep halving until the first stable case is found (Fs) Keep doubling SRF until first unstable case is found (Fu) Steps to test for stability or instability 1.  Do up to N steps and record unbalanced force ratio      (Ru) 2.  If Ru falls below 0.001 during stepping, exit as STABLE 3.  If |avg (Ru) - avg(Ru_old)| / avg (Ru) < 0.1, exit as      UNSTABLE 4.  If total iterations (steps 1-3) > 6, exit as UNSTABLE 5.  Set Ru = Ru_old and goto 1 Set F = (Fu + Fs) / 2 If stable Fs = FIf unstable Fu = F EXIT IF Fu - Fs < 0.005 CHECK CHECK CHECK Figure 3.9: Flowchart for UDEC-SSR showing stability criteria and bracketing search, modified from (Itasca, 2004c, 2008) 40 Chapter 3. Methodology 3.6.1 Displacement-Based SSR An important consideration for application of SSR to open pit stability problems is to incorporate components that are familiar and allow confidence in the result to be gained. For example, given the heavy use of geodetic monitoring of displacements using prisms mounted at strategic locations on a pit slope of interest, it makes sense to define failure in a SSR analysis using the same early warning trigger used in pit-slope monitoring, that is, sharp accelerations in slope displacements. Monitoring points are set in the model and displacements are tracked through- out the strength reduction loop. The point of failure is selected by the user from interpreting results such as those seen in figure 3.11. The point of failure (and cor- responding SRF) is chosen where the slope of the line rapidly changes from little displacement magnitudes to large displacements. This point is typically where lo- calization occurs and a full shear band has developed and significant deformations begin to occur with a rapid increase in displacement values with small increments in strength reduction factor. In this thesis, a monitoring point at the crest was used to track displacements. This strength reduction procedure has been implemented in UDEC using the in- ternal FISH programing language. An example of the strength reduction code can be seen in Appendix A. However, monitoring points are case specific and must be selected by the user at strategic locations (e.g. crest and toe) of the modelled slope. The reduction method always returns to an initial stable state before applying each successive strength reduction step such that the accumulated deformations of the previous reduction step are discounted (i.e. reset to zero). Thus, the strength reduction sequence looks like: initial→ ssr1+x → initial→ ssr1+2x → initial→ ...→ initial→ ssr1+(n)x Within this process, there are two different ways to cycle the model: automatically until a defined force ratio is reached (10−5 or some other value) or manually to a set number of steps. Using the force ratio option, each point is ultimately cycled a different amount with increasing number of cycles needed to reach the unbalanced force threshold the further the solution gets from SRF = 1. When cycling the model using a set number of cycles, the user must select a value and then ensure that each strength reduction step has reached equilibrium before moving on to 41 Chapter 3. Methodology the next reduction step. A conceptual code overview can be seen in Figure 3.10, detailing the displacement-based (D-SSR) strength reduction in UDEC. Stage 1 Stage 2 Stage 3 Initial model geometry Elastic gravity loading Plastic properties Removal of benches SSR initial state SSR 1+1(x) SSR 1+2(x) SSR 1+... SSR 1 + (n-1)(x) SSR 1 + (n)(x) Generate displacement charts (Results) - set initial strengths, SRF interval, #   of stages - define function: SRF_LOOP - write external loop to text file based   off of # of stages (n) and interval (x) - save int.sav with the above   parameters defined such that they do   not get deleted during restoring   previous stages of the model SRF_LOOP SRF_LOOP SRF_LOOP SRF_LOOP SRF_LOOP For each SRF_LOOP: - restores initial state (int.sav) - monitoring points defined (constant throughout model) - gets current SRF stage 1 + (n_current)(x) from external   text file - calculate new strength values based off of SRF stage for                    > different intact block materials                    > different joint materials - check exclusion flags - apply new strengths to the model (except exclusions) - step the model based on solution mode selected - save model state - read previous stages displacement values into memory   from file (for every stage but the first) - write displacement values + new displacement values to   file (this preserves data between initial state restoration) - exit SRF_LOOP Figure 3.10: Flowchart for general displacement-based SSR method This method requires interpretation by the user to select the point of failure rather 42 Chapter 3. Methodology than relying on an automatic numerical threshold such as the method implemented in UDEC by Itasca. The selection criteria used throughout this thesis is to take the first instance of significant deformation as the point of failure (i.e. the con- servative limit of movement). The method does not include any pre-determined range selection and requires the user to manually select the upper and lower SSR boundaries as well as the increment increase for each reduction step, in addition to choosing whether to automatically cycle the model to some ratio limit or to require the user to determine equilibrium at each stage. Thus, the D-SSR method devised here generally takes a longer time to determine the factor of safety value than UDEC-SSR. However, this method returns more data to the user by providing a save of each individual strength reduction stage rather than just the final stage as provided with UDEC-SSR. Also, as previously noted, it returns this data based on a more intuitive definition of failure (figure 3.11). It should be noted that extreme care must be taken when using any default and/or automated procedure, as certain assumptions are implied that may be valid for most cases but may not necessarily for all. 43 Chapter 3. Methodology Strength Reduction Factor D i s p l a c e m e n t Point of failure selected where the displacements begin to rapidly increase 1 1 + x Figure 3.11: Conceptual figure showing displacement response and selection of strength reduction factor indicating instability 44 Chapter 4 Parametric Study Given the relative lack of experience in applying the SSR method to discontin- uum problems, coupled with a potentially better SSR stability criteria within the commercial code UDEC, a parametric study was carried out to gain confidence in the method and the different controls. This was done to further refine the D-SSR technique and compare it to the UDEC-SSR criteria and to see how the method behaves conceptually on three basic slope configurations. General sensitivities of various input parameters and model setup options such as mesh density were carried out on one or more of the three distinct joint patterns defined in Table 4.1. Table 4.1: Slope geometries used in parametric study Geometry Dip Joint Set 1 Dip Joint Set 2 Daylighting joints -35o 45o Discontinuous daylighting joints -40o 45o Slope parallel joints -45o 45o All three of these base models share the same common geometry but differ with their respective joint networks. The slope is 250 m high and is dipping at 45o. Each model is also zoned in terms of jointing and meshing. This height was chosen as a typical depth of a small open pit as well as being constraining the problem in size for computational efficiency. A typical jointing and mesh setup, using triangular elements, can be seen in figure 4.1, in this case with 10 m spaced daylighting joints. This decrease in joint density and increase in mesh size away 45 Chapter 4. Parametric Study from the area of interest is for numerical efficiency during the computationally intensive SSR routines. The three different joint configurations represent three common kinematic fail- ure modes: the daylighting joint case is the most structurally controlled case as all the deformation is expected to occur as sliding along a critical discontinuity with rock block parameters playing little role; the discontinuous daylighting case represents a more complex rockmass failure that includes yielding of intact rock bridges as well as movement along key discontinuities; the slope parallel case lies between the other two extremes with expected crushing and breakout at the toe combined with movement along critical discontinuities parallel to the slope. The three examples of the different slopes can be seen in figure 4.2. Both the discontinuous daylighting and slope parallel joint cases share the material parameters found in table 4.2. This is a fictitious rockmass based on a low strength mudstone–sandstone with a GSI of 50 and UCS of 25 MPa, with a confinement valued calculated from a slope height of 250 m and uniform density of 2.6 kg/m3. Cohesion and friction values have been modified for the daylighting joint case in order to have a slope that is stable at the point of initialization. These changed parameters can be found in table 4.3. Joint shear stiffnesses are based off of repre- sentative values in Kulhawy (1975). Varying of initial parameters, covariance of parameters, different assumptions regarding Hoek-Brown equivalent continuum strength reductions, and ranges of Mohr-Coulomb parameters from Hoek-Brown parameters are beyond the scope of this thesis. Material factors other than the selection of initial parameters are the focus of this thesis. Joint tension has been included in order to include the SSR process on that variable, such that all 6 prop- erties (friction, cohesion, and tension for both intact rock and discontinuities) are factored equally. A range of different sensitivities are investigated in this chapter. Some of the sen- sitivities involve material parameters while others are related to numerical mod- elling setup options or D-SSR options. Many of the numerical modelling setup cases are analyzed using the 10 m discontinuous joint setup as it is the most gen- eral rockmass failure type (i.e. including shear through intact rock and slip along discontinuities) of the three cases being investigated. There are some cases where the setup options or D-SSR options are also run in conjunction with a material parameter sensitivity analysis. In these cases, the results are shown in the material parameter section rather than the general model setup or D-SSR option sensitivity sections, for clarity. 46 Chapter 4. Parametric Study    UDEC (Version 4.00) LEGEND    23-May-08  17:07  cycle    71470  time    1.935E+01 sec block plot -3.000 -1.000  1.000  3.000  5.000  7.000 (*10^2)  0.100  0.300  0.500  0.700  0.900  1.100 (*10^3) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada          (a) Joints    UDEC (Version 4.00) LEGEND    23-May-08  17:07  cycle    71470  time    1.935E+01 sec zones in fdef blocks -3.000 -1.000  1.000  3.000  5.000  7.000 (*1 ^2)  0.100  0.300  0.500  0.700  0.900  1.100 (*10^3) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada          (b) Mesh Figure 4.1: Example of scoping of joints (a) and meshing (b) within slope for computational efficiency 47 Chapter 4. Parametric Study    UDEC (Version 4.00) LEGEND    23-May-08  17:07  cycle    71470  time    1.935E+01 sec block plot  1.750  2.250  2.750  3.250  3.750  4.250  4.750  5.250  5.750 (*10^2)  3.000  3.500  4.000  4.500  5.000  5.500  6.000  6.500  7.000 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada (a) Continuous daylighting    UDEC (Version 4.00) LEGEND    20-May-08  13:44  cycle   137100  time    6.002E+01 sec block plot  1.750  2.250  2.750  3.250  3.750  4.250  4.750  5.250  5.750 (*10^2)  3.000  3.500  4.000  4.500  5.000  5.500  6.000  6.500  7.000 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada (b) Discontinuous daylighting    UDEC (Version 4.00) LEGEND    22-Apr-08  15:08  cycle    58750  time    4.229E+01 sec block plot  1.750  2.250  2.750  3.250  3.750  4.250  4.750  5.250  5.750 (*10^2)  3.000  3.500  4.000  4.500  5.000  5.500  6.000  6.500  7.000 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada (c) Slope parallel Figure 4.2: An example of the three joint networks used in the parametric study. Note that the view shown is zoomed to the area of interest compared to the entire model. 48 Chapter 4. Parametric Study Table 4.2: Discontinuous and slope parallel joint model parameters Parameter Value Unit Density 2600 kg/m3 E 7 x 109 Pa ν 0.33 - φrock 33o - crock 1.0 x 106 Pa trock 1.0 x 105 Pa φjoint 25o - cjoint 1.0 x 105 Pa tjoint 1.0 x 103 Pa Joint normal stiffness (kn) 1.0 x 1010 Pa Joint shear stiffness (ks) 8 x 109 Pa Table 4.3: Daylighting joint model parameters Parameter Value Unit Density 2600 kg/m3 E 7 x 109 Pa ν 0.33 - φrock 40o - crock 3.50 x 105 Pa trock 1.0 x 105 Pa φjoint 30o - cjoint 1.0 x 105 Pa tjoint 1.0 x 103 Pa Joint normal stiffness (kn) 10 x 109 Pa Joint shear stiffness (ks) 8 x 109 Pa 49 Chapter 4. Parametric Study In this study, the factor of safety is reported to three significant figures even though this level of accuracy is typically not used in industry as the range of model uncer- tainty and simplifying assumptions causes this level of confidence (i.e. precision) to be inadmissible. However, because the focus is on how the parameters affect the SSR factor of safety, three significant figures are quoted for the answer. Cases which cause large variations in factor of safety are highlighted. Those that do not affect the answer at the accuracy level used in industry (typically two significant figures) are noted as not likely being important to an analysis. All D-SSR based chart outputs are created from a monitoring point at the crest unless otherwise noted. An example of the input files used in UDEC can be seen in Appendix B 4.1 Cycling Mode When using the D-SSR method, there are two options for cycling the model at each strength reduction step: cycling to a specified force ratio or timestepping a specific number of steps. The force ratio is the ratio of the average unbalanced force magnitude to the average applied force magnitude for all block gridpoints in the model. The default value for this in UDEC corresponds to solving for a force ratio limit of 10−5. The implications of this is that each SRF stage in the routine is cycled a different amount (in terms of time steps). Typically, the further from SRF = 1, the more cycles it takes for the force ratio to reach the specified limit as deformations increase at each stage. The second option is to cycle the model a specified number of timesteps which results in each SRF strength stage being cycled the same amount. In both cases, the important consideration is that the model be cycled until equilibrium is assured for each strength reduction step. The force ratio route is an attempt to automate this process; the question is whether it does this adequately and at what values does this typically occur at, be it 10−4, 10−5, 10−6, or some other value. If the user bypasses the automated cycling and instead cycles for a specific number of timesteps, then the user must ensure that equilibrium conditions are being met before moving on to the next strength re- duction step. Typically, the user will select a number of timesteps that is quite large in order to ensure that equilibrium is likely being met at all stages with the obvious disadvantage of excessive cycling resulting in increased computational time. 50 Chapter 4. Parametric Study UDEC-SSR differs from the D-SSR method in that it is uses an automatic brack- eting solution method. As noted in section 3.6 a representative number of steps (Nr) is found for UDEC-SSR. There are several conditions which represent a state of equilibrium being met and moving on to the next stage which would require the model being cycled more than the Nr number of steps. However, on a basic level the UDEC-SSR routine is based on a force ratio limit that is checked against, after Nr number of steps is run (or some multiple of), rather than a specific number of cycles to determine stability at each stage. A 10 m discontinuous daylighting joint model was used to test a range of cycled values and compare against cycling for force ratio limits. The different cases run were 5000, 15000, 30000, 45000, 60000, and 100000 cycles in addition to running the D-SSR method for force ratio limits of 10−4, 10−5, and 10−6. D-SSR results for crest measurements can be seen below in figure 4.3 for cycle limits and figure 4.4 for force ratio limits. For cycle comparison purposes, table 4.4 gives the number of cycles until the 10−5 force ratio limit is met for each strength reduction stage. Table 4.4: Number of cycles for 10−5 force ratio limit condition in 10 m discon- tinuous daylighting joint case SRF stage Cycles for each stage 1.62 11470 1.64 13240 1.66 14970 1.68 15580 1.70 27150 1.72 25880 1.74 33120 1.76 32290 1.78 34000 1.80 34700 It is clear that there is a necessary minimum threshold of cycles that needs to occur in order for significant displacements to develop, as seen by the low dis- 51 Chapter 4. Parametric Study 18 00 16.00 . 12.00 14.00 e m e n t   ( m ) 10.00 t a l   d i s p l a c e 5000 15000 6.00 8.00 l e d   h o r i z o n 30000 45000 60000 2 00 4.00 M o d e l l 100000 0.00 . 1 60 1 65 1 70 1 75 1 80 1 85 1 90. . . . . . . Strength reduction factor Figure 4.3: Results of displacement-based SSR for varying number of cycles placements in each SRF stage in the 5,000 and 15,000 cycle cases in figure 4.3. Similarly, there is a distinct upper limit when solving for a set number of cycles where additional cycling does not significantly change the displacements that de- velop and the critical SRF value. In this case, at a SRF of 1.68 with 30,000 or more cycles. The primary difference between the 30,000, 45,000, 60,000, and 100,000 cycle/srf-interval cases is that displacements are correspondingly larger as the number of cycles increases until some cycling threshold is reached where the behaviour becomes the same. This can be seen by the marginally different values of displacement for the 60,000 and 100,000 cases despite the 40,000 dif- ference in cycles per SRF interval. This threshold appears to be the same as the maximum number of cycles required to reach the ratio limit of 10−5 ( 34700, for this particular model) seen in table 4.4. Similarly, from figure 4.4, there is a minimum ratio limit before significant dis- placements develop within the slope. The ratio limits of 10−3 and 10−4 are too large to allow the model to fully come to equilibrium before the next SRF stage begins. It is important to note that the actual value of displacements is not im- portant; rather the shift from low displacement to high displacement is important. 52 Chapter 4. Parametric Study 10.00 12.00 14.00 16.00 18.00 t a l   d i s p l a c e m e n t   ( m ) Ratio = 1E‐3 0.00 2.00 4.00 6.00 8.00 1.44 1.54 1.64 1.74 1.84 1.94 M o d e l l e d   h o r i z o n Strength reduction factor Ratio = 1E‐4 Ratio = 1E‐5 Ratio = 1E‐6 Figure 4.4: Results of displacement-based SSR for varying force ratio limits per strength reduction step 53 Chapter 4. Parametric Study There is agreement between the results shown in figure 4.3 and figure 4.4 with the 100,000 cycles case matching the response of the 10−6 ratio limit case. How- ever, the force ratio limit case took significantly less time to compute due to the lower number of cycles needed for each SRF interval, especially those at the early stages. UDEC-SSR returns the SRF as 1.74 for this joint network case when run with the default limit of 10−5. When looking at the plot in figure 4.3 it is possible for the D-SSR method to also return a result of 1.74 (in addition to the 1.70 indicating the limit of significant displacement) as the method is interpretive and requires the user to select a point at which the slope accelerates to failure. Looking at both the 10−5 force ratio case and 30,000 cycle case for D-SSR in the above results, the slope initially substantially deforms from previous stages at an SRF of 1.70 but then does not significantly deform between 1.70 and 1.74 before significantly deforming again at an SRF of 1.76. Remember that an increased SRF represents an increased reduction in strength properties in the model and therefore increased instability, even though a higher SRF also means a higher factor of safety as the difference between the initial parameters and those required for failure increases. This could be due to the SRF stage = 1.72 reaching the ratio limit earlier (less deformation) than the SRF = 1.70 stage due to slight variations in intact block deformations and movements along joints creating a slightly more stable force configuration at a SRF of 1.72 than at a SRF of 1.70. Though it is not clear from viewing other model information where the deformation difference could be occurring in the complicated block and joint network. The models that are cycled longer or use a lower force ratio limit do not pose a difficulty in subjective selection of the point of significant displacement. However, this subjective effect on critical SRF selection is small (less than 1%). In all of these cases, the model was investigated and the same failure surface was observed provided there was a significant amount of displacement in the model. For example, figure 4.5 shows the velocity plots for the 15,000 and 60,000 case have the same portion of the slope actively moving. Plots of velocity for other cases are similar. Occasionally a model may show anomalous displacements where a later reduction step displaces less than an earlier step. This phenomena can be seen in figure 4.6 which shows the results of a 15 m slope parallel jointing case solved for a ratio limit of 10−6. Note how the stronger materials in SRF = 1.40 displaces more than the weaker materials in the SRF = 1.42 case. This type of behaviour occasionally 54 Chapter 4. Parametric Study    UDEC (Version 4.00) LEGEND    22-May-08  14:44  cycle    61030  time    3.690E+01 sec velocity vectors      maximum =    7.713E-02 0   2E -1 block plot  2.000  2.500  3.000  3.500  4.000  4.500  5.000  5.500 (*10^2)  3.500  4.000  4.500  5.000  5.500  6.000  6.500  7.000 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada          (a) 15,000 cycles    UDEC (Version 4.00) LEGEND    23-May-08  20:44  cycle   146030  time    8.710E+01 sec velocity vectors      maximum =    8.271E-05 0   5E -4 block plot  2.000  2.500  3.000  3.500  4.000  4.500  5.000  5.500 (*10^2)  3.500  4.000  4.500  5.000  5.500  6.000  6.500  7.000 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada          (b) 1 ,000 cycles Figure 4.5: Velocity plots showing the same activated portions of slope (seen by the presence of large velocity vectors) for 10 m and 15 m discontinuous joint cases, at the onset of instability 55 Chapter 4. Parametric Study occurs for both the manual cycled and automated force ratio solving modes. Given the complex interactions between intact block deformation and movement along joints, it is believed that this arises from a more stable configuration of deformed blocks that does not occur in an earlier stage however this is difficult to investi- gate within the models. This behaviour does not occur with lower displacements modelled for two or more successive steps, though an upward trending sawtooth pattern as seen in figure 4.6 is common. Most importantly, this occurs after the selected stability criterion is exceeded and does not affect the selection of the ap- propriate critical strength reduction stage indicating instability. Similarly, UDEC-SSR can be run using a different force ratio limit such as 10−6. This was done for the slope parallel cases for the 10 m and 15 m spacing, along with UDEC-SSR for the 10−5 ratio, with results summarized in table 4.5. There is a 0.1 (∼7%) difference in factor of safety between the two methods but no difference in critical failure surface or failed volume. This is due to the higher threshold for stability in the 10−6 case leading to more cycles in the model per SRF stage which allow for greater plastic deformation to occur earlier leading to instability at an earlier SRF stage. This instability at an early stage using a ratio limit of 10−6 is similar to the discontinuous daylighting case using D-SSR and reinforces the choice of 10−6 force ratio limit to provide a more conservative factor of safety. Because of this, using a ratio limit of 10−6 for both D-SSR and UDEC-SSR is recommended. Table 4.5: UDEC-SSR results for two different ratio limits Case SRF (10−5) SRF(10−6) 10 m spacing slope parallel 1.53 1.43 15 m spacing slope parallel 1.53 1.43 56 Chapter 4. Parametric Study 2.50 3.00 3.50 4.00 4.50 t a l   d i s p l a c e m e n t   ( m ) 0.00 0.50 1.00 1.50 2.00 1.20 1.25 1.30 1.35 1.40 1.45 1.50 1.55 1.60 M o d e l l e d   h o r i z o n Strength reduction factor Figure 4.6: Example of sawtooth horizontal displacement results during D-SSR analysis 57 Chapter 4. Parametric Study 4.2 Element Size Six models were run using the discontinuous daylighting joint geometry with a spacing of 10 m for each of the joints (figure 4.2b). The purpose of these models is to test the sensitivity of the SSR method to element size as it is an important numerical feature that often controls model behaviour. The model geometry and input parameters are the same between all six models except for the generation of the mesh density throughout the slope. All models were run using the D-SSR method with a ratio limit of 10−5 unless otherwise noted, in order to reduce com- putation time. Six different models were created with mesh sizes of 1.0 , 2.5 , 5.0 , 7.5 , 10, and 15 m triangular elements. All mesh sizes are smaller than the intact block size with the exception of the 15 m mesh case. The intention of the 15 m case is to investigate behaviour when the element size is selected larger than the smallest block size. Figure 4.7 and table 4.6 show that the results of D-SSR are highly dependent on the element size used. These results show a general trend where the smaller the element size, the more conservative the factor of safety. In fact, an up to 25% difference in factor of safety was obtained for the range of element sizes tested. Between the different element sizes there is a slight difference in failed volumes (location of critical slip surfaces) with slightly deeper seated failures occurring with cases that fail at a later SRF stage (i.e. larger elements). This is due to the influence element size has on block deformations with a greater number of elements in a block (i.e. finer mesh) generally allowing for finer detail in shear surface localization and block deformations. Because of this, the blocks can behave differently between the grid sizes, causing small differences in deformation patterns which ultimately leads to a slightly different critical failure surface within the slope. Non-conservative, and difficult to interpret, results occur when the grid size is greater than the block size, as is seen by the 15 m case in figure 4.7. This is likely due to the large approximations of modelled block deformation and the inability of shear band localization when a large element size is used. Here there is a gradual increase in displacements as SRF increases and the strengths are decreased but there is no sharp acceleration of displacements as seen for the smaller element size models. Because of this, it is difficult to pick the exact stage at which failure occurs other than indicating that it may be occurring between an SRF of 1.88 and 2.0. It may be useful to restore models at these SRF intervals and look at 58 Chapter 4. Parametric Study diagnostic indicators such as gridpoint velocity, plastic states of the elements, and strain contours; however, following the D-SSR method as originally stated, an accelerated state is not discernible. 8.00 10.00 12.00 14.00 16.00 t a l   d i s p l a c e m e n t   ( m ) 1 m 2.5 m 0.00 2.00 4.00 6.00 1.60 1.70 1.80 1.90 2.00 2.10 M o d e l l e d   h o r i z o n Strength reduction factor 5 m 7.5 m 10 m 15 m Figure 4.7: Results of displacement based SSR for element size sensitivity Table 4.6: Summary of D-SSR FOS values for element size variation Element size D-SSR SRF 1 m NA 2.5 m 1.66 5 m 1.70 7.5 m 1.76 10 m 1.82 15 m 1.88 - 2.0 The 1.0 m element size is an interesting case as it does not follow the trend given by the larger element sizes. Following from the general trend, it is expected that 59 Chapter 4. Parametric Study the finer the element density the more detailed the computation becomes and thus the better approximation it is to the real system, especially with respect to tracking localization and the development of a critical shear failure surface. In this case, the 1.0 m element size has markedly different behaviour than the next larger el- ement size in this series at 2.5 m. This could be due to numerical instability but appears to be related to ultra localization where numerous localized shear bands are developing at the same time such that a single, dominant surface does not de- velop easily. This can be seen in figure 4.8 which is a zoomed in view of the same block showing the plastic state of each element in the grid for both the 1 m and 5 m element cases. The 1.0 m case (figure 4.8a) displays a clear shear band running through the block; however, the majority of the block has not yet reached the yield criteria of the material. Contrasting this is the 5 m case in figure 4.8b where all the zones in the block have failed in shear. In the first case, stresses have been able to redistribute on a small scale such that the intact blocks are still generally competent as stress concentrations occupy a smaller percentage of the block. In the second case, the blocks are more easily deformed as each particular zone is large and stress redistribution is not able to concentrate in such a small area as in the 1 m case, thus causing more of the block to be in a yielded state. Extrapolat- ing this over the entire slope, the 5 m case is able to deform much more readily as the blocks along the sliding surface are loaded beyond their strength and any additional stresses above this leads to large deformations. It is impossible to say at this stage which of these is the more correct SRF, as there is no stability basis to compare to. Back analysis of a case history in the next chapter is expected to help clarify this important element size issue as there is a significant computer cycling time involved (more computationally expensive) in finer mesh cases compared to those that are slightly coarser. For example, the 2.5 m grid case was computed in approximately 1 day while the 1.0 m grid case took 12 days to complete on a 2.40 GHz Intel Core 2 Quad CPU with 2 gigabytes of RAM in single precision mode. All the above cases were run for a ratio limit of 10−5. A test case was selected using the 7.5 m element grid size to see if the number of cycles could be contribut- ing to the large difference in results. Figure 4.9 shows the results for the 7.5 m case with data for 50,000 cycles/interval and the 10−5 force ratio cases. There are some minor differences between the two displacement responses; however, they both fail at an SRF of 1.78. Thus, mesh sensitivity appears to be independent of the solving method chosen. 60 Chapter 4. Parametric Study    UDEC (Version 4.00) LEGEND    17-May-08   2:09  cycle   260380  time    3.391E+01 sec zones in fdef blocks no. zones : total      7046 at yield surface (*)   2488 yielded in past  (X)      0 tensile failure  (o)      4  2.550  2.600  2.650  2.700  2.750  2.800  2.850 (*10^2)  6.050  6.100  6.150  6.200  6.250  6.300  6.350 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada (a) Grid size = 1 m    UDEC (Version 4.00) LEGEND     5-May-08  11:53  cycle    61000  time    3.689E+01 sec zones in fdef blocks no. zones : total       338 at yield surface (*)    264 yielded in past  (X)      0 tensile failure  (o)      4  2.550  2.600  2.650  2.700  2.750  2.800  2.850 (*10^2)  6.050  6.100  6.150  6.200  6.250  6.300  6.350 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada (b) Grid size = 5 m Figure 4.8: Plastic state of zones in a single block showing detail of shear band localization defined by elements yielding in shear (x yield in shear, o yield in tension) 61 Chapter 4. Parametric Study 6.00 8.00 10.00 12.00 t a l   d i s p l a c e m e n t   ( m ) 0.00 2.00 4.00 1.60 1.70 1.80 1.90 2.00 2.10 M o d e l l e d   h o r i z o n Strength reduction factor 7.5 m force equilibrium 7.5 m cycled  Figure 4.9: Comparison between cycled and force-equilibrium results for 7.5 m grid size 62 Chapter 4. Parametric Study 4.3 Element Types There are two different types of meshing possible within UDEC: quadrilateral elements and triangular elements. Quadrilateral elements provide a better solution for plasticity problems, provided that they can be generated for the different block shapes involved. A limitation of quadrilateral elements is that they can only be discretized within a quadrilateral block. Models were run for the 5 m and 7.5 m element size cases for both triangular and quadrilateral element types with the D-SSR method for a force ration of 10−6 and discontinuous joint spacing of 10 m with results provided in figure 4.10. The 5 m grid case predicts failure at SRF = 1.70 for both the triangular and quadri- lateral cases, with only slight differences in the displacements. The 7.5 m element case has triangular elements failing at SRF = 1.76 with quadrilateral elements failing at SRF = 1.74. Figures 4.11 and 4.12 show the horizontal displacement contours and plastic state of the 7.5 m triangular and quadrilateral elements at SRF = 1.76 and SRF = 1.74 respectively. A comparison of these figures shows that the failure surface is slightly bigger in the 7.5 m quadrilateral case. Also, the triangular elements predict a less continuous failure surface compared to the quadrilateral elements. UDEC-SSR was run for the 7.5 m cases giving a SRF of 1.78 with no difference in failure surface, unlike the D-SSR result. It is rec- ommended that quadrilateral elements are used whenever possible, particularly when using D-SSR. In reality this is not critical as a slope will likely consist of a mixture of elements as blocks which fail quadrilateral element discretization (non-quadrilateral blocks) will use triangular elements instead. 63 Chapter 4. Parametric Study 8.00 10.00 12.00 14.00 16.00 t a l   d i s p l a c e m e n t   ( m ) 5 m triangular 0.00 2.00 4.00 6.00 1.60 1.65 1.70 1.75 1.80 1.85 1.90 M o d e l l e d   h o r i z o n t Strength reduction factor 5 m quadrilateral 7.5 m triangular 7.5 m quadrilateral Figure 4.10: D-SSR results for element type sensitivity for discontinuous day- lighting joints 64 Chapter 4. Parametric Study    UDEC (Version 4.00) LEGEND     9-Jul-08  13:24  cycle    56290  time    5.200E+01 sec X displacement contours contour interval= 2.000E-01  2.000E-01 to  1.400E+00                   2.000E-01                   4.000E-01                   6.000E-01                   8.000E-01                   1.000E+00                   1.200E+00                   1.400E+00 block plot  1.500  2.000  2.500  3.000  3.500  4.000  4.500  5.000  5.500 (*10^2)  3.000  3.500  4.000  4.500  5.000  5.500  6.000  6.500  7.000 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada (a) Horizontal displacement contours (m)    UDEC (Version 4.00) LEGEND     9-Jul-08  13:24  cycle    56290  time    5.200E+01 sec block plot no. zones : total     12490 at yield surface (*)    565 yielded in past  (X)      0 tensile failure  (o)    164  1.500  2.000  2.500  3.000  3.500  4.000  4.500  5.000  5.500 (*10^2)  3.000  3.500  4.000  4.500  5.000  5.500  6.000  6.500  7.000 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada (b) Plastic states of elements (x yield in shear, o yield in tension) Figure 4.11: Failed volume indicators for 7.5 m triangular elements at SRF stage = 1.76 65 Chapter 4. Parametric Study    UDEC (Version 4.00) LEGEND     9-Jul-08  13:25  cycle    58450  time    4.946E+01 sec X displacement contours contour interval= 2.000E-01  2.000E-01 to  1.200E+00                   2.000E-01                   4.000E-01                   6.000E-01                   8.000E-01                   1.000E+00                   1.200E+00 block plot  1.500  2.000  2.500  3.000  3.500  4.000  4.500  5.000  5.500 (*10^2)  3.000  3.500  4.000  4.500  5.000  5.500  6.000  6.500  7.000 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada (a) Horizontal displacement contours (m)    UDEC (Version 4.00) LEGEND     9-Jul-08  13:25  cycle    58450  time    4.946E+01 sec block plot no. zones : total     16912 at yield surface (*)    747 yielded in past  (X)      0 tensile failure  (o)    246  1.500  2.000  2.500  3.000  3.500  4.000  4.500  5.000  5.500 (*10^2)  3.000  3.500  4.000  4.500  5.000  5.500  6.000  6.500  7.000 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada (b) Plastic states of elements (x yield in shear, o yield in tension) Figure 4.12: Failed volume indicators for 7.5 m quadrilateral elements at SRF stage = 1.74 66 Chapter 4. Parametric Study 4.4 Insitu Stress Three models were run using the discontinuous daylighting joint geometry with a spacing of 10 m for each of the joints (figure 4.2b). As before, this model was chosen as it is the least structurally controlled geometry and thus represents the most generalized rockmass failure case. The purpose of this series is to check if the SRF returned is sensitive to the insitu stress ratio in the model. This is im- portant as the insitu stresses in slopes are commonly unknown and assumed, thus any differences in reported model results would make the selection of this param- eter important. The model geometry and all parameters are identical between the three models except for the insitu stresses. The three insitu stresses chosen were a horizontal to vertical stress ratio (ko) of 0.5, 1, and 2. The horizontal and vertical stress distributions, after removal of the benches, for each of the three cases can be found in figures 4.13 and 4.14 respectively. Results for the D-SSR method monitoring the crest can be seen in figure 4.15, which returns a critical SRF of 1.70. UDEC-SSR returns a result of 1.71 for all three cases. Even though the stress conditions going into the strength reduction method are different between the models, the result are the same. The magnitude of the displacements varies slightly; however, the point of stability occurs at the same reduction step. D-SSR matches the behaviour of UDEC-SSR, both of which are unaffected by the insitu stress regime. Small changes in the initial distribution of stresses, on the cases tested, do not have much of an effect on the ultimate point of failure. This is likely due to the relatively large drop in material properties at the point of failure (a 59% reduction in strength parameters) causing significant stress concentrations in the toe of the slope compared to the rather small differences in initial horizontal and vertical stresses. This agrees with results by Hammah et al. (2005a) and Diederichs et al. (2007) for continuum models. 4.5 Elastic Properties A series of models were run using the discontinuous daylighting joint geometry with a spacing of 10 m for each of the joint sets (figure 4.2b). The purpose of this series is to test the sensitivity of the Poisson’s ratio, and by association the bulk and shear modulus, of the Mohr-Coulomb constitutive model on the SSR method. 67 Chapter 4. Parametric Study (a) Legend (Pa) (b) K = 0.5 (c) K = 1.0 (d) K = 2.0 Figure 4.13: Horizontal stress contours for varying modelled insitu stress ratios 68 Chapter 4. Parametric Study (a) Legend (Pa) (b) K = 0.5 (c) K = 1.0 (d) K = 2.0 Figure 4.14: Vertical stress contours for varying modelled insitu stress ratios 69 Chapter 4. Parametric Study 18 00 16.00 . 12.00 14.00 e m e n t   ( m )   10.00 t a l   d i s p l a c e k=0.5Point of failure 6.00 8.00 e d     h o r i z o n t k=1  k=2  at SRF= 1.7 4.00M o d e l l e 0.00 2.00 1.60 1.65 1.70 1.75 1.80 1.85 1.90 Strength Reduction Factor Figure 4.15: D-SSR results reduction results for varying horizontal to vertical stress ratios 70 Chapter 4. Parametric Study The models are set to take Young’s modulus (E) and Poisson’s ratio (ν) as inputs to calculate the bulk (K) and shear modulus (G) from the following equations, as recommended in (Itasca, 2004c): K = E 3(1− 2ν) (4.1) G = E 2(1 + ν) (4.2) Table 4.7 gives the model setup with associated rock modulus parameters for each case. D-SSR results can be seen in figure 4.16 which shows a similar point of fail- ure at a stress reduction factor of 1.8 for all cases. Displacements are not sensitive to various moduli (E, K, G) compared to the strength drop in parameters at a SRF of 1.8. This is because the elastic parameters have little influence on the shear fail- ure (plastic deformation) of the rockmass. This may not be true for the case where there are multiple stiffnesses within a model which vary significantly from each other. Dawson et al. (1999) states that to properly duplicate LEM results with the SSR method, the same value of E and ν should be used. However, duplication of the LEM result may not be the desired outcome of an SSR analysis. Table 4.7: Input parameters for Poisson’s Ratio sensitivity models ν E (GPa) K (GPa) G (GPa) 0.10 7.00 2.50 2.73 0.15 7.00 2.86 2.61 0.20 7.00 3.33 2.50 0.25 7.00 4.00 2.40 0.30 7.00 5.00 2.31 71 Chapter 4. Parametric Study 10 00 9.00 . 7.00 8.00 e m e n t   ( m ) 5.00 6.00 t a l   d i s p l a c e v = 0.10 v = 0.15 3 00 4.00 l e d   h o r i z o n v = 0.20 v = 0.25 v = 0 30User selected point 2.00 . M o d e l l     .    of failure 0.00 1.00 1 60 1 65 1 70 1 75 1 80 1 85 1 90. . . . . . . Strength reduction factor Figure 4.16: Results of displacement based SSR for Poisson’s Ratio sensitivity 4.6 Joint Shear Stiffness Four models were run to test the sensitivity of the D-SSR method to the joint shear stiffness (ks) property. The models were run on the 10 m discontinuous daylighting joint model as this model represents the most complicated, intact rock and joint failure model (figure 4.2b). Joint shear stiffness values tested were 8 x 107, 8 x 108, 8 x 109 , and 8 x 1010 Pa/m with the six Mohr-Coulomb parameters (friction, cohesion, tension for both intact rock and joints) factored equally. The model was initially run for a ratio limit of 10−5, giving results for D-SSR as seen in figure 4.17. These results indicate a relatively large range of factor of safety reported back with no clear general trend in terms of model response to stiffness change. The joint shear stiffness has a profound effect on the rate at which the model comes to equilibrium, as can be seen in figure 4.18 showing the number of cycles required for the ratio limit in the 8 x 108 and 8 x 1010 Pa/m cases. There is a significant difference in the number of cycles in each model and it has been shown in section 4.1 that the number of cycles can have a significant impact on the SRF stage indicating instability. 72 Chapter 4. Parametric Study The 8 x 108 and 8 x 1010 Pa/m cases were run with a ratio limit of 10−6 with results seen in figures 4.19 and 4.20 respectively. This was done to see if the differences in factor of safety persisted between the various values of joint shear stiffness when a more stringent equilibrium criterion is introduced. There is no difference between the ratio limits for the 8 x 1010 cases, as seen in figure 4.20. Further investigations reveal that the failure surface for the two ratio limits is the same. There is, however, a large difference between the ratio limits of the 8 x 108 ra- tio cases: a significant difference in SRF of 0.14 (∼9%) between the two force ratio limit cases. Figure 4.21 shows velocity vectors indicating actively moving portions of the slope at the SRF stage where failure is first indicated for the ra- tio limits in the 8 x 108 case. The difference in SRF reported is due to a deeper seated failure being predicted in the ks = 8 x 108 Pa/m case, which is seen in the velocity plots. Because the slope is sliding on a different surface, with a different percentage of rock bridges involved, the SRF reported is different. It is clear that the joint shear stiffness can have a profound effect on the SRF de- termined as well as the position of the critical sliding surface due to the influence the parameter has on modelled displacements and that this effect may or may not depend on the ratio limit selected for equilibrium at each SRF stage. 73 Chapter 4. Parametric Study 18 00 16.00 . 12.00 14.00 e m e n t   ( m ) 10.00 t a l   d i s p l a c e 8.00E+07 (Pa) 6.00 8.00 l e d   h o r i z o n 8.00E+08 (Pa) 8.00E+09 (Pa) 8.00E+10 (Pa) 2 00 4.00 M o d e l l 0.00 . 1 48 1 58 1 68 1 78 1 88 1 98. . . . . . Strength reduction factor Figure 4.17: Results of D-SSR analysis for joint shear stiffness sensitivity 74 Chapter 4. Parametric Study 60000 80000 100000 120000 r   o f   c y c l e s 0 20000 40000 1.50 1.55 1.60 1.65 1.70 1.75 1.80 1.85 N u m b e r Strength reduction factor 8.00E+08 (Pa) 8.00E+10 (Pa) Figure 4.18: Number of cycles for 10−5 ratio limit for ks = 8 x 108 and 8 x 1010 Pa/m cases 75 Chapter 4. Parametric Study 8.00 10.00 12.00 14.00 16.00 t a l   d i s p l a c e m e n t   ( m ) 0.00 2.00 4.00 6.00 1.48 1.58 1.68 1.78 1.88 1.98 M o d e l l e d   h o r i z o n Strength reduction factor Ratio 1E‐5 Ratio 1E‐6 Figure 4.19: Results of displacement-based SSR for varying ratio limits when ks = 8 x 108 Pa/m 76 Chapter 4. Parametric Study 18 00 16.00 . 12.00 14.00 e m e n t   ( m ) 10.00 t a l   d i s p l a c e 6.00 8.00 l e d   h o r i z o n Ratio 1E‐5 Ratio 1E‐6 2 00 4.00 M o d e l l 0.00 . 1 48 1 58 1 68 1 78 1 88 1 98. . . . . . Strength reduction factor Figure 4.20: Results of displacement-based SSR for varying ratio limits when ks = 8 x 1010 Pa/m 77 Chapter 4. Parametric Study    UDEC (Version 4.00) LEGEND    26-May-08  16:31  cycle    86660  time    5.629E+01 sec velocity vectors      maximum =    3.339E-01 0   2E  0 block plot  1.500  2.000  2.500  3.000  3.500  4.000  4.500  5.000  5.500  6.000 (*10^2)  3.000  3.500  4.000  4.500  5.000  5.500  6.000  6.500  7.000 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada (a) Velocity showing active portions of the slope for ratio limit of 10−5 at failure SRF of 1.8    UDEC (Version 4.00) LEGEND    30-Jun-08  15:43  cycle   109680  time    7.099E+01 sec velocity vectors      maximum =    3.990E-04 0   2E -3 block plot  1.500  2.000  2.500  3.000  3.500  4.000  4.500  5.000  5.500  6.000 (*10^2)  3.000  3.500  4.000  4.500  5.000  5.500  6.000  6.500  7.000 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada (b) Velocity showing active portions of the slope for ratio limit of 10−6 at failure SRF of 1.6 Figure 4.21: Comparison between active portions of the slope for different ratio limits when ks = 8 x 108 Pa/m 78 Chapter 4. Parametric Study 4.7 Strain-Softening Models It is possible to apply the SSR method to a constitutive model that allows for a strength drop after the peak strength has been reached. Hammah et al. (2005a) in- vestigates this briefly and comes to the conclusion that an elastic-perfectly plastic model best agrees with the LEM analysis. However, discontinuum models are typ- ically used where LEM methods are not applicable, or do not accurately capture the expected failure mechanics. This poses an interesting question: in some cases, why would one want the model to agree with a LEM analysis that is not able to properly represent rockmass behaviour? Considering that one of the strengths of the discontinuum model is that it allows for a better idealization of the rockmass and associated mechanical behaviour, the strain-softening model using the SSR method should be investigated. UDEC-SSR currently does not allow the use of softening materials to be included in a SSR analysis. A version of D-SSR was created in order to set residual parameters such that residual intact friction and cohesion were set between the properties for the peak intact parameters and joint parameters. Both peak and residual parameters are reduced in this series of models. The point at which these parameters are ap- plied depend on a linear relationship between residual parameters and modelled strain such that a specific strain limit is defined where strengths are set to be fully residual. Joints remained elastic perfectly-plastic, and are included in the strength reduction. It is obvious that the lower the residual strengths, the earlier the point of failure is; however, it is less obvious how the strain limit will affect the shape of the failure surface and the SRF reported. D-SSR was run for the 10 m discontinuous joint spacing case with parameters seen in table 4.8 for full residual parameters at strains of 0.1%, 1%, 5%, 50%, 100% with results seen in figure 4.22. From these results it is seen that the strain limit for post-peak behaviour is extremely sensitive, with failure at an earlier SRF stage as lower strain limits are set for mobilizing the full residual parameters. Also, there is a slight change in sliding surface location with slightly shallower failures occurring with smaller peak to residual strain limits. 79 Chapter 4. Parametric Study Table 4.8: Peak-residual parameters used for strain-softening sensitivity Property Peak Residual φrock 33 28 φjoint 25 - crock 1.01 x 106 Pa 5 x 105 Pa cjoint 1 x 105 Pa - 8.00 10.00 12.00 14.00 16.00 t a l   d i s p l a c e m e n t   ( m ) 100% 10% 0.00 2.00 4.00 6.00 1.10 1.20 1.30 1.40 1.50 1.60 1.70 1.80 1.90 2.00 M o d e l l e d   h o r i z o n t Strength reduction factor 5% 1% 0.1% Perfectly plastic Figure 4.22: D-SSR for strain limits for strain-softening cases 80 Chapter 4. Parametric Study 4.8 Joint Network Cases This section presents the reasoning and results for a sequence of models that tests the sensitivity of the SSR method to the joint spacing. As noted in section 3.2, the joints in the rockmass provide significant control over the mechanical response of the stressed rockmass. It is expected that the higher the fracture density, the lower the factor of safety returned due to the fact that the joint properties are lower than the properties of the continuum blocks between the discontinuities. For all of the models with varying joint spacing or other properties (persistence, etc.), the intact properties are kept as the same ‘equivalent continuum’ material as there is no clear relationship between the modelled fracture density and the amount of ‘equivalent continuum’ strength reduction needed to account for the different amount of discrete features being modelled. 4.8.1 Daylighting Joints This set of models represents an extreme case where there are long, continuous daylighting joints with a continuous, back-dipping cross joint. This model is typ- ical of a case that may be modelled using simpler analytical kinematic stability methods and probabilistic planar sliding routines. The purpose of this model is to examine the influence of joint spacing on the factor of safety with particular attention paid to ultimate failure volumes in the model. A set of models were run where the joint spacing for both sets (daylighting and orthogonal), were varied by the same amount. A typical joint distribution for this joint network case can be seen previously in figure 4.2a in section 4, showing a 10 m spacing on both joint sets. Figure 4.23 shows the results of D-SSR for varying joint spacings run for a force ratio of 10−5. UDEC-SSR method was also run for these cases and the results can be seen in table 4.9 along with values of SRF for D-SSR. The results show that there is some sensitivity to the joint spacing; the more heav- ily jointed the rockmass, the lower the factor of safety. An exception however is the 10 and 15 m cases using UDEC-SSR method that are actually less stable than the 5.0 m case. A closer look at the failed states, shown as horizontal dis- placement contours with shear along the major slipping joint in figure 4.24 reveals why: the joints are slipping on a similar critical joint whose placement may differ 81 Chapter 4. Parametric Study 6.00 5.00 e m e n t   ( m ) 4.00 t a l   d i s p l a c e 2.5 m 5 m 2.00 3.00 l e d   h o r i z o n 10 m 15 m 20 m 1.00 M o d e l l   0.00 1 00 1 05 1 10 1 15 1 20 1 25 1 30 1 35 1 40. . . . . . . . . Strength reduction factor Figure 4.23: D-SSR results for different daylighting joint spacings slightly in each spacing case as the joint generator within UDEC does not neces- sarily create joints at the exact same place within the model for different spacing realizations. All models slip on the longest continuously daylighting joint above the toe, re- gardless of how many joints there are in the slope and thus share a similar ultimate failure volume. Some D-SSR models, namely 10, 15 and 20 m spacings, show no clear trend for where displacements begin to change with increasing SRF stage. These models initially deform but then plateau to a series of stable stages before a second large increase in displacements/SRF stage occurs. The difference between the possible critical SRF stages for these models can be as great as 10%. This has potential implications for design. Within the model, blocks could be deforming and locking giving these quasi-stable states despite the slope being weaker (higher SRF stage). However, this is similar to behaviour that is seen in real rock slopes where movements begin and the slope appears to be at the onset of instability but ultimately hangs up on critical blocks and asperities. UDEC-SSR does not show this behaviour and agrees with the conservative D-SSR factor of safety se- lected from when significant movement first occurs (the lower bound of the range 82 Chapter 4. Parametric Study (a) 2.5m (b) 5 m (c) 10 m Figure 4.24: Horizontal displacement and shear along joints for varying daylight- ing joint spacings (m) 83 Chapter 4. Parametric Study (d) 15 m (e) 20 m Figure 4.24: Horizontal displacement and shear along joints for varying daylight- ing joint spacings (m) (cont.) returned). 84 Chapter 4. Parametric Study Table 4.9: SSR results for different daylighting joint spacings Spacing (m) UDEC-SSR SRF Displacement-based SSR SRF 2.5 1.06 1.08 5.0 1.10 1.10 10 1.08 1.12 - 1.20 15 1.09 1.14 - 1.26 20 1.12 1.14 - 1.26 4.8.2 Discontinuous Daylighting Joints This set of models represents a typical rockmass containing daylighting joints with a back-dipping cross joint, both of which are non-continuous. This repre- sents a more complicated rockmass failure that involves both failure of intact rock and movement along discontinuities. This type of failure is not captured in other methods such as LEM or continuum modelling methods. A set of models was run where the joint spacing was varied. The joints in UDEC were generated using the settings seen in table 4.10 to create a slope as previously seen in figure 4.2. An explanation of each of the terms used to generate the joint set can be found back in section 3.5. These models physically represent the rock bridges as intact blocks using a Mohr-Coulomb constitutive model, as opposed to using continuous joints and modifying joint cohesion properties to simulate the presence of rock bridges, such as in Kim et al. (2007). Scaling the gap length and trace length of the joint by a factor of 4 compared to the spacing was determined to give a better match (number and position of rock bridges) than just changing the spacing. Given the way the terminating joints within a block are removed in UDEC, it is difficult to keep a consistent joint network between the different spacings. D-SSR results for force ratio of 10−5 can be seen in figure 4.25. There is a scatter in the displacement results with critical SRF values ranging from 1.62 – 1.74. It was expected that the slope would become more unstable as the fracture frequency increased as there would be more lower strength joints combined with smaller in- 85 Chapter 4. Parametric Study Table 4.10: Discontinuous joint generation settings Spacing (m) Angle trace length(m) gap length (m) µ, σ µ, σ 7 -40 28, 0 7, 7 7 45 28, 0 0, 0 8 -40 32, 0 8, 8 8 45 32, 0 0, 0 9 -40 36,0 9, 9 9 45 36, 0 0, 0 10 -40 40, 0 10, 10 10 45 40, 0 0, 0 11 -40 44, 0 11, 11 11 45 44, 0 0, 0 12 -40 48, 0 12, 12 12 45 48, 0 0, 0 13 -40 52, 0 13, 13 13 45 52, 0 0, 0 86 Chapter 4. Parametric Study tact rock bridges, reducing the overall rockmass strength. The results in figure 4.25 do not show this consistently. For example, the 8 m case is less stable than the 7 m case due to the complex joint network generated. Given the complex- ity involved and because only fractures that completely pass through a block in UDEC are included, there is no way to guarantee that the same number of rock bridges will occur in the slope for each different joint spacing. Small changes in gap lengths, trace lengths, or spacing can cause significant differences in the joint network created. Because of this it is possible that a more fractured rockmass will have rock bridges and interlocking blocks providing strength in critical places for one spacing that may not be present in models with a smaller joint spacing. This can be clearly seen by looking at the joint network in the toe for the 7 m and 8 m spacing cases (figures 4.26a and 4.26b). Note the difference between the size and distribution of intact rock bridges (blocks without joints dipping into the slope) even though the difference in joint spacing is only one meter. Despite these differences, it is also interesting to note that the spread in the SRF values for the range of joint network configurations tested is less than 8%. This narrow range would provide a degree of confidence in the factor of safety calcu- lated for a forward analysis, due to the limited sensitivity. 87 Chapter 4. Parametric Study 2 00 1.80 . 1.40 1.60 e m e n t   ( m ) 1.00 1.20 t a l   d i s p l a c e 7m 8m 9m 0 60 0.80 l e d   h o r i z o n 10m 11m 12m 0.40 . M o d e l l 13m 0.00 0.20 1 55 1 60 1 65 1 70 1 75 1 80 1 85. . . . . . . Strength reduction factor Figure 4.25: D-SSR results for discontinuous daylighting joint cases with varying joint network properties 88 Chapter 4. Parametric Study    UDEC (Version 4.00) LEGEND     9-Apr-08  16:21  cycle    45130  time    2.868E+01 sec block plot  2.300  2.500  2.700  2.900  3.100  3.300  3.500 (*10^2)  5.700  5.900  6.100  6.300  6.500  6.700  6.900  7.100 (*10^2) JOB TITLE :  UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada          (a) 7 m spacing    UDEC (Version 4.00) LEGEND    13-May-08  10:40  cycle    45560  time    3.110E+01 sec block plot  2.300  2.500  2.700  2.900  3.100  3.300 (*10^2)  5.700  5.900  6.100  6.300  6.500  6.700  6.900  7.100 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada          (b) 8 m spacing Figure 4.26: Joint network differences between 7 m and 8 m spacing cases 89 Chapter 4. Parametric Study 4.8.3 Slope Parallel Joints This model represents a case where the rockmass is comprised of dipping, lay- ered sedimentary sequences with a continuous, slope parallel joint set and a back- dipping cross joint, shown previously in figure 4.2. This represents a difficult case to assess in traditional LEM methods as the expected failure mechanism is shear- ing of intact rock at the toe followed by slip along the joint sets, which cannot be captured fully using LEM or continuum models due to the important influence of structure. Joints are expected to play a significant role. A set of models were run where the joint spacing was varied to see the influence of joint spacing on the SSR method. D-SSR was run for 5, 6, 7, 10, 11, 12, 13, 15, 20 m cases in the cycled mode for 20,000 cycles. This cycle limit was investigated and it was determined that equilibrium was being met for each model stage tested. UDEC-SSR routine was also run for 5, 10, 15, and 20 m cases. D-SSR results can be seen in figure 4.27 while UDEC-SSR routine results can be seen for several select cases in table 4.11. 2.50 3.00 3.50 4.00 4.50 5.00 t a l   d i s p l a c e m e n t   ( m ) 5 m 6 m 7 m 10 m 0.00 0.50 1.00 1.50 2.00 1.30 1.35 1.40 1.45 1.50 M o d e l l e d   h o r i z o n Strength reduction factor 11 m 12 m 13 m 15 m 20 m Figure 4.27: Results of displacement based SSR for slope parallel joint spacing cases 90 Chapter 4. Parametric Study Table 4.11: Selected results from UDEC-SSR for slope parallel joint case Spacing (m) SRF(10−5) SRF(10−6) 5 1.53 1.43 10 1.53 1.43 15 1.53 1.43 20 1.53 1.43 From figure 4.27, the D-SSR routine predicts failure between an SRF of 1.34 and 1.40 for all cases. Each model was investigated and it was determined that each model satisfies equilibrium with enough cycles having been performed for a force ratio of at least 10−6 to be satisfied. Selecting the critical SRF stage is clear for some cases, such as the 15 m case with an SRF of 1.36. Others, such as the 5 m case, show a general increase in displacements before critical displacements occur resulting in a range of possible instability points between 1.36 – 1.40. Regardless of the difference in joint spacing, the range of results is tightly clustered between 1.32 – 1.40 for the critical SRF stage. UDEC-SSR results for the select cases in table 4.11 give a consistent answer regardless of joint spacing. The answer of 1.43 for the 10−6 UDEC-SSR ratio case is similar to the result for D-SSR (which has a ratio limit of at least 10−6). This result indicates that the selection of critical point of failure in UDEC-SSR may not be as rigorous as the D-SSR criteria, as evidenced by the insensitivity to joint spacing in this case. Generally, the spacing of the joint does not have much of an affect on factor of safety returned and it is not possible to predict whether a higher or lower spacing is going to necessarily be more or less stable. For example, the 7 m spacing has a higher SRF than the 5 m spacing case even though the modelled strength of the rockmass as a whole for the 5 m case is lower due to having more joints at a lower strength. Each case was investigated to find the joint surface that represented the failure surface by looking at different diagnostic plots over the range of the factor of safety reported. From these plots, the distance from the slope crest to the critical sliding joint (table 4.12), can be found. From this table, it is seen that the critical joint distance does not change much with different joint spacings. Thus, there is a structure located about 55 m back from the slope that the majority of 91 Chapter 4. Parametric Study the failure slides on, regardless of how many other joints are between it and the slope surface. Small differences in critical SRF values are due to the slight change in critical joint position between the different spacings which affects the amount of intact rock at the toe that must yield for failure to occur. The failure mode for each case is identical and involves crushing at the toe, a tension crack opening up near a critical joint, and sliding along a critical joint setback from the crest as given in table 4.12. A representative plasticity plot, for example for the 15 m case, can be seen in 4.28 along with shear along the controlling joint showing the typical failure mechanism and sliding surface. The UDEC-SSR cases show the same sliding surfaces as those for the D-SSR; that is, the shear surface and failure volume are the same but the SRF reported is different due to the differences in criteria for stability. However, the UDEC-SSR cases using a force ratio of 10−6 return a closer result to the D-SSR results than the default 10−5 cases. Table 4.12: Distance of critical sliding surface from crest for Slope parallel cases Spacing (m) Distance from crest (m) 5 55 6 53 7 52 10 56 11 47 12 50 13 50 15 63 20 58 92 Chapter 4. Parametric Study Figure 4.28: Plasticity plot showing crushing at the toe, tension cracks at the crest, and location of critical discontinuity for 15 m case (x yield in shear, o yield in tension) 93 Chapter 4. Parametric Study 4.9 Summary The analysis presented in this chapter looks at the sensitivity of various parame- ters, including both model setup options and material properties, and the resulting affect on the SSR-predicted factor of safety and failure surface. Results show that there is generally little difference in the results for changes in material properties or site-specific conditions and that the parameters that do have a bigger influence on the factor of safety, and slip surface are the model settings. The exception to this is complex discontinuous joint networks where the spacing can have a signif- icant impact on both the factor of safety reported and the location of the failure surface due to limitations in the way discontinuities are handled within UDEC. Ta- ble 4.13 below summarizes the results of the sensitivity analysis with comments on their effect on the SSR method. 94 Chapter 4. Parametric Study Ta bl e 4. 13 :S um m ar y of pa ra m et ri c st ud y re su lts Se tt in g /P ar am et er SR F Se ns iti vi ty E ff ec t on Po si tio n of Fa ilu re su rf ac e C om m en ts N um be ro fc yc le s or ra - tio lim it us ed fo r ea ch st ag e L ar ge Sm al l T he hi gh er th e nu m be ro fc yc le s or lo w er th e ra tio lim it us ed , th e m or e co ns er va tiv e th e fa ct or of sa fe ty re tu rn ed .F ai lu re su rf ac es ar e oc ca si on al ly m or e sh al lo w th e m or e cy cl es th at oc cu r. E le m en tS iz e L ar ge Sm al l T he el em en t si ze ha s a m aj or in flu en ce on th e fa ct or of sa fe ty .E xt re m el y sm al le le m en t si ze s ar e to be av oi de d. E le m en ts iz e af fe ct s fa ilu re su rf ac e po si tio n du e to in flu en ce on lo ca liz at io n an d sh ea rb an di ng . E le m en tT yp e N eg lig ib le Sm al l Q ua dr ila te ra l el em en ts ar e re po rt ed ly be tte r fo r pl as tic m at er ia ls ;h ow ev er m os tp ro bl em do m ai ns w ill re qu ir e bo th qu ad ri la te ra l an d tr ia ng ul ar el em en tt yp es du e to lim ita tio ns in sp at ia ld is cr et iz at io n. In si tu St re ss N eg lig ib le N eg lig ib le A bs ol ut e di sp la ce m en t va lu es ar e af fe ct ed bu tt he on se to fi ns ta bi lit y oc cu rs at th e sa m e SR F st ag e. E la st ic Pr op er tie s N eg lig ib le N eg lig ib le A bs ol ut e di sp la ce m en t va lu es ar e af fe ct ed bu tt he on se to fi ns ta bi lit y oc cu rs at th e sa m e SR F st ag e w ith in a sm al l ra ng e. L ite ra tu re ca ut io ns ag ai ns t us in g m at er ia ls w ith dr as - tic al ly di ff er en t st iff ne ss pr op er tie s w he n co m pa ri ng to L E M so lu tio ns . 95 Chapter 4. Parametric Study Ta bl e 4. 13 :S um m ar y of pa ra m et ri c st ud y re su lts (c on t.) Se tt in g /P ar am et er SR F Se ns iti vi ty E ff ec t on Po si tio n of Fa ilu re su rf ac e C om m en ts Jo in tS he ar St iff ne ss V ar ia bl e Sm al l St iff ne ss af fe ct s th e fa ct or of sa fe ty , ho w - ev er ,t he re is no sy st em at ic re la tio ns hi p be - tw ee n th e tw o. T hi s pr op er ty ha s a la rg e af - fe ct on th e tim e re qu ir ed fo ra m od el to co m e to eq ui lib ri um an d oc ca si on al ly sh ow s di f- fe re nt be ha vi ou rd ep en di ng on th e fo rc e ra tio lim it us ed fo re ac h SR F st ep . St ra in -s of te ni ng m od - el s (i nt ac tb lo ck s) L ar ge Sm al l T he tr an si tio n to re si du al st re ng th pa ra m e- te rs re su lts in a lo w er fa ct or of sa fe ty ov er - al l. L ow er st ra in va lu es be fo re on se to f fu ll re si du al st re ng th s re su lt in a lo w er fa ct or of sa fe ty . T he re is a sl ig ht ef fe ct on po si tio n of fa ilu re su rf ac es w ith sl ig ht ly sh al lo w er fa il- ur e su rf ac es oc cu rr in g at lo w er st ra in lim it va lu es . Jo in tN et w or k V ar ia bl e Sm al l Fo r an y co nt in uo us jo in ts ,t he fr ac tu re de n- si ty ha s an al m os t ne gl ig ib le af fe ct on th e fa ct or of sa fe ty an d po si tio n of cr iti ca ls lid - in g su rf ac e. Fo r di sc on tin uo us jo in ts , fr ac - tu re de ns ity co nt ro ls th e si ze an d di st ri bu - tio n of ro ck br id ge s w ith in th e m od el an d ca n ha ve a la rg e af fe ct on th e fa ct or of sa fe ty re tu rn ed . H ea vi ly fr ac tu re d ro ck m as se s ar e no tn ec es sa ri ly le ss st ab le th an le ss fr ac tu re d on es du e to di st ri bu tio n of ro ck br id ge s an d lim ita tio ns of jo in te d ro ck w ith in U D E C . 96 Chapter 5 Case History 5.1 Introduction The case history presented will remain unnamed due to a confidentiality agree- ment, but involves an open pit failure that occurred in the south-western United States. The purpose of using the case history is to apply the SSR method to a complex rock slope problem using material parameters and knowledge from site investigations of an operational mine. The slope failure case provides a known failure surface that can be used for model validation of the SSR method and al- lows for comparisons between SSR and LEM solutions. It also allows a check of some of the results from the parametric study, as well as the study of a few new options due to the increased complexity of the slope and failure mechanism. Thus, this series of models and discussions is to test whether SSR is a valuable tool for analysis of large slopes as well as helping to further refine the guidelines for using SSR presented in the next chapter. 5.1.1 Geology The regional geology of the case history site is comprised of early and middle paleozoic siliclastic sedimentary and volcanic rocks that have been thrust upon shallow water, carbonate-rich rocks of the continental platform during the late Devonian. There are two main rock units that comprise the failure area, given the monikers upper and lower units based on their locations in the stratigraphic section. The upper unit is a locally carbonaceous, bedded, siliclastic mudstone with chert, 180 m in thickness, and can be broken into two different units: 97 Chapter 5. Case History • Basal unit - gray carbonate-cemented quartz sandstone, 0 – 21 m in thick- ness. This unit is absent in the section of the failed slope. • Upper unit - interbedded, dark grey, sillaceous mudstone and light grey, laminated silty, limestone with sub-meter bedding spacing The lower unit, a 180 m thick limestone, is similarly broken into two different units: • Upper limestone - dark grey to black, carbonaceous, planar laminated, thinly bedded to massive, micritic, silty limestone with sub-meter bedding spacing • Lower limestone - dark grey, carbonaceous, micritic, muddy limestone with local green-gray marble alteration with sub-meter bedding spacing The deformed transitional contact zone between the upper and lower units is struc- turally important for the failure case being investigated. The contact zone varies in thickness (2 – 8 m) and is composed of many different rock types: bedded siltstones, sandstones, cherts, argillite, slatty mudstone, phyllite rock, massive muddy limestones, limy mudstones, and red–knobby breccias. The contact is also described as containing tectonically-derived breccias formed from shear zones, faults, and folding . Geological reports mention ripped up breccia in this contact area with 0.5 - 2.0 cm angular clasts of argillite and shell fragments in a sandy to silty matrix of similar material, and may be locally silicified. Variations in the contact zone are not captured in the geotechnical data. 5.1.2 Slope Failure The failure transpired as a one hour event and the area had been cleared of per- sonnel and equipment due to the early warnings provided by monitored elevated movements in the weeks prior to the failure event. The failure is approximately 120 m high, 250 m wide and has an estimated depth of 25 - 30 m. A cross-section, seen in figure 5.1, shows the geology of the pre-failure slope along a section ori- ented through the thickest part of the failure, with a the lower boundary extended for modelling purposes. 98 Chapter 5. Case History Waste Upper Unit Contact Zone Lower Unit H = 120 m Figure 5.1: Cross-section of pre-failure pit, showing major geological units Slope failure is thought to have initiated as significant toe heave (upwards of 12 m) in the contact zone, that was documented on May 18th, 2005. The failed zone has some structural controls and appears to be: • bounded at the crest by some high-angle persistent structure perpendicular to the dip direction of the pit slope • bounded to the east by a large fault • bounded to the west by two large faults A photo taken on July 13th, 2005 showing the failure and major structures can be seen in figure 5.2. The failure is attributed to a range of conditions that adversely affected this portion of the pit wall: 1. Daylighting of the low strength carbonaceous, black, silt/clay contact that is adversely dipping (averagely) 16o into the pit at the toe of the failure. 2. 115 mm of rainfall over an 11-day period prior to failure, with surface water infiltrating tensions cracks open at the crest. During a peak rainfall event of 53 mm over 24 hours on May 16th, an estimated 400 l/min of water was entering these tension cracks. 3. Pore pressure build-up due to water influx. 99 Chapter 5. Case History At the time of failure, there were no piezometers monitoring water levels within the slope, which resulted in limited hydrological inputs and constraints for the models. Figure 5.2: View of the south wall showing the failure, taken July 13th 2005. 5.1.3 Structure The structural mapping data has been provided by Piteau Associates through var- ious internal reports. The failure occurred within a single structural domain. The rockmass contains a large amount of structures, summarized in the following stereonets: • Figure 5.3 showing fault orientations • Figure 5.4 showing bedding orientations • Figure 5.5 showing joint orientations There is no detailed information on persistence, aperture, infilling, trace length, and spacing of the various discontinuities in the failure area. Discussions with engineers that have been on site indicate that structure is non-continuous, and there are many small structures with a small spacing in addition to the larger structures that have been mapped in the field. Only large scale features will be modelled in the analysis because: 100 Chapter 5. Case History BY: DATE: APPROVED: DWG: P I T E A U   A S S O C I A T E S GEOTECHNICAL & HYDROGEOLOGICAL CONSULTANTS VANCOUVER                                                         CALGARY PREPARED SOLELY FOR THE USE OF OUR CLIENT AND NO REPRESENTATION OF ANY KIND IS MADE TO OTHER PARTIES WITH WHICH PITEAU ASSOCIATES ENGINEERING LTD. HAS NOT ENTERED INTO A CONTRACT. LIMA fil e:  h :/p ro je ct /8 9- 14 8a /2 00 5/ m ay 05 /s tru ct /D ip s- C om bi ne d/ B U ZZ -F .g rf IRS JUL 05STRUCTURAL ANALYSIS OF FAULTS IN BUZZARD DOMAIN (RODEO CREEK FORMATION) BARRICK GOLDSTRIKE MINES INC. OPEN PIT GEOTECHNICAL ASSESSMENTS 11TH WEST LAYBACK EQUAL  AREA LOWER HEMISPHERE N E S W CONTOUR LEGEND SCHMIDT POLE CONCENTRATIONS % of total per   1.0 % area Minimum Contour =   0.5 Contour Interval =   0.5 Max.Concentration =  10.5 842 Poles Plotted 842 Data  Entries MAJOR   PLANES  ORIENTATIONS SET   DIP/DIR.   %        # FA1 FA2 FA3 FA4 FC1 FC2 FD1 FB1 FB2 66/245 67/231 64/260 64/280 55/296 79/299 77/137 64/108 59/072 10.5 6 6 4.5 3.5 2 2 2 1.5 88 51 51 38 29 17 17 17 13 FA1 FA2 FA3 FA4 FC1 FC2 FB2 FB1 FD1 B20 Figure 5.3: Stereonet showing faults in the failure area 101 Chapter 5. Case History BY: DATE: APPROVED: DWG: P I T E A U   A S S O C I A T E S GEOTECHNICAL & HYDROGEOLOGICAL CONSULTANTS VANCOUVER                                                         CALGARY PREPARED SOLELY FOR THE USE OF OUR CLIENT AND NO REPRESENTATION OF ANY KIND IS MADE TO OTHER PARTIES WITH WHICH PITEAU ASSOCIATES ENGINEERING LTD. HAS NOT ENTERED INTO A CONTRACT. LIMA fil e:  h :/p ro je ct /8 9- 14 8a /2 00 5/ m ay 05 /s tru ct /D ip s- C om bi ne d/ B U ZZ -B .g rf IRS JULY 05STRUCTURAL ANALYSIS OF BEDDING IN BUZZARD DOMAIN (RODEO CREEK FORMATION) BARRICK GOLDSTRIKE MINES INC. OPEN PIT GEOTECHNICAL ASSESSMENTS 10TH WEST LAYBACK EQUAL  AREA LOWER HEMISPHERE N E S W CONTOUR LEGEND SCHMIDT POLE CONCENTRATIONS % of total per   1.0 % area Minimum Contour =     2 Contour Interval =     2 Max.Concentration =  23.2 512 Poles Plotted 512 Data  Entries MAJOR   PLANES  ORIENTATIONS SET   DIP/DIR.   %        # BG1 BG2 BG3 15/353 19/058 23/324 23.2 12 12 119 61 61 BG1 BG2 BG3 FIG: B22 Figure 5.4: Stereonet showing bedding in the failure area 102 Chapter 5. Case History BY: DATE: APPROVED: DWG: P I T E A U   A S S O C I A T E S GEOTECHNICAL & HYDROGEOLOGICAL CONSULTANTS VANCOUVER                                                         CALGARY PREPARED SOLELY FOR THE USE OF OUR CLIENT AND NO REPRESENTATION OF ANY KIND IS MADE TO OTHER PARTIES WITH WHICH PITEAU ASSOCIATES ENGINEERING LTD. HAS NOT ENTERED INTO A CONTRACT. LIMA fil e:  h :/p ro je ct /8 9- 14 8a /2 00 5/ m ay 05 /s tru ct /D ip s- C om bi ne d/ B U ZZ -J .g rf IRS JUL 05STRUCTURAL ANALYSIS OF JOINTS IN BUZZARD DOMAIN (RODEO CREEK FORMATION) BARRICK GOLDSTRIKE MINES INC. OPEN PIT GEOTECHNICAL ASSESSMENTS 11TH WEST LAYBACK EQUAL  AREA LOWER HEMISPHERE N E S W CONTOUR LEGEND SCHMIDT POLE CONCENTRATIONS % of total per   1.0 % area Minimum Contour =   0.5 Contour Interval =   0.5 Max.Concentration =  12.4 485 Poles Plotted 485 Data  Entries MAJOR   PLANES  ORIENTATIONS SET   DIP/DIR.   %        # JA1 JD1 JA2 JA3 JD2 JD3 JC1 JB1 JC2 JC3 JB2 75/240 74/135 73/223 77/248 78/156 80/173 76/298 84/020 84/338 75/326 69/055 12.4 10 6 6 6 2.5 1.5 1.5 1.5 1.5 1 60 49 29 29 29 12 7 7 7 7 5 JD1 JD2 JD3 JA2 JA1 JA3 JB2 JB1 JC2 JC3 JC1 B21 Figure 5.5: Stereonet showing joints in the failure area 103 Chapter 5. Case History • It is impossible to model every possible discontinuity in the slope due to memory and computational time limits • The parametric study shows that only key discontinuities are important • Key discontinuities are likely to be large scale, multi-bench features Also, major structures from field mapping done on site were plotted on pit plans showing contours and then used to give a starting point for what joint spacing and persistence should be used to generate the fracture network in the UDEC model. The modelled section through the failure is orientated at 010o. Up to three struc- tural sets were chosen to be modelled: bedding ranging from a dip of 15 to 20 degrees, cross-jointing dipping into the pit at orientations of 70-85 degrees, and a second set of cross-joints dipping away from the pit at 70-85 degrees. These ranges are given as being representative of the variation and uncertainty present in the actual slope and were incorporated into the UDEC models as such. 5.1.4 Geotechnical Properties A series of site investigations have been done over time at the site, most recently in 2008 in response to movements in the pit wall. Triaxial, uniaxial, and brazillian tests were preformed on core samples from the pit. Piteau Associates have pro- vided a summary of their testing and mapping data to generate the Hoek-Brown properties for the pertinent rockmasses in the failure area. Some limited data on the strength of discontinuities was provided along with recommendations for strength ranges that have been used in previous modelling throughout the pit. Ta- ble 5.1 summarizes the available strength data used in the modelling. A single value of Young’s Modulus (E) and Poisson’s ratio (ν) was chosen from lab testing results for each unit even though they can vary significantly within units. Models were run to test the influence of this parameter and it was found to be negligible over the range of possible E and ν values. This non-sensitivity was also shown in the parametric analysis. The linear Mohr-Coulomb strengths found in table 5.1 are calculated from the 2002 Hoek-Brown strength criterion and empirical corrections, using a distur- bance factor of 1.0 and confinement of 1.0 MPa. These two parameters were selected after discussions with Piteau Associates about previous modelling work 104 Chapter 5. Case History Table 5.1: Summary of geotechnical parameters used in the case history Unit Density ν E UCS RMR mi mb s a Cohesion φ Tensile Strength (kg/m3) (GPa) (MPa) (kPa) (kPa) Upper unit 2400 0.112 2.39 24 46 12 0.255 1.25 x 10−4 0.507 185 29.7 8.3 Upper Limestone 2560 0.276 8.00 38 50 17 0.499 2.66 x 10−4 0.506 278 39.3 11.8 Contact Zone 2240 - - - - - - - - 48 25 8.3 Waste 1922 0.2 3 - - - - - - 0 39 0 Joints & Bedding Planes - - - - - - - - - 0 - 70 23-26 0 done in the pit and the use of Mohr-Coulomb approximations of Hoek-Brown properties. The same conservative ‘equivalent continuum’ strength reduction to account for unmodelled joints is used regardless of the joint network modelled as there is no clear relationship between modelled joint density and intact strength reduction in the Hoek-Brown criteria. There is limited data for the contact zone properties. A starting point of φ = 25o and cohesion = 48 kPa was selected. These values coincide with the strength of bedding discontinuities within the Upper unit. Discontinuity strengths vary, particularly with the cohesion component. There is limited data on the strength of discontinuities in the pit, with most testing focused on fault gouge strengths. The distribution of these faults is not known within the slope section, requiring the use of the representative values of φ = 25o and cohesion = 48 kPa to be used. Geometric properties of the joint network are not well documented in terms of persistence, continuity, rock bridging, etc which when combined with modelling limitations (i.e. inability to model every joint present in the field) makes it possible for a range of plausible joint networks. Joint normal stiffness (kn) and joint shear stiffness (ks) are selected from similar units in reference data provided in Kulhawy (1975). A ks value of 8 x 108 Pa/m and a kn value of 10 x 108 Pa/m are used as a starting point unless otherwise noted. A range of shear stiffnesses is tested within the models due to the influence of ks on failure surface and factor of safety as seen previously in the parametric 105 Chapter 5. Case History study, as well as the range of ks values reported in Kulhawy (1975) for similar rock types. 5.1.5 Modelling Approach The case history was undertaken here with the following three objectives: 1) to compare the SSR method to the LEM method, 2) to predict the failure with no a priori surface assumptions, and 3) to provide a factor of safety. These were then used to help further refine SSR guidelines (as presented in Chapter 6). The UDEC models build from a simple unjointed configuration to more complicated joint models and ends with the introduction of a water table and strain-softening effects, while keeping the three objectives in mind. The failure occurred within the Contact zone and the Upper unit only. Because of this, only those units undergo shear strength reduction in the analysis in order to constrain the development of the shear surface. The waste pile on top of the pit does not undergo SSR to inhibit small scale failures that are not important to the problem. This assumption has been checked and shown to be correct. It is used to constrain failure and for numerical efficiency purposes. Similar to the parametric study, the problem is zoned in terms of element size and jointing. An element size of 5 m is chosen for the Upper and Contact zone areas. A coarser 10 m element size is used for the waste and Upper unit limestone. A 15 m mesh is used for the benches that are removed. Finally, a 20 m mesh size is used at the base of the problem far away from the area of interest. Joints away from the area of interest are typically twice the spacing of the joints in the Upper unit and Contact zone areas. This setup, for computational efficiency, can be seen in figure 5.6 for the configurations that included jointing. The problem follows the SSR outline and modelling procedure given in section 3.6.1: geometry, jointing, elastic gravity loading, plastic loading, bench removal, and finally SSR. A Mohr-Coulomb constitutive model is used for the intact blocks with a Coulomb-slip model for joints. Strength parameters are given in table 5.1. Boundary conditions were set to zero velocities for the left and right boundaries of the model in the x-direction, and zero velocities for the bottom of the model in the y-direction. This allows for movements in the y-direction on the left and right edges, and movements in the x-direction at the base of the model. An automatic solution limit of 1 x 106 for force ratio was used for both D-SSR and UDEC-SSR 106 Chapter 5. Case History    UDEC (Version 4.00) LEGEND     7-Aug-08  11:28  cycle    86920  time    4.387E+01 sec block plot -1.500 -0.500  0.500  1.500  2.500  3.500  4.500  5.500 (*10^2)  0.500  1.500  2.500  3.500  4.500  5.500  6.500  7.500 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada          (a) Joints    UDEC (Version 4.00) LEGEND     7-Aug-08  11:28  cycle    86920  time    4.387E+01 sec zones in fdef blocks -1.500 -0.500  0.500  1.500  2.500  3.500  4.500  5.500 (*10^2)  0.500  1.500  2.500  3.500  4.500  5.500  6.500  7.500 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada          (b) Mesh Figure 5.6: Example of telescoping of joints (a) and meshing (b) within slope for the case history 107 Chapter 5. Case History unless otherwise stated. 108 Chapter 5. Case History 5.1.6 Selection of Critical Sliding Surface Selection of the sliding surface in UDEC requires explanation as there are many different indicators that can be used but none that do not require some degree of interpretation. . UDEC provides two convenient output plots that can be used to infer a sliding surface: yielded elements and velocity vectors. The criteria used in this thesis for the selection of the final sliding surface is agreement between the sliding surfaces predicted by both the yielded elements and the velocity vec- tors. Each individual element can exist in one of three states based on the elasto-plastic constitutive criteria used: unyielded, yielded in shear, yielded in tension. A com- mon indicator of failure used for the case history analysis was the development of a tension crack opening at the crest that intersects at depth with a large shear band through the Upper unit and Contact zone. However, there are many intermediate stages of yielding that may develop apart from the ultimate shear zone. Table 5.2 below lists a series of model states and associated figures showing yielded ele- ments, which can be used to track the development of the shear band at different cycle stages of a particular model. Table 5.2: Model and cycle stages of yielded elements Model SRF State Figure UDEC-SSR 1.12 10−6 limit 5.7 D-SSR 1.11 5000 cycles 5.8 D-SSR 1.11 10−6 ratio limit 5.9 D-SSR 1.11 10−6 ratio limit plus 5,000 cycles 5.10 D-SSR 1.11 10−6 ratio limit plus 15,000 cycles 5.11 D-SSR 1.11 10−6 ratio limit plus 30,000 cycles 5.12 D-SSR 1.11 10−6 ratio limit plus 50,000 cycles 5.13 UDEC-SSR and D-SSR at the SRF stage of failure cycled for 5000 time steps gives the same plastic shear zone with a second defined internal shear zone, seen in figures 5.7 and 5.8. In these plots the surface is quite wide and it is not possible to define a distinct, localized surface that the failure is occurring along. Looking at 109 Chapter 5. Case History    UDEC (Version 4.00) LEGEND     6-Aug-08  11:07  cycle    86278  time    5.375E+01 sec block plot no. zones : total     13530 at yield surface (*)   2984 yielded in past  (X)      0 tensile failure  (o)    795  1.750  2.250  2.750  3.250  3.750  4.250 (*10^2)  2.500  3.000  3.500  4.000  4.500  5.000 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.7: Yielded elements at 10−6 ratio limit for UDEC-SSR giving FS = 1.12 (x yield in shear, o yield in tension)    UDEC (Version 4.00) LEGEND     7-Aug-08  13:23  cycle    87062  time    5.424E+01 sec block plot no. zones : total     13530 at yield surface (*)   3101 yielded in past  (X)      0 tensile failure  (o)    843  1.750  2.250  2.750  3.250  3.750  4.250 (*10^2)  2.500  3.000  3.500  4.000  4.500  5.000 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.8: Yielded elements at 5,000 cycles for D-SSR at SRF = 1.11 (x yield in shear, o yield in tension) 110 Chapter 5. Case History    UDEC (Version 4.00) LEGEND     7-Aug-08  13:03  cycle   113880  time    7.121E+01 sec block plot no. zones : total     13530 at yield surface (*)    925 yielded in past  (X)      0 tensile failure  (o)    140  1.750  2.250  2.750  3.250  3.750  4.250 (*10^2)  2.500  3.000  3.500  4.000  4.500  5.000 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.9: Yielded elements at 10−6 ratio limit for D-SSR at SRF = 1.11 (x yield in shear, o yield in tension)    UDEC (Version 4.00) LEGEND     7-Aug-08  16:16  cycle   118880  time    7.437E+01 sec block plot no. zones : total     13530 at yield surface (*)   2484 yielded in past  (X)      0 tensile failure  (o)    724  1.750  2.250  2.750  3.250  3.750  4.250 (*10^2)  2.500  3.000  3.500  4.000  4.500  5.000 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.10: Yielded elements at 10−6 ratio limit plus 5,000 cycles for D-SSR at SRF = 1.11 (x yield in shear, o yield in tension) 111 Chapter 5. Case History    UDEC (Version 4.00) LEGEND     7-Aug-08  12:34  cycle   128880  time    8.070E+01 sec block plot no. zones : total     13530 at yield surface (*)   2544 yielded in past  (X)      0 tensile failure  (o)    609  1.750  2.250  2.750  3.250  3.750  4.250 (*10^2)  2.500  3.000  3.500  4.000  4.500  5.000 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.11: Yielded elements at 10−6 ratio limit plus 15,000 cycles for D-SSR at SRF = 1.11 (x yield in shear, o yield in tension)    UDEC (Version 4.00) LEGEND     7-Aug-08  12:27  cycle   143880  time    9.018E+01 sec block plot no. zones : total     13530 at yield surface (*)   1615 yielded in past  (X)      0 tensile failure  (o)    209  1.750  2.250  2.750  3.250  3.750  4.250 (*10^2)  2.500  3.000  3.500  4.000  4.500  5.000 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.12: Yielded elements at 10−6 ratio limit plus 30,000 cycles for D-SSR at SRF = 1.11, showing the ultimate failure surface (x yield in shear, o yield in tension) 112 Chapter 5. Case History    UDEC (Version 4.00) LEGEND     7-Aug-08  13:48  cycle   163880  time    1.028E+02 sec block plot no. zones : total     13530 at yield surface (*)   1643 yielded in past  (X)      0 tensile failure  (o)    221  1.750  2.250  2.750  3.250  3.750  4.250 (*10^2)  2.500  3.000  3.500  4.000  4.500  5.000 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.13: Yielded elements at 10−6 ratio limit plus 50,000 cycles for D-SSR at SRF = 1.11, showing the ultimate failure surface (x yield in shear, o yield in tension) the output from the D-SSR analysis in figure 5.9, the plastic indicators show that a shallower failure is occurring; however, cycling the model beyond this point (fig- ures 5.10–5.13) reveals that this is an intermediate shear zone representing internal shearing of the failing mass. As cycling continues, a clearly defined ultimate fail- ure surface eventually forms as in figure 5.12 and this surface does not change with additional timesteps as evidenced in figure 5.13. This shows that the 10−6 ratio limit is a good criteria to use to determine the factor of safety, but additional timesteps may be required in order for an ultimate shear surface to be clearly seen in the yielded element output. Velocity plots can be used to look at portions of the slope that are actively moving by comparing the relative size of velocity vectors within the model. There is a clear break between the activated slope (large vectors) and inactive slope (small vectors) that can be used to define a failure in the analysis output. Figure 5.14 provides an example. These can then be externally plotted for different models and model states and compared, as in figure 5.15. This plot shows that using velocity vectors to determine the location and shape of the sliding surface results in a similar result regardless of the state (current number of cycles, current force ratio) of the model. If the yielded elements are superimposed on this plot, as in figure 5.16, it is clear that the surface from the velocity plot is the same as the 113 Chapter 5. Case History ultimate shear surface predicted by yielded elements. The intermediate surfaces, which may be important, are not seen in the velocity plots. As these outputs represent model states that are at the limit of stability, the velocity vectors in the D-SSR results have a tendency to be difficult to interpret with no clear break between vector sizes. However, figure 5.16 shows that time stepping beyond the model state which reaches the ratio limit provides a clearly defined failure surface. Throughout the following sections, SRF values are often included with the differ- ent plots presented. It is assumed that these SRF values are for the critical SRF stage (i.e. the factor of safety), unless otherwise noted.    UDEC (Version 4.00) LEGEND     6-Aug-08  11:07  cycle    86278  time    5.375E+01 sec velocity vectors      maximum =    9.559E-02 0   5E -1 block plot  1.750  2.250  2.750  3.250  3.750  4.250 (*10^2)  2.500  3.000  3.500  4.000  4.500  5.000 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Slip surface dened by velocity vectors representing active portion of the slope Figure 5.14: Example of selecting sliding surface from UDEC velocity vector output 114 Chapter 5. Case History 275 325 375 UDEC-SSR [1E-6] D-SSR SRF = 1.11 [1E-6 + 15,000 cycles] D-SSR SRF = 1.11 [1E-6 + 30,000 cycles] D-SSR SRF = 1.11 [1E-6 + 50,000 cycles] D-SSR SRF=1.11 [5,000 cycles] 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters Figure 5.15: Velocity plot derived surfaces from models states listed in table 5.2 115 Chapter 5. Case History 275 325 375 Surfaces from velocity Intermediate shear surface [at SRF 1.11 (1E‐6)] 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters Ultimate shear surface [at SRF 1.11 (1E‐6 + 15,000 ‐ 50,000 cycles)] Figure 5.16: Velocity plot derived surfaces from models states listed in table 5.2 with yielded elements overlain 116 Chapter 5. Case History 5.2 Results This section outlines the results of models of increasing geologic complexity. The LEM solution allows for a comparison to be made with the SSR method in UDEC for both the position of the failure surface and the value of FOS reported. The depth of the actual failure surface is largely unknown due to the remaining volume; however, the headscarp is known and the deepest possible sliding surface can be assumed by projecting the headscarp into the Contact zone area. Thus, the sliding surfaces predicted in UDEC can also be compared to the location of the headscarp and whether or not the rest of the sliding surface lies deeper than the largest possible sliding surface. 5.2.1 Limit Equilibrium A LEM model was provided by Piteau Associates for reference purposes to test the the SSR method. Limits were set to search for the critical surface in order to constrain it more readily to the actual failure surface. The linear Mohr-Coulomb strengths found in table 5.1 were used. A tension crack is explicitly defined in some models to force the failure headscarp to match the actual headscarp. The model is run both dry, and with pore pressures via the use of an Ru value of 0.15 for the Upper unit and contact zone. The latter is to account for the influx of water that occurred in the week prior to the failure. The model was run for a variety of cases, with results in table 5.3 and lowest fac- tor of safety sliding surfaces in figure 5.17 the non-circular block search method. There are two possible sliding surfaces predicted by SLIDE. Surface One, which matches the head scarp found in the field is only possible if a tension crack is explicitly modelled in that area and the Janbu Corrected solution is used. If there is no head scarp, or the General Limit Equilibrium (GLE) method is used, the predicted failure surface is the deeper-seated Surface Two. The inclusion of pore pressures results in the same failure surfaces but with lower FOS values. Phys- ically modelling a water surface rather than using the Ru parameter results in a wide range of factor of safeties due to the unknown position of the water table within the slope. Failure surfaces for these cluster around LEM Surface Two re- gardless of solution mode. This makes the Ru = 0.15 condition questionable as it is not one that can be reproduced using a modelled water table In these modelled 117 Chapter 5. Case History water table results, factor of safeties range from 1.30 - 1.50 for surfaces matching the actual failure surface. Table 5.3: Summary of LEM results, surfaces can be seen in figure 5.17 Model Tension Crack Water Table Method FOS Surface Identifier A Yes Yes Janbu Corrected 0.99 1 A Yes Yes GLE 0.98 2 B Yes No Janbu Corrected 1.20 1 B Yes No GLE 1.18 2 C No Yes Janbu Corrected 1.00 2 C No Yes GLE 0.98 2 D No No Janbu Corrected 1.20 2 D No No GLE 1.20 2 118 Chapter 5. Case History 275 325 375 Approximate failure surface  based on field observations and  geologic projects (exact depth  unknown) referred to as the  failure surface LEM Surface One: ‐ Tension crack ‐ Dry or Ru = 0.15 ‐ Janbu Corrected LEM Surface Two: ‐ Tension crack 125 175 225 150 200 250 300 350 400 450 500 550 ‐ GLE/MP ‐ Dry or Ru = 0.15 OR ‐ No tension crack ‐ Janbu Corrected, GLE/MP ‐ Dry or Ru = 0.15 m e t e r s meters Figure 5.17: Critical LEM failure surfaces predicted by SLIDE 119 Chapter 5. Case History 5.2.2 UDEC: Unjointed Slope Models A set of UDEC models was created to compare the results to those from the LEM analysis, by omitting all joints or other discontinuities in the model other than the contacts between the geologic models themselves. This was taken as being as close as possible in terms of geotechnical representation, for a UDEC and LEM model to be, and is also the most basic geometry assessed for the case history. The simplified block geometry can be seen in figure 5.18. Intact properties and contact properties (joints) are as per table 5.1. D-SSR was run for a force ratio limit of 10−6 and UDEC-SSR was run for a ratio limit of 10−6 and 10−5. The following model series were run and compared to the LEM solution with results seen in table 5.4: • U1-U3: 2.5, 5, and 10 m element sizes in the area of interest to determine which element size best represents the failure. • U4: Insitu stress (σh : σv) = 0.5. As factor of safety is near 1 these properties may be important. These models are run to double-check the findings from section 4.4 of the parametric study. • U5: All materials have the same E and ν of 2.4GPa and 0.11 respectively in order to give the same bulk and shear modulus. The guidelines in Hammah et al. (2005a) state that elastic properties should be the same for each unit in order to compare a SSR solution to an LEM solution. These models check if this is true or not for the case history in question. Figure 5.19 shows the results of D-SSR used for the cases listed in table 5.4. In all cases, the D-SSR returns a slightly lower value for the SRF than UDEC-SSR; however, the difference is small. Similarly, UDEC-SSR run at a lower force ratio (10−6) gives a lower SRF which agrees with results from the parametric study. The results of the SSR method in UDEC are more conservative, but closely resembles those obtained using LEM (1.18-1.20 depending on the solution method). The 10 m element size case gives a larger FOS value for D-SSR and UDEC-SSR at a ratio limit of 10−6 and a significantly larger FOS for UDEC-SSR at a ratio limit of 10−5. This is likely due to the poor shear band localization when using larger elements, though the UDEC-SSR result of 1.42 is anomalous and indicators as to why are not present within UDEC. Figure 5.19 shows that the model is not particularly sensitive to the insitu stress ratio, even for cases where the factor of safety is close to 1, which reinforces the 120 Chapter 5. Case History    UDEC (Version 4.00) LEGEND    28-Jul-08   5:59  cycle    77440  time    5.696E+01 sec block plot -1.500 -0.500  0.500  1.500  2.500  3.500  4.500  5.500 (*10^2)  0.500  1.500  2.500  3.500  4.500  5.500  6.500  7.500 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.18: Unjointed model geometry after bench removal Table 5.4: Summary of unjointed SSR results Model Element σh : σv Bulk/Shear D-SSR UDEC-SSR UDEC-SSR size (m) Modulus 10−5 10−6 U1 2.5 2.0 Variable 1.12–1.14 1.13 1.12 U2 5.0 2.0 Variable 1.13–1.15 1.17 1.14 U3 10.0 2.0 Variable 1.17 1.42 1.21 U4 5.0 0.5 Variable 1.13–1.15 1.16 1.14 U5 5.0 2.0 Consistent 1.13–1.15 1.16 1.16 121 Chapter 5. Case History 8 00 7.00 . 6.00 m e n t   ( m ) 4.00 5.00 n t a l   d i s p l a c e U1 U2 3.00 d e l l e d   h o r i z o n U3 U4 U5 1.00 2.00M o d 0.00 1.08 1.10 1.12 1.14 1.16 1.18 1.20 1.22 1.24 Strength reduction factor Figure 5.19: Results of displacement-based SSR for models U1-U5 showing a critical SRF range of 1.12 – 1.17 122 Chapter 5. Case History findings from the parametric study. Also, there is negligible difference between having different elastic constants for each geologic unit compared to those with the same values used throughout (U5). There are slight differences in the displace- ment values at each SRF stage between U5 and U2; however, the point of failure is identical. Lastly, the smaller the element size, the lower the strength reduction factor at the point of failure for both D-SSR and UDEC-SSR. This agrees with findings from the parametric study. The size and location of the failed surface is an equally important output of the SSR analysis as the absolute value of factor of safety. Using the methodology set out in section 5.1.6, failure surfaces can be compared for the unjointed models as follows:: • U1, using plasticity indicators, for displacement-based SSR at SRF = 1.13 in figure 5.20 • U1, using velocity vectors, for UDEC-SSR velocity surfaces for ratio limits of 10−5 and 10−6 in figure 5.21 • U2, comparing the use of plasticity indicators and velocity vectors, for displacement-based SSR at SRF = 1.14 in figure 5.22 • U2, using velocity vectors, for UDEC-SSR velocity surfaces for ratio limits of 10−5 and 10−6 in figure 5.23 The failure mechanism is the same in all cases, showing the same general failure as in the previous plasticity plots, i.e., tension cracks at some distance behind the crest with a shear band extending from the toe through the contact zone before breaking up through the Upper unit. However, looking at figures 5.20 – 5.23 it can be seen that none of the estimated critical slip surfaces match the failure sur- face or the LEM solution which includes the explicitly modelled tension crack (i.e. Surface One in figure 5.17). This indicates that the complete failure mecha- nism is not being captured in these simple unjointed models beyond capturing the influence of the weak contact zone. Some models do match LEM Surface Two; however they are dependent on ele- ment size and force ratio selection. In both UDEC-SSR models in U1/U2, using a force ratio of 10−5, the failure surface is significantly deeper than that of LEM Surface Two. Similarly, U2 with an element size of 5 m using the D-SSR method gives a poor match to Surface Two by predicting a deeper failure. However, both the D-SSR and UDEC-SSR methods when cycled to a force ratio limit of 10−6 123 Chapter 5. Case History using a 2.5 m element size (figures 5.20 and 5.21 respectively) predict failure surfaces identical to LEM Surface Two. The D-SSR solution shows a band of yielded elements which brackets the surface while the velocity-derived surface from UDEC-SSR follow the LEM solution almost exactly. The simple unjointed models shows that some amount of structure needs to be defined in the model in order for the failure surface to match the failure, which likely lies between LEM Surface One and the projected deepest possible failure surface. However, the simple model (in terms of geology and problem geometry) does verify that both SSR methods used in the UDEC modelling produce results that match those found using a LEM solution, both in terms of a similar factor of safety and failure surface predicted, provided a ratio limit of 10−6 and element size of 2.5 m is used for the SSR analysis. 275 325 375 Actual failure surface LEM Surface One LEM Surface Two Yielded elements (U1) SRF = 1 13 Internal shearing of  failed mass 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters     . Figure 5.20: Comparison of LEM surfaces and D-SSR SRF = 1.13 yielded ele- ments for model U1 124 Chapter 5. Case History 275 325 375 Actual failure surface LEM Surface One LEM Surface Two UDEC‐SSR 1E‐5 (U1) (velocity) UDEC‐SSR 1E‐6 (U1) 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters     (velocity) Figure 5.21: Comparison of LEM surfaces and UDEC-SRF velocity-derived sur- faces for model U1 125 Chapter 5. Case History 275 325 375 Actual failure surface LEM Surface One LEM Surface Two Yielded elements (U2) SRF = 1 14 Internal shearing of  failed mass 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters     . Velocity (U2) SRF = 1.14 Figure 5.22: Comparison of LEM surfaces and D-SSR SRF = 1.14 yielded ele- ments for model U2 126 Chapter 5. Case History 275 325 375 Actual failure surface LEM Surface One LEM Surface Two UDEC‐SSR 1E‐5 (U2) (velocity) UDEC‐SSR 1E‐6 (U2) 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters     (velocity) Figure 5.23: Comparison of LEM surfaces and UDEC-SRF velocity-derived sur- faces for model U2 127 Chapter 5. Case History 5.2.3 UDEC: Continuous Bedding Models Complexity was next added to the model in the form of continuous bedding in order to see how it affects both the value of factor of safety returned and the failure surface predicted. The models were compared to LEM Surface One and the actual projected failure surface in order to see if the models would predict a tension crack opening at the crest of the slope together with a slip surface that lies somewhere between the LEM solution and the actual failure surface. From the above results, a 2.5 m element size in the Upper unit and Contact zone was based on previous experience and because a 2.5 m element size predicted a shallower failure than other values investigated. UDEC-SSR is run using a ratio limit of 10−6 for the same reasons. Four models are investigated, with properties defined in table 5.1: • C1: 5 m spaced continuous joints dipping at 15o into the pit • C2: 10 m spaced continuous joints dipping at 15o into the pit • C3: 10 m spaced continuous joints dipping at 20o into the pit • C4: 20 m spaced continuous joints dipping at 15o into the pit Table 5.5 summarizes the results of the above model runs. The critical SRF re- turned is 1.12–1.13 which is identical to the result for the 2.5 m element size unjointed model (U1). Failure surfaces are similar to that of the U1 model, that is, a match to the LEM Surface Two without a tension crack but not a match to the actual failure surface. The bedding, surprisingly, does not affect the failure surface, as models C1 – C4 return the same result. Insights into this are partly provided through figure 5.24, which shows the shear displacement (in red) on joints in the C4 model at the point of failure for UDEC-SSR, where joint displacement is occurring only at a bench scale in the model. The maximum value of shear displacement is 7.82 x 10−2 m; however the value is not as important as the location. There is very little movement along the other bedding joints aside from some isolated bench scale movements in the upper slope. The majority of the movement localizes and concentrates along the Contact zone. There is no large scale joint displacement anywhere else in the model even though joint strengths are lower than that of the intact rock blocks. This indicates that a more complicated joint network is needed in order to predict the failure correctly and that bedding is not the only factor playing a significant 128 Chapter 5. Case History role in the predicted failure mechanism. Table 5.5: Summary of continuous bedding SSR results Model Joint Dip towards D-SSR UDEC-SSR spacing (m) pit 10−6 C1 5 15o 1.12 1.13 C2 10 15o 1.13 1.13 C3 10 20o 1.13 1.13 C4 20 15o 1.11 – 1.12 1.12    UDEC (Version 4.00) LEGEND    29-Jul-08  17:24  cycle    99202  time    6.048E+01 sec block plot shear displacement on joint max shear disp =  7.825E-02 each line thick = 2.000E-03  2.400  2.800  3.200  3.600  4.000  4.400 (*10^2)  2.800  3.200  3.600  4.000  4.400  4.800 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.24: Plot showing shear displacements, in red, on bedding surfaces at the point of failure in UDEC-SSR for model C3 5.2.4 UDEC: Discontinuous Joint Models These sets of models represent the more complicated models developed for the case history. The cross-cutting joints are expected to play a significant role due to 129 Chapter 5. Case History the presence of a complex joint network. This is where the SSR method is tested for it’s viability –both in terms of FOS reported and slip surface returned. Several different joint networks were generated and tested to offset limitations in the site investigation data where details of the joint network (e.g. persistence) were not fully captured. Also, limitations in the modelling code UDEC make it impractical for a 1:1 representation of the joint frequency in the model coupled with the inability to form true non-persistent joints terminating within a block. Continuous and discontinuous bedding are investigated with one or two different discontinuous high angle joints using orientations taken from stereonets presented in section 5.1.3 and discussions with site engineers. The limitations in represent- ing the joint network in UDEC makes it possible for a near-infinite number of different joint networks to be created from the combination of persistence and gap lengths with the above two to three joint sets modelled. Because of this, a range of models with different spacings, joint persistence, and other variations of the joint network are run to see how changes in this joint network affect the FOS and slip surface returned. This would be how the analysis would be carried out in a forward analysis; however, as this is a back analysis there is a failure surface to compare to, allowing the verification of the SSR method through the prediction of the correct failure surface, which LEM analysis does not. Table 5.6 shows the input parameters used to create the different joint networks in this series of models. The joint network generated for each model can be found in the figures in Appendix C. Displacement-SSR results with some select UDEC- SSR results can be found in Table 5.7, which summarizes the results in figure 5.25. None of the results produce a FOS below 1.0; however, pore-pressures are not included at this stage. A narrow range of FOS (1.06 – 1.14) is reported even considering the significant difference in possible joint works which is between 4 to 9% of the LEM FOS. Velocity-derived failure surfaces from UDEC-SSR can be separated into three categories, as was done for the above table 5.7: good matches to the actual fail- ure surface (figure 5.26), matches to LEM Surface Two (figure 5.27), and deeper seated failures (figure 5.28). From these figures, and table 5.7, we see that the joint network has a large control on the factor of safety and the position of the critical slip surface. Some results (J1, J6, J7) provide a good match to LEM Surface Two with a factor of safety ranging from 1.10 - 1.14. Many models predict a deep seated failure surface with factors of safety ranging from 1.08 - 1.14. Only model J8 predicts 130 Chapter 5. Case History Table 5.6: Joint generation parameters for jointed case history model Model Joint dip Trace length (m) Gap length (m) Spacing (m) m std m std m std m std J1 -15 0 100 0 0 0 10 0 -80 7 30 11 12 15 10 5 J2 -15 0 100 0 0 0 10 0 -74 7 27 11 12 15 10 5 J3 -15 0 100 0 0 0 5 0 75 5 17 10 15 15 5 0 -80 7 23 11 12 15 5 0 J4 -20 0 120 0 0.5 1 15 1 75 7 30 5 15 10 12 3 -80 7 23 7 12 15 13 5 J5 -20 0 120 0 0.5 1 3 1 75 7 8 5 15 10 7 0 -80 7 8 7 12 15 5 0 J6 -15 0 100 0 0 0 10 0 80 5 25 10 15 15 10 4 -80 5 20 10 10 15 10 5 J7 -18 3 80 30 0.5 1 3 1 72 13 12 5 15 10 7 1 -71 11 8 7 12 15 5 0 J8 -15 0 100 0 0 0 7.5 0 75 7 35 10 15 15 10 4 -75 6 40 14 10 15 14 5 J9 -15 0 100 0 0 0 4 0 75 7 35 10 15 15 10 4 -75 6 28 14 5 15 17 5 J10 -16 0 60 20 0.5 0.5 3.5 1.5 70 5 25 5 5 15 15 4 -78 6 25 38 10 10 14.5 4 J11 -15 0 100 0 0 0 7.5 0 75 5 35 10 15 15 10 5 -75 6 30 10 10 6 14 5 J12 -15 0 100 0 0 0 2.75 0.5 75 5 10 5 5 0 10 0 -77 7 4 5 4 7 4 5 -80 3 6 3 3 5 5 3 J13 -20 0 120 0 0.5 1 3 1 75 7 8 5 15 10 7 0 -80 6 8 9 6 15 4.5 1.2 131 Chapter 5. Case History Table 5.7: Summary of discontinuous SSR results with a wide range of joint net- works, revealing small range in factor of safety values Model D-SSR UDEC-SSR Failure surface location J1 1.12 1.13 Matches LEM Surface Two J2 1.14 1.14 Deeper seated failure J3 1.09 1.09 Deeper seated failure J4 1.12 1.11 Deeper seated failure J5 1.13 1.13 Deeper seated failure J6 1.10 1.10 Matches LEM Surface Two J7 1.11 1.10 Matches LEM Surface Two J8 1.06 1.06 Good match to actual failure surface J9 1.07 1.08 Deeper seated failure J10 1.07 1.08 Deeper seated failure J11 1.08 1.07 Fair match to actual failure surface J12 1.10 1.10 Matches LEM Surface Two J13 1.10 1.11 Deeper seated failure 132 Chapter 5. Case History 5.00 6.00 7.00 8.00 9.00 n t a l   d i s p l a c e m e n t   ( m ) J1 J2 J3 J4 J5 J6 Large joint influence on  stability Moderate joint influence on  stability 0.00 1.00 2.00 3.00 4.00 1.00 1.05 1.10 1.15 1.20 1.25 M o d e l l e d   h o r i z o Strength reduction factor J7 J8 J9 J10 J11 J12 J13Small joint influence on  stability   Figure 5.25: Results of D-SSR for jointed case history models showing range of critical SRF values from 1.06 – 1.15 133 Chapter 5. Case History 275 325 375 J8 SRF = 1.06 J11 SRF = 1.07 LEM Surface Two LEM Surface One Actual failure surface 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters Figure 5.26: Velocity surfaces representing critical sliding surfaces, from UDEC- SSR, that provide a good match to the actual failure surface 134 Chapter 5. Case History 275 325 375 J1 SRF = 1.12 J7 SRF = 1.11 J12 SRF = 1.10 J6 SRF = 1.10 LEM Surface Two LEM Surface One Actual failure surface 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters Figure 5.27: Velocity surfaces representing critical sliding surfaces, from UDEC- SSR, that provide a good match to LEM Surface Two 135 Chapter 5. Case History 275 325 375 J2 SRF = 1.14 J3 SRF = 1.09 J4 SRF = 1.12 J5 SRF = 1.13 J9 SRF = 1.07 LEM Surface Two LEM Surface One Actual failure surface 125 175 225 150 200 250 300 350 400 450 500 550 J10 SRF = 1.07 J13 SRF = 1.10m e t e r s meters Figure 5.28: Velocity surfaces representing critical sliding surfaces, from UDEC- SSR, that do not provide a good match and predict a deeper seated failure 136 Chapter 5. Case History a failure surface that is close to the actual failure surface; however it requires continuous steeply-dipping structures placed directly behind the crest of the pit in order for this to occur. Other models, such as J11, lie somewhere between the LEM Surface Two and the actual failure surface. None of the models tested match the LEM Surface One. Even the most heavily jointed slope (J12) prefers to break back to a deeper failure surface. These results show that UDEC predicts a much deeper final failure surface than that in the actual slope. Figure 5.29 shows the yielded elements at the point just prior to failure in model J4 superimposed with the actual failure surface. At this timestep, the model indicates that there is a shear band that closely matches the projected failure surface. This model, however, eventually develops a deeper ulti- mate slip surface (seen by velocity vectors in figure 5.28 with the surface shown in figure 5.29 instead forming an intermediate slip surface representing internal shearing which never fully develops as a failure surface. Material inputs and other behaviours can be varied, however, within a sensible range to attempt to influence the position of the final failure surface to make it shallower. This includes us- ing weaker material properties in the Upper unit and moving to a strain softening model, in addition to including the effects of water.    UDEC (Version 4.00) LEGEND     7-Aug-08  13:03  cycle   113880  time    7.121E+01 sec block plot no. zones : total     13530 at yield surface (*)    925 yielded in past  (X)      0 tensile failure  (o)    140  1.750  2.250  2.750  3.250  3.750  4.250 (*10^2)  2.500  3.000  3.500  4.000  4.500  5.000 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.29: Yielded elements showing internal shearing of the failing mass in the J4 joint network for a SRF stage just prior to catastrophic failure 137 Chapter 5. Case History 5.2.4.1 Strain-Softening Behaviour Select discontinuous models were run using a strain softening model to attempt to account for the generation and propagation of new fractures within the rockmass due to tensile failure behind the crest of the slope. Residual rockmass tension is set to zero to mimic the behaviour of a large discontinuity when tensile failure occurs within a block. Models J7, J10, J11, J12 representing a range of responses from the different joint networks tested, were modified to include the following parameters defining strain-softening in the upper unit with joints using a simple elastic perfectly-plastic constitutive model: 1. Strain limit for full residual strength: 1% strain 2. Relationship between strain and residual parameters: Linear 3. Residual friction: φres = φpeak 4. Residual cohesion: 25% reduction (to 141 kPa) 5. Residual tension: 0 kPa Both peak properties and residual properties are reduced in the SSR routine. UDEC-SSR does not allow for SSR of strain softening constitutive models and is thus not included in this analysis. Results of D-SSR for models J7, J10, J11, J12 can be seen in figure 5.30. The factor of safety returned is lower than running an elastic perfectly-plastic model, but not significantly. It was hypothesized that the strain-softening model would allow for some of the internal shearing occurring during failure to initiate softening of the material. This would then preferentially develop into a shallow failure surface instead of the deeper seated failures predicted using a perfectly plastic materially. However, none of the failure surfaces predicted in UDEC match the projected failure sur- face or the behaviour of the model as closely as that with a steep persistent joint included behind the crest. The position of the failure surface can be modified slightly, as is the case in the strain-softening model using the J12 joint network, as evidenced in figure 5.31 comparing the strain softening and perfectly plastic failure surfaces, but again not significantly. This J12 joint network model is used for the following experiments on the position of the failure surface predicted using D-SSR, to attempt to mimic the effects of a continuous joint propagating behind the crest through a continuum (i.e. the intact 138 Chapter 5. Case History 8 00 7.00 . 6.00 e m e n t   ( m ) 4.00 5.00 n t a l   d i s p l a c e J7 strain softening SRF = 1.09 3.00 d e l l e d   h o r i z o J10 strain softening SRF = 1.05 J11 strain softening SRF = 1.10 J12 strain softening SRF = 1.08 1.00 2.00M o d 0.00 1.00 1.02 1.04 1.06 1.08 1.10 1.12 1.14 Strength reduction factor Figure 5.30: Results of displacement-based SSR for jointed strain-softening case history models 139 Chapter 5. Case History 275 325 375 J12 strain softening [SRF = 1.08] J12 perfectly plastic [SRF = 1.10] Actual failure surface 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters Figure 5.31: Velocity surfaces representing critical sliding surfaces for the J12 joint network for elastic perfectly-plastic behaviour and strain softening be- haviour, at the point of instability as defined by D-SSR 140 Chapter 5. Case History blocks): 1. R1 – Base case presented above with 1% strain limit for full residual pa- rameters listed above (tension = 0) 2. R2 – 0.5% strain limit for full residual parameters with residual cohesion equal to joint cohesion (48 kPa) with zero tensile strength 3. R3 – 1% strain limit for full residual parameters with residual cohesion equal to joint cohesion (48 kPa) with zero tensile strength Figure 5.32 shows the results of the D-SSR method, and figure 5.33 shows velocity- derived failure surfaces for strain-softening R1–R3 cases. Factors of safety range from 1.05–1.08, as seen in figure 5.32, with no clear link between strain limit for full residual parameters and effect on factor of safety. The range returned, how- ever, is small. The more brittle behaviour, R2 with residual parameters occurring at lower values of strain, return a higher factor of safety than the less brittle be- haviour, R3 with a 1% strain limit. This could be due to slight variations in the stress path and deformation of intact blocks that occurs in one model, that does not occur in the other. From figure 5.33 we see shallower failures with lower residual cohesion properties; however, the difference is slight. Regardless, none of the strain-softening models used to incorporate the effect of fracture generation through intact rock are able to reproduce a match to the failure surface mapped in the field, as seen in the J8 model presented earlier. 141 Chapter 5. Case History 5 00 4.50 . 3.50 4.00 e m e n t   ( m ) 2.50 3.00 n t a l   d i s p l a c e R1 SRF = 1.08 1.50 2.00 d e l l e d   h o r i z o R2 SRF = 1.07 R3 SRF = 1.05 0 50 1.00 M o d 0.00 . 1.00 1.02 1.04 1.06 1.08 1.10 1.12 1.14 Strength reduction factor Figure 5.32: Results of displacement-based SSR for models jointed strain soft- ening models R1-R3 with low residual properties to capture fracture generation within a continuum 142 Chapter 5. Case History 275 325 375 R1 SRF = 1.08 R2 SRF = 1.07 R3 SRF = 1.05 Actual failure surface 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters Figure 5.33: Velocity surfaces representing critical sliding surfaces for R1–R3 strain-softening models 143 Chapter 5. Case History 5.2.4.2 Material Property Variations In this section initial parameters are varied slightly to see the effect on the pre- dicted failure surface. Because all factor of safety calculations depend on initial values to compare against, modifications to the initial parameters change the fac- tor of safety, regardless of the method (SSR or LEM) used. This set of models investigates more the effect of varied values on the position of the failure surface rather than the ultimate value of factor of safety due to the modification of the starting parameters. Three models M1, M2, M3 were created that had initial strength parameters 10% lower than the base case for either the Upper unit, joint strengths, or Contact zone respectively. These three models shared the same joint network as model J12, which consists of a very dense discontinuous fracture network without a large persistence steeply dipping structure behind the crest. Figure 5.34 shows the velocity-derived failure surfaces at the point of instability. The failure sur- face locations are relatively insensitive to the upper unit and joint strengths, but deeper seated failures are predicted when the contact zone strength is reduced (in the absence of a continuous steeply dipping structure). The factors of safety, as expected, are lower with the initial parameters reduced, with the joint strengths (M2) having the smallest influence on the factor of safety. None of the models were able to capture the actual failure surface indicating that a large controlling structure must be present in order to obtain the actual failure surface. The model that most closely matches the failure surface, J8 with very persistent back dipping structures, was investigated for sensitivity to joint shear stiffness (ks) and initial joint strength properties. Joint shear stiffness models were varied to 8 x 108, 8 x 109, 8 x 1010 Pa/m respectively, with D-SSR results seen in figure 5.35. This figure shows that the joint shear stiffness has a mild impact on the strength reduction factor at failure, despite greater differences in absolute values of displacement at each SRF stage. The failure surfaces for the three ks cases tested were investigated and are identical, indicating that this failure mechanism would have been the same for any reasonable value of kstaken from the literature such as Kulhawy (1975). Joint strength properties were varied for two cases, also for the J8 fracture net- work: • Case one: φjoint = 27, cjoint = 53 kPa 144 Chapter 5. Case History 275 325 375 J12 M1 M2 M3 Base case FS = 1.10 10% lower initial upper unit parameters FS = 1.05 10% lower initial joint parameters FS = 1.08 10% lower initial contact zone parameters FS = 1.04 Failure surface position  insensitive to joint and upper  unit strengths Actual failure surface 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters Failure surface position sensitive to  contact zone strengths Figure 5.34: Failure surface and factor of safety sensitivities with respect to initial parameters before strength reduction 145 Chapter 5. Case History • Case two: φjoint = 27, cjoint = 80 kPa The base case returns an SRF at failure of 1.06 while cases one and two return values of 1.07 and 1.08 respectively. This is not a large increase in factor of safety despite the relatively large increases in friction and cohesion, indicating that joint strengths alone are not the most significant properties in the model. The velocity-derived failure surfaces, seen in figure 5.36, indicate that the joint strengths can play a significant role in the shape and location of the failure sur- face. If the strengths of the discontinuities are high, but still lower than those of the upper layer, the failure does not propagate upwards through the large steeply dipping structure. Instead a more general, deeper seated rockmass failure is pre- dicted. 1.50 2.00 2.50 3.00 n t a l   d i s p l a c e m e n t   ( m ) 8E8 (Pa) All models fail at the same  critical SRF stage despite  differences in absolute values of  0.00 0.50 1.00 1.00 1.02 1.04 1.06 1.08 1.10 1.12 M o d e l l e d   h o r i z o Strength reduction factor   8E9 (Pa) 8E10 (Pa) displacement Figure 5.35: Results of displacement-based SSR for models ks sensitivity for joint network model J8 146 Chapter 5. Case History 275 325 375 J8 Base SRF = 1.06 J8 Case one SRF = 1.07 J8 Case two SRF = 1.08 Actual failure surface 125 175 225 150 200 250 300 350 400 450 500 550 m e t e r s meters Figure 5.36: Variations in failure surface location with initial joint strength pa- rameters for joint network model J8 147 Chapter 5. Case History 5.2.4.3 Water Influence Water played a significant role in the actual failure with 115 mm of rainfall having occurred over the eleven days prior to failure and an estimated 400 l/min of water entering tension cracks behind the crest at one point. Unfortunately, there was no pore pressure monitoring in this section of the pit during this time, making the water table analysis more an investigation of potential effects rather than one that examines the absolute effect on the factor of safety. Also, simulating the direct effects of precipitation runoff entering the tension cracks at the surface in the model is beyond the scope of this thesis. A simple static water table was implemented in UDEC, after discussion with Piteau Associates regarding site conditions and anecdotal evidence (standing wa- ter in blast holes, etc.), and can be seen outlined in blue in figure 5.37. Model J12, with finely spaced bedding and cross-joints, was selected along with the model J8, which best predicted the actual failure to be run with a static water table and are labeled W12 and W8 respectively. D-SSR results can be seen in figure 5.38 showing a large amount of displacement in the first SRF stage. This, coupled with the state of the model after the removal of the last bench, indicate that the factor of safety with these water conditions is at, or just less than 1.0. Evidence of this can be seen in the plasticity plots after removal of the last bench for the W12 model (figure 5.39), which shows a clear ultimate failure surface defined by sheared elements, and in the velocity plots, which shows movement along the critical discontinuity behind the crest in W8 at the same SRF stage in the model before strength reduction (figure 5.40). Updated LEM results using this static water table in SLIDE return a factor of safety of 1.04 and 1.09 for GLE and Janbu Corrected respectively. Both models, with a tension crack explicitly defined, follow LEM Surface Two. Trial surfaces of the actual failure zone return results of 1.26 and 1.47 for the same two methods. This indicates that the J8/W8 model is a failure that is not properly predicted in the LEM analysis. 148 Chapter 5. Case History Figure 5.37: Estimated static water table location, in blue, for case history 149 Chapter 5. Case History 4.00 5.00 6.00 7.00 8.00 n t a l   d i s p l a c e m e n t   ( m ) Water table models W8 and  W12 have factors of safety  at, or below, 1.0  0.00 1.00 2.00 3.00 1.00 1.02 1.04 1.06 1.08 1.10 1.12 1.14 M o d e l l e d   h o r i z o Strength reduction factor W8 W12 Figure 5.38: Results of displacement-based SSR for models W8 and W12 with static water tables   UDEC (Version 4.00) LEGEND    22-Aug-08  10:08  cycle    81067  time    3.901E+01 sec block plot no. zones : total     17814 at yield surface (*)   1643 yielded in past  (X)      0 tensile failure  (o)    464  2.000  2.500  3.000  3.500  4.000 (*10^2)  2.500  3.000  3.500  4.000  4.500 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.39: Plasticity indicators showing clear shear band after removal of last bench in W12 model 150 Chapter 5. Case History    UDEC (Version 4.00) LEGEND    19-Aug-08  13:57  cycle   113190  time    6.715E+01 sec block plot velocity vectors      maximum =    8.090E-02 0   2E -1  2.000  2.500  3.000  3.500  4.000 (*10^2)  2.500  3.000  3.500  4.000  4.500  5.000 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada Figure 5.40: Velocity plot showing slope failure along critical discontinuity be- hind the crest after removal of last bench in W8 model 151 Chapter 5. Case History 5.3 Continuous Cycling Solution Mode It is possible modify the D-SSR routine to continuously drop strengths from SRF = 1 rather than returning to an initial state between each SRF stage. Running the routine in this method means that the response modelled in each stage is directly influenced by the stages before it. Yielding and yielded elements are carried for- ward between each SRF stage. This code is included as an option in the D-SSR UDEC code developed (found in Appendix A). It is important that the SRF code is always run from an SRF stage of 1, even if the ultimate factor of safety of the slope is greater than 1, in order to capture the path dependency of the model or at least to keep it the same as much as possible between models. The benefit of running the routine this way is computational efficiency as each stage introduces only small changes to the strengths of the materials which allows the force ratio limits during cycling to be satisfied quicker. The models in table 5.8 were investigated in order to assess whether or not the path dependency of the continuous cycling displacement-based SSR (CD-SSR) routine affects the failure surface and factor of safety returned, with results summarized in the same table from D-SSR output in figure 5.41 Table 5.8: Summary of continuous displacement-based SSR for jointed model cases Model Joint Network D-SSR CD-SSR CD1 J1 1.12 1.12 CD3 J3 1.09 1.09 CD5 J5 1.13 1.13 CD8 J8 1.06 1.06 From this analysis, it is seen that the continuously cycled mode gives the same result as the method involving an initial state between SRF stages, despite path dependancy issues. Velocity surfaces and plots of yielded elements return the same failure surface between both methods as well. This is a useful result as the computational efficiency gained in terms of solution runtimes can be significant 152 Chapter 5. Case History over the normal D-SSR method. 4.00 5.00 6.00 7.00 8.00 n t a l   d i s p l a c e m e n t   ( m ) J1 CD1 J3 0.00 1.00 2.00 3.00 1.00 1.05 1.10 1.15 1.20 1.25 M o d e l l e d   h o r i z o Strength reduction factor CD3 J5 CD5 J8 CD8 Figure 5.41: Results of displacement-based SSR for models CD1–CD8 showing identical behaviour between continuously-cycled and initial state displacement- based SSR 5.4 Discussion The above series of analyses, incorporating many different models and parame- ter variations, show that the SSR method returns a range of safety factors from 1.06 to 1.14 with several possible failure surfaces depending on the specifics of the model parameters, joint network, and constitutive model used. As the anal- ysis moves away from a traditional LEM analysis, the amount of model inputs increases requiring additional information from a site investigation. The materials being modelled are natural geomaterials with known and unknown variations in properties like material strengths, joint distributions and connectivity, etc. Mod- elling assumptions required in response to limitations in the constitutive models used, inability to generate joints terminating within a block, and inability to model 153 Chapter 5. Case History fracture generation and propagation in the model results, in a large number of cases that would be required (prudent) for a forward analysis in order to ensure capturing a representative failure mode for the slope For this reason, the design factor of safety is never taken to be exactly 1.0 due to variations and unknowns within the analysis. The back analysis presented here can be equated to a forward analysis of the slope in question as it is based on known rockmass strength properties (albeit, without an analysis to quantify the uncertainty of initial parameter choice) together with that information pertaining to the joint network in the slope that was available prior to failure. While only one model, which includes a large persistent steeply dipping structure at the crest of the slope, predicts failure, the other models do reveal that the slope in question has the potential for deep seated failure given the low factors of safety returned (1.07 – 1.14). The result of this modelling shows that there are two distinct, related failure mechanisms occurring within the slope. The first mechanism is the large deep seated failure caused by the weak carbona- ceous black silt layer between the upper and lower units in the slope. This zone provides a mechanism of crushing at the toe and allows for significant sliding of the upper mass along it, particularly as the dip of this zone increases towards the pit at depth. Regardless of the specifics of the joint network, material properties, water table, or constitutive strength model of the upper slope, this failure propa- gates along the weaker unit before some limiting sliding mass volume is reached which allows localization and internal shearing through the Upper unit to surface. The weak contact unit represents a known geological problem identified during the site investigation. The second mechanism pertains to the J8 model with a large continuous struc- ture inserted directly behind the crest. The failure mechanism still relies on the movement of the upper mass along the weak contact unit facilitated by crushing in the toe; however, the structure controls the breakout of the failure to the sur- face and is responsible for a lower factor of safety. Figure 5.42, is a photograph showing the large continuous fault which acted as a release surface in the case history failure. This suggests that the J8 model is actually representative of the failure mechanism that likely occurred with additional kinematic control on the failure being provided by the structure behind the crest. This photograph (and thus, information) was not included in the previous analysis so as not to bias the results of the forward analysis mindset. It should be noted that this feature was absent in the structural data collected before the failure. Determining the persis- 154 Chapter 5. Case History tence of features, at depth using borehole investigation, can be difficult if there are no clear defining characteristics for the discontinuity. The feature does not ap- pear on surface behind the crest, nor does it appear in bench face mapping. This represents an unknown in the forward analysis, and was missed entirely except for one model (J8) which likely would not have been considered in the absence of additional data identifying the presence of this feature. Regardless, a large failure would have been predicted for this slope given the range of FOS from 1.07 to 1.14 determined for the various dry scenarios, even when the major structure behind the crest is omitted leaving only the weak contact zone for kinematic control on the failure. Figure 5.42: Photograph roughly perpendicular to failure slope showing presence of a steep, continuous structure The key advantage of carrying out more complex discontinuum modelling in con- junction with the SSR method is that this failure had the potential to be predicted. The failure mechanism is similar to the deeper seated failure, except that there is a more convenient release surface forming the back scarp that requires only a small volume of rock in the upper unit to shear. This formed a weaker path than that predicted in most of the models that involved more extensive internal shearing. It is possible that this case could have been run in a typical forward analysis, such as the one presented in this chapter. Attempts using strain softening behaviour to simulate the fracturing of intact rock and joint propagation by setting residual 155 Chapter 5. Case History cohesion and tension properties in the continuum portion of the models does not fully capture this failure mechanism. At the same time, no physically modelled water tables in LEM analysis, even with a defined tension crack, would have pre- dicted the actual failure surface. Only the curious condition of Ru = 0.15 with the Janbu Corrected method defines this as a low factor of safety surface, which could be seen as an unreliable answer within a typical analysis. This highlights the fact that only a discontinuum model would have fully realized this failure, and that the SSR method provides a suitable means to quantify the factor of safety and potential volume of failure, for similar discontinuous rock slope problems, pro- vided well defined joint geometry information from a site investigation is avail- able. 156 Chapter 6 Discussion and Conclusions 6.1 Framework Geological engineering increasingly depends on risk-based assessments including both the consequence of a slope failure and the probability of the failure, which prompts the use of a variety of tools in order to achieve the overall goal. These range from simple kinematic stereonet investigations, to deterministic limit equi- librium methods, numerical continuum and discontinuum (e.g. distinct element) failure mechanism investigations and the SSR method, in addition to a probabilis- tic analysis involving the use of geostatistics and incorporation of GIS databases for input. There is no single all-encompassing method for slope stability; instead a range of analyses are undertaken and the relative merits of each are used. The SSR method is one such tool, to be used in combination with other investigative and analytical tools. A general framework, shown in figure 6.1, depicts a typical approach to risk as- sessment for a rock slope using a variety of techniques, and what insights are typ- ically gathered at each stage. Ultimately, probabilistic and limit equilibrium anal- yses form the basis for most of the design calculations with information gathered from more complex numerical modelling methods being used for failure mode investigations or estimates of volumes of unstable material. The SSR method en- ters this framework as a means to incorporate the advanced rockmass behaviours possible in finite-element/finite-difference continuum and discontinuum codes to create models with more realistic idealized conditions that can return a factor of safety, in addition to capturing complex failure mechanisms that are not possible using more basic methods. 157 Chapter 6. Discussion and Conclusions Data Collection -Geology - Groundwater -Geometry - Dynamic loads -Strength Kinematic Analysis - Discontinuity sets governing the slope stability Numerical Modelling Continuum - Potential failure mechanism - Factor of Safety (via SSR) - Depth of slide Discontinuum - Potential failure mechanism - Factor of Safety (via SSR) - Depth of slide Sensitivity Analysis - material parameters - modelling parameters - Failure surface(s) Uncertainty Analysis - Quantification of uncertainties in sensitive paramteters Probabilistic Modelling - Limit Equilibrium Analysis Evaluation of slide magnitude Vulnerability Analysis - Magnitude (M) - Scale (S) - Elements at Risk (E) Risk Analysis - Computation of risk Risk Evaluation - Assessment of acceptability/tolerability of risk Figure 6.1: Typical risk assessment framework for rock slopes, modified from Duzgun and Lacasse (2005) 158 Chapter 6. Discussion and Conclusions 6.2 Discussion This thesis focused on UDEC, a tool typically used to compare hypotheses of fail- ure mechanisms, material parameters, and joint networks, to assist in determining the factor of safety for a given rock slope. Chapter 4 involved a rigorous investi- gation of some of the major modelling parameters and material variation ranges and their effects on the calculated factor of safety. The major impacts on the SSR solution, assuming that the starting conditions were kept consistent were found to be related to model settings and complications arising from simple idealized behaviours of a complex natural system. It was found that the results were not that sensitive to variations in the rockmass parameters within the range of preci- sion typically applied in practice. At first glance, the results of Chapter 4 show that as long as there is a reasonable amount of field mapping and lab testing data available, then the SSR method can be confidently applied to find potential critical shearing surfaces and associated factors of safety, as the input parameters do not affect the answer significantly in these simpler models. Still, the difficulty in capturing the proper failure surface in the case history pre- sented in Chapter 5 highlights that the method requires a minimum threshold of information in order for it to become useful as a stability analysis tool. The major- ity of these models presented in sections 5.2.2–5.2.4 result in a surface predicted that is no closer to the actual failure surface than those from LEM solutions, and took significantly longer to obtain. This set of results, over quite a diverse range of different possible joint networks, indicates that the jointed rockmass can be- have as a continuum over a larger scale when a deeper seated failure mechanism is involved. The advantage of a numerical discontinuum analysis, such as that pre- sented, is that it can capture failure mechanisms and failure propagation that other methods cannot through explicit inclusion of the discontinuity network. These joints enable the extra degrees of freedom present when the failure evolves, giv- ing the model ability to capture the correct failure surface. If this network is not well defined for the modeller, then the Shear Strength Reduction method is not as reliable a tool for stability analysis as critical discontinuities may not be accounted for which invalidates the key benefit of the discontinuum analysis. The back analysis reveals that there must be a steep persistent structure behind the crest in order for the predicted surface to match the actual failure surface. A blind forward analysis, done using the discontinuity data available beforehand, that includes many sampled structural orientations but limited geometric data such 159 Chapter 6. Discussion and Conclusions as persistence, spacing and statistical variations, does not provide an exact match to the actual failure surface or additional insights over quicker, simpler methods. However, the analysis would have provided a prediction of failure, although more conservative in terms of the volume of the slope mobilized, and together with the LEM results would have provided added confidence in the conclusion that the slope was critically unstable. Based on these results, it was found that the Shear Strength Reduction method using discontinuities is applicable when: 1. A factor of safety value (or range) is important to the design of the project. For example to meet some design safety criteria, or use in some risk assess- ment framework, or to quantitatively compare design options. 2. The rockmass is expected to behave as a complex interaction between joints and intact rock rather than dominated by either one. 3. Some amount of lab testing is available, combined with field observations and GSI/RMR parameters, in order to give reasonable initial rockmass and joint properties for the relevant constitutive models. 4. Good discontinuity data exists beyond simple orientations. Fully being able to describe the continuity and spacing of the fracture network is integral to the SSR method producing meaningful results. Examining these, there must be a reason why an analysis is being carried out, which is covered with point one above. The Shear Strength Seduction method does not add any additional information over a traditional discontinuum analysis if only potential failure mechanisms are being investigated. The factor of safety value must be of interest, and thus the SSR method is not a tool to be used in an initial analysis with little data, outside of roughly scoping the problem. Point two means that the SSR method is one of many tools, and that it is only ap- plicable at certain times in the overall framework. The expected behaviour of the slope (or potential slope) being investigated should be neither purely structurally controlled nor heavily fractured/disturbed, as other methods exist that are more efficient (and can possibly include probabilistic techniques), that will give similar answers. The SSR models would help provide verification of the simpler models as the results from Chapter 5 show that the SSR method can capture continuum type failure modes without the assumption of a predefined failure surface. With all factor of safety methods, the starting point is important (point three). 160 Chapter 6. Discussion and Conclusions The more constrained the information in terms of statistical variability, zones of different strengths, etc., the greater the confidence in the values of likely factors of safety returned. With the SSR method, due to model sensitivities related to the constitutive models and other modelling properties (e.g. element size), a range of possible safety factors is returned. The better the material properties are defined, the smaller and more meaningful this range will be. Using GSI or other criteria to define an equivalent continuum is a common practice; however care should be taken in choosing assumed values of confining stresses and disturbance values as these can have a large influence on the starting parameters. Point four highlights the fact that the discontinuity network is essential in using the distinct element based SSR method. As previously noted, most other material factors do not play a significant role on the predicted failure surface in the models evaluated in this thesis. Material properties and initial conditions affect mainly the absolute value of the factor of safety. The result returned for the particular case history investigated was, at least in part, that of a continuum, with a reduced factor of safety compared to LEM results due the presence of in the SSR analysis of lower strength joints within the rockmass. To capture the actual failure surface, good discontinuity information is needed with particular respect to the continuity of the joints. In the analysis presented, the forward analysis was unable to capture the failure; only with the incorporation of additional structural data discovered post-failure was it possible for the correct failure surface to be evaluated using the SSR method. Transforming the discontinuity network information into UDEC, given the programs inability to terminate a joint within a block, can result in some difficulties in capturing the joint network completely within the model. However, with care a reasonable joint model of the actual network can be made. There has been no defined relationship between the modelled joint density and the equivalent continuum strength reduction via GSI/Hoek-Brown; however, as this analysis has shown only critically located discontinuities as dominating the modelled response. This relationship may not be as important with respect to the prediction of the failure surface. More advanced constitutive model involving strain softening should be used with caution. Lab data and field observations are needed in order to constrain the resid- ual properties and define how they change with respect to strain within the rock- mass, which was not present in the data used for the case history analysis. There are some conceptual problems with using a strain softening material in order to find the limit of instability after a material has softened and failed in a traditional sense. Technically, the limiting equilibrium is still being found because both the 161 Chapter 6. Discussion and Conclusions peak and residual parameters are being reduced, but the validity or meaning of the factor of safety returned is debatable. Using softened parameters, particularly cohesion and tension, in order to capture the behaviour of fracture generation through a continuum (i.e. the intact block) was attempted, rather than a true anal- ysis using softened materials to determine the absolute value of factor of safety. The factor of safety returned was similar to that of an elastic-perfectly plastic material, however the method failed to capture the expected behaviour of joint propagation behind the crest in a continuum. The method is applicable to both natural and engineered slopes. In order to be used for a natural slope, there must be sufficient discontinuity data available from outcrop exposures on a sufficient scale to provide insights into the features con- trolling stability. Similarly, in the early stages of a site investigation for an engineered slope (e.g. a preliminary pit design), orientated borehole data will not return the requisite in- formation on persistence and spacing of key discontinuities. In a modelled system with the limitation of not being able to model every discontinuity, the inclusion of major features is of prime importance. Correlating multi-bench scale features using limited borehole data is not possible. Thus, the SSR method is a specific tool that can be useful for certain site con- ditions and phases in design, including ongoing operation (or safety assurance) of the slope in question. Confidence in the results, for practical applications, re- quires mine experience to be gained in applying the method to design cases and in back analysis situations. However, it has been shown that it needs to be used correctly otherwise the method can give faulty results, resulting in reduced confi- dence in the method within industry. Further validation of the SSR method could benefit from additional treatments using monitoring data to calibrate and verify the models –something that was not done in this thesis due to time constraints and data availability. Situations that meet the above four criteria, but that also include slope displacement or other monitoring data is where the SSR method can really add value to the engineering process. If the model is calibrated, then future design decisions, changes, and new instabilities can be modelled and the factor of safety can be returned, with confidence, in order to compare the merits of particular strategies or mitigation methods. Thus, the SSR method is not only useful for a preliminary analysis, but may be useful over the operational life of the slope involved to find the factor of safety of various design and operational decisions provided the model is calibrated to monitoring data and updated with 162 Chapter 6. Discussion and Conclusions new material constraints and discontinuity information as it becomes available. Also, the method would be useful in cases where the failure mechanism changes over time as only discontinuum models can capture the shift from block yielding to movement along joints. A prime example of where this method could be used as a predictive tool for fail- ure is with the ongoing operation of an open pit. In this case, the model is updated as the pit is constructed. Actual displacement monitoring points are modelled, adopting key concepts from the Terzaghi and Peck (1948) observational approach to design. This is where other parameters that affect displacements but not the actual SRF method become important as these would also need to be tested and refined to better correlate modelled with the actual displacements. The joint net- work could be redefined and updated as better discontinuity information becomes available from bench face mapping and/or the incorporation of LIDAR (Light Intensity Detection and Ranging) and photogrammetry remote sensing methods, especially with respect to better defining multi-bench scale features that have an impact on overall slope stability. In essence, the model evolves with the pit instead of being looked at once and then discarded. The strength of the SSR method, in this situation, is that the pit design and operational stability strategies can be con- stantly updated, and problems not foreseen in the original design can be uncovered as new information becomes available. Again, this links to the use of numerical modelling to Terzaghi and Peck’s observational approach. Because the method is deterministic and feeds directly into a risk analysis via the factor of safety con- cept, decisions and options can be quantitatively evaluated. The monitoring data is essential for this in order to gain confidence in the results, as the forward model predictions using an idealized behaviour will be firmly rooted in a model that matches the controlling conditions of the slope. 6.3 Conclusion This thesis reports efforts to investigate the applicability of new developments in the application of Shear Strength Reduction (SSR) techniques for stability as- sessments of jointed rock slopes. For this, a wide variety of model and material parameters were investigated for several conceptual slope geometries followed by applying the method to a case history with more complex site parameters. There is the common misconception that because discontinuum models require 163 Chapter 6. Discussion and Conclusions more inputs than a typical LEM analysis, that the discontinuum model is unre- liable, or more sensitive to a larger number of parameters. The analysis in this thesis has shown that this is not always the case. Although it is true that there are more model and material parameters in these advanced models, it was found that the majority of them had little effect on the value of factor of safety when the initial shear strength parameters are kept constant before strength reduction. Sim- ilarly, the failure surface that is predicted through the generation of strains, block displacements, and both shear and tensile failure was insensitive to the majority of input parameters. Furthermore, using the displacement-based criteria to define the point of instability agrees with the numerical unbalanced force instability cri- terion for the majority of the models. This allows for the SSR method to distance the analysis from problems with numerical instability criterions, and has the po- tential to be used in finite-element methods which traditionally apply SSR using convergence criteria for defining failure that are seen to be imprecise by Krahn (2007). The rockmasses investigated in this work assume predominantly low to medium strengths. These findings may not hold true for brittle rock failure involving high strength rocks, particularly in light of work done by Martin and Chandler (1994) and Carter et al. (2007) showing that frictional strengths are not mobilized un- til significant cohesion loss occurs. In high strength rockmass cases, the joint network and joint properties may play a more significant role in the failure mech- anism than in the yielding of intact rock dominated shear failure mechanism ana- lyzed in this thesis. 6.3.1 Guidelines for Use The results from this thesis allow for several key recommendations to made for the practical application of the SSR method, involving limiting the range of factor of safety values returned as well as optimizing computational efficiency. Because of the element density, joint network assumptions, and constituitive models, the method will always result in a range of modelled responses. As with any slope sta- bility method, the better the input data, the narrower the range of results. The three most significant influences on the value of factor of safety, for the low strength rockmass cases tested, were found to be: • Element size – The smaller the element size, typically, the lower the FOS. The case history gave good surface matches to the LEM surfaces for simple 164 Chapter 6. Discussion and Conclusions unjointed models using an element size to slope height ratio of 0.02. Using larger ratios predicted deeper seated failures with higher factors of safety. Evidence from the case history indicates that models using element size to slope height ratios as small as 0.01 give good results however smaller ratio sizes may give anomalous results at greatly increased computation time. • Joint network (if discontinuous jointing is modelled) – A range of joint network realizations needs to be modelled in order to constrain the factor of safety within a possible range based on the mapped joint network properties. • Ratio or cycle limit for each stage – A lower ratio limit (or larger number of cycles) results in a lower factor of safety and a shallower failure. A ratio limit of 10−6 rather than the default 10−5 is recommended as the failure surface predicted in simple unjointed models are verified with LEM models. The other parameters are important, but these are the three that require the most at- tention. Element size has a controlling influence on the development of localized shear, which in turn, defines the rupture surface that the slope fails along. In failure mechanisms that involve yielding of the weak intact rock, this is the most impor- tant control on the value of FOS. Similarly, a discontinuous joint network imposes kinematic controls on those failure mechanisms involving significant movement along discontinuities. In the case history presented, the effect of jointing was mostly muted and only slightly affected the factor of safety. Mainly, the failure mechanism for this case involved localization of the shear band through intact blocks representing weak, disturbed rock. The joint network is expected to play a larger role in failures with more structural control, such as failure mechanisms involving flexural toppling. The distinct element SSR analysis may take a long time to complete, particularly if many different scenarios are being investigated. Thus, it is useful to detail the controls on computation time and advantages gained in order to optimize the time budgeted and to make the method more accessible to the engineering practitioner. There are five main controls on computation time: 1. Element size plays a large role; however, element size is required to be small due to the influence it has on shear band localization. This increases computation time. 2. The joint network can also be important; however, if the failure mecha- nism involves a large amount of intact rock yielding then a widely spaced joint network may be accepted. This would help reduce computing time but 165 Chapter 6. Discussion and Conclusions should include a limited number of model runs with closely spaced joints to ensure the assumption is valid. 3. The 10−6 cycle limit is recommended for the SSR method, both for D-SSR and UDEC-SSR, even though it requires more computational cycles (time). The work in this thesis has shown that is required to ensure that failure is not occurring for a specific SSR stage before erroneously moving on to the next stage as this was something that was seen to occur for the default ratio limit of 10−5. 4. Stiffer joints increase the computational speed due to their effect on the numerical timestep within UDEC. However, the shear stiffness of the joints can, but not necessarily always, have an affect on the failure mechanism, FOS returned, and position of the failure surface. Caution should be used here. 5. The continuously cycled displacement-based SSR (CD-SSR) method has shown, for limited cases with critical SRF stages near 1, to return the same answer as the D-SSR method with significant time gains. Also, UDEC-SSR is computationally efficient; however, it does not return as much information as the displacement-based methods. 6.3.2 Recommendations for Further Work Although a large variety of model scenarios were tested in this thesis covering most of the major parameters and other considerations a typical modeller would need to investigate, there is room for further studies into the applicability of SSR in distinct element modelled, particularly with respect to different rockmass strengths and different failure mechanisms. These include: 1. An SSR sensitivity analysis on stronger rockmasses with greater slope heights (i.e greater stresses) that involve brittle rockmass failure rather than the lower strength rockmasses analyzed in this thesis. 2. A better understanding of the relationship between the factor of safety re- turned by SSR and the element size. This thesis recommends a reasonably small element size, but not too small, and work needs to be done to con- strain how small is too small, possibly with respect to the minimum block size in the model. 166 Chapter 6. Discussion and Conclusions 3. The relationship between the ‘equivalent continuum’ strength reduction via Hoek-Brown and the modelled joint density, particularly as the block size gets small and approaches the properties of intact rock. 4. As with any method, additional case histories involving a wider scope of geological and site conditions (particularly differing rockmass strengths, failure mechanisms, structural influence, and stress conditions (controlled by the slope height)). The results obtained in this thesis are by no means conclusive and focus on a small range of low to medium strength rockmasses failing in shear and the limitations of the analysis should be accounted for in future works. Some important model op- tions, particularly in regards to choice of constitutive mode that take into account strength related to confinement, were left out of this research because either the application of SSR to more complicated constitutive models are not well-defined or due to time constraints. These involve: 1. Using the state of stress within the model to set confinement-based intact rock properties using a Hoek-Brown failure criterion (currently not avail- able within the UDEC model, but could be implemented using the internal programming language FISH). 2. Influence of dilation on the SSR method. 3. Covariant parameters (e.g. joint shear stiffness with the joint cohesion and friction angle) and their inclusion in the strength reduction process 4. Using the SSR method with a non-linear failure criterion such as Hoek- Brown, for intact rock representation in the model. 5. Getting around the limitations in the way UDEC handles joints (cannot ter- minate within a block) by creating fully continuos joints and modifying the select joint contacts directly with increased properties (φ, cohesion, tension) to simulate rock bridges, rather than physically defining the rock bridges as was done in this thesis. This would require a drastic change to the SSR code presented as each contact would have to be modified individually. Other nu- merical options such as hybrid finite-element / discrete particle codes (e.g. ELFEN), do allow for true non-persistent joints and joint propagation and SSR should be investigated in these more advanced, but uncommon mod- elling tools. 6. Application of the displacement-based stability criteria to finite-element 167 Chapter 6. Discussion and Conclusions based methods in order to compare to the non-convergence criteria. This could later be used for comparison purposes between continuum and dis- continuum SSR models. 168 Bibliography Barton, N. R. (1973). Review of a new shear strength criterion for rock joints. Engineering Geology, 7:287–322. Bawden, W. (1993). The use of rock mechanics in Canadian underground hard rock mine design, volume 2 of Comprehensive Rock Engineering, pages 247– 289. Oxford-Pergamon. Bieniawski, Z. (1976). Rock mass classifiction in rock engineering. Exploration for rock engineering, 1:97–106. Bieniawski, Z. (1989). Engineering rock mass classification. Wiley, New York. Brady, B. and Brown, E. (1993). Rock Mechanics for Underground Mining. Chapman and Hall, London, second edition. Cai, M., Kaiser, P., Uno, H., Tasaka, Y., and Minami, M. (2004). Estimation of rock mass strength an deformation modulus of jointed hard rock masses using the gsi system. International Journal of Rock Mechanics and Mining Sciences, 41(1):3–19. Carter, T., Diederichs, M., and Cavalho, J. (2007). A unified procedure for hoek- brown prediction of strength and post yield behaviour for rockmasses at the ex- treme ends of the rock competency scale. In Luis, S., Olalla, C., and Grossman, N., editors, 11th Congress of the international society for rock mechanics, vol- ume 1, pages 161–164, Lisbon. Taylor & Francis. Carvalo, J. L., Carter, T., and Diederichs, M. (2007). An approach for prediction of strength and post yield behaviour for rock masses of low intact strength. In Eberhardt, E., Stead, D., and Morrison, T., editors, 1ST Canada-US Rock Me- chanics Symposium, volume 1, pages 277 – 285, Vancouver. 169 Bibliography Coulomb, C. A. (1773). Sur une application des regles de maximus et minimis a quelques probleme de statique relatifs a l’architecture. Acad. Roy. des Sciences Memoires de math. et de physique par diversavans, 7, 343-82. Dawson, E., Roth, W., and Drescher, A. (1999). Slope stability analysis by strength reduction. Geotechnique, 49:835–840. Diederichs, M., Lato, M., Hammah, R., and Quinn, P. (2007). Shear strength reduction for slope stability. In Eberhardt, Stead, and Morrison, editors, 1ST Canada-US Rock Mechanics Symposium, volume 1, pages 319 – 325, Vancouver. Duzgun, H. S. B. and Lacasse, S. (2005). Vulnerability and acceptable risk in intergrated risk assessment framework. In Hungr, O., Fell, R., Couture, R., and Eberhardt, E., editors, Landslide Risk Management, pages 505–515, Vancouver. Taylor & Francis. Elishakoff, I. (2004). Safety factors and reliability: friends or foes? Springer, 1st edition. Franklin, J. and Dusseault, M. (1989). Rock Engineering. McGraw-Hill, New York. Fredlund, D. and Krahn, J. (1977). Comparison of slope stability methods of analysis. Canadian Geotech Journal, 14:429–439. Griffiths, D. and Lane, P. (1999). Slope stability analysis by finite elements. Geotechnique, 49(3):387–403. Hammah, R. and Yacoub, T. (2006). Investigating the performance of the shear strength reduction (ssr) method on the analysis of reinforced slopes. In 59th Canadian Geotechnical and 7th Joint IAH–CNC and CGS Groundwater Spe- ciality Conferences, Vancouver. Hammah, R., Yacoub, T., Corkum, B., and Curran, J. (2005a). A comparison of finite-element slope stability analysis with cnventional limit-equilibrium in- vestigation. In 58th Canadian Geotechnical and 6th Joint IAH–CNC and CGS Groundwater Speciality Conferences, Saskatoon, Canada. GeoSask. Hammah, R., Yacoub, T., Corkum, B., and Curran, J. (2005b). The shear strength reduction method for the generalized hoek-brown criterion. In 40th U.S. Sym- posium on Rock Mechanics (USRMS), Anchorage, Alaska. American Rock Me- chanics Association. 170 Bibliography Hart, R. (1993). An introduction to distinct element modeling for rock engineer- ing. In Hutchinson, J., editor, Comprehensive Rock Engineering: Principles, Practice, and Projects, volume 2, pages 245–261. Pergamon Press, Oxford. Havenith, H., Strom, A., Calvetti, F., and Jongmans, D. (2003). Seismic trig- gering of landslides. part b: Simulation of dynamic failure processes. Natural Hazards and Earth System Sciences, 3:663–682. Hoek, E. (2007). Practical rock engineering. http://www.rocscience. com/hoek/PracticalRockEngineering.asp. Hoek, E., Carranza-Torres, C., and Corkum, B. (2002). Hoek-brown failure criterion 2002 edition. In North American Rock Mechanics Society Meetings, pages 267–73, Toronto, Canada. Hoek, E., Kaiser, P., and Bawden, W. (1995). Support of Underground Excava- tions in Hard Rock. Balkema, Rotterdam. Hoek, E., Marinos, P., and Benissi, M. (1998). Applicability of the geological strength index (gsi) classification for very weak and sheared rock masses. the case of the athens schist formation. Bulletin of Engineering Geology and the Environment, pages 151–160. Hudson, J. and Harrison, J. (2005). Engineering rock mechanics. An introduction to principles. Pergamon. Itasca (2004a). UDEC 4.0: Command Reference, volume 2. Itasca Consulting Group Inc, Minneapolis, Minnesota. Itasca (2004b). UDEC 4.0: Theory and Background, volume 2. Itasca Consulting Group Inc, Minneapolis, Minnesota. Itasca (2004c). UDEC 4.0: Universal Distinct Element Code User’s Guide, volume 2. Itasca Consulting Group Inc, Minneapolis, Minnesota. Itasca (2008). Udec training course: Basic concepts and recommended proceu- dres for geotechnical numerical analysis. Short course. Itasca Inc. (2007). Udec - universal distinct element code. http://www. itascacg.com/udec.html. Kim, B., Cai, M., and Kaiser, P. (2007). Rock mass strength with non-persistent joints. In Eberhardt, E., Stead, D., and Morrison, T., editors, 1ST Canada-US 171 Bibliography Rock Mechanics Symposium, volume 1, pages 241–248, Vancouver. Taylor & Francis. Krahn, J. (2003). The 2001 r.m. hardy lecture: The limits of limit equilibrium analysis. Canadian Geotech Journal, 40:643–660. Krahn, J. (2007). Limit equilibrium, strength summation and strength reduction methods for assessing slope stability. In Eberhardt, Stead, and Morrison, edi- tors, 1ST Canada-US Rock Mechanics Symposium, volume 1, pages 311 – 318, Vancouver. Kulhawy, F. (1975). Stress deformation properties of rock and rock discontinu- ities. Engineering Geology, 9:327–350. Marinos, P., Marinos, V., and Hoek, E. (2005). The geological strength index: applications and limitations. Bulletin of Engineering Geology and the Environ- ment, 64:55–65. Martin, C. and Chandler, N. (1994). The progressive fracture of lac du bon- net granite. International Journal of Rock Mechanics and Mining Sciences, 31(6):643–659. Mostyn, G. and Douglas, K. (2000). Strength of intact rock and rock masses. In Geological Engineering, volume 1, Melbourne, Australia. International confer- ence on Geotechnical and Geological Engineering, Technomic Publishing. Nichol, S., Hungr, O., and Evans, S. (2002). Large-scale brittle and ductile toppling of rock slopes. Canadian Geotech Journal, 39:773–788. Pande, G., Beer, G., and Williams, J. (1990). Numerical methods in rock me- chanics. John Wiley & Sons. Scott, R. (1987). Failure. Geotechnique, 37:423–466. Terzaghi, K. and Peck, R. (1948). Soil Mechanics in Engineering Practice. John Wiley & Sons, New York. Wyllie, D. and Mah, C. (2004). Rock Slope Engineering. Taylor & Francis, 4 edition. Zienkiewicz, O., Humpheson, C., and Lewis, R. W. (1975). Associated and non-associated visco-plasticity and plasticity in soil mechanics. Geotechnique, 25(4):671–689. 172 Appendix A Displacement-based SSR Code This section contains sample source code for the displacement-based SSR method implemented in UDEC using the FISH programming language. Options are given in the source for returning to an initial state between SRF stages or to continuously drop strengths without returning to the initial state. ; STRENGTH REDUCTION FACTOR ; TO THE LETTER APPROACH ; AND CONTINUOUS ; X−DISPLACEMENT VALUES ; LAST UPDATED: JULY 2nd 2008 ; CREATED BY: MATHEW FOURNIER ; ATTACHED TO THESIS : INVESTIGATIONS INTO THE SHEAR STRENGTH REDUCTION METHOD ; USING DISTINCT ELEMENT MODELS : UNIVERSITY OF BRITISH COLUMBIA ; 2008 ; ; USE AT YOUR OWN RISK ; ; ; ; AUTHOR ASSUMES NO LIABILITY FOR ANYTHING THAT MAY GO WRONG ; ; COLLAPSE, LOSS OF LIFE , MONEY, ETC . ; ; ; STRENGTHS PREVIOUSLY DEFINED ; eg . r o d f r i c t i o n d e f i n e d i n an e a r l i e r s t a g e by ; d e f b l a h ; r o d f r i c t i o n = 29 ; end ; b l a h ; prop jma t = 2 f r i c t i o n = r o d f r i c t i o n ; r e q u i r e s STABLE STAGE t o run from , c a l l e d 3 . sav ; a t bot tom of b e n c h i n g r e 3 . sav ; s e t max c y c l e l i m i t h i g h e r a s d e f a u l t i s r e a c h e d f o r s o l v i n g f o r ; f o r c e r a t i o f o r some l e v e l s ( q u i t s b e f o r e f o r c e r a t i o c r i t e r i a i s met ) s o l v e c y c l e 160000 ; s e t p o i n t s t o check f o r d i s p l a c e m e n t e q u i l i b r i u m ; t h r e e p o i n t s t h a t match t h o s e g i v e n i n t h e 173 Appendix A. Displacement-based SSR Code ; d e f i n e d f u n c t i o n s below ; t o p o f s l o p e , midd le o f s l o p e , t o e h i s t n c y c l e 100 x d i s 3 2 2 . 0 , 352 .14 h i s t n c y c l e 100 y d i s 3 2 2 . 0 , 352 .14 h i s t n c y c l e 100 x d i s 3 6 8 . 0 , 316 .938 h i s t n c y c l e 100 y d i s 3 6 8 . 0 , 316 .938 h i s t n c y c l e 100 x d i s 4 1 5 . 0 , 280 .362 h i s t n c y c l e 100 y d i s 4 1 5 . 0 , 280 .362 h i s t n c y c l e 100 x d i s 4 6 2 . 0 , 243 .786 h i s t n c y c l e 100 y d i s 4 6 2 . 0 , 243 .786 ; r e s e t t r a c k i n g d a t a r e s e t d i s p v e l j d i s p ; ; ; ; ; DEFINED FUNCTIONS ; ; ; ; ; ; ; ; ; ; ; GENERATE LOOP CREATES A FILENAME TO CONTROL THE LOOPING SEQUENCE ; ; ; SET n i t e r a t i o n s = x f o r t h e amount o f i n t e r a t i o n s ; ; ; SET n i n c r e m e n t = x f o r t h e amount t o i n c r e m e n t each t ime ; ; ; SET n c y c l e s = x f o r how many c y c l e s t o s t e p each model ; ; ; SET n c o n t = 0 / 1 t o t u r n o f f / on c o n t i n u o u s c y c l i n g ( no i n i t i a l s t a t e ) ; ; ; G e n e r a t e s a f i l e s r f . t x t which c o n t r o l s t h e m u l t i f o s f u n c t i o n ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; d e f g e n e r a t e l o o p a r r a y makeloop ( 1 0 0 ) ; open t h e f i l e f o r t h e SRF loop s t a t u s = open ( ’ s r f . t x t ’ , 1 , 1 ) ; THESE VARIABLES CONTROL THE LOOP SETTINGS ; AND VARIOUS OTHER THINGS ABOUT THE SRF RUN n i t e r a t i o n s = 13 ; how many l o o p s t o run n i n c r e m e n t = 0 . 0 1 ; how much t o i n c r e m e n t each loop n c y c l e s = 20000 ; how many c y c l e s f o r each model n r a t i o = 1e−6 ; i f s o l v i n g f o r f o r c e r a t i o n f o r c e s o l v e = 1 ; f l a g f o r s o l v i n g f o r f o r c e r a t i o : 1 on 0 o f f ; w i l l be i s s u e d ( warn ing ) n s t a r t = 0 . 9 9 ; s t a r t i n g p o s i t i o n o f t h e loop n c o n t = 1 ; t u r n on c o n t i n u o u s c y c l i n g ; no r e t u r n i n g t o an i n i t i a l s t a t e be tween SRF s t e p s ; MUST have n s t a r t s e t t o 1 ; working on a z on in g method i f e x p e c t e d SRF i s f a r from 1 n c o n t s a v e = 1 . 0 2 ; a t which SRF s t a g e t o t u r n on s a v i n g f i l e s ; d i s p l a c e m e n t d a t a w i l l be w r i t t e n f o r e v e r y s t a g e ; however s a v i n g f i l e s f o r e v e r y s t a g e i s a was t e 174 Appendix A. Displacement-based SSR Code ; o f s p a c e i f FOS i s f a r from 1 when do ing ; c o n t i n u o u s c y c l i n g . ; i f u s i n g c o n t i n u o u s c y c l i n g , f o r c e n s t a r t t o have f i r s t s t a g e a t ; SRF = 1 i f n c o n t = 1 n s t a r t = 1 − n i n c r e m e n t e n d i f ; GLOBAL FOR HOW MANY VARIABLES ARE BEING TRACKED n u m v a r i a b l e s = 8 ; how many p o i n t s a r e t r a c k e d i n t h e model ; THESE VARIABLES MAKE THE LOOP BEHAVE SIMILAR TO THE ’EXCLUDE’ ; KEYWORD IN UDEC’ S NORMAL FOS SOLVER ; 1 − t u r n e d on ; 0 − t u r n e d o f f ; d e f i n e p o i n t s i n model t o be t r a c k e d ; t h e s e s h o u l d be t h e same as t h e h i s t o r y p o i n t s p1 = g p n e a r ( 3 2 2 . 0 , 3 5 2 . 1 4 ) p2 = g p n e a r ( 3 6 8 . 0 , 3 1 6 . 9 3 8 ) p3 = g p n e a r ( 4 1 5 . 0 , 2 8 0 . 3 6 2 ) p4 = g p n e a r ( 4 6 2 . 0 , 2 4 3 . 7 8 6 ) ; ; Upper u n i t t o g g l e s ; ; 1 / 0 = on / o f f t o be i n c l u d e d i n t h e r e d u c t i o n r o d r f f r i c t i o n = 1 ; rockmass f r i c t i o n r o d r f c o h e s i o n = 1 ; rockmass c o h e s i o n r o d r f t e n s i o n = 1 ; rockmass t e n s i o n ; ; f r i c t i o n f o r ALL j o i n t s r o d j f f r i c t i o n = 1 ; j o i n t f r i c t i o n r o d j f c o h e s i o n = 1 ; j o i n t c o h e s i o n r o d j f t e n s i o n = 0 ; j o i n t t e n s i o n ; t u r n e d OFF as j o i n t s have z e r o t e n s i o n v a l u e ; ; c o n t a c t zone t o g g l e s d p r f f r i c t i o n = 1 ; rockmass f r i c t i o n d p r f c o h e s i o n = 1 ; rockmass c o h e s i o n d p r f t e n s i o n = 1 ; rockmass t e n s i o n ; ; no j o i n t s i n c o n t a c t zone so t u r n e d OFF ; d p j f f r i c t i o n = 1 ; j o i n t f r i c t i o n ; d p j f c o h e s i o n = 1 ; j o i n t c o h e s i o n ; d p j f t e n s i o n = 1 ; j o i n t t e n s i o n ; C r e a t e t h e loop v i a b a t c h UDEC r u n s ; two d i f f e r e n t v e r s i o n s 175 Appendix A. Displacement-based SSR Code ; − i n i t i a l s t a t e ; − c o n t i n u o u s i f n c o n t = 0 t h e n ; i n i t i a l s t a t e loop n ( 1 , n i t e r a t i o n s ) makeloop ( 1 ) = ’ r e i n t . sav ’ makeloop ( 2 ) = ’ d e f c o u n t l o o p ’ makeloop ( 3 ) = ’ n= ’+ s t r i n g ( i n t ( n ) ) ; n c o n t r o l s c u r r e n t s t a t e o f t h e loop makeloop ( 4 ) = ’ end ’ makeloop ( 5 ) = ’ c o u n t l o o p ’ makeloop ( 6 ) = ’ m u l t i f o s ’ ; c a l l s t h e FOS f u n c t i o n , p a s s i n g t h e p a r a m e t e r makeloop ( 7 ) = ’ ’ ; n s t a t u s = w r i t e ( makeloop , 7 ) ; w r i t e f i l e end loop e n d i f i f n c o n t = 1 t h e n ; c o n t i n u o u s loop n ( 1 , n i t e r a t i o n s ) makeloop ( 1 ) = ’ d e f c o u n t l o o p ’ makeloop ( 2 ) = ’ n= ’+ s t r i n g ( i n t ( n ) ) ; n c o n t r o l s c u r r e n t s t a t e o f t h e loop makeloop ( 3 ) = ’ end ’ makeloop ( 4 ) = ’ c o u n t l o o p ’ makeloop ( 5 ) = ’ m u l t i f o s ’ ; c a l l s t h e FOS f u n c t i o n , p a s s i n g t h e p a r a m e t e r makeloop ( 6 ) = ’ ’ ; n s t a t u s = w r i t e ( makeloop , 6 ) ; w r i t e f i l e end loop e n d i f s t a t u s = c l o s e end ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; MULTI FOS IS THE FACTOR OF SAFETY REDUCTION FUNCTION ; ; ; − t h e f u n c t i o n g e n e r a t e l o o p must be run f i r s t ; ; ; − n i n c r e m e n t must be s e t ; ; ; − t r a c k e d v a r i a b l e s a r e s e t and h a n d l e d ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; d e f m u l t i f o s ; i n i t i a l FOS s e t t i n g s d rad = 0.0174532952 ; d e g r e e t o r a d i a n c o n v e r s i o n rdeg = 57.29577951 ; r a d i a n t o d e g r e e c o n v e r s i o n ; ; ; ; i n t a c t rock s t r e n g t h r e d u c t i o n ; ; ; ; s t r e n g t h r e d u c e d by n i n c r e m e n t each loop f o s c u r r e n t = n s t a r t + n i n c r e m e n t ∗ n ; c a l c u l a t e new v a l u e s f o r i n t a c t uppe r u n i t rod ph inew = a t a n ( t a n ( r o d r f r i c t i o n ∗ d rad ) / ( f o s c u r r e n t ) )∗ r deg rod cohnew = r o d r c o h e s i o n / ( f o s c u r r e n t ) 176 Appendix A. Displacement-based SSR Code r o d t e n n e w = r o d r t e n s i o n / ( f o s c u r r e n t ) ; w r i t e new v a l u e s t o memory i f p r o p e r f l a g s a r e s e t ; i n g e n e r a t e l o o p f u n c t i o n i f r o d r f f r i c t i o n = 1 t h e n command prop mat=1 f r i c = rod ph inew end command e n d i f i f r o d r f c o h e s i o n = 1 t h e n command prop mat=1 coh= rod cohnew end command e n d i f i f r o d r f t e n s i o n = 1 t h e n command prop mat=1 t e n = r o d t e n n e w end command e n d i f ; c a l c u l a t e new v a l u e s f o r c o n t a c t zone dp phinew = a t a n ( t a n ( d p r f r i c t i o n ∗ d rad ) / ( f o s c u r r e n t ) )∗ r deg dp cohnew = d p r c o h e s i o n / ( f o s c u r r e n t ) dp tennew = d p r t e n s i o n / ( f o s c u r r e n t ) ; w r i t e new v a l u e s t o memory i f p r o p e r f l a g s a r e s e t ; i n g e n e r a t e l o o p f u n c t i o n i f d p r f f r i c t i o n = 1 t h e n command prop mat=3 f r i c = dp phinew end command e n d i f i f d p r f c o h e s i o n = 1 t h e n command prop mat=3 coh=dp cohnew end command e n d i f i f d p r f t e n s i o n = 1 t h e n command prop mat=3 t e n = dp tennew end command e n d i f ; j o i n t s t r e n g t h r e d u c t i o n ; c a l c u l a t e new v a r i a b l e s f o r uppe r u n i t j 1 177 Appendix A. Displacement-based SSR Code r o d j 1 p h i n e w = a t a n ( t a n ( r o d j 1 f r i c t i o n ∗ d rad ) / ( f o s c u r r e n t ) )∗ r deg r o d j 1 c o h n e w = r o d j 1 c o h e s i o n / ( f o s c u r r e n t ) r o d j 1 t e n n e w = r o d j 1 t e n s i o n / ( f o s c u r r e n t ) ; c a l c u l a t e new v a r i a b l e s f o r uppe r u n i t j 2 r o d j 2 p h i n e w = a t a n ( t a n ( r o d j 2 f r i c t i o n ∗ d rad ) / ( f o s c u r r e n t ) )∗ r deg r o d j 2 c o h n e w = r o d j 2 c o h e s i o n / ( f o s c u r r e n t ) r o d j 2 t e n n e w = r o d j 2 t e n s i o n / ( f o s c u r r e n t ) ; w r i t e new v a r i a b l e s t o memory i f r o d j f f r i c t i o n = 1 t h e n command prop jma t =1 j f r i c = r o d j 1 p h i n e w prop jma t =2 j f r i c = r o d j 2 p h i n e w end command e n d i f i f r o d j f c o h e s i o n = 1 t h e n command prop jma t =1 j c o h = r o d j 1 c o h n e w prop jma t =2 j c o h = r o d j 2 c o h n e w end command e n d i f i f r o d j f t e n s i o n = 1 t h e n command prop jma t =1 j t e n = r o d j 1 t e n n e w prop jma t =2 j t e n = r o d j 2 t e n n e w end command e n d i f ; s t e p t h e model by n c y c l e s o r s o l v e f o r ; f o r c e e q u i l i b r i u m ; c o n t r o l l i n g f l a g s ; n f o r c e s o l v e = 1 t u r n s on f o r c e r a t i o s o l v i n g ; = 2 j u s t i s s u e s a s o l v e command ( may hang ) ; = 0 s t e p s n c y c l e s command p r i n t n f o s c u r r e n t end command i f n f o r c e s o l v e = 1 t h e n ; s o l v e f o r a f o r c e r a t i o command p r i n t n r a t i o s o l v e r a t i o n r a t i o end command e n d i f i f n f o r c e s o l v e = 0 t h e n ; s o l v e f o r a s e t number o f c y c l e s command p r i n t n c y c l e s s t e p n c y c l e s end command 178 Appendix A. Displacement-based SSR Code e n d i f ; pause ; w r i t e save f i l e s based o f f n p a r a m e t e r ; f o s c i f c o n t i n u o u s ; f o s n i f i n i t i a l s t a t e i f n c o n t = 0 t h e n sav name = ’ f o s n ’+ s t r i n g ( i n t (1000∗ f o s c u r r e n t ) ) + ’ . sav ’ command save sav name end command e n d i f i f n c o n t = 1 t h e n i f f o s c u r r e n t > n c o n t s a v e t h e n sav name = ’ f o s c ’+ s t r i n g ( i n t (1000∗ f o s c u r r e n t ) ) + ’ . sav ’ command save sav name end command e n d i f e n d i f ; g e n e r a t e t a b l e s f o r c h a r t s t o be c r e a t e d ; w r i t e t r a c k i n g v a r i a b l e s t o f i l e ; i f new v a r i a b l e s added , u p d a t e n u m v a r i a b l e s ; a p p r o p r i a t e l y and f o l l o w a r r a y c o n v e n t i o n ; f o r x d i s p l a c e m e n t s ; a r r a y t o c a r r y t r a c k e d v a r i a b l e s a r r a y x d i s p l a c e m e n t s ( 4 0 0 0 ) ; c a n n o t append t o a f i l e u s i n g FISH f i l e i / o ; i n s t e a d w r i t e s a f i l e and r e p l a c e s t h e f i l e on n e x t w r i t e ; so i n s t e a d r e a d t h e f i l e i n t o memory each loop and t h e n ; w r i t e o l d p a r a m e t e r s a s w e l l a s t h e new p a r a m e t e r s ; t h i s i s n o t n e c e s s a r y on t h e f i r s t i t e r a t i o n i f n > 1 s t a t u s = open ( ’ x d i s . da t ’ , 0 , 0 ) s t a t u s = r e a d ( x d i s p l a c e m e n t s , n u m v a r i a b l e s ∗ ( n−1)) s t a t u s = c l o s e e n d i f ; u p d a t e t h i s a c c o r d i n g l y wi th how many v a r i a b l e s i n n u m v a r i a b l e s ; w r i t e s t h e d a t a b e i n g t r a c k e d x d i s p l a c e m e n t s ( ( n u m v a r i a b l e s ∗ ( n−1))+1 )= s t r i n g ( abs ( g p x d i s ( p1 ) ) ) x d i s p l a c e m e n t s ( ( n u m v a r i a b l e s ∗ ( n−1))+2 )= s t r i n g ( abs ( g p y d i s ( p1 ) ) ) x d i s p l a c e m e n t s ( ( n u m v a r i a b l e s ∗ ( n−1))+3 )= s t r i n g ( abs ( g p x d i s ( p2 ) ) ) x d i s p l a c e m e n t s ( ( n u m v a r i a b l e s ∗ ( n−1))+4 )= s t r i n g ( abs ( g p y d i s ( p2 ) ) ) x d i s p l a c e m e n t s ( ( n u m v a r i a b l e s ∗ ( n−1))+5 )= s t r i n g ( abs ( g p x d i s ( p3 ) ) ) 179 Appendix A. Displacement-based SSR Code x d i s p l a c e m e n t s ( ( n u m v a r i a b l e s ∗ ( n−1))+6 )= s t r i n g ( abs ( g p y d i s ( p3 ) ) ) x d i s p l a c e m e n t s ( ( n u m v a r i a b l e s ∗ ( n−1))+7 )= s t r i n g ( abs ( g p x d i s ( p4 ) ) ) x d i s p l a c e m e n t s ( ( n u m v a r i a b l e s ∗ ( n−1))+8 )= s t r i n g ( abs ( g p y d i s ( p4 ) ) ) ; w r i t e s bo th t h e o l d d a t a and t h e new d a t a ( hence n u m v a r i a b l e s ∗n ) s t a t u s = open ( ’ x d i s . da t ’ , 1 , 0 ) s t a t u s = w r i t e ( x d i s p l a c e m e n t s , n u m v a r i a b l e s ∗n ) s t a t u s = c l o s e ; s a v e s t h e f i l e based o f f o f p a r a m e t e r n end ; ; ; ; ; ; FUNCTION CONVERT ARRAY ; ; ; ; ; ; Th i s f u n c t i o n r e a d s d a t a w r i t t e n t o f i l e i n FISH b i n a r y ; i n t o memory . U n f o r t u n e t l y FISH h a n d l e s t h i s a s a s t r i n g ; d e s p i t e i t c l e a r l y b e i n g a f l o a t once i t i s b r o u g h t i n . ; t h e i n n e r f u n c t i o n CONVINT w r i t e s a f i l e t h a t e x p l i c i t y ; b r i n g s t h e v a l u e s i n t o FISH as f l o a t s r a t h e r t h a n s t r i n g s ; i n o r d e r t o p r o c e s s them p r o p e r l y ; r e q u i r e d : n u m v a r i a b l e s − number o f t a b l e s b e i n g t r a c k e d ; n − number o f SRF l o o p s ; d i r e c t l y w r i t i n g t h e i n f o r m a t i o n a t t h e end of t h e loop ; c r a s h e s UDEC f o r an unknown r e a s o n . G e t t i n g around t h i s ; happens v i a t h e f u n c t i o n c r e a t e t a b l e s which , f o r some ; r ea son , doesn ’ t c r a s h d e s p i t e c o n t a i n i n g t h e e x a c t ; same code t h a t would c r a s h UDEC. d e f c o n v e r t a r r a y ; s e t up a r r a y s a r r a y x r e a d ( 4 0 0 0 ) a r r a y xconv ( 4 0 0 0 ) a r r a y xtemp ( 1 0 ) ; open t h e f i l e f o r r e a d i n g i n FISH b i n a r y mode s t a t u s = open ( ’ x d i s . da t ’ , 0 , 0 ) ; r e a d e v e r y t h i n g i n t o memory s t a t u s = r e a d ( xread , n u m v a r i a b l e s ∗ ( n ) ) s t a t u s = c l o s e ; g e n e r a t e t a b l e s ; c o n v e r t s t r i n g a r r a y ( x r e a d ) i n t o some th ing I can a c t u a l l y use ; open f i l e t o w r i t i n g ( a s c i i ) and g e n e r a t e h e a d e r 180 Appendix A. Displacement-based SSR Code s t a t u s = open ( ’ c o n v e r t . da t ’ , 1 , 1 ) xtemp ( 1 ) = ’ d e f c o n v i n t ’ s t a t u s = w r i t e ( xtemp , 1 ) ; w r i t e t h e d a t a t o t h e f i l e l i n e by l i n e loop i ( 1 , n u m v a r i a b l e s ∗ ( n ) ) xtemp ( 1 ) = ’ xconv ( ’ + s t r i n g ( i ) + ’ ) = ’+ s t r i n g ( x r e a d ( i ) ) command p r i n t xtemp ( 1 ) end command s t a t u s = w r i t e ( xtemp , 1 ) e n d l o o p xtemp ( 1 ) = ’ end ’ xtemp ( 2 ) = ’ c o n v i n t ’ ; f i n i s h t h e f i l e s t a t u s = w r i t e ( xtemp , 2 ) s t a t u s = c l o s e ; b r i n g t h e f i l e i n t o re−w r i t e v a l u e s as f l o a t s command c a l l c o n v e r t . d a t end command end ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; FUNCTION CREATE TABLES ; ; ; r e a d s d a t a from memory s t o r e d i n xconv a r r a y i n t o s e p a r a t e t a b l e s ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; d e f c r e a t e t a b l e s ; goes t h r o u g h loop f o r t h e amount o f d a t a w r i t t e n loop i ( 1 , n ) ; f i r s t m o n i t o r i n g p o i n t x t a b l e ( 1 0 0 , i ) = n s t a r t + n i n c r e m e n t ∗ i y t a b l e ( 1 0 0 , i ) =xconv ( 1 + ( ( n u m v a r i a b l e s )∗ ( i −1) ) ) x t a b l e ( 1 0 1 , i ) = n s t a r t + n i n c r e m e n t ∗ i y t a b l e ( 1 0 1 , i ) = xconv ( 2 + ( ( n u m v a r i a b l e s )∗ ( i −1) ) ) ; second m o n i t o r i n g p o i n t x t a b l e ( 2 0 0 , i ) = n s t a r t + n i n c r e m e n t ∗ i y t a b l e ( 2 0 0 , i ) =xconv ( 3 + ( ( n u m v a r i a b l e s )∗ ( i −1) ) ) 181 Appendix A. Displacement-based SSR Code x t a b l e ( 2 0 1 , i ) = n s t a r t + n i n c r e m e n t ∗ i y t a b l e ( 2 0 1 , i ) =xconv ( 4 + ( ( n u m v a r i a b l e s )∗ ( i −1) ) ) ; t h i r d m o n i t o r i n g p o i n t x t a b l e ( 3 0 0 , i ) = n s t a r t + n i n c r e m e n t ∗ i y t a b l e ( 3 0 0 , i ) =xconv ( 5 + ( ( n u m v a r i a b l e s )∗ ( i −1) ) ) x t a b l e ( 3 0 1 , i ) = n s t a r t + n i n c r e m e n t ∗ i y t a b l e ( 3 0 1 , i ) =xconv ( 6 + ( ( n u m v a r i a b l e s )∗ ( i −1) ) ) ; f o u r t h m o n i t o r i n g p o i n t x t a b l e ( 4 0 0 , i ) = n s t a r t + n i n c r e m e n t ∗ i y t a b l e ( 4 0 0 , i ) =xconv ( 7 + ( ( n u m v a r i a b l e s )∗ ( i −1) ) ) x t a b l e ( 4 0 1 , i ) = n s t a r t + n i n c r e m e n t ∗ i y t a b l e ( 4 0 1 , i ) =xconv ( 8 + ( ( n u m v a r i a b l e s )∗ ( i −1) ) ) e n d l o o p end ; end of f u n c t i o n ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ;−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; END FUNCTION DEFINITIONS ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ;−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−; ; c a l l g e n e r a t e l o o p f u n c t i o n t o c r e a t e b a t c h f i l e f o r l o o p i n g g e n e r a t e l o o p ; s a v e s t h e i n i t i a l s t a t e , p l a c e m e n t i s i m p o r t a n t , a s i t ; keeps a l l t h e above s e t t i n g s i n memory so t h a t whent ; i n t . sav i s r e s t o r e d a l l t h e i m p o r t a n t v a r i a b l e s such as ; n i n c r e m e n t e t c a r e t r a c k e d save i n t . sav ; c a l l s t h e b a t c h f i l e f o r e x t e r n a l l o o p i n g c a l l s r f . t x t ; t o g e n e r a t e t a b l e s , a t any t ime r e s t o r e l a s t s ave f i l e ; even i f t h e f i l e c r a s h e s o u t due t o c o n t a c t o v e r l a p t o o g r e a t ; e . g . r e f o s n 1 0 7 9 . sav ; and c a l l t h e f u n c t i o n s : ; c o n v e r t a r r a y −−−> f o l l o w e d by ; c r e a t e t a b l e s ; c a l l o u t p u t . t x t t o dump t a b l e s c r e a t e d by c r e a t e t a b l e s 182 Appendix A. Displacement-based SSR Code ; t o t e x t f i l e s f o r e x t e r n a l use ( e x c e l ) 183 Appendix A. Displacement-based SSR Code The next listing is a bit of code to dump the contents of the table to file: ; BEGIN OUTPUT . TXT ; r e s t o r e t h e l a s t s ave f i l e i n t h e s r f l oop ; t h e n c a l l o u t p u t . t x t c o n v e r t a r r a y c r e a t e t a b l e s t a b l e 100 w r i t e 0 d a t a t 1 0 0 . t x t t a b l e 101 w r i t e 0 d a t a t 1 0 1 . t x t t a b l e 200 w r i t e 0 d a t a t 2 0 0 . t x t t a b l e 201 w r i t e 0 d a t a t 2 0 1 . t x t t a b l e 300 w r i t e 0 d a t a t 3 0 0 . t x t t a b l e 301 w r i t e 0 d a t a t 3 0 1 . t x t t a b l e 400 w r i t e 0 d a t a t 4 0 0 . t x t t a b l e 401 w r i t e 0 d a t a t 4 0 1 . t x t 184 Appendix B Example Input Files for Parametric Study This section contains sample source code for the parametric study. Many different joint networks and other parameters were run and this is merely one example of the 10 m slope parallel jointing case. Contains three sets of code, 1.txt, 2.txt, and 3.txt, which are different sequences in the model generation. The first file: ; BEGIN 1 . t x t ; c r e a t e s b l o c k and s e t s geomet ry ; c r e a t e s j o i n t s and e l e m e n t meshing ; 10 M d i p p a r a l l e l c a s e new t i t l e ’UDEC p a r a m e t r i c s t u d i e s − 10 m d i p p a r a l l e l ’ ; ; d e f i n e geomet ry round 0 . 1 ; d e f a u l t =0 .5 ; ; s p e c i f y a minimum v a l u e o f b l o c k edge l e n g t h s e t edge 1 ; d e f a u l t ; ; c r e a t e a s i n g l e r i g i d b l o c k d e f i n i n g t h e o r i g i n a l boundary ; c r e a t e u s i n g block , l e f t hand r u l e b l 0 ,−250 0 ,500 1350 ,500 1350 ,−250 ; c r 438 ,500 688 ,250 c r 682 ,250 1350 ,250 ; c r 6 1 1 . 7 , 326 .28 1350 , 326 .28 ; c r 5 2 8 . 9 9 7 , 409 .01 1350 , 409 .01 c r 5 9 7 . 6 9 4 , 335 .680 1350 , 335 .680 c r 5 2 8 . 5 1 6 , 404 .764 1350 , 404 .764 185 Appendix B. Example Input Files for Parametric Study j r e g i o n i d 7 0 , 500 175 , 500 175 , −250 0 , −250 j r e g i o n i d 6 175 , 500 438 ,500 788 , 150 175 , 150 j r e g i o n i d 8 682 , 250 920 , 250 880 , 150 788 , 150 j r e g i o n i d 9 175 , 150 1350 , 150 1350 , −250 175 , −250 j r e g i o n i d 10 880 , 250 1350 , 250 1350 , 150 880 , 150 ; ; g e n e r a t e j o i n t p a t t e r n j s e t −45,0 100 ,0 0 ,0 10 ,0 r a n g e j r e g i o n 6 j s e t 45 ,0 100 ,0 0 ,0 10 ,0 r a n g e j r e g i o n 6 j s e t −45,0 100 ,0 0 ,0 10 ,0 r a n g e j r e g i o n 8 j s e t 45 ,0 100 ,0 0 ,0 10 ,0 r a n g e j r e g i o n 8 j s e t −45,0 100 ,0 0 ,0 40 ,0 r a n g e j r e g i o n 7 j s e t 45 ,0 100 ,0 0 ,0 40 ,0 r a n g e j r e g i o n 7 j s e t −45,0 100 ,0 0 ,0 40 ,0 r a n g e j r e g i o n 9 j s e t 45 ,0 100 ,0 0 ,0 40 ,0 r a n g e j r e g i o n 9 j s e t −45,0 100 ,0 0 ,0 40 ,0 r a n g e j r e g i o n 10 j s e t 45 ,0 100 ,0 0 ,0 40 ,0 r a n g e j r e g i o n 10 j d e l ; main j o i n t s p a c i n g ; 5 ) Mesh G e n e r a t i o n ch mat 2 ch mat 3 r a n g e b l 5923 ch mat 3 r a n g e b l 6369 ch mat 3 r a n g e b l 6848 ; c o a r s e mesh away gen edge 25 r a n g e j r e g i o n 7 gen edge 25 r a n g e j r e g i o n 9 gen edge 25 r a n g e j r e g i o n 10 gen edge 25 r a n g e mat 3 ; f i n e mesh on s l o p e gen edge 5 r a n g e j r e g i o n 6 gen edge 5 r a n g e j r e g i o n 8 ch mat 2 save 1 . sav ca 2 . t x t 186 Appendix B. Example Input Files for Parametric Study The second file: ; BEGIN 2 . t x t ; s e t s s t r e n g t h s ; e l a s t i c g r a v i t y l o a d i n g r e s 1 . sav ; s e t c o n s t i t u t i v e models f o r i n t a c t r o c k s and j o i n t s ; i n t a c t rock models : ; 2 ( cons =1) − e l a s t i c , i s o t r o p i c ; 3 ( cons =3) − Mohr−Coulomb p l a s t i c i t y ; j o i n t models : ; 2 ( j c o n s = 2) − j o i n t a r e a c o n t a c t − coulomb s l i p ; 1 ( j c o n s = 1) − p o i n t c o n t a c t − coulomb s l i p ; 3 ( j c o n s = 5) − coloumb s l i p wi th r e s . s t r e n g t h ; s e t i n t a c t rock t o MC p l a s i c i t y change cons 3 ; s e t j o i n t s t o j o i n t a r e a c o n t a c t , coulomb s l i p change j c o n s 2 ;−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− ; DEFAULT PROPERITES AND CONST . RELATION FOR NE CONTACTS s e t j m a t d f =1 j c o n d f =2 ;−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− ; ; ; ; ; ; ; ; INPUT OF ROCK MASS MATERIAL PROPERTIES ; ; ; ; ; ; ; ; ; ; ; ; ; ; ;∗∗ l i s t o f p r o p e r t i e s used ; mat = 1 − d e f a u l t m a t e r i a l ; mat = 2 − j o i n t e d rock ; mat = 3 − u n j o i n t e d rock o f t y p e 2 ;−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− ; s i g c i = 25 ; g s i = 50 ; mi = 10 ; MR = 275 ; s l o p e s ; s l o p e h e i g h t 250m ; ; ; D e f a u l t m a t e r i a l d e f Emat1 b l k =m1E/(3∗(1−2∗m1pr ) ) ; bu lk modulus she =m1E/ ( 2∗ ( 1 + m1pr ) ) ; s h e a r modulus command prop mat=1 d=2600 bu lk = b l k s h e a r = she f r i c =56 coh =1.771 e6 t e n =0 zone dens =2600 bu lk = b l k s h e a r = she f r i c =56 coh =1.771 e6 t e n =0 r a n g e mat 1 end command end s e t m1E=7 e9 m1pr =0 .33 Emat1 187 Appendix B. Example Input Files for Parametric Study ; ; ; J o i n t e d rock d e f Emat2 b l k =m2E/(3∗(1−2∗m1pr ) ) ; bu lk modulus she =m2E/ ( 2∗ ( 1 + m1pr ) ) ; s h e a r modulus command prop mat=2 d=2600 bu lk = b l k s h e a r = she f r i c =33 coh =1.01 e6 t e n =100 e3 end command end s e t m2E=7 e9 m1pr =0 .27 Emat2 ; ; ; U n j o i n t e d rock d e f Emat3 b l k =m3E/(3∗(1−2∗m1pr ) ) ; bu lk modulus she =m3E/ ( 2∗ ( 1 + m1pr ) ) ; s h e a r modulus command prop mat=3 d=2600 bu lk = b l k s h e a r = she f r i c =33 coh =1.01 e6 t e n =100 e3 end command end s e t m3E=7 e9 m1pr =0 .27 Emat3 ; ; ; ; ; ; ; ; ; ; ; INPUT OF JOINT MATERIAL PROPERTIES ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ;∗∗ l i s t o f p r o p e r t i e s used ; jma t = 1 − a r t i f i c i a l f r a c t u r e s ; jma t = 5 − JS # 1 ; jma t = 6 − JS # 2 ; jma t = 7 − JS # 3 ;−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− ; ; ; S e t a l l a r t i f i c a l j o i n t s t o jma t =1 change jma t =1 ; a r t i f i c i a l j o i n t s prop jma t =1 j k n =10 e9 j k s =8 e9 j f r i c =65 j c o h =1 e6 j t e n =100 e3 ; ; ; JS #1 prop jma t =2 j k n =10 e9 j k s =8 e9 j f r i c =25 .0 j c o h =100 e3 j t e n =0 ; ; ; JS #2 prop jma t =3 j k n =10 e9 j k s =8 e9 j f r i c =25 .0 j c o h =100 e3 j t e n =0 ; ; ; ; ; ; ; ; ; ; ; ; BOUNDARY CONDITIONS ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; NEEDS TO BE CHANGED EVERY TIME GEOMETRY IS CHANGED ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; Def Bound Buf fe r ; f o r s p e c i f y i n g boundary c o n d i t i o n s x y b u f f e r = 2 ; 1 0 . 0 x l e f t m = 0.00− x y b u f f e r x l e f t p = 0 .00+ x y b u f f e r x r i g h t m = 1350− x y b u f f e r x r i g h t p = 1350+ x y b u f f e r y bot tom m =−250−x y b u f f e r y b o t t o m p =−250+ x y b u f f e r y top m = 500− x y b u f f e r 188 Appendix B. Example Input Files for Parametric Study y t o p p = 500+ x y b u f f e r End Bound Buf fe r ; ; ; R o l l e r Boundary on s i d e s , bot tom bound x v e l = 0 . 0 r a n g e x l e f t m , x l e f t p y bottom m , y t o p p ; l e f t s i d e ; bound x v e l = 0 . 0 r a n g e x r i g h t m , x r i g h t p y bottom m , y t o p p ; r i g h t s i d e ; bound y v e l = 0 . 0 r a n g e x l e f t m , x r i g h t p y bottom m , y b o t t o m p ; a t t h e bot tom s e t g r av 0 −9.81 ; ; ; ; ; ; ; ; ; ; ; INITIAL STRESS STATE ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; a v e r a g e d e n s i t y assumed t o be : 2500 kg / m3 Def C a l I n s i t u S t r e s s a d e n s i t y = 2100 .0 ; a v e r a g e d e n s i t y a g r = −9.81 ; g r a v i t y a r a t i o = 1 . 0 ; r a t i o n o f h o r i z o n t a l t o v e r t i c a l s t r e s s a t o p = 500 ; t o p h e i g h t o f model , needs t o be changed a ne g = −1 ; e v e r y t i m e model geomet ry i s changed syy 0 = a d e n s i t y ∗ a g r ∗ a t o p ; syy a t z e r o e l e v a t i o n sxy 0 = 0 . 0 ; sxy a t z e r o e l e v a t i o n sxx 0 = sy y 0 ∗ a r a t i o ; sxx a t z e r o e l e v a t i o n s z z 0 = syy 0 ∗ a r a t i o ; s z z a t z e r o e l e v a t i o n ; y s t r e s s g r a d i e n t s y s y y g = a d e n s i t y ∗ a g r ∗ a ne g y s x x g = a d e n s i t y ∗ a g r ∗ a r a t i o ∗ a ne g y s z z g = a d e n s i t y ∗ a g r ∗ a r a t i o ∗ a ne g End C a l I n s i t u S t r e s s ; s e t i n s i t u s t r e s s reg ime i n s i t u s t r sxx 0 , sxy 0 , syy 0 yg y sxx g , 0 , y s y y g ; ywt t a b 100 ; w a t e r t a b l e NOT used i n s i t u s z z s z z 0 zg 0 , y s z z g ;∗ A l l m a t e r i a l s a r e e l a s t i c t o r e a c h i n i t i a l e q u i l i b r i u m ch cons =1 Emat1 ; d e f a u l t Emat2 ; j o i n t e d Emat3 ; u n j o i n t e d damp a u t o h i s n c y c l e 100 u n b a l ; 1 h i s n c y c l e 100 damp ; 2 189 Appendix B. Example Input Files for Parametric Study h i s n c y c l e 100 vmax ; 3 s o l v e ; e l a s t i c s t e p 10000 sav 2 . sav c a l l 3 . t x t 190 Appendix B. Example Input Files for Parametric Study The third file: ; BEGIN 3 . t x t ; s e t s r e a l s t r e n g t h s ; benches down ; change rockmass t o MC b e h a v i o u r r e 2 . sav ch cons =3 ; zone model mohr ch j c o n s =2 s e t j m a t d f =2 j c o n d f =2 ch jma t =2 Emat1 ; d e f a u l t Emat2 ; j o i n t e d rock Emat3 ; u n j o i n t e d rock s o l v e ; S e t j o i n t s t r e n g t h s , k e e p i n g a r t i f i c i a l j o i n t s a s jma t1 ch j c o n s =2 prop jma t =1 j k n =10 e9 j k s =8 e9 j f r i c =65 j c o h =1 e6 j t e n =0 ; ; ; JS #1 prop jma t =2 j k n =10 e9 j k s =8 e9 j f r i c =25 .0 j c o h =100 e3 j t e n =0 ; ; ; JS #2 prop jma t =3 j k n =10 e9 j k s =8 e9 j f r i c =25 .0 j c o h =100 e3 j t e n =0 ch jma t =2 s o l v e s t e p 5000 d e l b l 5923 s o l v e s t e p 8000 d e l b l 6369 d e l b l 170591 s o l v e s t e p 8000 d e l b l 6848 s o l v e s t e p 15000 sav 3 . sav ; Ready f o r SSR from now ; e i t h e r s o l v e r a t i o =1e−6 f o s ; o r c a l l d i s p l a c e m e n t based SSR method 191 Appendix C Case History Joint Networks    UDEC (Version 4.00) LEGEND    30-Jul-08   9:54  cycle    75801  time    7.263E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J1 joint network used in case history modelling 192 Appendix C. Case History Joint Networks    UDEC (Version 4.00) LEGEND     5-Aug-08  22:34  cycle    78690  time    6.158E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J2 joint network used in case history modelling    UDEC (Version 4.00) LEGEND     5-Aug-08  22:58  cycle    86580  time    4.557E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J3 joint network used in case history modelling 193 Appendix C. Case History Joint Networks    UDEC (Version 4.00) LEGEND     5-Aug-08  22:46  cycle    82060  time    5.108E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J4 joint network used in case history modelling    UDEC (Version 4.00) LEGEND     6-Aug-08  11:13  cycle    83330  time    5.576E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J5 joint network used in case history modelling 194 Appendix C. Case History Joint Networks    UDEC (Version 4.00) LEGEND     6-Aug-08  17:06  cycle    83150  time    5.002E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J6 joint network used in case history modelling    UDEC (Version 4.00) LEGEND     6-Aug-08  16:49  cycle    86880  time    4.492E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J7 joint network used in case history modelling 195 Appendix C. Case History Joint Networks    UDEC (Version 4.00) LEGEND     8-Aug-08  14:11  cycle    84360  time    4.962E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J8 joint network used in case history modelling    UDEC (Version 4.00) LEGEND     9-Aug-08  16:03  cycle    84530  time    4.764E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J9 joint network used in case history modelling 196 Appendix C. Case History Joint Networks    UDEC (Version 4.00) LEGEND     9-Aug-08  11:03  cycle    85790  time    4.752E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J10 joint network used in case history modelling    UDEC (Version 4.00) LEGEND     8-Aug-08  17:08  cycle    84320  time    4.971E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J11 joint network used in case history modelling 197 Appendix C. Case History Joint Networks    UDEC (Version 4.00) LEGEND    11-Aug-08  12:52  cycle    87980  time    4.122E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J12 joint network used in case history modelling    UDEC (Version 4.00) LEGEND     7-Aug-08  11:28  cycle    86920  time    4.387E+01 sec block plot  2.000  2.400  2.800  3.200  3.600  4.000 (*10^2)  2.600  3.000  3.400  3.800  4.200  4.600 (*10^2) JOB TITLE : University of Brish Columbia Vancouver, BC Canada J13 joint network used in case history modelling 198

Cite

Citation Scheme:

    

Usage Statistics

Country Views Downloads
China 13 81
United States 9 1
Japan 7 0
Germany 3 0
Canada 2 0
Iran 2 2
United Kingdom 2 0
Bosnia and Herzegovina 1 0
Chile 1 0
Australia 1 0
Malaysia 1 1
City Views Downloads
Beijing 7 0
Tokyo 7 0
Shenyang 6 81
Ashburn 4 0
Sunnyvale 3 0
Aachen 3 0
Washington 2 0
Addlestone 2 0
Khalilabad 2 2
Unknown 2 1
Brossard 2 0
Santiago 1 0
Sarajevo 1 0

{[{ mDataHeader[type] }]} {[{ month[type] }]} {[{ tData[type] }]}

Share

Share to:

Comment

Related Items