Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Investigations into the Shear Strength Reduction method using distinct element models Fournier, Mathew 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 [ 36.29MB ]
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
Original Record: 1.0052558 +original-record.json
Full Text
1.0052558.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 Reduction (SSR) method to determine factor of safety values in discontinuum models 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 displacementmonitoring based method. A parametric study was first undertaken, using a simple homogeneous rock slope, with three different joint networks representing common 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 constitutive models that better idealize the rockmass than simpler methods such as limit equilibrium (e.g. either method of slices or wedge solutions) and even numerical continuum models (e.g. finite difference, finite element). Joints are explicitly 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. Strainsoftening models are also discussed with respect to the SSR method. The results presented highlight several important relationships to consider related to both numerical 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 factor 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. Softening 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 Problem Statement . . . . . 1.2 Research Objectives . . . . 1.2.1 Sensitivity Analysis 1.2.2 Case History . . . . 1.3 Thesis Structure . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  2  Background . . . . . . . . . . . . . . . . . . . . . . . . 2.1 Definition of Failure . . . . . . . . . . . . . . . . . 2.2 Use of Factor of Safety . . . . . . . . . . . . . . . . 2.2.1 Method of Slices . . . . . . . . . . . . . . . 2.2.2 Numerical Models . . . . . . . . . . . . . . 2.2.2.1 Continuum Modeling . . . . . . . 2.2.2.2 Distinct-Element Modeling . . . . 2.2.2.3 Shear Strength Reduction Method  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  . 5 . 5 . 6 . 9 . 11 . 13 . 13 . 13  3  Methodology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1 Geomechanical Behavior of Rockmasses and Intact Rock . . . 3.1.1 Elasticity . . . . . . . . . . . . . . . . . . . . . . . . . 3.1.2 Plasticity . . . . . . . . . . . . . . . . . . . . . . . . .  . 15 . 15 . 17 . 19  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  1 1 2 2 3 3  iii  Table of Contents 3.1.3  3.2 3.3 3.4 3.5  3.6  Yield Criterion . . . . . . . . . . . 3.1.3.1 Mohr-Coulomb . . . . . 3.1.3.2 Hoek-Brown . . . . . . 3.1.4 Rockmass Classification Schemes . 3.1.4.1 Rock Mass Rating . . . . 3.1.4.2 Geological Strength Index Behaviour of Discontinuities . . . . . . . . Discontinuity Network . . . . . . . . . . . Limit Equilibrium . . . . . . . . . . . . . Universal Distinct Element Code . . . . . 3.5.1 Constitutive Models . . . . . . . . 3.5.1.1 Intact Rocks . . . . . . . 3.5.1.2 Joints . . . . . . . . . . Shear Strength Reduction Using UDEC . . 3.6.1 Displacement-Based SSR . . . . .  . . . . . . . . . . . . . .  . . . . . . . . . . . . . . .  . . . . . . . . . . . . . . .  . . . . . . . . . . . . . . .  . . . . . . . . . . . . . . .  . . . . . . . . . . . . . . .  . . . . . . . . . . . . . . .  . . . . . . . . . . . . . . .  . . . . . . . . . . . . . . .  . . . . . . . . . . . . . . .  . . . . . . . . . . . . . . .  . . . . . . . . . . . . . . .  20 20 22 23 26 26 27 28 30 30 35 36 37 38 41  4  Parametric Study . . . . . . . . . . . . . . . 4.1 Cycling Mode . . . . . . . . . . . . . . 4.2 Element Size . . . . . . . . . . . . . . . 4.3 Element Types . . . . . . . . . . . . . . 4.4 Insitu Stress . . . . . . . . . . . . . . . . 4.5 Elastic Properties . . . . . . . . . . . . . 4.6 Joint Shear Stiffness . . . . . . . . . . . 4.7 Strain-Softening Models . . . . . . . . . 4.8 Joint Network Cases . . . . . . . . . . . 4.8.1 Daylighting Joints . . . . . . . . 4.8.2 Discontinuous Daylighting Joints 4.8.3 Slope Parallel Joints . . . . . . . 4.9 Summary . . . . . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  . . . . . . . . . . . . .  45 50 58 63 67 67 72 79 81 81 85 90 94  5  Case History . . . . . . . . . . . . . . . . . . . 5.1 Introduction . . . . . . . . . . . . . . . . . 5.1.1 Geology . . . . . . . . . . . . . . 5.1.2 Slope Failure . . . . . . . . . . . . 5.1.3 Structure . . . . . . . . . . . . . . 5.1.4 Geotechnical Properties . . . . . . 5.1.5 Modelling Approach . . . . . . . . 5.1.6 Selection of Critical Sliding Surface  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  . . . . . . . .  97 97 97 98 100 104 106 109 iv  Table of Contents 5.2  5.3 5.4 6  Results . . . . . . . . . . . . . . . . . . . . . 5.2.1 Limit Equilibrium . . . . . . . . . . . 5.2.2 UDEC: Unjointed Slope Models . . . 5.2.3 UDEC: Continuous Bedding Models . 5.2.4 UDEC: Discontinuous Joint Models . . 5.2.4.1 Strain-Softening Behaviour . 5.2.4.2 Material Property Variations 5.2.4.3 Water Influence . . . . . . . Continuous Cycling Solution Mode . . . . . . Discussion . . . . . . . . . . . . . . . . . . .  Discussion and Conclusions . . . . . . . . . . 6.1 Framework . . . . . . . . . . . . . . . . . 6.2 Discussion . . . . . . . . . . . . . . . . . 6.3 Conclusion . . . . . . . . . . . . . . . . . 6.3.1 Guidelines for Use . . . . . . . . . 6.3.2 Recommendations for Further Work  . . . . .  . . . . . .  . . . . . . . . . .  . . . . . . . . . .  . . . . . . . . . .  . . . . . . . . . .  . . . . . . . . . .  . . . . . . . . . .  . . . . . . . . . .  . . . . . . . . . .  . . . . . . . . . .  . . . . . . . . . .  117 117 120 128 129 138 144 148 152 153  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  . . . . . .  157 157 159 163 164 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 4.1 4.2 4.3 4.4 4.5 4.6 4.7 4.8 4.9 4.10 4.11 4.12 4.13 4.13 5.1 5.2 5.3 5.4 5.5 5.6 5.7 5.8  Statics satisfied and interslice forces in various methods modified from (Krahn, 2003) . . . . . . . . . . . . . . . . . . . . . . . . . Slope geometries used in parametric study . . . . . . . . . . . . Discontinuous and slope parallel joint model parameters . . . . Daylighting joint model parameters . . . . . . . . . . . . . . . Number of cycles for 10−5 force ratio limit condition in 10 m discontinuous daylighting joint case . . . . . . . . . . . . . . . UDEC-SSR results for two different ratio limits . . . . . . . . . Summary of D-SSR FOS values for element size variation . . . Input parameters for Poisson’s Ratio sensitivity models . . . . . Peak-residual parameters used for strain-softening sensitivity . . SSR results for different daylighting joint spacings . . . . . . . Discontinuous joint generation settings . . . . . . . . . . . . . . Selected results from UDEC-SSR for slope parallel joint case . . Distance of critical sliding surface from crest for Slope parallel cases . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Summary of parametric study results . . . . . . . . . . . . . . . Summary of parametric study results (cont.) . . . . . . . . . . .  10  . 45 . 49 . 49 . . . . . . . .  51 56 59 71 80 85 86 91  . 92 . 95 . 96  Summary of geotechnical parameters used in the case history . . . Model and cycle stages of yielded elements . . . . . . . . . . . . Summary of LEM results, surfaces can be seen in figure 5.17 . . . Summary of unjointed SSR results . . . . . . . . . . . . . . . . . Summary of continuous bedding SSR results . . . . . . . . . . . Joint generation parameters for jointed case history model . . . . Summary of discontinuous SSR results with a wide range of joint networks, revealing small range in factor of safety values . . . . . Summary of continuous displacement-based SSR for jointed model cases . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .  105 109 118 121 129 131 132 152 vi  List of Figures 2.1  Continuum and discontinuum model examples . . . . . . . . . . .  Scale effects of rockmass failures and properties, modified from (Hoek, 2007) . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2 Stress strain relationships of varying rockmass quality, from Hoek (2007) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.3 Graphical representation of the Mohr-Coulomb failure envelope, modified from (Hudson and Harrison, 2005) . . . . . . . . . . . 3.4 Approximate M-C failure envelope from H-B parameters, modified from (Hoek et al., 2002) . . . . . . . . . . . . . . . . . . . 3.5 Disturbance factor guidelines for slopes, modified from (Hoek et al., 2002) . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.6 Computational cycle used in UDEC, from (Itasca, 2004c) . . . . 3.7 Difference between actual jointing and modelled joints . . . . . 3.8 Nomenclature for JSET joint generator in UDEC, from (Itasca, 2004a) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.9 Flowchart for UDEC-SSR showing stability criteria and bracketing search, modified from (Itasca, 2004c, 2008) . . . . . . . . . 3.10 Flowchart for general displacement-based SSR method . . . . . 3.11 Conceptual figure showing displacement response and selection of strength reduction factor indicating instability . . . . . . . . .  12  3.1  4.1 4.2  4.3 4.4  .  16  .  18  .  21  .  24  . 25 . 31 . 33 .  34  . 40 . 42 .  44  Example of scoping of joints (a) and meshing (b) within slope for computational efficiency . . . . . . . . . . . . . . . . . . . . . . 47 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 Results of displacement-based SSR for varying number of cycles . 52 Results of displacement-based SSR for varying force ratio limits per strength reduction step . . . . . . . . . . . . . . . . . . . . . 53 vii  List of Figures 4.5  4.6 4.7 4.8  4.9 4.10 4.11 4.12 4.13 4.14 4.15 4.16 4.17 4.18 4.19 4.20 4.21 4.22 4.23 4.24 4.24  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 . . . . . . . . . . Example of sawtooth horizontal displacement results during DSSR analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . Results of displacement based SSR for element size sensitivity . . 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) . . . . . . . . . . . . . . . . . . . . . . Comparison between cycled and force-equilibrium results for 7.5 m grid size . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . D-SSR results for element type sensitivity for discontinuous daylighting joints . . . . . . . . . . . . . . . . . . . . . . . . . . . . Failed volume indicators for 7.5 m triangular elements at SRF stage = 1.76 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Failed volume indicators for 7.5 m quadrilateral elements at SRF stage = 1.74 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Horizontal stress contours for varying modelled insitu stress ratios Vertical stress contours for varying modelled insitu stress ratios . . D-SSR results reduction results for varying horizontal to vertical stress ratios . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Results of displacement based SSR for Poisson’s Ratio sensitivity Results of D-SSR analysis for joint shear stiffness sensitivity . . . Number of cycles for 10−5 ratio limit for ks = 8 x 108 and 8 x 1010 Pa/m cases . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Results of displacement-based SSR for varying ratio limits when ks = 8 x 108 Pa/m . . . . . . . . . . . . . . . . . . . . . . . . . . Results of displacement-based SSR for varying ratio limits when ks = 8 x 1010 Pa/m . . . . . . . . . . . . . . . . . . . . . . . . . Comparison between active portions of the slope for different ratio limits when ks = 8 x 108 Pa/m . . . . . . . . . . . . . . . . . . . D-SSR for strain limits for strain-softening cases . . . . . . . . . D-SSR results for different daylighting joint spacings . . . . . . . Horizontal displacement and shear along joints for varying daylighting joint spacings (m) . . . . . . . . . . . . . . . . . . . . . Horizontal displacement and shear along joints for varying daylighting joint spacings (m) (cont.) . . . . . . . . . . . . . . . . .  55 57 59  61 62 64 65 66 68 69 70 72 74 75 76 77 78 80 82 83 84  viii  List of Figures 4.25 D-SSR results for discontinuous daylighting joint cases with varying 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 5.2 5.3 5.4 5.5 5.6 5.7 5.8 5.9 5.10 5.11 5.12  5.13  5.14 5.15 5.16 5.17  Cross-section of pre-failure pit, showing major geological units . . View of the south wall showing the failure, taken July 13th 2005. . Stereonet showing faults in the failure area . . . . . . . . . . . . . Stereonet showing bedding in the failure area . . . . . . . . . . . Stereonet showing joints in the failure area . . . . . . . . . . . . . Example of telescoping of joints (a) and meshing (b) within slope for the case history . . . . . . . . . . . . . . . . . . . . . . . . . Yielded elements at 10−6 ratio limit for UDEC-SSR giving FS = 1.12 (x yield in shear, o yield in tension) . . . . . . . . . . . . . . Yielded elements at 5,000 cycles for D-SSR at SRF = 1.11 (x yield in shear, o yield in tension) . . . . . . . . . . . . . . . . . . . . . Yielded elements at 10−6 ratio limit for D-SSR at SRF = 1.11 (x yield in shear, o yield in tension) . . . . . . . . . . . . . . . . . . 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) . . . . . . . . . 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) . . . . . . . . . Yielded elements at 10−6 ratio limit plus 30,000 cycles for DSSR at SRF = 1.11, showing the ultimate failure surface (x yield in shear, o yield in tension) . . . . . . . . . . . . . . . . . . . . . Yielded elements at 10−6 ratio limit plus 50,000 cycles for DSSR at SRF = 1.11, showing the ultimate failure surface (x yield in shear, o yield in tension) . . . . . . . . . . . . . . . . . . . . . Example of selecting sliding surface from UDEC velocity vector output . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Velocity plot derived surfaces from models states listed in table 5.2 Velocity plot derived surfaces from models states listed in table 5.2 with yielded elements overlain . . . . . . . . . . . . . . . . . Critical LEM failure surfaces predicted by SLIDE . . . . . . . . .  99 100 101 102 103 107 110 110 111 111 112  112  113 114 115 116 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 elements 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 elements 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 UDECSSR, that provide a good match to the actual failure surface . . . . 134 5.27 Velocity surfaces representing critical sliding surfaces, from UDECSSR, that provide a good match to LEM Surface Two . . . . . . . 135 5.28 Velocity surfaces representing critical sliding surfaces, from UDECSSR, 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 softening behaviour, at the point of instability as defined by D-SSR . . 140 5.32 Results of displacement-based SSR for models jointed strain softening models R1-R3 with low residual properties to capture fracture 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 parameters for joint network model J8 . . . . . . . . . . . . . . . . 5.37 Estimated static water table location, in blue, for case history . . . 5.38 Results of displacement-based SSR for models W8 and W12 with static water tables . . . . . . . . . . . . . . . . . . . . . . . . . . 5.39 Plasticity indicators showing clear shear band after removal of last bench in W12 model . . . . . . . . . . . . . . . . . . . . . . . . 5.40 Velocity plot showing slope failure along critical discontinuity behind the crest after removal of last bench in W8 model . . . . . . 5.41 Results of displacement-based SSR for models CD1–CD8 showing identical behaviour between continuously-cycled and initial state displacement-based SSR . . . . . . . . . . . . . . . . . . . 5.42 Photograph roughly perpendicular to failure slope showing presence of a steep, continuous structure . . . . . . . . . . . . . . . . 6.1  147 149 150 150 151  153 155  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 particular, needs to be thanked for supplying the case history data and model guidance. Many thanks to Doug Stead for time and effort on my behalf as my external examiner. 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 increase 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 establish 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 numerical 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 excellent reference against which these advanced numerical models can be tested. Slope stability is of prime concern to geological engineering from both an economic standpoint and with respect to the risk large slopes can pose to mine personnel 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 limitations 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 locations of interslice forces, and the inability to include deformation in the definition of failure. Complex rock deformation, including both intact rock failure and structurally controlled 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 proposed 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 constitutive 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 parametric 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 configurations 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, behaviour 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 common 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 conditions. 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 mechanism with time dependence. However, the sequence of events is particularly essential. 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 potential 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 mechanism are involved. Key information has been destroyed by the failure, resulting 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, analysis 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 measures (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 predicting the load, variations in material properties, uncertainties in material properties, and differences between modelled idealized behaviour and true behaviour (El6  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. Essentially, 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 geomaterials, it is possible to apply probability theory for some methods. Other techniques 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 confidence 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 complicated assumptions and methodologies whose implications are yet to be fully understood. Modeling of any system requires good boundary conditions, reasonable knowledge of material input parameters, and selection of proper constitutive 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 insitu measurements are highly variable, include some sort of testing bias, and are difficult to constrain. Some models include parameters that are difficult to measure 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 weaknesses: 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, prediction 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 complicated rockmass problems to be modeled. This has introduced the temptation to model very detailed rock behavior that may not be useful for general rock engineering 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 oversimplified. 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 important to understand how the fundamental behaviour of the modelling code and constitutive models of materials affects the single value of factor of safety returned. 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 factor 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  Σ(c βR + (N − uβ)R tan φ ) for moment equilibrium (2.1) Σ(W x − N f ) ± Dd) Σ(c β cos α + (N − uβ) tan φ cos α) F SF = for force equilibrium (2.2) Σ(N sin α) − D cos w F Sm =  where the slice base normal force is: N=  W + (Xr − Xl ) − cos α +  c β sin α+uβ sin α tan φ FS 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 equilibrium equilibrium normal (E) Ordinary Yes No No Bishop’s simplified Yes No Yes Janbu’s simplified No Yes Yes Spencer Yes Yes Yes Morgenstern-Price Yes Yes Yes Corps of Engineers -1 No Yes Yes  Interslice shear (X) No No No Yes Yes Yes  Corps of Engineers -2 No  Yes  Yes  Yes  Lowe-Karafiath  Yes  Yes  Yes  No  Inclination of X/E No force Horizontal Horizontal Constant Variable Inclination of a line from crest to toe Slice top ground surface inclination Average of ground surface slope and slice base inclination  It is important to note that the entire analysis is done using static forces on assumed 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 reasonable factors of safety.  10  Chapter 2. Background  2.2.2  Numerical Models  Numerical models are computational methods that attempt to capture the mechanical response (e.g. deformation) of a rockmass using a set of material parameters, 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 definitions 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 located 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  JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2)  UDEC (Version 4.00)  5.150  Chapter 2. Background LEGEND  5.050  28-Aug-08 14:25 cycle 0 time 0.000E+00 sec  4.950  zones in fdef blocks 4.850  4.750  4.650  4.550  4.450  4.350 (*10^2)  JOB TITLE : UDEC parametric studies - dip daylighting  UDEC (Version 4.00) University of Brish Columbia Vancouver, BC Canada LEGEND  5.150 4.250  4.000  (a) Continuum zoning4.400 example connected 4.200 4.300 4.500 with4.600 4.700 zones 4.800  4.100  4.900  (*10^2)  5.050  28-Aug-08 14:24 cycle 0 time 0.000E+00 sec  4.950  block plot 4.850  4.750  4.650  4.550  4.450  4.350  4.250  University of Brish Columbia Vancouver, BC Canada 4.000  (b) Discontinuum example showing4.600 discrete 4.200 4.300 4.400 4.500 4.700blocks 4.800  4.100  4.900  (*10^2)  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 arguably 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 models also rely on using ‘equivalent continuum’ properties for the material, such as those derived from rockmass characterization procedures, to account for the unmodelled discontinuities. There are no explicitly modeled discontinuities; however, 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 discrete 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 disturbance 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 motion 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 reduce the strength parameters of a material in a numerical model until some stabil13  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, maximum 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):  1 c SRF measured 1 φtrial = tan−1 tan φ SRF i τmax = ctrial + σn tan φtrial ctrial =  (2.4) (2.5) for the ith material  (2.6)  where the factor of safety is the SRF when the stability threshold has been overcome. This method can be applied to both continuum and discontinuum models. It has also been applied to more advanced constitutive models such as the HoekBrown 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 explicitly 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 combined 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 controls. The third category of failure mode exists between these two end members and includes stress controlled failure of intact rock bridges between joints to create 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 include 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 discontinuities. Structurally controlled failure is a bench-scale problem in the open pit environment with larger, deeper-seated failure involving intact rock only occurring 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 occurring 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 important 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 below 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:         1 −v −v 0 0 0 σxx  yy  −v 1 −v   0 0 0  σyy       zz  −v −v 1   σzz  1 0 0 0  =     (3.2) γxy  E  0   0 0 2(1 + v) 0 0      σxy  γyz  0   σyz  0 0 0 2(1 + v) 0 γzx 0 0 0 0 0 2(1 + v) σzx xx  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 material behaviour. Pande et al. (1990) lays out the four basic tenets of plastic theory: 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 coincide 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 criterion 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 representation 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 materials 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. However, 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 + σn tanφ  (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, modified 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. Further 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 + s σci  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 = mi exp  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 a=  GSI − 100 9 − 3D  −20 1 1 − GSI + (e 15 − e −3 ) 2 6  (3.9) (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: a  σc = σci s  (3.11)  Similarly, the tensile strength (σt ) can be found by setting σ1 = σ3 = σt in the same equation giving: sσci (3.12) σt = − mb 22  Chapter 3. Methodology Geotechnical software is generally written in terms of Mohr-Coulomb failure criterion, necessitating the need to be able to determine equivalent Mohr-Coulomb parameters from Hoek-Brown parameters. The process involves fitting an average 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  c =  6amb (s + mb σ3n )a−1 2(1 + a)(2 + a) + 6amb (s + mb σ3n )a−1  σ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.13)  (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 deformability. 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 structural 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 RMR 89 = GSI − 5  (3.15) (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 orientation. 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 rockmass behaves as an isotropic mass. The GSI system should not be applied to rockmasses with a clear dominant structural control. Also, the GSI system requires 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 structurally controlled. Conditions for movement along these discontinuities are controlled 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 asperities 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 dependent attribute. Also, testing of field scale joints (those controlling stability within the rockmass) is exceptionally difficult. True shear performance of the discontinuity 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 dilatency 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 discon28  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 discontinuities 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, however, 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 directions 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 methods 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 distribution 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 forcedisplacement 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 become 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, corners, 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 important 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 + 43 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  JOB TITLE : UDEC parametric studies - 10 m dip parallel  (*10^2)  UDEC (Version 4.00) 4.690  LEGEND 3-Sep-08 10:48 cycle 0 time 0.000E+00 sec  4.670  4.650  block plot  4.630  4.610  4.590  4.570  4.550  University of Brish Columbia JOB TITLE : UDEC parametric studies - 10 m dip parallel Vancouver, BC Canada  UDEC (Version 4.00)  4.530 (*10^2)  (a) Real joint network, with joints terminating within a block 4.250  4.270  4.290  4.310  4.330 (*10^2)  4.350  4.370  4.390  4.410  4.690  LEGEND 3-Sep-08 10:49 cycle 0 time 0.000E+00 sec  4.670  4.650  block plot  4.630  4.610  4.590  4.570  4.550  University of Brish Columbia Vancouver, BC Canada  4.530  (b) Same jointing as modelled in UDEC with incomplete joints joints removed 4.250  4.270  4.290  4.310  4.330 (*10^2)  4.350  4.370  4.390  4.410  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. Quadrilateral elements provide a better solution to plastic deformations (Itasca Inc., 2007). The block constitutive model, such as elastic or elasto-plastic, then controls 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 proportional to the mass velocity but in the opposite direction. This, typically acts on lower frequency vibrations and is described by the following equa34  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 α + βwn 2 wn  (3.21)  Alternatively, the λ value is adjusted as a fraction of the rate of energy displacement 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 constitutive models that include joint closure, dilation, continuously yielding materials, etc., only simple constitutive models were used that require simple strength parameters 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 − σ3 Nφ + 2c Nφ  (3.22)  with a tensile failure envelope: f t = σ t − σ3  (3.23)  1 + sinφ 1 − sinφ  (3.24)  c tanφ  (3.25)  where: Nφ = and: t σmax =  The non-associated flow rule uses a shear potential function defined as: g s = σ1 − σ3 Nψ  (3.26)  1 + sinψ 1 − sinψ  (3.27)  where: Nψ =  An associated flow rule for the tensile potential function such that: g t = −σ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 controlled 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)  ∆τs = Ks ∆ues  (3.30)  |τs | ≥ τmax  (3.31)  τs = sign(∆us )τmax  (3.32)  then  else if  then  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 displacement 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 tensile 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 friction angle selected depend on joint confinement which can vary widely within the rockmass, making the selection of a single value for these two parameters difficult. 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 discontinuum 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 thesis, 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 discontinuities) are all included in the SSR routine, and all factored equally. The built-in UDEC factor of safety calculator (UDEC-SSR) only allows for solving 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 reduction 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 brackets 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 system 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  Initial model geometry  Stage 2  Elastic gravity loading  Stage 3  Plastic properties Removal of benches Steps to test for stability or instability  SSR initial state  1. Do up to N steps and record unbalanced force ratio (Ru) Determine representative number of steps (N) for the characteristic response time of the system  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 CHECK  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)  Set F = (Fu + Fs) / 2  5. Set Ru = Ru_old and goto 1  CHECK  CHECK  If stable Fs = F If unstable Fu = F  EXIT IF Fu - Fs < 0.005  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 throughout 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 corresponding SRF) is chosen where the slope of the line rapidly changes from little displacement magnitudes to large displacements. This point is typically where localization 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 internal 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  Initial model geometry  Stage 2  Elastic gravity loading  Stage 3  Plastic properties Removal of benches  - 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  SSR initial state SRF_LOOP  For each SRF_LOOP:  SSR 1+1(x)  - 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  SRF_LOOP SSR 1+2(x) SRF_LOOP SSR 1+... SRF_LOOP SSR 1 + (n-1)(x) SRF_LOOP SSR 1 + (n)(x)  Generate displacement charts (Results)  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 conservative 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  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  Strength Reduction Factor  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 discontinuum 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 Daylighting joints Discontinuous daylighting joints Slope parallel joints  Dip Joint Set 1  Dip Joint Set 2  -35o -40o -45o  45o 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 failure 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 representative 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 properties (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 sensitivities involve material parameters while others are related to numerical modelling 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 general 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  JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2)  UDEC (Version 4.00)  Chapter 4. Parametric Study  7.000  LEGEND 23-May-08 17:07 cycle 71470 time 1.935E+01 sec  5.000  block plot  3.000  1.000  -1.000  JOB TITLE : UDEC parametric studies - dip daylighting  -3.000 (*10^2)  UDEC (Version 4.00) University of Brish Columbia Vancouver, BC Canada LEGEND  0.100  0.300  (a) Joints0.700  0.500  7.000 0.900  1.100  (*10^3)  23-May-08 17:07 cycle 71470 time 1.935E+01 sec  5.000  zones in fdef blocks  3.000  1.000  -1.000  -3.000  University of Brish Columbia Vancouver, BC Canada 0.100  0.300  (b) Mesh0.700  0.500  0.900  1.100  (*10^3)  Figure 4.1: Example of scoping of joints (a) and meshing (b) within slope for computational efficiency 47  Chapter 4. Parametric Study  JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2) 5.750  UDEC (Version 4.00) LEGEND  5.250  23-May-08 17:07 cycle 71470 time 1.935E+01 sec  4.750  block plot 4.250  3.750  3.250  2.750  2.250  JOB TITLE : UDEC parametric studies - dip daylighting University of Brish Columbia Vancouver, BC Canada  (*10^2) JOB TITLE : UDEC parametric studies - dip daylighting  UDEC (Version 4.00)  UDEC (Version 4.00) 3.000  3.500  4.000  4.500  5.000 (*10^2)  LEGEND  LEGEND 22-Apr-08 15:08 cycle 58750 time 4.229E+01 sec  block plot  block plot  (*10^2) 5.750  5.500  6.000  (a) Continuous daylighting  20-May-08 13:44 cycle 137100 time 6.002E+01 sec  University of Brish Columbia Vancouver, BC Canada  1.750  5.750 6.500  7.000  5.250  5.250  4.750  4.750  4.250  4.250  3.750  3.750  3.250  3.250  2.750  2.750  2.250  2.250  1.750  1.750  University of Brish Columbia Vancouver, BC Canada 3.000  3.500  4.000  4.500  5.000 (*10^2)  5.500  6.000  6.500  (b) Discontinuous daylighting  7.000  3.000  3.500  4.000  4.500  5.000 (*10^2)  5.500  (c) Slope parallel  6.000  6.500  7.000  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 Density E ν φrock crock trock φjoint cjoint tjoint Joint normal stiffness (kn ) Joint shear stiffness (ks )  Value  Unit  2600 7 x 109 0.33 33o 1.0 x 106 1.0 x 105 25o 1.0 x 105 1.0 x 103 1.0 x 1010 8 x 109  kg/m3 Pa Pa Pa Pa Pa Pa Pa  Table 4.3: Daylighting joint model parameters Parameter  Value  Density 2600 E 7 x 109 ν 0.33 φrock 40o crock 3.50 x 105 trock 1.0 x 105 φjoint 30o cjoint 1.0 x 105 tjoint 1.0 x 103 Joint normal stiffness (kn ) 10 x 109 Joint shear stiffness (ks ) 8 x 109  Unit kg/m3 Pa Pa Pa Pa Pa Pa 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 uncertainty 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 reduction 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 bracketing 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 discontinuous daylighting joint case SRF stage 1.62 1.64 1.66 1.68 1.70 1.72 1.74 1.76 1.78 1.80  Cycles for each stage 11470 13240 14970 15580 27150 25880 33120 32290 34000 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 dis51  Chapter 4. Parametric Study 18 00 18.00  Modellled horizontal displace ement (m)  16.00 14.00 12.00 5000  10.00  15000 30000  8.00  45000 60000  6.00  100000 4.00 2 00 2.00 0.00 1 60 1.60  1 65 1.65  1 70 1.70  1 75 1.75  1 80 1.80  1 85 1.85  1 90 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 develop 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 difference 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 displacements 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 important; rather the shift from low displacement to high displacement is important. 52  Chapter 4. Parametric Study  18.00  Modelled horizontal displacement (m)  16.00 14.00 12.00 10.00  Ratio = 1E‐3 Ratio = 1E‐4  8.00  Ratio = 1E‐5 Ratio = 1E‐6  6.00 4.00 2.00 0.00 1.44  1.54  1.64  1.74  1.84  1.94  Strength reduction factor  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. However, 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  JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2)  UDEC (Version 4.00)  5.500  Chapter 4. Parametric Study LEGEND 5.000  22-May-08 14:44 cycle 61030 time 3.690E+01 sec  4.500  velocity vectors maximum = 7.713E-02 4.000  0  2E -1  block plot 3.500  3.000  2.500  JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2)  UDEC (Version 4.00) University of Brish Columbia Vancouver, BC Canada LEGEND  2.000 5.500 3.500  4.000  4.500  (a)5.000 15,000 cycles 5.500  6.000  6.500  7.000  (*10^2)  23-May-08 20:44 cycle 146030 time 8.710E+01 sec  5.000  velocity vectors maximum = 8.271E-05  0  4.500  5E -4  4.000  block plot 3.500  3.000  2.500  University of Brish Columbia Vancouver, BC Canada  2.000 3.500  4.000  4.500  (b) 100,000 cycles 5.000 5.500  6.000  6.500  7.000  (*10^2)  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 investigate 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 appropriate 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 10 m spacing slope parallel 15 m spacing slope parallel  SRF (10−5 )  SRF(10−6 )  1.53 1.53  1.43 1.43  56  Chapter 4. Parametric Study  4.50  Modelled horizontal displacement (m)  4.00 3.50 3.00 2.50 2.00 1.50 1.00 0.50 0.00 1.20  1.25  1.30  1.35  1.40  1.45  1.50  1.55  1.60  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 computation 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. 16.00  Modelled horizontal displacement (m)  14.00 12.00 10.00 1 m 2.5 m  8.00  5 m 7.5 m  6.00  10 m 15 m  4.00 2.00 0.00 1.60  1.70  1.80  1.90  2.00  2.10  Strength reduction factor  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 1m 2.5 m 5m 7.5 m 10 m 15 m  D-SSR SRF NA 1.66 1.70 1.76 1.82 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 element 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 develop 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. Extrapolating 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 contributing 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  JOB TITLE : UDEC parametric studies - dip daylighting (*10^2)  JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2)  UDEC (Version 4.00)  UDEC (Version 4.00)  2.850  2.850  2.800  2.800  2.750  2.750  2.700  2.700  2.650  2.650  2.600  2.600  LEGEND  LEGEND 17-May-08 2:09 cycle 260380 time 3.391E+01 sec  5-May-08 11:53 cycle 61000 time 3.689E+01 sec  zones in fdef blocks no. zones : total 7046 at yield surface (*) 2488 yielded in past (X) 0 tensile failure (o) 4  zones in fdef blocks no. zones : total 338 at yield surface (*) 264 yielded in past (X) 0 tensile failure (o) 4  2.550  2.550  (b) Grid size = 5 m  (a) Grid size University = 1 mof Brish Columbia Vancouver, BC Canada  University of Brish Columbia Vancouver, BC Canada 6.050  6.100  6.150  6.200 (*10^2)  6.250  6.300  6.350  6.050  6.100  6.150  6.200 (*10^2)  6.250  6.300  6.350  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  12.00  Modelled horizontal displacement (m)  10.00  8.00  6.00 7.5 m force equilibrium 7.5 m cycled  4.00  2.00  0.00 1.60  1.70  1.80  1.90  2.00  2.10  Strength reduction factor  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 quadrilateral 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 recommended 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  16.00  Modelled horizonttal displacement (m)  14.00  12.00  10.00  5 m triangular  8.00  5 m quadrilateral 7.5 m triangular 6.00  7.5 m quadrilateral  4.00  2.00  0.00 1.60  1.65  1.70  1.75  1.80  1.85  1.90  Strength reduction factor  Figure 4.10: D-SSR results for element type sensitivity for discontinuous daylighting joints  64  Chapter 4. Parametric Study JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2)  UDEC (Version 4.00) 5.500  LEGEND 5.000  9-Jul-08 13:24 cycle 56290 time 5.200E+01 sec  4.500  X displacement contours contour interval= 2.000E-01 2.000E-01 to 1.400E+00  4.000  2.000E-01 4.000E-01 6.000E-01 8.000E-01 1.000E+00 1.200E+00 1.400E+00  3.500  3.000  block plot 2.500  2.000  JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2)  UDEC (Version 4.00) University of Brish Columbia Vancouver, BC Canada  1.500  (a) Horizontal displacement contours (m) 3.000  LEGEND  3.500  4.000  4.500  5.000 (*10^2)  5.500  6.000  5.500 6.500  7.000  5.000  9-Jul-08 13:24 cycle 56290 time 5.200E+01 sec  4.500  block plot no. zones : total 12490 at yield surface (*) 565 yielded in past (X) 0 tensile failure (o) 164  4.000  3.500  3.000  2.500  2.000  University of Brish Columbia Vancouver, BC Canada  1.500 3.500 4.000 4.500 5.000 5.500 7.000 (b)3.000 Plastic states of elements (x yield in6.000shear,6.500o yield in (*10^2) tension)  Figure 4.11: Failed volume indicators for 7.5 m triangular elements at SRF stage = 1.76  65  Chapter 4. Parametric Study JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2)  UDEC (Version 4.00) 5.500  LEGEND 5.000  9-Jul-08 13:25 cycle 58450 time 4.946E+01 sec  4.500  X displacement contours contour interval= 2.000E-01 2.000E-01 to 1.200E+00  4.000  2.000E-01 4.000E-01 6.000E-01 8.000E-01 1.000E+00 1.200E+00  3.500  3.000  block plot  2.500  2.000  JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2)  UDEC (Version 4.00) University of Brish Columbia Vancouver, BC Canada  1.500  (a) Horizontal displacement contours (m) 3.000  LEGEND  3.500  4.000  4.500  5.000 (*10^2)  5.500  6.000  5.500 6.500  7.000  5.000  9-Jul-08 13:25 cycle 58450 time 4.946E+01 sec  4.500  block plot no. zones : total 16912 at yield surface (*) 747 yielded in past (X) 0 tensile failure (o) 246  4.000  3.500  3.000  2.500  2.000  University of Brish Columbia Vancouver, BC Canada  1.500 3.500 4.000 4.500 5.000 5.500 7.000 (b)3.000 Plastic states of elements (x yield in6.000shear,6.500o yield in (*10^2) 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 important as the insitu stresses in slopes are commonly unknown and assumed, thus any differences in reported model results would make the selection of this parameter 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 18.00  Modelle ed  horizonttal displace ement (m)   16.00 14.00 12.00 .00 10.00 Point of failure at SRF= 1.7  8.00  k=0.5 k=1  k=2   6.00 4.00 2.00 0.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):  E 3(1 − 2ν) E G= 2(1 + ν)  K=  (4.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 failure 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 failure (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 ν 0.10 0.15 0.20 0.25 0.30  E (GPa) K (GPa) 7.00 7.00 7.00 7.00 7.00  2.50 2.86 3.33 4.00 5.00  G (GPa) 2.73 2.61 2.50 2.40 2.31  71  Chapter 4. Parametric Study 10 00 10.00  Modellled horizontal displace ement (m)  9.00 8.00 7.00 6.00 v = 0.10 5.00  v = 0.15 v = 0.20  4.00  v = 0.25  3 00 3.00  v = 0 30 v = 0.30  User selected point  selected point of failure  2.00 1.00 0.00 1 60 1.60  1 65 1.65  1 70 1.70  1 75 1.75  1 80 1.80  1 85 1.85  1 90 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 ratio 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 ratio 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 determined 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 18.00  Modellled horizontal displace ement (m)  16.00 14.00 12.00 10.00  8.00E+07 (Pa) 8.00E+08 (Pa)  8.00  8.00E+09 (Pa) 8.00E+10 (Pa)  6.00 4.00 2 00 2.00 0.00 1 48 1.48  1 58 1.58  1 68 1.68  1 78 1.78  1 88 1.88  1 98 1.98  Strength reduction factor  Figure 4.17: Results of D-SSR analysis for joint shear stiffness sensitivity  74  Chapter 4. Parametric Study  120000  100000  Numberr of cycles  80000  60000 8.00E+08 (Pa) 8.00E+10 (Pa)  40000  20000  0 1.50  1.55  1.60  1.65  1.70  1.75  1.80  1.85  Strength reduction factor  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  16.00  Modelled horizontal displacement (m)  14.00 12.00 10.00 8.00 Ratio 1E‐5 Ratio 1E‐6  6.00 4.00 2.00 0.00 1.48  1.58  1.68  1.78  1.88  1.98  Strength reduction factor  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 18.00  Modellled horizontal displace ement (m)  16.00 14.00 12.00 10.00 Ratio 1E‐5  8.00  Ratio 1E‐6 6.00 4.00 2 00 2.00 0.00 1 48 1.48  1 58 1.58  1 68 1.68  1 78 1.78  1 88 1.88  1 98 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  JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2) 6.000  UDEC (Version 4.00)  Chapter 4. Parametric Study 5.500  LEGEND 26-May-08 16:31 cycle 86660 time 5.629E+01 sec  5.000  velocity vectors maximum = 3.339E-01  0  4.500  4.000  2E 0  block plot 3.500  3.000  2.500  2.000  JOB TITLE : UDEC parametric studies - dip daylighting University of Brish Columbia UDEC (Version 4.00) Vancouver, BC Canada LEGEND  6.000 1.500 3.500showing 4.000 4.500 portions 5.000 6.000 (a) 3.000 Velocity active of5.500the slope for6.500ratio 7.000 limit of (*10^2) 10−5 at failure SRF of 1.8  30-Jun-08 15:43 cycle 109680 time 7.099E+01 sec  5.500  5.000  velocity vectors maximum = 3.990E-04  0  (*10^2)  4.500  4.000  2E -3  block plot 3.500  3.000  2.500  2.000  University of Brish Columbia Vancouver, BC Canada  1.500 3.500 showing 4.000 4.500 5.000 5.500the slope 6.000 (b)3.000 Velocity active portions of for6.500ratio 7.000 limit of (*10^2) 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) investigates this briefly and comes to the conclusion that an elastic-perfectly plastic model best agrees with the LEM analysis. However, discontinuum models are typically 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 applied 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 φr ock φj oint cr ock cj oint  Peak  Residual  33 25 1.01 x 106 Pa 1 x 105 Pa  28 5 x 105 Pa -  16.00  Modelled horizonttal displacement (m)  14.00  12.00  10.00  100% 10%  8.00  5% 1% 6.00  0.1% Perfectly plastic  4.00  2.00  0.00 1.10  1.20  1.30  1.40  1.50  1.60  1.70  1.80  1.90  2.00  Strength reduction factor  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 typical 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 heavily 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 displacement 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  Modellled horizontal displace ement (m)  6.00  5.00  4.00 2.5 m 5 m  3.00  10 m 15 m  2.00  20 m 20 m  1.00  0.00 1 00 1.00  1 05 1.05  1 10 1.10  1 15 1.15  1 20 1.20  1 25 1.25  1 30 1.30  1 35 1.35  1 40 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 necessarily 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, regardless 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 selected 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 daylighting 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 daylighting joint spacings (m) (cont.) returned).  84  Chapter 4. Parametric Study  Table 4.9: SSR results for different daylighting joint spacings Spacing (m) 2.5 5.0 10 15 20  4.8.2  UDEC-SSR SRF  Displacement-based SSR SRF  1.06 1.10 1.08 1.09 1.12  1.08 1.10 1.12 - 1.20 1.14 - 1.26 1.14 - 1.26  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 represents 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 in85  Chapter 4. Parametric Study  Table 4.10: Discontinuous joint generation settings Spacing (m)  7 7 8 8 9 9 10 10 11 11 12 12 13 13  Angle  trace length(m) µ, σ  gap length (m) µ, σ  -40 45 -40 45 -40 45 -40 45 -40 45 -40 45 -40 45  28, 0 28, 0 32, 0 32, 0 36,0 36, 0 40, 0 40, 0 44, 0 44, 0 48, 0 48, 0 52, 0 52, 0  7, 7 0, 0 8, 8 0, 0 9, 9 0, 0 10, 10 0, 0 11, 11 0, 0 12, 12 0, 0 13, 13 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 complexity 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 calculated for a forward analysis, due to the limited sensitivity.  87  Chapter 4. Parametric Study  2 00 2.00  Modellled horizontal displace ement (m)  1.80 1.60 1.40 7m  1.20  8m 1.00  9m 10m  0.80  11m 12m  0 60 0.60  13m 0.40 0.20 0.00 1 55 1.55  1 60 1.60  1 65 1.65  1 70 1.70  1 75 1.75  1 80 1.80  1 85 1.85  Strength reduction factor  Figure 4.25: D-SSR results for discontinuous daylighting joint cases with varying joint network properties  88  JOB TITLE : UDEC parametric studies - dip daylighting  (*10^2)  UDEC (Version 4.00)  Chapter 4. Parametric Study 3.500  LEGEND 9-Apr-08 16:21 cycle 45130 time 2.868E+01 sec  3.300  block plot 3.100  2.900  2.700  2.500  JOB TITLE :  (*10^2)  UDEC (Version 4.00) University of Brish Columbia Vancouver, BC Canada LEGEND  2.300 3.300 5.700  5.900  6.100  (a)6.300 7 m spacing 6.500  6.700  6.900  7.100  (*10^2)  13-May-08 10:40 cycle 45560 time 3.110E+01 sec  3.100  block plot 2.900  2.700  2.500  2.300  University of Brish Columbia Vancouver, BC Canada 5.700  5.900  6.100  (b)6.300 8 m spacing 6.500  6.700  6.900  7.100  (*10^2)  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, layered sedimentary sequences with a continuous, slope parallel joint set and a backdipping 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 shearing 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. 5.00  Modelled horizontal displacement (m)  4.50 4.00 3.50  5 m 6 m  3.00  7 m 2.50  10 m 11 m  2.00  12 m 13 m  1.50  15 m 1.00  20 m  0.50 0.00 1.30  1.35  1.40  1.45  1.50  Strength reduction factor  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) 5 10 15 20  SRF(10−5 )  SRF(10−6 )  1.53 1.53 1.53 1.53  1.43 1.43 1.43 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) 5 6 7 10 11 12 13 15 20  Distance from crest (m) 55 53 52 56 47 50 50 63 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 parameters, 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 significant 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. Table 4.13 below summarizes the results of the sensitivity analysis with comments on their effect on the SSR method.  94  SRF Sensitivity  Large  Large  Negligible  Negligible  Negligible  Setting / Parameter  Number of cycles or ratio limit used for each stage  Element Size  Element Type  Insitu Stress  Elastic Properties  Negligible  Negligible  Small  Small  Small  Effect on Position of Failure surface The higher the number of cycles or lower the ratio limit used, the more conservative the factor of safety returned. Failure surfaces are occasionally more shallow the more cycles that occur. The element size has a major influence on the factor of safety. Extremely small element sizes are to be avoided. Element size affects failure surface position due to influence on localization and shear banding. Quadrilateral elements are reportedly better for plastic materials; however most problem domains will require both quadrilateral and triangular element types due to limitations in spatial discretization. Absolute displacement values are affected but the onset of instability occurs at the same SRF stage. Absolute displacement values are affected but the onset of instability occurs at the same SRF stage within a small range. Literature cautions against using materials with drastically different stiffness properties when comparing to LEM solutions.  Comments  Table 4.13: Summary of parametric study results  Chapter 4. Parametric Study  95  SRF Sensitivity  Variable  Large  Variable  Setting / Parameter  Joint Shear Stiffness  Strain-softening models (intact blocks)  Joint Network  Small  Small  Small  Effect on Position of Failure surface Stiffness affects the factor of safety , however, there is no systematic relationship between the two. This property has a large affect on the time required for a model to come to equilibrium and occasionally shows different behaviour depending on the force ratio limit used for each SRF step. The transition to residual strength parameters results in a lower factor of safety overall. Lower strain values before onset of full residual strengths result in a lower factor of safety. There is a slight effect on position of failure surfaces with slightly shallower failure surfaces occurring at lower strain limit values. For any continuous joints, the fracture density has an almost negligible affect on the factor of safety and position of critical sliding surface. For discontinuous joints, fracture density controls the size and distribution of rock bridges within the model and can have a large affect on the factor of safety returned. Heavily fractured rockmasses are not necessarily less stable than less fractured ones due to distribution of rock bridges and limitations of jointed rock within UDEC.  Comments  Table 4.13: Summary of parametric study results (cont.)  Chapter 4. Parametric Study  96  Chapter 5 Case History 5.1  Introduction  The case history presented will remain unnamed due to a confidentiality agreement, 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 allows 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 thickness. 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 structurally 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 personnel 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 oriented 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 various 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  file: h:/project/89-148a/2005/may05/struct/Dips-Combined/BUZZ-F.grf  Chapter 5. Case History  EQUAL AREA LOWER HEMISPHERE  842 842  N  Poles Plotted Data Entries  FD1 FA2  FA1 FB1 FA3  E  W FA4 FB2 FC1 FC2  MAJOR PLANES  CONTOUR LEGEND SCHMIDT POLE CONCENTRATIONS % of total per 1.0 % area Minimum Contour Contour Interval Max.Concentration  ORIENTATIONS  S  SET DIP/DIR. %  = 0.5 = 0.5 = 10.5  Figure 5.3: Stereonet showing faults in the  FA1 66/245 FA2 67/231 FA3 64/260 FA4 64/280 FC1 55/296 FC2 79/299 FD1 77/137 FB1 64/108 59/072 failureFB2area  #  10.5 6 6 4.5 3.5 2 2 2 1.5  88 51 51 38 29 17 17 17 13  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.  BARRICK GOLDSTRIKE MINES INC. OPEN PIT GEOTECHNICAL ASSESSMENTS 11TH WEST LAYBACK  STRUCTURAL ANALYSIS OF FAULTS IN BUZZARD DOMAIN (RODEO CREEK FORMATION)  PITEAU ASSOCIATES GEOTECHNICAL & HYDROGEOLOGICAL CONSULTANTS  LIMA CALGARY  VANCOUVER BY:  DATE:  IRS APPROVED:  JUL 05 DWG:  101  B20  file: h:/project/89-148a/2005/may05/struct/Dips-Combined/BUZZ-B.grf  Chapter 5. Case History  EQUAL AREA LOWER HEMISPHERE  512 512  N  Poles Plotted Data Entries  E  W BG2 BG1 BG3  MAJOR PLANES  CONTOUR LEGEND SCHMIDT POLE CONCENTRATIONS % of total per 1.0 % area Minimum Contour Contour Interval Max.Concentration  S  ORIENTATIONS SET DIP/DIR. %  BG1 15/353 BG2 19/058 BG3 23/324  = 2 = 2 = 23.2  #  23.2 119 12 61 61 12  Figure 5.4: Stereonet showing bedding in the failure area 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.  BARRICK GOLDSTRIKE MINES INC. OPEN PIT GEOTECHNICAL ASSESSMENTS 10TH WEST LAYBACK  STRUCTURAL ANALYSIS OF BEDDING IN BUZZARD DOMAIN (RODEO CREEK FORMATION)  PITEAU ASSOCIATES GEOTECHNICAL & HYDROGEOLOGICAL CONSULTANTS  LIMA CALGARY  VANCOUVER BY:  DATE:  IRS APPROVED:  JULY 05 FIG: DWG:  102  B22  file: h:/project/89-148a/2005/may05/struct/Dips-Combined/BUZZ-J.grf  Chapter 5. Case History  EQUAL AREA LOWER HEMISPHERE  485 485  N  Poles Plotted Data Entries  JD3  JD2  JA2 JD1 JA1 JA3  E  W  JB2  JC1  JC3 JB1 JC2  MAJOR PLANES ORIENTATIONS SET DIP/DIR. %  CONTOUR LEGEND  S  SCHMIDT POLE CONCENTRATIONS % of total per 1.0 % area Minimum Contour Contour Interval Max.Concentration  = 0.5 = 0.5 = 12.4  Figure 5.5: Stereonet showing joints in the 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.  BARRICK GOLDSTRIKE MINES INC. OPEN PIT GEOTECHNICAL ASSESSMENTS 11TH WEST LAYBACK  STRUCTURAL ANALYSIS OF JOINTS IN BUZZARD DOMAIN (RODEO CREEK FORMATION)  JA1 75/240 JD1 74/135 JA2 73/223 JA3 77/248 JD2 78/156 JD3 80/173 JC1 76/298 JB1 84/020 JC2 84/338 75/326 failureJC3 area JB2 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  PITEAU ASSOCIATES GEOTECHNICAL & HYDROGEOLOGICAL CONSULTANTS  LIMA CALGARY  VANCOUVER BY:  DATE:  IRS APPROVED:  JUL 05 DWG:  103  B21  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 structural 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 provided 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. Table 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 disturbance 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 Density (kg/m3 )  ν  Upper unit  2400  0.112  2.39  Upper Limestone  2560  0.276  Contact Zone  2240  Waste  Unit  Joints & Bedding Planes  E UCS (GPa) (MPa)  RMR  mi  mb  s  a  Cohesion (kPa)  φ  24  46  12  0.255  1.25 x 10−4  0.507  185  29.7  8.3  8.00  38  50  17  0.499  2.66 x 10−4  0.506  278  39.3  11.8  -  -  -  -  -  -  -  -  48  25  8.3  1922  0.2  3  -  -  -  -  -  -  0  39  0  -  -  -  -  -  -  -  -  -  0 - 70  23-26  0  Tensile Strength (kPa)  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  JOB TITLE :  (*10^2)  UDEC (Version 4.00) 5.500  Chapter 5. Case History LEGEND 4.500  7-Aug-08 11:28 cycle 86920 time 4.387E+01 sec  3.500  block plot  2.500  1.500  0.500  -0.500  JOB TITLE :  (*10^2) -1.500  UDEC (Version 4.00) University of Brish Columbia Vancouver, BC Canada LEGEND  5.500  0.500  1.500  2.500  (a) Joints 4.500  3.500  5.500  6.500  7.500  (*10^2) 4.500  7-Aug-08 11:28 cycle 86920 time 4.387E+01 sec  3.500  zones in fdef blocks  2.500  1.500  0.500  -0.500  -1.500  University of Brish Columbia Vancouver, BC Canada 0.500  1.500  2.500  (b) Mesh 4.500  3.500  5.500  6.500  7.500  (*10^2)  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 vectors. 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 common 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 elements, 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  UDEC-SSR D-SSR D-SSR D-SSR D-SSR D-SSR D-SSR  1.12 1.11 1.11 1.11 1.11 1.11 1.11  10−6 limit 5000 cycles 10−6 ratio limit 10−6 ratio limit plus 5,000 cycles 10−6 ratio limit plus 15,000 cycles 10−6 ratio limit plus 30,000 cycles 10−6 ratio limit plus 50,000 cycles  Figure 5.7 5.8 5.9 5.10 5.11 5.12 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  JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.250  LEGEND 6-Aug-08 11:07 cycle 86278 time 5.375E+01 sec  3.750  block plot no. zones : total 13530 at yield surface (*) 2984 yielded in past (X) 0 tensile failure (o) 795  3.250  2.750  2.250  1.750  University of Brish Columbia Vancouver, BC Canada 2.500  3.000  3.500 (*10^2)  4.000  4.500  5.000  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) JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.250  LEGEND 7-Aug-08 13:23 cycle 87062 time 5.424E+01 sec  3.750  block plot no. zones : total 13530 at yield surface (*) 3101 yielded in past (X) 0 tensile failure (o) 843  3.250  2.750  2.250  1.750  University of Brish Columbia Vancouver, BC Canada 2.500  3.000  3.500 (*10^2)  4.000  4.500  5.000  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  JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.250  LEGEND 7-Aug-08 13:03 cycle 113880 time 7.121E+01 sec  3.750  block plot no. zones : total 13530 at yield surface (*) 925 yielded in past (X) 0 tensile failure (o) 140  3.250  2.750  2.250  1.750  University of Brish Columbia Vancouver, BC Canada 2.500  3.000  3.500 (*10^2)  4.000  4.500  5.000  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) JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.250  LEGEND 7-Aug-08 16:16 cycle 118880 time 7.437E+01 sec  3.750  block plot no. zones : total 13530 at yield surface (*) 2484 yielded in past (X) 0 tensile failure (o) 724  3.250  2.750  2.250  1.750  University of Brish Columbia Vancouver, BC Canada 2.500  3.000  3.500 (*10^2)  4.000  4.500  5.000  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  JOB TITLE :  (*10^2)  Chapter 5. Case History  UDEC (Version 4.00) 4.250  LEGEND 7-Aug-08 12:34 cycle 128880 time 8.070E+01 sec  3.750  block plot no. zones : total 13530 at yield surface (*) 2544 yielded in past (X) 0 tensile failure (o) 609  3.250  2.750  2.250  1.750  University of Brish Columbia Vancouver, BC Canada 2.500  3.000  3.500 (*10^2)  4.000  4.500  5.000  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) JOB TITLE : (*10^2)  UDEC (Version 4.00) 4.250  LEGEND 7-Aug-08 12:27 cycle 143880 time 9.018E+01 sec  3.750  block plot no. zones : total 13530 at yield surface (*) 1615 yielded in past (X) 0 tensile failure (o) 209  3.250  2.750  2.250  1.750  University of Brish Columbia Vancouver, BC Canada 2.500  3.000  3.500 (*10^2)  4.000  4.500  5.000  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  JOB TITLE :  (*10^2)  UDEC (Version 4.00)  Chapter 5. Case History  4.250  LEGEND 7-Aug-08 13:48 cycle 163880 time 1.028E+02 sec  3.750  block plot no. zones : total 13530 at yield surface (*) 1643 yielded in past (X) 0 tensile failure (o) 221  3.250  2.750  2.250  1.750  University of Brish Columbia Vancouver, BC Canada 2.500  3.000  3.500 (*10^2)  4.000  4.500  5.000  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 (figures 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 failure 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 different plots presented. It is assumed that these SRF values are for the critical SRF stage (i.e. the factor of safety), unless otherwise noted. JOB TITLE :  (*10^2)  UDEC (Version 4.00)  Slip surface defined by velocity vectors representing active portion of the slope  LEGEND 6-Aug-08 11:07 cycle 86278 time 5.375E+01 sec  4.250  3.750  velocity vectors maximum = 9.559E-02 3.250  0  5E -1  block plot 2.750  2.250  Figure 5.14: Example of selecting sliding surface from UDEC velocity vector output 1.750  University of Brish Columbia Vancouver, BC Canada 2.500  3.000  3.500 (*10^2)  4.000  4.500  5.000  114  Chapter 5. Case History  UDEC-SSR [1E-6] D-SSR SRF = 1.11 [1E-6 + 15,000 cycles]  375  D-SSR SRF = 1.11 [1E-6 + 30,000 cycles] D-SSR SRF = 1.11 [1E-6 + 50,000 cycles] 325  D-SSR SRF=1.11 [5,000 cycles]  meters  275  225  150 175  200  250  300  350  400  450  500  550  meters  125  Figure 5.15: Velocity plot derived surfaces from models states listed in table 5.2  115  Chapter 5. Case History  375  Surfaces from velocity  Intermediate shear surface [at SRF 1.11 (1E‐6)]  325  meters  275  225  150 175  Ultimate shear surface [at SRF 1.11 (1E‐6 + 15,000 ‐ 50,000 cycles)]  200  250  300  350  400  450  500  550  meters  125  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 factor 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. Physically 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 regardless 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 Identifier A A B B C C D D  Yes Yes Yes Yes No No No No  Water Table  Method  FOS  Surface  Yes Yes No No Yes Yes No No  Janbu Corrected GLE Janbu Corrected GLE Janbu Corrected GLE Janbu Corrected GLE  0.99 0.98 1.20 1.18 1.00 0.98 1.20 1.20  1 2 1 2 2 2 2 2  118  Chapter 5. Case History  Approximate failure surface  based on field observations and  geologic projects (exact depth  unknown) referred to as the  failure surface  375  325  meters  275  225  150 175  LEM Surface One: ‐ Tension crack ‐ Dry or Ru = 0.15 ‐ Janbu Corrected LEM Surface Two: ‐ Tension crack ‐ GLE/MP ‐ Dry or Ru = 0.15 OR ‐ No tension crack ‐ Janbu Corrected, GLE/MP ‐ Dry or Ru = 0.15 200  250  300  350  400  450  500  meters  125  Figure 5.17: Critical LEM failure surfaces predicted by SLIDE  119  550  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  JOB TITLE :  (*10^2)  Chapter 5. Case History  UDEC (Version 4.00)  5.500  LEGEND 28-Jul-08 5:59 cycle 77440 time 5.696E+01 sec  4.500  3.500  block plot  2.500  1.500  0.500  -0.500  -1.500  University of Brish Columbia Vancouver, BC Canada 0.500  1.500  2.500  3.500  4.500  5.500  6.500  7.500  (*10^2)geometry after bench removal Figure 5.18: Unjointed model  Table 5.4: Summary of unjointed SSR results Model  U1 U2 U3 U4 U5  Element σh : σv size (m) 2.5 5.0 10.0 5.0 5.0  2.0 2.0 2.0 0.5 2.0  Bulk/Shear Modulus  D-SSR  UDEC-SSR 10−5  UDEC-SSR 10−6  Variable Variable Variable Variable Consistent  1.12–1.14 1.13–1.15 1.17 1.13–1.15 1.13–1.15  1.13 1.17 1.42 1.16 1.16  1.12 1.14 1.21 1.14 1.16  121  Chapter 5. Case History  8 00 8.00  Mod delled horizon ntal displacement (m)  7.00  6.00  5.00  U1 4.00  U2 U3  3.00  U4 U5  2.00  1.00  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 displacement 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 surface 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 mechanism 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 element 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. LEM Surface Two 375  Actual failure surface  325  LEM Surface One Internal shearing of  failed mass  meters  275  Yielded elements (U1) SRF = 1 13 SRF = 1.13  225  150 175  200  250  300  350  400  450  500  550  meters  125  Figure 5.20: Comparison of LEM surfaces and D-SSR SRF = 1.13 yielded elements for model U1  124  Chapter 5. Case History  LEM Surface Two 375  Actual failure surface  325  LEM Surface One  UDEC‐SSR 1E‐5 (U1) (velocity) 275  meters  UDEC SSR 1E UDEC‐SSR 1E‐6 6 (U1) (U1) (velocity)  225  150 175  200  250  300  350  400  450  500  550  meters  125  Figure 5.21: Comparison of LEM surfaces and UDEC-SRF velocity-derived surfaces for model U1  125  Chapter 5. Case History  LEM Surface Two 375  Actual failure surface LEM Surface One  325  meters  275  Internal shearing of  failed mass  Yielded elements (U2) SRF = 1 14 SRF = 1.14  225  150 175  Velocity (U2) SRF = 1.14 200  250  300  350  400  450  500  550  meters  125  Figure 5.22: Comparison of LEM surfaces and D-SSR SRF = 1.14 yielded elements for model U2  126  Chapter 5. Case History  LEM Surface Two 375  Actual failure surface  325  LEM Surface One  UDEC‐SSR 1E‐5 (U2) (velocity) 275  meters  UDEC SSR 1E UDEC‐SSR 1E‐6 6 (U2) (U2) (velocity)  225  150 175  200  250  300  350  400  450  500  550  meters  125  Figure 5.23: Comparison of LEM surfaces and UDEC-SRF velocity-derived surfaces 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 returned 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  C1 C2 C3 C4  Joint spacing (m)  Dip towards pit  D-SSR  UDEC-SSR 10−6  5 10 10 20  15o 15o 20o 15o  1.12 1.13 1.13 1.11 – 1.12  1.13 1.13 1.13 1.12  JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.400  LEGEND 29-Jul-08 17:24 cycle 99202 time 6.048E+01 sec  4.000  block plot shear displacement on joint max shear disp = 7.825E-02 each line thick = 2.000E-03  3.600  3.200  2.800  2.400  University of Brish Columbia Vancouver, BC Canada 2.800  3.200  3.600 (*10^2)  4.000  4.400  4.800  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 representing 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 UDECSSR 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 failure 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 J1 J2 J3  J4  J5  J6  J7  J8  J9  J10  J11  J12  J13  Joint dip m std -15 0 -80 7 -15 0 -74 7 -15 0 75 5 -80 7 -20 0 75 7 -80 7 -20 0 75 7 -80 7 -15 0 80 5 -80 5 -18 3 72 13 -71 11 -15 0 75 7 -75 6 -15 0 75 7 -75 6 -16 0 70 5 -78 6 -15 0 75 5 -75 6 -15 0 75 5 -77 7 -80 3 -20 0 75 7 -80 6  Trace length (m) m std 100 0 30 11 100 0 27 11 100 0 17 10 23 11 120 0 30 5 23 7 120 0 8 5 8 7 100 0 25 10 20 10 80 30 12 5 8 7 100 0 35 10 40 14 100 0 35 10 28 14 60 20 25 5 25 38 100 0 35 10 30 10 100 0 10 5 4 5 6 3 120 0 8 5 8 9  Gap length (m) Spacing (m) m std m std 0 0 10 0 12 15 10 5 0 0 10 0 12 15 10 5 0 0 5 0 15 15 5 0 12 15 5 0 0.5 1 15 1 15 10 12 3 12 15 13 5 0.5 1 3 1 15 10 7 0 12 15 5 0 0 0 10 0 15 15 10 4 10 15 10 5 0.5 1 3 1 15 10 7 1 12 15 5 0 0 0 7.5 0 15 15 10 4 10 15 14 5 0 0 4 0 15 15 10 4 5 15 17 5 0.5 0.5 3.5 1.5 5 15 15 4 10 10 14.5 4 0 0 7.5 0 15 15 10 5 10 6 14 5 0 0 2.75 0.5 5 0 10 0 4 7 4 5 3 5 5 3 0.5 1 3 1 15 10 7 0 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 networks, revealing small range in factor of safety values Model J1 J2 J3 J4 J5 J6 J7 J8 J9 J10 J11 J12 J13  D-SSR  UDEC-SSR  1.12 1.14 1.09 1.12 1.13 1.10 1.11 1.06 1.07 1.07 1.08 1.10 1.10  1.13 1.14 1.09 1.11 1.13 1.10 1.10 1.06 1.08 1.08 1.07 1.10 1.11  Failure surface location Matches LEM Surface Two Deeper seated failure Deeper seated failure Deeper seated failure Deeper seated failure Matches LEM Surface Two Matches LEM Surface Two Good match to actual failure surface Deeper seated failure Deeper seated failure Fair match to actual failure surface Matches LEM Surface Two Deeper seated failure  132  Chapter 5. Case History  9.00  Modelled horizontal displacement (m)  8.00 Moderate joint influence on  stability  7.00  J1 J2  6.00  J3 J4  5.00  J5  Large joint influence on  stability stability   4.00  J6 J7 J8  3.00  J9 J10  2.00  J11 J12  1.00 Small joint influence on  stability  0.00 1.00  1.05  1.10  1.15  1.20  J13  1.25  Strength reduction factor  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  LEM Surface Two  375  Actual failure surface  LEM Surface One  325  J8 SRF = 1.06 J11 SRF = 1.07  meters  275  225  150 175  200  250  300  350  400  450  500  550  meters  125  Figure 5.26: Velocity surfaces representing critical sliding surfaces, from UDECSSR, that provide a good match to the actual failure surface  134  Chapter 5. Case History  LEM Surface Two  375  Actual failure surface  LEM Surface One  325  J1 SRF = 1.12 J7 SRF = 1.11 J12 SRF = 1.10  275  meters  J6 SRF = 1.10  225  150 175  200  250  300  350  400  450  500  550  meters  125  Figure 5.27: Velocity surfaces representing critical sliding surfaces, from UDECSSR, that provide a good match to LEM Surface Two  135  Chapter 5. Case History  LEM Surface Two 375  Actual failure surface  J2 SRF = 1.14 LEM Surface One  325  J3 SRF = 1.09 J4 SRF = 1.12 J5 SRF = 1.13 J9 SRF = 1.07  275  meters  J10 SRF = 1.07 J13 SRF = 1.10 225  150 175  200  250  300  350  400  450  500  550  meters  125  Figure 5.28: Velocity surfaces representing critical sliding surfaces, from UDECSSR, 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 ultimate 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 JOB TITLE the: position of the final failure surface to make it shallower. This includes usUDECing (Version 4.00) material properties in the Upper unit and moving to a strain softening weaker model, in addition to including the effects of water. (*10^2)  4.250  LEGEND  7-Aug-08 13:03 cycle 113880 time 7.121E+01 sec  3.750  block plot no. zones : total 13530 at yield surface (*) 925 yielded in past (X) 0 tensile failure (o) 140  3.250  2.750  2.250  1.750  University of Brish Columbia Vancouver, BC Canada 2.500  3.000  3.500 (*10^2)  4.000  4.500  5.000  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 surface 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 8.00  Mod delled horizontal displace ement (m)  7.00  6.00  5.00  4.00  J7 strain softening SRF = 1.09 J10 strain softening SRF = 1.05 J11 strain softening SRF = 1.10  3.00  J12 strain softening SRF = 1.08 2.00  1.00  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  Actual failure surface  375  J12 strain softening [SRF = 1.08] J12 perfectly plastic [SRF = 1.10]  325  meters  275  225  150 175  200  250  300  350  400  450  500  550  meters  125  Figure 5.31: Velocity surfaces representing critical sliding surfaces for the J12 joint network for elastic perfectly-plastic behaviour and strain softening behaviour, 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 parameters 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 velocityderived 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, however, 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 behaviour, 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 5.00 4.50  Mod delled horizontal displace ement (m)  4.00 3.50 3.00 2.50  R1 SRF = 1.08 R2 SRF = 1.07  2.00  R3 SRF = 1.05 1.50 1.00 0 50 0.50 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 softening models R1-R3 with low residual properties to capture fracture generation within a continuum  142  Chapter 5. Case History  Actual failure surface  375  R1 SRF = 1.08 R2 SRF = 1.07  325  R3 SRF = 1.05  meters  275  225  150 175  200  250  300  350  400  450  500  550  meters  125  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 predicted failure surface. Because all factor of safety calculations depend on initial values to compare against, modifications to the initial parameters change the factor 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 surface 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 ks taken from the literature such as Kulhawy (1975). Joint strength properties were varied for two cases, also for the J8 fracture network: • Case one: φjoint = 27, cjoint = 53 kPa  144  Chapter 5. Case History  J12Base case FS = 1.10  Actual failure surface  375  M1 10% lower initial upper unit parameters FS = 1.05 M2 10% lower initial joint parameters FS = 1.08 M3 10% lower initial contact zone parameters FS = 1.04  325  Failure surface position  insensitive to joint and upper  unit strengths  meters  275  Failure surface position sensitive to  contact zone strengths 225  150 175  200  250  300  350  400  450  500  550  meters  125  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 surface. 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 predicted. 3.00  Modelled horizontal displacement (m)  2.50  2.00  All models fail at the same  critical SRF stage despite  differences in absolute values of  displacement  1.50  8E8 (Pa) 8E8 (Pa) 8E9 (Pa) 8E10 (Pa)  1.00  0.50  0.00 1.00  1.02  1.04  1.06  1.08  1.10  1.12  Strength reduction factor  Figure 5.35: Results of displacement-based SSR for models ks sensitivity for joint network model J8  146  Chapter 5. Case History  375  Actual failure surface J8 Base SRF = 1.06 325  J8 Case one SRF = 1.07 J8 Case two SRF = 1.08  meters  275  225  150 175  200  250  300  350  400  450  500  550  meters  125  Figure 5.36: Variations in failure surface location with initial joint strength parameters 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 water 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 8.00  Modelled horizontal displacement (m)  7.00  Water table models W8 and  W12 have factors of safety  at, or below, 1.0   6.00  5.00  4.00 W8 W12  3.00  2.00  1.00  0.00 1.00  1.02  1.04  1.06  1.08  1.10  1.12  1.14  Strength reduction factor  JOB TITLE :  Figure 5.38: Results of displacement-based SSR for models W8 and W12 with UDEC (Version 4.00) static water tables (*10^2)  LEGEND 4.000  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  3.500  3.000  2.500  2.000  University of Brish Columbia Vancouver, BC Canada 2.500  3.000  3.500 (*10^2)  4.000  4.500  Figure 5.39: Plasticity indicators showing clear shear band after removal of last bench in W12 model  150  Chapter 5. Case History  JOB TITLE :  (*10^2)  UDEC (Version 4.00) LEGEND 4.000  19-Aug-08 13:57 cycle 113190 time 6.715E+01 sec block plot velocity vectors maximum = 8.090E-02  0  3.500  2E -1 3.000  2.500  2.000  University of Brish Columbia Vancouver, BC Canada 2.500  3.000  3.500 (*10^2)  4.000  4.500  5.000  Figure 5.40: Velocity plot showing slope failure along critical discontinuity behind 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 forward 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 CD1 CD3 CD5 CD8  Joint Network  D-SSR  CD-SSR  J1 J3 J5 J8  1.12 1.09 1.13 1.06  1.12 1.09 1.13 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. 8.00  Modelled horizontal displacement (m)  7.00  6.00  5.00  J1 CD1  4.00  J3 CD3 J5  3.00  CD5 J8  2.00  CD8 1.00  0.00 1.00  1.05  1.10  1.15  1.20  1.25  Strength reduction factor  Figure 5.41: Results of displacement-based SSR for models CD1–CD8 showing identical behaviour between continuously-cycled and initial state displacementbased SSR  5.4  Discussion  The above series of analyses, incorporating many different models and parameter 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 analysis 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. Modelling 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 carbonaceous 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 propagates 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 structure 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 surface 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 persis154  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 appear 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 conjunction 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 predicted 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, provided well defined joint geometry information from a site investigation is available.  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 equilibrium methods, numerical continuum and discontinuum (e.g. distinct element) failure mechanism investigations and the SSR method, in addition to a probabilistic 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 assessment for a rock slope using a variety of techniques, and what insights are typically gathered at each stage. Ultimately, probabilistic and limit equilibrium analyses 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 enters 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 failure mechanisms, material parameters, and joint networks, to assist in determining the factor of safety for a given rock slope. Chapter 4 involved a rigorous investigation 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 precision 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 presented 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 majority 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 behave 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 presented, 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, giving 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 assessment 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 applicable 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 residual properties and define how they change with respect to strain within the rockmass, 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 analysis 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 controlling 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 information 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 conditions and phases in design, including ongoing operation (or safety assurance) of the slope in question. Confidence in the results, for practical applications, requires 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 confidence 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 failure 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 network 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 constantly 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 concept, 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 assessments 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 unreliable, 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. Similarly, 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 criterion 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 potential 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 until significant cohesion loss occurs. In high strength rockmass cases, the joint network and joint properties may play a more significant role in the failure mechanism than in the yielding of intact rock dominated shear failure mechanism analyzed 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 stability 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 attention. 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 important 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 mechanism 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 returned 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 constrain 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 options, 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 available 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 HoekBrown, for intact rock representation in the model. 5. Getting around the limitations in the way UDEC handles joints (cannot terminate 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 numerical 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 modelling 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 discontinuum 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 hoekbrown prediction of strength and post yield behaviour for rockmasses at the extreme ends of the rock competency scale. In Luis, S., Olalla, C., and Grossman, N., editors, 11th Congress of the international society for rock mechanics, volume 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 Mechanics 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 Speciality Conferences, Vancouver. Hammah, R., Yacoub, T., Corkum, B., and Curran, J. (2005a). A comparison of finite-element slope stability analysis with cnventional limit-equilibrium investigation. 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. Symposium on Rock Mechanics (USRMS), Anchorage, Alaska. American Rock Mechanics Association.  170  Bibliography Hart, R. (1993). An introduction to distinct element modeling for rock engineering. 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 triggering 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 Excavations 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 proceudres 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, editors, 1ST Canada-US Rock Mechanics Symposium, volume 1, pages 311 – 318, Vancouver. Kulhawy, F. (1975). Stress deformation properties of rock and rock discontinuities. 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 Environment, 64:55–65. Martin, C. and Chandler, N. (1994). The progressive fracture of lac du bonnet 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 conference 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 mechanics. 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 2 nd 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 def blah r o d f r i c t i o n = 29 end blah prop jmat = 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 r u n from , c a l l e d 3 . s a v ; a t bottom of benching re 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 to 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 ; top of slope , middle of slope , toe h i s t n c y c l e 100 x d i s 3 2 2 . 0 , 3 5 2 . 1 4 h i s t n c y c l e 100 y d i s 3 2 2 . 0 , 3 5 2 . 1 4 h i s t n c y c l e 100 x d i s 3 6 8 . 0 , 3 1 6 . 9 3 8 h i s t n c y c l e 100 y d i s 3 6 8 . 0 , 3 1 6 . 9 3 8 h i s t n c y c l e 100 x d i s 4 1 5 . 0 , 2 8 0 . 3 6 2 h i s t n c y c l e 100 y d i s 4 1 5 . 0 , 2 8 0 . 3 6 2 h i s t n c y c l e 100 x d i s 4 6 2 . 0 , 2 4 3 . 7 8 6 h i s t n c y c l e 100 y d i s 4 6 2 . 0 , 2 4 3 . 7 8 6 ; reset tracking data reset disp vel jdisp ; ; ; ; ; 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 e a c h t i m e ; ; ; 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 e a c h 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 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; def generate loop a r r a y makeloop ( 1 0 0 ) ; open t h e f i l e f o r t h e SRF l o o p 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 n n n n  i t e r a t i o n s = 13 increment = 0.01 cycles = 20000 ratio = 1 e−6 force solve = 1  n s t a r t = 0.99 n cont = 1  ; ; ; ;  ; ; ; ; ;  how many l o o p s t o r u n how much t o i n c r e m e n t e a c h l o o p how many c y c l e s f o r e a c h model if solving for force ratio 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 ( w a r n i n g ) ; s t a r t i n g p o s i t i o n of the loop  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 b e t w e e n SRF s t e p s MUST h a v e n s t a r t s e t t o 1 w o r k i n g on a z o n i n 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 w a s 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 d o i n g ; continuous cycling .  ; 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 cont = 1 n s t a r t = 1 − n increment end if ; GLOBAL FOR HOW MANY VARIABLES ARE BEING TRACKED num variables ; ; ; ;  = 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 − turned off  ; 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 a s t h e h i s t o r y p o i n t s p1 p2 p3 p4  = = = =  gp gp gp gp  near (322.0 near (368.0 near (415.0 near (462.0  , , , ,  352.14) 316.938) 280.362) 243.786)  ; ; 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 cohesion 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  rod jf friction = 1 ; rod jf cohesion = 1 ; rod jf tension = 0 ; ;  joint friction j o i n t cohesion joint tension t u r n e d OFF a s j o i n t s h a v e 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 cohesion 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 z o n e s o t u r n e d OFF ; dp jf friction = 1 ; joint friction ; dp jf cohesion = 1 ; j o i n t cohesion ; dp jf tension = 1 ; joint tension  ; ;  C r e a t e t h e l o o p 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 ; ;  − initial state − continuous i f n cont = 0 then ; 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 l o o p 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 ) ; write f i l e endloop end if i f n cont = 1 then ; continuous 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 l o o p 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 ) ; write f i l e endloop end if status = close  end  ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; ; ; MULTI FOS I S 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 r u n 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 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; def m ul ti f os ; i n i t i a l FOS s e t t i n g s drad = 0.0174532952 ; degree to 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 to degree 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 e a c h l o o p  fos current = n s t a r t + n increment ∗ 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 u p p e r u n i t rod phinew = atan ( tan ( r o d r f r i c t i o n ∗drad ) / ( f o s c u r r e n t ))∗ rdeg 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 rod tennew = 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 ; in generate loop function if  r o d r f f r i c t i o n = 1 then command p r o p mat =1 f r i c = r o d p h i n e w end command  end if i f rod rf cohesion = 1 then command p r o p mat =1 coh = r o d c o h n e w end command end if i f r o d r f t e n s i o n = 1 then command p r o p mat =1 t e n = r o d t e n n e w end command end if ; 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 z o n e dp phinew = a t a n ( t a n ( d p r f r i c t i o n ∗ drad ) / ( f o s c u r r e n t ) ) ∗ rdeg 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 ; in generate loop function if  d p r f f r i c t i o n = 1 then command p r o p mat =3 f r i c = d p p h i n e w end command  end if i f dp rf cohesion = 1 then command p r o p mat =3 coh = dp cohnew end command end if i f d p r f t e n s i o n = 1 then command p r o p mat =3 t e n = d p t e n n e w end command end if  ; joint strength reduction ; c a l c u l a t e new v a r i a b l e s f o r u p p e r u n i t j 1  177  Appendix A. Displacement-based SSR Code rod j1 phinew = atan ( tan ( r o d j 1 f r i c t i o n ∗drad ) / ( f o s c u r r e n t ))∗ rdeg rod j1 cohnew = rod j1 cohesion / ( f o s c u r r e n t ) rod j1 tennew = rod j1 tension / ( fos current ) ; c a l c u l a t e new v a r i a b l e s f o r u p p e r u n i t j 2 rod j2 phinew = atan ( tan ( r o d j 2 f r i c t i o n ∗drad ) / ( f o s c u r r e n t ))∗ rdeg rod j2 cohnew = rod j2 cohesion / ( f o s c u r r e n t ) rod j2 tennew = rod j2 tension / ( fos current ) ; w r i t e new v a r i a b l e s t o memory if  r o d j f f r i c t i o n = 1 then command p r o p j m a t =1 j f r i c = r o d j 1 p h i n e w p r o p j m a t =2 j f r i c = r o d j 2 p h i n e w end command end if i f rod jf cohesion = 1 then command p r o p j m a t =1 j c o h = r o d j 1 c o h n e w p r o p j m a t =2 j c o h = r o d j 2 c o h n e w end command end if i f r o d j f t e n s i o n = 1 then command p r o p j m a t =1 j t e n = r o d j 1 t e n n e w p r o p j m a t =2 j t e n = r o d j 2 t e n n e w end command end if  ; 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 ; force equilibrium ; controlling flags ; n force solve = 1 ; = 2 ; = 0 command print n fos end command  t u r n s on f o r c e r a t i o s o l v i n g j u s t i s s u e s a s o l v e command ( may hang ) steps n cycles current  i f n force solve = 1 then ; solve for a force r a t i o command print n ratio solve r a t i o n ratio end command end if 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 print n cycles step n cycles end command  178  Appendix A. Displacement-based SSR Code end if ; pause ; w r i t e save f i l e s based o f f n parameter ; fos c i f continuous ; fos n if i n i t i a l state i f n cont = 0 then s a v n a m e = ’ f o s n ’ + s t r i n g ( i n t ( 1 0 0 0 ∗ f o s c u r r e n t ) ) + ’ . sav ’ command save sav name end command end if i f n cont = 1 then i f f o s c u r r e n t > n cont save then s a v n a m e = ’ f o s c ’ + s t r i n g ( i n t ( 1 0 0 0 ∗ f o s c u r r e n t ) ) + ’ . sav ’ command save sav name end command end if end if ; ; ; ; ;  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 write tracking variables to 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 for x displacements  ; array to carry tracked variables array x displacements (4000) ; ; ; ; ;  c a n n o t a p p e n d 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 s o i n s t e a d r e a d t h e f i l e i n t o memory e a c h l o o p 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  if n > 1 s t a t u s = open ( ’ x d i s . d a 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)) status = close end if ; u p d a t e t h i s a c c o r d i n g l y w i t h 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 ; writes the data being tracked 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 ( a b s ( 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 ( a b s ( 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 ( a b s ( 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 ( a b s ( 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 ( a b s ( 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 ( a b s ( 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 ( a b s ( 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 ( a b s ( g p y d i s ( p4 ) ) ) ; w r i t e s b o t h t h e o l d d a t a and t h e new d a t a ( h e n c e n u m v a r i a b l e s ∗n ) s t a t u s = open ( ’ x d i s . d a 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 ) status = close ; saves the f i l e based o f f of parameter n  end ; ; ; ; ; ; FUNCTION CONVERT ARRAY ; ; ; ; ; ; T h 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 being 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 a s 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 o f t h e l o o p c r a s h e s UDEC f o r an unknown r e a s o n . Getting around t h i s h a p p e n s 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 e a s o n , 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 c o d e t h a t would c r a s h UDEC .  def co nv er t ar r ay ; s e t up a r r a y s array xread (4000) 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 . d a 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 = read ( xread , n u m v a r i a b l e s ∗( n ) ) status = close ; generate tables ; c o n v e r t s t r i n g a r r a y ( xread ) i n t o something 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 . d a 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 ) end loop xtemp ( 1 ) = ’ end ’ xtemp ( 2 ) = ’ c o n v i n t ’ ; f i n i s h the f i l e s t a t u s = w r i t e ( xtemp , 2 ) status = close ; b r i n g t h e f i l e i n t o r e−w r i t e v a l u e s a s f l o a t s command c a l l convert . dat 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 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;  def c r e a t e t a b l e s ; g o e s t h r o u g h l o o p 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 monitoring point xtable (100 , i ) = n s t a r t + n increment ∗ 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 ) ) )  xtable (101 , i ) = n s t a r t + n increment ∗ 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 monitoring p o i n t xtable (200 , i ) = n s t a r t + n increment ∗ 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  xtable (201 , i ) = n s t a r t + n increment ∗ 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 monitoring point xtable (300 , i ) = n s t a r t + n increment ∗ 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 ) ) ) xtable (301 , i ) = n s t a r t + n increment ∗ 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 ) ) ) ; fourth monitoring point xtable (400 , i ) = n s t a r t + n increment ∗ 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 ) ) ) xtable (401 , i ) = n s t a r t + n increment ∗ 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 ) ) ) end loop end ; end o f f u n c t i o n  ;;;;;;;;;;;;;;;;;;;;  ;−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; ; ; ; ; ; ; END FUNCTION DEFINITIONS ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−; ; c a l l generate loop function to create batch f i l e for looping generate loop ; ; ; ;  saves the i n i t i a l s t a t e , placement i s important , as i t k e e p s a l l t h e a b o v e s e t t i n g s i n memory s o t h a t whent i n t . sav i s r e s t o r e d a l l the i m p o r t a n t v a r i a b l e s such as n increment etc are tracked  save i n t . sav ; c a l l s the batch f i l e for external looping call srf . txt ; ; ; ; ; ;  t o g e n e r a t e t a b l e s , a t any t i m e r e s t o r e l a s t s a v e f i l e e v e n 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 . re fos n 1079 . sav and c a l l t h e f u n c t i o n s : convert array −−−> f o l l o w e d by create tables  ; 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 ; to t e x t f i l e s f o r e x t e r n a l use ( excel )  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 the l a s t save f i l e in the s r f loop ; then c a l l output . t x t  convert array create tables table table table table table table table table  100 101 200 201 300 301 400 401  write write write write write write write write  0 0 0 0 0 0 0 0  data data data data data data data data  t100 t101 t200 t201 t300 t301 t400 t401  . . . . . . . .  txt txt txt txt txt txt txt txt  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 g e o m e t r y ; c r e a t e s j o i n t s and e l e m e n t m e s h i n g ; 10 M d i p p a r a l l e l c a s e  new title  ’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 geometry 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 e d g e 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 block d e f i n i n g the o r i g i n a l boundary ; c r e a t e u s i n g b l o c k , l e f t hand r u l e b l 0 , −250 0 , 5 0 0 1 3 5 0 , 5 0 0 1350 , −250 ; c r 438 ,500 688 ,250 c r 682 ,250 1350 ,250 ; c r 611.7 , 326.28 1350 , 326.28 ; c r 528.997 , 409.01 1350 , 409.01 c r 597.694 , 335.680 1350 , 335.680 c r 528.516 , 404.764 1350 , 404.764  185  Appendix B. Example Input Files for Parametric Study jregion jregion jregion jregion jregion  id id id id id  7 0 , 500 1 7 5 , 500 1 7 5 , −250 0 , −250 6 1 7 5 , 500 4 3 8 , 5 0 0 7 8 8 , 150 1 7 5 , 150 8 6 8 2 , 250 9 2 0 , 250 8 8 0 , 150 7 8 8 , 150 9 1 7 5 , 150 1 3 5 0 , 150 1 3 5 0 , −250 1 7 5 , −250 10 8 8 0 , 250 1 3 5 0 , 250 1 3 5 0 , 150 8 8 0 , 150  ; ; generate joint pattern  j s e t −45 ,0 1 0 0 , 0 jset 45 ,0 100 ,0 j s e t −45 ,0 1 0 0 , 0 jset 45 ,0 100 ,0  0 ,0 0 ,0 0 ,0 0 ,0  10 ,0 10 ,0 10 ,0 10 ,0  range range range range  jregion jregion jregion jregion  6 6 8 8  j s e t −45 ,0 1 0 0 , 0 0 , 0 4 0 , 0 r a n g e j r e g i o n jset 45 ,0 100 ,0 0 ,0 40 ,0 range j r e g i o n  7 7  j s e t −45 ,0 1 0 0 , 0 0 , 0 4 0 , 0 r a n g e j r e g i o n jset 45 ,0 100 ,0 0 ,0 40 ,0 range j r e g i o n  9 9  j s e t −45 ,0 1 0 0 , 0 0 , 0 4 0 , 0 r a n g e j r e g i o n jset 45 ,0 100 ,0 0 ,0 40 ,0 range j r e g i o n  10 10  jdel ; 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 ch ch ch  mat mat mat mat  2 3 r a n g e b l 5923 3 r a n g e b l 6369 3 r a n g e b l 6848  ; c o a r s e mesh away gen gen gen gen  edge edge edge edge  25 25 25 25  range range range range  jregion jregion jregion mat 3  7 9 10  ; f i n e mesh on s l o p e gen e d g e 5 r a n g e j r e g i o n 6 gen e d g e 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 ; sets strengths ; e l a s t i c gravity loading r e s 1. sav ; s e t c o n s t i t u t i v e m o d e l s 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 ( c o n s = 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 w i t h r e s . s t r e n g t h ; s e t i n t a c t r o c k 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 of p r o p e r t i e s used ; mat = 1 − default material ; mat = 2 − j o i n t e d rock ; mat = 3 − u n j o i n t e d rock of type 2 ;−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−  ; ; ; ; ; ;  s i g c i = 25 g s i = 50 mi = 10 MR = 275 slopes s l o p e h e i g h t 250m  ; ; ; Default material d e f Emat1 b l k =m1E/(3∗(1 −2∗ m1pr ) ) ; b u l k modulus s h e =m1E / ( 2 ∗ ( 1 + m1pr ) ) ; s h e a r modulus command p r o p mat =1 d =2600 b u l k = b l k s h e a r = s h e f r i c =56 coh = 1 . 7 7 1 e6 t e n =0 z o n e d e n s =2600 b u l k = b l k s h e a r = s h e f r i c =56 coh = 1 . 7 7 1 e6 t e n =0 r a n g e mat 1 end command end s e t m1E=7 e9 m1pr = 0 . 3 3 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 ) ) ; b u l k modulus s h e =m2E / ( 2 ∗ ( 1 + m1pr ) ) ; s h e a r modulus command p r o p mat =2 d =2600 b u l k = b l k s h e a r = s h e f r i c =33 coh = 1 . 0 1 e6 t e n =100 e3 end command end s e t m2E=7 e9 m1pr = 0 . 2 7 Emat2 ; ; ; Unjointed rock d e f Emat3 b l k =m3E/(3∗(1 −2∗ m1pr ) ) ; b u l k modulus s h e =m3E / ( 2 ∗ ( 1 + m1pr ) ) ; s h e a r modulus command p r o p mat =3 d =2600 b u l k = b l k s h e a r = s h e f r i c =33 coh = 1 . 0 1 e6 t e n =100 e3 end command end s e t m3E=7 e9 m1pr = 0 . 2 7 Emat3 ; ; ; ; ; ; ; ; ; ; ; INPUT OF JOINT MATERIAL PROPERTIES ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ∗ ∗ l i s t of p r o p e r t i e s used ; jmat = 1 − a r t i f i c i a l f r a c t u r e s ; j m a t = 5 − JS # 1 ; j m a t = 6 − JS # 2 ; j m a 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 j m a t =1 c h a n g e j m a t =1 ; artificial joints p r o p j m a 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 p r o p j m a t =2 j k n =10 e9 j k s =8 e9 j f r i c = 2 5 . 0 j c o h =100 e3 j t e n =0 ; ; ; JS #2 p r o p j m a t =3 j k n =10 e9 j k s =8 e9 j f r i c = 2 5 . 0 j c o h =100 e3 j t e n =0  ; ; ; ; ; ; ; ; ; ; ; ; BOUNDARY CONDITIONS ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; ; NEEDS TO BE CHANGED EVERY TIME GEOMETRY I S CHANGED ; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; Def B o u n d B u f f e 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 xy buffer =2;10.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 b o t t o m 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 top p  = 500+ x y b u f f e r  End Bound Buffer ; ; ; R o l l e r Boundary on s i d e s , b o t t o m 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 b o t t o m s e t g r a v 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 g r = −9.81 a r a t i o = 1.0 a t o p = 500 a n e g = −1 syy sxy sxx szz  0 0 0 0  = = = =  ; average density ; gravity ; r a t i o n of h o r i z o n t a l to v e r t i c a l s t r e s s ; t o p h e i g h t o f model , n e e d s t o be c h a n g e d ; e v e r y t i m e model g e o m e t r y i s c h a n g e d  a density ∗ a gr ∗ a top 0.0 syy 0 ∗ a r a t i o syy 0 ∗ a r a t i o  ; y stress y syy g = y sxx g = y szz g =  ; syy ; sxy ; sxx ; szz  at at at at  zero zero zero zero  elevation elevation elevation elevation  gradients a d e n s i t y ∗ a gr ∗ a neg a d e n s i t y ∗ a gr ∗ a r a t i o ∗ a neg a d e n s i t y ∗ a gr ∗ a r a t i o ∗ a neg  End Cal Insitu Stress ; s e t i n s i t u s t r e s s regime i n s i t u s t r s x x 0 , s x y 0 , s y y 0 yg y s x x g , 0 , y s y y g i n s i t u s z z s z z 0 zg 0 , y s z z g  ; ywt t a b 100 ; w a t e r t a b l e NOT u s e d  ;∗ All m a t e r i a l s are e l a s t i c to reach i n i t i a l equilibrium ch c o n s =1 Emat1 Emat2 Emat3  ; default ; jointed ; unjointed  damp a u t o h i s n c y c l e 100 u n b a l h i s n c y c l e 100 damp  ;1 ;2  189  Appendix B. Example Input Files for Parametric Study h i s n c y c l e 100 vmax solve ; e l a s t i c s t e p 10000 sav 2. sav call 3. txt  ;3  190  Appendix B. Example Input Files for Parametric Study The third file: ; BEGIN 3 . t x t ; sets real strengths ; b e n c h e s down  ; c h a n g e r o c k m a s s t o MC b e h a v i o u r re 2. sav ch c o n s =3 ; z o n e 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 j m a t =2 Emat1 Emat2 Emat3  ; default ; j o i n t e d rock ; u n j o i n t e d rock  solve ; Set j o i n t strengths , keeping a r t i f i c i a l  j o i n t s as jmat1  ch j c o n s =2 p r o p j m a 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 p r o p j m a t =2 j k n =10 e9 j k s =8 e9 j f r i c = 2 5 . 0 j c o h =100 e3 j t e n =0 ; ; ; JS #2 p r o p j m a t =3 j k n =10 e9 j k s =8 e9 j f r i c = 2 5 . 0 j c o h =100 e3 j t e n =0 ch j m a t =2 solve s t e p 5000 d e l b l 5923 solve s t e p 8000 d e l b l 6369 d e l b l 170591 solve s t e p 8000 d e l b l 6848 solve 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 b a s e d SSR method  191  Appendix C Case History Joint Networks  JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 30-Jul-08 9:54 cycle 75801 time 7.263E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J1 joint network used in case history modelling (*10^2)  192  JOB TITLE :  Appendix C. Case History Joint Networks  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 5-Aug-08 22:34 cycle 78690 time 6.158E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J2 joint network used in case history modelling (*10^2)  JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 5-Aug-08 22:58 cycle 86580 time 4.557E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J3 joint network used in case history modelling (*10^2)  193  JOB TITLE :  Appendix C. Case History Joint Networks  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 5-Aug-08 22:46 cycle 82060 time 5.108E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J4 joint network used in case history modelling (*10^2)  JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 6-Aug-08 11:13 cycle 83330 time 5.576E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J5 joint network used in case history modelling (*10^2)  194  JOB TITLE :  Appendix C. Case History Joint Networks  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 6-Aug-08 17:06 cycle 83150 time 5.002E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J6 joint network used in case history modelling (*10^2)  JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 6-Aug-08 16:49 cycle 86880 time 4.492E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J7 joint network used in case history modelling (*10^2)  195  JOB TITLE :  Appendix C. Case History Joint Networks  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 8-Aug-08 14:11 cycle 84360 time 4.962E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J8 joint network used in case history modelling (*10^2)  JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 9-Aug-08 16:03 cycle 84530 time 4.764E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J9 joint network used in case history modelling (*10^2)  196  JOB TITLE :  Appendix C. Case History Joint Networks  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 9-Aug-08 11:03 cycle 85790 time 4.752E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J10 joint network used in case history modelling (*10^2)  JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 8-Aug-08 17:08 cycle 84320 time 4.971E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J11 joint network used in case history modelling (*10^2)  197  JOB TITLE :  Appendix C. Case History Joint Networks  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 11-Aug-08 12:52 cycle 87980 time 4.122E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J12 joint network used in case history modelling (*10^2)  JOB TITLE :  (*10^2)  UDEC (Version 4.00) 4.000  LEGEND 7-Aug-08 11:28 cycle 86920 time 4.387E+01 sec  3.600  block plot  3.200  2.800  2.400  2.000  University of Brish Columbia Vancouver, BC Canada 2.600  3.000  3.400  3.800  4.200  4.600  J13 joint network used in case history modelling (*10^2)  198  

Cite

Citation Scheme:

    

Usage Statistics

Country Views Downloads
China 25 86
United States 20 1
Canada 10 0
Japan 8 0
South Africa 5 0
France 3 0
Australia 3 1
Germany 3 0
Brazil 2 0
Bolivia 2 0
Chile 2 0
Republic of Lithuania 2 0
Indonesia 2 1
City Views Downloads
Unknown 32 3
Shenyang 9 84
Beijing 7 0
Tokyo 7 0
Ashburn 6 0
Chengdu 5 2
Aachen 3 0
Burnaby 3 0
Sunnyvale 3 0
Shanghai 3 0
Santiago 2 0
Durban 2 0
La Paz 2 0

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

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.24.1-0052558/manifest

Comment

Related Items