A STRESS-STRAIN MODEL FOR THE UNDRAINED RESPONSE OF OIL SAND B.Sc, University of Manchester, 1983 A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF APPLIED SCIENCE in FACULTY OF GRADUATE STUDIES Department of C i v i l Engineering We accept this thesis as conforming by KA FAI HENRY/CHEUNG to the required standard THE UNIVERSITY OF BRITISH COLUMBIA April 1985 © KA FAI HENRY CHEUNG » ?985 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at THE UNIVERSITY OF BRITISH COLUMBIA, I agree that the Library shall make i t freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the Head of my Department or by his or her representative. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department of C i v i l Engineering THE UNIVERSITY OF BRITISH COLUMBIA 2075 Wesbrook Place Vancouver, B.C. CANADA, V6T 1W5 Date: April, 1985 ABSTRACT An efficient undrained model for the deformations analyses of o i l sand masses upon undrained loading i s presented i n this thesis. An analysis which couples the s o i l skeleton and pore fluids i s used. The s o i l skeleton i s modelled as a non-linear elastic-plastic isotropic material. In undrained conditions, the constitutive relationships for the pore fluids are formulated based on the ideal gas laws. The coupling between the s o i l skeleton and the pore fluids i s based upon volume compatibility. The undrained model was verified with the experimental results and one dimensional expansion of s o i l sand cores. Comparisons between computed and measured responses are i n good agreement and suggest that this model may prove useful as a tool in evaluating undrained response of o i l sand. The response of a wellbore i n o i l sand upon unloading was analysed using the developed model. Such analyses are important in the rational design of o i l recovery systems i n o i l sand. i i i TABLE OF CONTENTS Page Abstract i i Tabele of Contents i i i List of Tables v i i Li s t of Figures v i i i L ist of Symbols x Acknowledgements x i i i CHAPTER 1 Introduction 1.1 Introduction 1 1.2 Behaviour of Oil Sand 2 1.3 The Scope 5 1.4 The Organization of Thesis 6 CHAPTER 2 Review of Previous Work 2.1 Introduction 8 2.2 Mathematical Model for Undrained Behaviour of O i l Sand 2.2.1 Harris and Sobkowicz 8 2.2.2 Dusseault 12 2.2.3 Byrne and Janzen 14 CHAPTER 3 Stress-Strain Model 3.1 Introduction 17 3.2 Development of the Undrained Model 17 3.3 Constitutive Relations 3.3.1 Incremental Non-linear Elastic Soils Model 22 3.3.2 Incremental Non-linear Elastic Fluid Model 25 i v Page 3.3.2.1 General 25 3.3.2.2 Partly Miscible Gas/Liquid Mixture 27 a. Air/Water Mixture 28 b. Carbon Dioxide, Air/Water 30 Mixture c. Gas/Bitumen and Water Mixture 33 CHAPTER 4 Finite Element Formulations 4.1 Introduction 36 4.2 The Plane Strain Formulation 4.2.1 The Constitutive Matrix [D] 36 4.2.2 The Strain Displacement Matrix [B] 40 4.2.3 The Stiffness Matrix [K] 44 4.3 The Axisymmetric Formulation 4.3.1 The Constitutive Matrix [D] 46 4.3.2 The Strain Displacement Matrix [B] 50 4.3.3 The Stiffness Matrix [K] 53 4.4 Load Shedding Formulation 55 CHAPTER 5 Comparisons with Existing Solutions 5.1 Introduction 60 5.2 Comparisons with Theoretical Solutions 5.2.1 Elastic Closed Form Solutions 60 5.2.2 Elastic-Plastic Closed Form Solutions 61 5.3 Comparisons with Observed Data 5.3.1 One Dimensional Unloading of O i l Sand 68 Page 5.3.2 Triaxial Tests on Gassy Soils 72 CHAPTER 6 Stresses Around a Wellbore or Shaft i n O i l Sand 6.1 Introduction 79 6.2 General Model Description 80 6.3 Theoretical Solutions for Stresses Around a Wellbore 80 6.3.1 Stresses 6.3.1.1 Stresses in Elastic Zone 83 6.3.1.2 Stresses in Plastic Zone 85 6.3.1.3 Radius of Plastic Zone 87 6.3.2 Stability 89 6.3.3 Pore Pressure Profile 90 6.4 Comparisons of Predicted Response and Closed Form Solution 93 6.5 Analysis of Response of Borehole in Oil Sand on Unloading 6.5.1 Undrained Response 93 6.5.2 Drained Response 100 6.5.3 Implications of Undrained and Drained Analysis 104 6.6 Application of o i l recovery 104 CHAPTER 7 Summary and Conclusions 108 Bibliography 110 Appendix A 114 Appendix B 116 Appendix C 119 VI Page Appendix D Appendix E Appendix F Appendix G Appendix H Appendix J . _, 136 LIST OF TABLES Table v i i Page A comparison of computed and measured results (Test 11, Sobkowicz) 74 V l l l LIST OF FIGURES Figure _ __ Page 1.1 Schematic spring analogy for o i l sand 3 3.1 Stress-strain curves for drained t r i a x i a l tests on loose sand 24 3.2 Phase diagram for gassy s o i l 32 3.3 Phase diagram for o i l sand 34 4.1 Stresses associated with load shedding 56 5.1 thick wall cylinder 62 5.2 Stresses and displacements around circular opening in an elastic material 63 5.3 Stresses and displacements around circular opening in an elastic-plastic material 64 5.4a Model for one dimensional unloading of o i l sand 70 5.4b Response of o i l sand to one dimensional unloading .... 71 5.5 Comparisons of predicted and observed pore pressure .. 76 5.6 Comparisons of predicted and observed strains 77 6.1 Outline of the wellbore problem 81 6.2a Finite element mesh for wellbore problem 82 6.3b Mohr Coulomb failure envelope 84 6.3 Idealised flow to a wellbore 92 LIST OF FIGURES - Continued Figure Page 6.4 Stresses around a wellbore i n an elastic-plastic material 94 6.5 Undrained response of wellbore i n o i l sand on unloading 96 6.6 Pore pressure profile around a borehole 101 6.7 Drained response of a wellbore in o i l sand on unloading 102 6.8 Comparisons of undrained and drained response of a wellbore in o i l sand at a support pressure of 3500 kPa 105 6.9 Comparisons of undrained and drained response of a wellbore in o i l sand at a support pressure of 3200 kPa 106 X LIST OF SYMBOLS The following is a l i s t of the commonly used symbols in this thesis. Multiple use of several symbols is unavoidable because of the complexity of the formulations. The symbol w i l l be defined immediately in the text where the use of symbols differs from those l i s t e d below. SYMBOL MEANING 3 compressibility A change a stress e strain v Poisson's ratio (j) f r i c t i o n angle a temperature solubility constants a radial stress r o"g tangential stress a vertical stress z B bulk modulus e . void ratio E Young's modulus G shear modulus H Henry's constant k permeability K apparent bulk f l u i d modulus x i bulk f l u i d modulus n porosity P pressure q flow rate r radius (variable) R radius (constant) R radius of plastic zone c r radius of wellbore w S saturation T temperature u pore f l u i d pressure V volume x i i Subscripts 1 f i n a l a a i r f pore f l u i d g gas i i n i t i a l (internal i n Chapter 6) 0 o i l (outer in Chapter 6) S s o i l skeleton v volumetric w water cb combined dg dissolved gas fg free gas Tg total gas Superscripts e elastic f f i n a l 1 i n i t i a l p plastic x i i i ACKNOWLEDGEMENTS I am greatly indebted to my supervisor, Professor P.M. Byrne, for his guidance, encouragement end enthusiastic interest throughout this research. I would also like to thank Professor Y.P.Vaid for reviewing the manuscript and making valuable suggestions. My colleagues, U. Atukorala, F. Salgado, J. She, H. Va z i r i and especially, C. Lum, shared a common active interest in s o i l mechanics. I thank them a l l for their helpful discussions and constructive criticisms. Appreciation is extended to Ms. S.N. Krunic for typing the manuscript and her patience during the preparation of this thesis. Support and assistance provided by the Natural Science and Engineering Research Council of Canada i s acknowledged with deep appreciation. A special thanks i s extended to my family and P r i s c i l l a , for their constant support and encouragement. 1 1.1 INTRODUCTION Many schemes for o i l recovery require open excavations, tunnels or wellbore i n o i l sand. As a result, an accurate and efficient analysis of stresses and deformations around these openings i s becoming increasingly important. Mathematical models have been developed by Byrne et a l (1980), Dusseault (1979), and Harris and Sobkowicz (1977) to investigate these problems. An efficient undrained model for analysing the stresses and deformations around open excavations, tunnels and wellbores i n o i l sand i s presented i n this thesis. Oil sand i s comprised of a dense sand skeleton with i t s pore spaces f i l l e d with bitumen, water and free or dissolved gas. The presence of bitumen reduces the effective permeability of o i l sand, hence undrained conditions occur on rapid unloading. Gas evolves from pore fluids during unloading when the pore fluid pressure i s below gas saturation pressure. Because gas exsolution takes some time to occur, two undrained conditions arise, (1) an immediate or short term condition i n which there i s no time for gas exsolution, and (2) a long term or equilibrium condition in which complete gas exsolution has occured. Both of these conditions are considered in this thesis. The rate of gas exsolution is not considered herein. The drained analysis with pore pressures under steady-state conditions i s also addressed. The sand skeleton is modelled as a non-linear elastic-plastic porous material. The fluids stress-strain relationships for undrained condition are formulated on the basis of ideal gas law. For the undrained condition, the pore pressure changes ;.re computed from the constant of volume compatibility between the sand skeleton and the pore fluids. It i s assumed that the pore f l u i d pressures are known in the drained analysis. 2 Validation of the stress strain model i s made by comparing responses with theoretical solutions and observed data. The observed data are the experimental results in Sobkowicz's doctorate thesis. The stresses and deformations around a wellbore in o i l sand upon unloading i s investigated using the new stress-strain model. 1.2 BEHAVIOUR OF OIL SAND Oil sand is comprised of a dense, highly incompressible, uncemented, interlocked skeleton with pore spaces f i l l e d by water, bitumen, and dissolved or free gas. The interpenetrative structure leads to the low in-situ void ratio and high shear strength. The response of o i l sand is mainly governed by the rate of loading. Undrained conditions occur as a result of: 1) low effective permeability to pore fluids 2) large amount of dissolved gas in pore fluids 3) rapid unloading The response of o i l sand may be physically modelled by a set of springs shown i n Figure 1.1. The o i l sand i s s p l i t into two load carrying components - s o i l skeleton and pore fluids. The s o i l skeleton compressibility, 3 g, characterizes the deformation of s o i l skeleton, which results i n a change i n effective stress. The pore f l u i d compressibility, 8^, characterizes the deformation of pore fluids (because of free gas) which results i n a change i n pore pressure. 3 a + Aa X X X X 4 4 X 4 6 f K=> O § Soil § Skeleton o § Water «o § Oil (bitumen) 5 \ § Stress = a'+ Aa § ) o § Gas Pore Fluid Pressure = u + Au Rigid Frame, No Lateral Yield a = a' + u A o = A a ' + A u A V = / ( A a , A u ) Fig.1.1 - Schematic spring analogy for o i l sand ( After Dusseault , 1979 ) 4 Strain compatibility between the pore fluids and s o i l skeleton controls the relative magnitudes of the changes in pore pressures and effective stress which together equal the change i n total stress. Stress changes (unloading) are shared between the sand skeleton and pore fluids according to their compressibilities. When the pore pressure is above the gas saturation pressure and the o i l sand is 100% saturated, the stress changes w i l l be accommodated by the pore fluids because their compressibilities are lower. Once the pore pressure drops below the gas saturation pressure, gas starts to evolve which increases the fluids' compressibilities. Therefore the sand skeleton becomes the less compressible phase and takes up the stresses rather than the pore fluids. When the effective stress drops to zero, the s o i l skeleton i s very compressible relative to the pore fluids; hence any further decrease in total stress is entirely accommodated by the pore pressure. The gas exsolution takes some time to occur, hence two undrained conditions (Sobkowicz, 1982) arise. The expressions 'short term' and 'long term' w i l l be applied to these processes exclusively. 1) 'Short term' undrained in which there i s no time for gas exsolution 2) 'Long term' undrained i n which equilibrium state has been reached, i.e. completion of gas exsolution. In the f i e l d , the unusal behaviour of o i l sand manifests in a number of ways. They include: 5 1) Volumetric expansion of 5 to 15% occurred when core samples were l e f t in an unconfined state, i.e. not retained by plastic core sleeves (Dusseault 1980; Byrne et a l 1980). 2) Core samples spontaneously split longitudinally and perpendicularly to the core axis, effervescence was observed on several freshly recovered cores (Hardy and Hemstock 1963). 3) Retrogression of slopes i n o i l sand on rapid excavation. 4) Oil sand at the base of excavation subjects to softening and heaving, followed by settlement on reloading. When the decrease in external stress occurs over a period of time, the evolved gas has time to drain off and effective stress does not go to zero which results in an undisturbed sand skeleton. Its high in - s i t u density and hence high shear strength i s retained. Such situations can be seen on the exposures of o i l sand deposits along the Valley of Athabasca River where erosion (unloading) by the river has occurred over thousands of years. These o i l sand deposits are standing on steep stable slopes with slope angles i n excess of 60° and heights up to 60 metres, exhibiting high strength (Harris and Sobkowicz 1977). 1.3 THE SCOPE The purpose of this study i s to present a stress-strain model (Vaziri, 1985) for the deformation analyses of o i l sand masses, i.e. to explain the behaviour of o i l sand as described in Section 1.2. Modelling the undrained response of o i l sand requires the pore fluids pressure to be numerically evaluated. For the undrained condition, pore f l u i d pressures are computed from the constraint of volumetric compatibility between the sand skeleton and pore f l u i d phases. Compressibilities of sand skeleton and pore fluids have to be evaluated in the new stress-strain model. A non-linear stress-strain relationship (Duncan et a l 1970) is adopted for the sand skeleton. The compressiblities of pore fluids are formulated on ths basis of ideal gas laws. The new stress-strain model i s incorporated into f i n i t e element programmes (INCOIL, MHANS). For the validation of the stress-strain relationship, the computed results are compared with the theoretical solutions and observed data. The validated programme was used to study the unloading response of stress and deformations around a deep wellbore in o i l sand. 1.4 ORGANIZATION OF THE THESIS This thesis consists of seven chapters. A review of previous work on undrained models for o i l sand i s given i n Chapter 2. This chapter concentrates on the examinations of their capabilities and shortcomings. A stress-strain model which was developed by Vaziri i s presented in Chapter 3. Appropriate s o i l skeleton and pore fluids constitutive relationships are also recommended in this chapter. Chapter 4 summarizes the f i n i t e element formulations used i n the development of the programme. Validation of the developed model by comparing predicted response with existing theoretical solutions, f i e l d and laboratory observations i s mentioned in Chapter 5. The response of a wellbore o i l sand upon unloading is investigated in Chapter 6 using the validated f i n i t e element programme. A summary of work, and major conclusions are presented i n Chapter 7. 8 CHAPTER 2 : REVIEW OF PREVIOUS WORK 2.1 INTRODUCTION Theoretical solutions for the undrained response of o i l sand were not available u n t i l 1977 because of i t s unusual behaviour. Due to the increase in demand for construction in o i l sand formation, such as open pit mining, tunnels and deep shafts, a considerable amount of research has been done on this topic since 1977. These developed theoretical relations share the same basic approach of coupling the so i l skeleton and pore fluids together. Pore pressure changes are computed from the constraint of volumetric compatibility between the s o i l skeleton and pore fluids. Harris and Sobkowicz's model i s capable of evaluating pore pressure change upon a stress or/and temperature change. This model was then extended by Byrne et a l and incorporated in a f i n i t e element programme. Dusseault's approach is restricted to one dimensional problems. A careful examination of these theoretical models and their applications is made and several shortcomings are also reviewed. 2.2 MATHEMATICAL MODELS FOR OIL SAND 2.2.1 Harris and Sobkowicz (1977) Harris and Sobkowicz presented a mathematical model to analyse the undrained response of o i l sand subjected to changes of stress and/or temperature. This is the f i r s t analytical model developed to explain the behaviour of o i l sand such as: 1) Movement and sta b i l i t y of slopes and tunnels formed in o i l sand. 2) Settlements or heave of structures placed in o i l sand (e.g. hot o i l tank). 3) Heave at the_base of excavations i n o i l sand. The model can be explained by the same spring analogy as shown in Figure 1.1. The response of o i l sand to changes i n stresses or/and to changes in temperature may be computed by the following equations: a) One-Dimensional Analysis Au can be obtained from a quadratic equation L * Au 2 + M * Au + N = 0 2.2 where L = 3 s Tl M = 3 (P. - Aoi) + n + Pi — (n H + n H ) s i 1 / g 1 T v w w o o ' a N = P ± [n (1 - - A 0 l 3 s ] " Pi AT ^ - ( n ^ + n ^ ) 10 b) Two-Dimensional Analysis from strain compatibility again Aa 3 + A (Ao-! - Ao 3) Au = and where and 1 + *s n6f 2.3 P Au z + Q Au + R = 0 2.4 P = B s T 1 Q = 8 (P. - X) + n + Pi — (n H + n H ) s v i g 1 T w w o o ' R = P. [n (1 - — ) - X 0 ] - P, AT — (n a + n a ) l L g T s J i T v w w 0 0 X = A02 + Ap (Aai - Ao"2) = volumetric strain B g = s o i l skeleton compressibility n = porosity of s o i l n Q = porosity of o i l n w = porosity of water ^ n = porosity of gas Aa = change in total stress 11 P = i n i t i a l pore f l u i d pressure Au = change in pore fluid pressure Ao"i = change i n major principal stress Ao"3 = change in minor principal stress T a = standard temperature (288 °K) = i n i t i a l pore f l u i d temperature °K Tl = f i n a l pore f l u i d temperature °K AT = change i n temperature = temperature solubility constant for gas in water = temperature solubility constant for gas in o i l H W = Pressure solubility constants for gas i n water H q = Pressure solubility constants for gas i n o i l Ap = s o i l skeleton dilatancy factor Harris and Sobkowicz examined the response of o i l sand during excavation and reloading of a square footing, and the behaviour of a tunnel excavation in the same material. Their results w i l l allow an assessment of the applicability and shortcoming of the theory. 1) This solution..incorporated a linear constitutive relationship for the s o i l . As o i l sand behaves like an elastic-plastic material, there w i l l be a plastic zone developed adjacent to the tunnel wall on unloading when the effective.stresses are such that failure (Mohr Coulomb Criteria) occurs. 2) The extent of the plastic zone i s only an approximation because redistribution of stresses has not been taken into account during the formation of plastic zone. 3) An iterative procedure is required to obtain Au 2.2.2 DUSSEAULT Dusseault (1979) extended the one dimensional Skempton's B equation to analyse the behaviour of cohesionless materials with large amount of free or dissolved gas in pore fluids. He presented a more rigorous derivation for the compressibility of pore fluids, and coupled the compressibility with that of s o i l skeleton ( e - a' relationship) in Skempton's B equation. This model shares the same basic idea as the previous one (Harris and Sobkowicz), which has two load carrying components - s o i l skeleton and pore fluids. The relative compressibilities of these components control the magnitude of changes in pore pressure and effective stress which together balance the change in total stress i n an undrained state. The traditional Skempton's B equation is Aa 1 + n-2.5 8^ = compressibility of pore fluid (water) 3 = compressibility of s o i l skeleton The extended one i s ^ = l / [ l + f(u,o)] 2.6 a+bAn(a-u)-e -e +H e -Hi e f(u,a) = [ e g + e 3 + J . o o w w J D U a, b = void ratio - effective stress relationship parameters e , e , e = void ratios : o i l , water, gas o' w' g » » e> H , H = Henry's constants : o i l , water o' w 3 Q, 3 w = compressibility : o i l , water u, a = current values of pore pressure and total stress Dusseault applied this model to investigate the response of an element of o i l sand upon unloading. He examined a shallow case (15 m) and a deeper one (500 m). His results w i l l allow an assessment of the applicability and shortcoming of the theory: 1) In order for the solution to be numerically stable, further reduction in total stress i s assumed to be entirely taken by the pore pressure once effective stress drops to zero. 2) This solution incorporated a non-linear constitutive relationship for the s o i l . 3) Accurate e - a' relationship or compressibility of o i l sand is extremely d i f f i c u l t to obtain since they are very sensitive to sample disturbances. No undisturbed o i l sand samples have been cored so far. 4) An iterative procedure i s required to get Au. 2.2.3. Byrne and Janzen Byrne et al (INCOIL, 1983) developed an incremental analytical f i n i t e element method for predicting stresses and deformations in excavations and around tunnels in o i l sand using a nonlinear elastic sand skeleton with shear dilation. An extension to Harris and Sobkowicz's model was used to evaluate the pore pressure change, Au. The f i n i t e element program, INCOIL, can handle both undrained and drained analyses. For the undrained condition, pore flu i d pressures are computed from the constraint of volumetric compatibility between the sand skeleton and the pore fluids. For the fu l l y drained condition, i t i s assumed that the pore f l u i d pressures are known. The general framework of the f i n i t e element model for o i l sand are as follows: [K] (fi) = (Af) - (K ) (Au) 2.7 w where u Ti | — * AT * (n a +n a ) - u fn *(1- — ) -T v w w o o ' « L a v T ' o L g A e 1 v J Au = u Ti a L 2.8 (n H +n H ) + n - Ae w w o o g v Derivation of equation 2.7 is presented in Appendix A 15 The effective stress may be evaluated from (Ae) = [c] (A6) 2.9 (Aa') = [D»] (Ae) 2.10 where [c] = the matrix which depends on element geometry [D'J = the matrix property matrix (affective stress) [Aa'] = the change in effective stress vector [Ae] = the change in strain vector [K] =» the element stiffness matrix (A6) = the incremental element nodal deflections vector (Af) = the incremental element nodal forces vector (K ) = the pore pressure load vector n w, n Q, n^ = porosity : water, o i l and gas phase a w > a Q = temperature solubility constant : water, o i l H W , H Q = pressure solubility constant : water, o i l u = reference (atmospheric) pressure T a = reference temperature (288°K) T Q = i n i t i a l temperature (°K) Tl - f i n a l temperature (°K) AT = change in temperature (Ti - T ) (°K) o U Q = i n i t i a l absolute pore pressure Ae = volumetric strain (compression positive) They examined the response of cylindrical shaft in o i l sand on reduction of support pressure and the response of an element of o i l sand to one dimensional unloading. The latter case i s to simulate core samples of o i l sand l e f t in an unconfined state. Their results w i l l allow an assessment of the applicability and shortcoming of the theory: 1) The predicted expansion of the core on unloading is small compared with those measured in the f i e l d which i s 5-15% when core samples of o i l sand are l e f t in an unconfined state. 2) An iterative procedure i s required to obtain Au. CHAPTER 3 : STRESS-STRAIN MODEL 3.1 INTRODUCTION It is noted that the mathematical models described in Chapter 2 have quite a few shortcomings. Therefore, a more sophisticated and efficient undrained model was developed (Naylor 1973, Vaziri 1985) and i s incorporated i n a f i n i t e element programme. A total stress approach coupling the s o i l skeleton and pore fluids i s used . Pore pressure changes are computed from the constraint of volume compatibility between the s o i l skeleton and pore fluids under undrained conditions. Separate constitutive relationships for s o i l skeleton and pore fluids are required in the new undrained model. An incremental non-linear elastic and isotropic stress-strain model as described by Duncan et al is adopted for the s o i l skeleton. Depending on the component of the pore fluids, different formulations for the non-linear elastic and isotropic stress-strain relationships of pore fluid are derived. 3.2 DEVELOPMENT OF UNDRAINED MODEL The effective stress concepts (Terzaghi) in conventional s o i l mechanics seem to be applicable to determine the shear strength of o i l sand (Hardy and Hemstock). However, the pore pressures of fluids in the o i l sand have to be numerically evaluated i n the effective stress approach. Hence, a quantitative relationship between the magnitudes of stress release and pore pressure i s required. When a saturated s o i l mass is subjected to undrained loading, 18 stress change must be shared between the s o i l skeleton and pore f l u i d (Bishop and Eldin 1950). A theoretical expression for the relationship between total stress change and resulting pore f l u i d pressure change was derived by Skempton (1954) which is the Skempton's B equation. Using the f i n i t e element method, Christian (1968) introduced an effective stress approach for soils subjected to undrained loading, enabling resulting pore f l u i d pressure to be evaluated. Programmes incorporating this alternative are relatively inefficient. Naylor (1973) developed a more elegant approach which allows excess pore pressure to be computed explicitly in terms of material skeleton stiffness parameters and an independently specified pore f l u i d stiffness. However, Naylor only considers soils that are two-phase system - solids and water. Due to the presence of bitumen, free and dissolved gas i n o i l sand, the above mentioned approaches are inadequate to describe the undrained behaviour of o i l sand. Not u n t i l 1977, Harris and Sobkowicz developed an analytical expression incorporating a linear constitutive relationship for the o i l sand, to relate the change i n pore pressure to the change in total stress. Dusseault (1979) extended Skempton's one dimensional B equation to model the equilibrium behaviour of o i l sand. Byrne et al (1980) studied the behaviour of o i l sand by using f i n i t e element method. They extended Harris and Sobkowicz's model and incorporated i t into Christian's f i n i t e element formulation. As there are shortcomings of the approaches developed by Harris and Sobkowicz, Dusseault and Byrne et a l , a more sophisticated numerical model i s required. Vaziri adopted Naylor's approach and extended i t to model the undrained behaviour of o i l sand using f i n i t e elements. The stress-strain model for o i l sand behaviour i s based on the following assumptions. 1) Volumetric change of s o i l skeleton i s governed by the effective stress. 2) Liquids and gas in the voids are at the same pressure, i.e. effect of surface tension between the pore flu i d phases are neglected. 3) Free gas in the pores behaves i n accordance with classic gas laws with respect to pressure. Gas comes out from solution i n accordance with Henry's law. Henry solubility constants are constant. 4) The compressibility of s o i l grains i s negligible, and has no contribution to the volume changes. 5) The gas Is in the form of occluded bubbles inside the pore f l u i d . The total stress constitutive law may be written as: (Aa) = [D] (Ae) 3.1 where (Aa), (Ae) are the incremental total stress vectors and strain vectors respectively, and [D] is the material property matrix (total stress). Computation of element stiffness matrix, assembly and solutions for displacements proceeds along the standard lines. The analysis yields a total stress f i e l d . Since the s o i l skeleton and the pore fluids deform together when conditions are undrained, strains - in a macroscopic sense - are the same in each phase. Thus in addition to Equation 3.1 (Aa') = [D»] (Ae) 3.2 (Au) = K ( A e n + Ae 22 + Ae 3 3) 3.3a = K (Ae ) 3.3b a v v' in which prime means effective, [D'] is the material property matrix (affective stress), K is the apparent pore fluid bulk modulus and (Au is the pore pressure change vector. The actual pore f l u i d modulus, K^ , i s related to the apparent one as K f K = — 3.4 a n where n is the s o i l porosity. Derivations of equations 3.3 and 3.4 are presented i n Appendix B. Equation 3.3 may be expressed in a form compatible with equation 3.1 and 3.2 as: (A° f) = [D f] (Ae) 3.5 where (Aa f) = [Au A U A U 0 0 0 ] T Since the pore f l u i d cannot transmit shear, [D^] can be expressed in terms of apparent bulk modulus, [D F] = K f a where I 3 is as 3 x 3 matrix with a l l the elements equal to 1, whereas O3 are 3 x 3 null matrices. The principle of effective stress may be used to relate the changes in effective stress and pore presusre caused by the applied loads to the corresponding change in total stress (Ao) = (Aa') + (Aa f) 3.7 Substituting from Equations 3.1, 3.2, and 3.3 into 3.7, yields: [D] = [D«] + [D f] 3.8 The elastic material is now considered to be two phase, with the stiffness defined by effective stress moduli, E, B and pore f l u i d apparent builk modulus K^. The material properties of [D'] and are read in separately. They are combined in the programme automatically using equations 3.6 and 3.8. Thereafter, computation of element stiffness, assembly and solution of displacements and hence strains proceed along the standard lines. The effective stress and pore pressure are obtained by Equations 3.2 and 3.5 respectively. r3 °3 0 3 0 3 3.6 This approach allows stresses, porepressure and deformation in s o i l mass, with nearly incompressible to highly compressible pore fluids, to be evaluated by f i n i t e elements for undrained conditions. Naylor studied the response of undrained t r i a x i a l test on clay using this model. The end platens were rigi d and rough. But the computed excess pore pressures near the centre of the sample were in good agreement with the theoretical solutions (Cam-Clay). In the case of drained analysis, Equations 2.7 to 2.10 are used instead, assuming a l l the pore f l u i d pressures are known. 3.3 CONSTITUTIVE RELATIONS 3.3.1 Incremental Non-linear Elastic Soil Skeleton An incremental non-linear elastic-plastic and isotropic stress-strain s o i l model as described by Duncan et a l (1970) i s employed in this thesis. In this approach, two independent elastic parameters are required to represent non-linear stress-strain and volume change behaviour. These are usually the Young's modulus, E, and the Poisson's ratio, v. The shear modulus, G, and the bulk modulus B, are the more fundamental parameters because they separate shear or distortion and volume components of strain and would be the most desirable ones to use. However, the shear modulus is d i f f i c u l t to obtain directly in laboratory testings and for this reason Duncan et a l (1980) used the Young's modulus and the bulk modulus as their two parameters. The Young's modulus i s very similar in character to the shear modulus as both are a measure of distortional response. Therefore, E and B are used i n this thesis. The stress-dependent E and B are usually obtained from laboratory tests. A typical example for sand is shown in Figure 3.1. The distortional response can be reasonably approximated by modified hyperbolas (Konder) and the volumetric response in exponential form. They are expressed as follows: The tangent Young's modulus and R (l-sin4.) ( a i - 0 3 ) 2 E = [ l - — =-j — 1 E, 3.9 t L 2oj sin<j> J i where E ± = Kg P a ( p - ) n 3.10 a <J. = <(>!- A* log (—) 3.11 The tangent bulk modulus B t = K B P a < F - > m 3 ' 1 2 a where E i n = i n i t i a l Young's modulus = Young's modulus number = Young's modulus exponent 24 16 Fig.3.1 - S t r e s s - S t r a i n curves f o r d r a i n e d t r i a x i a l t e s t s on loose sand ( A f t e r Byrne and E l d r i d g e , 1982 ) = bulk modulus number = bulk modulus exponent = atmospheric pressure = major and minor principal effective stress = failure ratio = f r i c t i o n angle at confining stress of 1 atm = decrease in f r i c t i o n angle for a tenfold increase in confining stress The procedures for evaluating these parameters from laboratory tests are described i n detail by Duncan et a l (1980) and Byrne and Eldridge (1982). 3.3.2 Incremental Non-linear Fluid Modulus 3.3.2.1 General An incremental non-linear elastic and isotropic stress-strain f l u i d model i s employed here. Since the f l u i d cannot transmit shear, the stress-strain relations are defined only by one elastic parameter, the bulk f l u i d modulus, K^ , which is a measure of volumetric response. Before the derivation of K^ , Henry's law and Boyle's law must be mentioned because these laws govern the derivations. Pore fluids may be immiscible, miscible, or a combination of the two. Examples of both miscible and immiscible fluids w i l l be considered in this thesis. That i s , water undersaturated with air and carbon dioxide, water and bitumen saturated with gas (methane, (X^)-The f i r s t combination i s to describe the pore fluids of the gassy soils which Sobkowicz (1982) used in his laboratory testings. The second combination i s to similate o i l sand pore fluids. \ m P a °{> R f • l A<{>. If both free gas and liquids are present i n the pore fluids, and the gas is soluble in the pore liquid in a certain extent, the pore f l u i d compressibility w i l l be both pressure dependent and influenced by the solubility relationship. Hence, Boyle's and Henry's laws are appropriate for describing these volume and pressure relationships. 1) Boyle's law (Laidler et a l ) : The volume of a free gas i s inversely proportional to the pressure applied to i t when the temperature is kept constant. Mathematically, V a j 3.13 where V is the volume, P i s the absolute pressure 2) Henry's law (Laidler et a l ) : The mass of gas m dissolved by a given volume of solvent at constant temperature, i s proportional to the pressure of the gas i n equilibrium with the solution. Mathematically, m (dissolved gas) = H * P 3.14 where H is Henry's solubility constant, P is the absolute pressure. In other words, the volume of dissolved gas is constant in a fixed volume of solvent at constant temperature when the volume i s measured at P 27 (dissolved gas) = H * V (solvent) 3.15 Most gases obey Henry's law when the temperature is not too low and the pressure i s moderate. If several gases from a mixture of gases dissolve in a solution, Henry's law applied to each gas independently, regardless of the pressure of the other gases present i n the mixture. H is both temperature and pressure dependent, particularly for natural gases i n hydrocarbons (Burcik, 1956). Since the variation of H on pressure is not very significant, i t is assumed that H is independent of pressure in this thesis. 3.3.2.2 Partly Mlsclble Gas/Liquid Mixture Definitions: 1) Pore fluids compressibility, 3 3.16a where V is the volume of pore fluids, P i s the absolute pressure, assuming surface tension effects are neglected, 3.16b 28 2) Compressibility of a gas g 1 d v 8 where V is the volume of gas. a) Air/water Mixture Let the i n i t i a l volume of free gas and water i n a s o i l element be: V* and V fg w Thus the total volume of free and dissolved gas v i = + H V 3.18 Tg fg w For a change in pressure, the volume of water is assumed constant as the compressibility of water i s insignificant compared with the pore gas. Applying Boyle's law to the .,1 total volume of gas (Fredlung 1973, Sobkowicz 1982) in the element, f i P i V = V * — 3.19 VTg VTg P f Then the new volume of free gas after change in pressure f f f v = V - V V f g VTg Vdg 29 P Tg P. Vdg ( V f g + Vdg ) * Pf " V d g 3.20 Change of free gas in the s o i l element fg fg fg p fg dg' P f dg v f g - ( V L + V,!) * ra» 3.21 fg dg' P +AP By definition g i AP fg (v* + V, 1) fg dg' # 1 v i P ± +AP fg Also, by definition 1 * d V f V f dP 3.22 30 , dV. dV V LdP dP J V f L P + AP V w J 1 -S + SH w P, + AP + S e - 3.23a As AP approaches zero, 3 f P + Sg w 3.23B and b) Carbon Dioxide, Air/water Mixture By definition, the compressibility of pore fluids, 3 f, 1_ ^ f*f 3 f " V * dP , dVr dV = - — f — & + — V LdP dP J in which there may have air and carbon dioxide as free gas, V* g. Consider the phase diagram in Figure 3.2, the volume of solids i s assumed to be 1 unit, the volume of void i s e units according to definition of void ratios Therefore V* e - e fg = w V f " e V* H e + H „ e dg _ a w co2 w V f 3.25a 3.25a Substituting equations 3.25 into 3.24a, yields , ( e - e ) + H e + H _ e , _ 1 r w a w co2 w L Q -i „ 0 / , !, = — I „ , r g e 3.24b f e L P + AP w wJ As AP approaches zero, 32 •+ • w Air & C02 Water Solids v 9 Fig.3 .2 - Phase diagram for gassy s o i l 33 . ( e - e ) + H e + H _ e a _ 1 r w a w co2 w , . -1 - . . 3 = _ ^ + 8 e J 3.24c f e L P w wJ and K - 2-1 s £ c) Gas (Methane)/Bitumen and Water Mixture Again, by definition the compressibility of pore fluids 1 d V f 8 = - — * -P f V f dP i d v ^ d vx V f L dP dP J (V* + V* ) rr [ p g , ° 8 - V 0 - V. B. ] 3.26a v ^ P^ + AP ww b b J Consider the phase diagrams in Figure 3.3, the volume of solids i s assumed to be 1 unit, the volume of voids i s e units according to definition of void ratios Therefore (e - e w - e b) 3.27a Fig. 3 . 3 - Phase diagram for o i l sand where Hg/W> ^g/b a r e Henry's solubility constants of gas in water and bitumen respectively. Substituting Equations 3.27 into 3.26a, and AP approaches zero, yields, (e - e - e, ) + H , e + H e, B = I [ 2 * - ^ ! _ w g/b_J> + , f e L P w w b b J 3.26b and CHAPTER 4 : FINITE ELEMENT FORMULATIONS 4.1 INTRODUCTION Two types of formulations are presented herein which are suitable for modelling a variety of problems encountered in practice. Depending on the nature of the problems, they w i l l f a l l into one of the following categories: 1) Plane strain - 2 Dimensional, e.g. tunnels, shaft, etc. 2) Axisymmetric - 3 Dimensional, e.g. t r i a x i a l test, wellbore, etc. The s o i l i s modelled by isoparametric quadrilateral or triangular elements. Stress distribution formulations are added to cope with problems where plastic zones are developed during loading or unloading. 4.2. THE PLANE SRAIN FORMULATIONS 4.2.1 The Constitutive Matrix [P] Plane strain problems are characterized by the following two properties: 1) no deflection i n z direction 2) f i r s t derivative of the x and y deflections with respect to z are zero. Therefore where u, v, w are the displacements in x, y and z directions respectively with corresponding subscripts, 1, 2 and 3 respectively. From the generalized Hooke's law for incremental elasticity Aei j = [Aon - v (A022 + Ac 3 3)]/E 4.2a Ae 22 = [AO"22 ~ v ^ a i l + Ao"33)]/E 4.2b Ae 3 3 = [Aa 3 3 - v (Ao^ + Aa 2 2)]/E 4.2c Ae 1 2 = Ao 1 2/G 4.2d Ae 2 3 = Aa 2 3/G 4.2e Ae 3 1 = Aa 3 1/G 4.2f After substituting the conditions from equation 4.1b into equations 4.2, i t follows that: Aoi 3 = Ao 3 1 = 0 4.3a A o 2 3 = Aa 32 = 0 4.3b where Aa^ .. is the incremental shear stress with the direction indicated by the subscripts. With the above eliminations, the incremental stress and strain vectors become (Aa) = [Aau Ao"22 Ao-12] 4.4a (Ae) = [ A e n Ae 2 2 Ae 1 2] 4.4b The constitutive relations for a plane strain problem in total stress analysis are written (Naylor, 1973) as: A o n Ao"2 Aa 1 2 (l+v)(l-2v) 1 1 0 + K 1 1 0 a 0 0 0 r l-v V l-v 0 1 A e l l Ae Ae 22 12 0 0 d-2v) 2 4.5a or (Aa) = ([D'l + [D f]) (Ae) 4.5b where E = tangent Young's modulus v = tangent Poisson's ratio K = apparent tangent bulk f l u i d modulus cL The constitutive relations can also be written in another equivalent form which i s adopted in INCOIL. * A a 2 2 < A o 1 2 + K B' + G' B' - G' 0 B1 - G B* - G' 0 0 0 G' < 1 1 0 A e n 1 1 0 • A £22 0 0 0 A e 1 2 or (Aa) - ([ n']) + [D f]) (Ae) Where •a i _ 3B 2(l+v) G' = E 2(l+v) B = tangent bulk modulus E = tangent Young's modulus v = tangent Poisson's ratio K = apparent tangent bulk f l u i d modulus Equation 4.6 can also be written as (Aa) = [D] ( A S ) 40 where [D] is the constitutive matrix 4.2.2 The Strain Displacement Matrix [B] For isoparametric elements, the geometry (x,y) and displacement (u,v) are both expressed by te same shape functions and are approximated as: (*) = [N] («) 4.10 and O - [M] («') 4.11 where Nj 0 N 2 0 N 3 0 Uk 0 W to Ni 0 N 2 0 N 3 0 NJ ( 6 ) - [x! y! x 2 y 2 x 3 y 3 X l + y j -( 6 ' ) = [ u j V! u 2 v 3 u 3 v 3 uit v 4] in which ( 6 ) is the nodal coordinate vector and ( 6 ' ) is the incremental nodal displacement vector and N x = (l-s) (l-t)/4 N 2 = (l-s) (l+t)/4 N3 = (1+s) (l+t)/4 - (1+s) (l-t)/4 x y = nodal coordinates i n x and y directions 1, i respectively u. v. = incremental nodal displacement i n x and y 2» 3 directions respectively. s,t are local coordinates The incremental strain vector can be expressed in terms of displacement as follows: * 3u/3x Ae 1 2 = 3v/8y A e13 3u/3y + 3v/3x Substitution of u and v from equation 4.11 into equation 4.12 yield Ae Ae 11 22 A e13 3Ni 0 W 3Nj 3Ni 3y 3x u l v l u 2 v 2 u3~ v3 3N2 3x~ 3N2 W 3N 2 3y~~ 3N 3 3x~ 0 3N3 3No 3Nq 3N. 3x 3y 3x 3N^ 3x~ 3N 4 3y~ 0 3N\ 3N1+: 3 7 " I 4.10 Equation 4.10 can also be written in matrix notation as (Ae) - [B] (6) 4.13 where [B] is the strain displacement matrix However, the shape functions for isoparametric elements are defined with respect to the local coordinates s and t and therefore cannot be differentiated directly with respect to the global x, y axes. In order to overcome this d i f f i c u l t y i t i s necessary to obtain the derivatives of the two sets of coordinates and this can be achieved through the chain rule of partial differentiation. For plain strain problems, the derivatives are related as 3_ 3s 3_ 3t 3_ 3x 3_ 3y 4.14 where M -3x 3y 3s 3s 3x 3t 3y 3t ) is called the Jacobian matrix Hence the derivatives w*r»t» x and y can be expressed as derivatives w»r«t s and t as follows 3_ 3x 3_ 3y [j] -1 3s 3t 4.15 The strain displacement matrix in Equation 4.13 are evaluated numerically, using Gaussian quadrature over quadrilateral regions. The quadrature rules are a l l of the form n r. // f (s,t) ds dt - E I K K. f(s , t.) i-1 j-1 1 2 1 J where K^ , K\ are weighting functions and s^, t^ are coordinate position within the element. A 2 x 2 Gauss quadrature i s used to evaluate the strain displacement. 4.2.3 The Stiffness Matrix [K] The stiffness matrix for the force displacement relationship i s obtained by the principle of virtual work. For a virtual nodal displacement vector (6) g, the external work done, W(ext), by the external force vector (f) caused by virtual displacements i s written The virtual strain vector caused by the virtual displacements vector i s written as Hence the internal work done, W(ext), caused by the virtual strain i s written as as W(ext) - {6)1 (f) 4.16 (^) = [B] ( 6 ) e 4.17 W(int) = J (Ae) 1 (Ao) t dA 4.18a where t = thickness of the element A = area of the element Substituting equation 4.17 into 4.18a yields W(int) = J A [«]* [ B ] T [D] [B] [6] t dA 4.18b Applying principle of virtual work W(ext) = W(int) 4.19 Therefore (f) = / A [B] T [D] [B] [6] t dA 4.20 (f) = [K] (6) ' 4.21 T where [K] = / [B] [D] [B] t dA is called the stiffness matrix for the element. The global stiffness matrix, [K ], i s obtained by assembling a l l the element stiffness matrices together. The procedure of assembling the element matrix i s based on the requirement of 'compatibility' at the element nodes. This means that at the nodes where elements are connected, the values of the unknown nodal degrees of freedom are the same for a l l the elements joining at that node. The global force vector, [F], i s assembled by adding nodal loads of each of the elements sharing the node. Displacements are calculated using standard procedure (e.g. Gaussian elimination) to solve the simultaneous equations,[K ] [s] = [F],represented by the global stiffness matrix and the force vector. Strains can be computed from equation 4.13 or 4.36 after knowing the elements nodal displacements. After solving for displacements and strains, the effective stress and pore pressure can be computed from (Aa') = [D'j (Ae) and (Au) = [D f] (Ae) 4.3 THE AXISYMMETRIC FORMULATIONS 4.3.1 The Constitutive Matrix [D| Axisymmetric problems are characterized by the following properties: 1) Symmetry of both geometry and loading 2) Stress components are independent of*the angular (0) coordinates. Hence v = 0 4.22a A e 1 3 " A e 3 1 = 0 Ae 2 3 = Ae 3 2 = 0 4.22b 4.22c where u, v, w are displacements in the r, z and 0 directions with corresponding cubscripts, 1, 2 and 3 respectively. From the generalized Hooke's law for incremental e l a s t i c i t y A e n = [ l o u - v (Ao22 + Ao 3 3)]/E 4.23a Ae 2 2 = [ACT 2 2 - v (A033 + A a n ) ] / E 4.23b A e 3 3 " [ A o"33 ~ v ( A a n + Aa 2 2)]/E 4.23c A e i 3 = Ao13/G 4.23d Ae 1 2 = Aa 1 2/G 4.23e Ae 2 3 = Aa 2 3/G 4.23f After substituting the conditions from Equations 4.22b and 4.22c into Equations 4.23, i t follows that Aai3 = Aa 3i = 0 4.14a Aa 2 3 = Aa 3 2 = 0 4.15a where A o i j is the incremental shear stress with the direction indicated by the subscripts. With the above eliminations, the incremental stress and strain vectors become (Ao) = [ A a n Aa 2 2 Aa 3 3 Aa 1 2] 4.25 48 (Ae) = [ A e n Ae 22 A e 3 3 A e 1 2 ] 4.26 The constitutive relations for an axisymmetric problem in total stress analysis are written (Naylor, 1973) as Also Aa 11 A a 2 2 A a 3 3 A a 1 2 + K (l+v)(l-2v) 1 1 1 0 1 1 1 0 1 1 1 0 0 0 0 0 1-v V V 1-v V V 1-v 0 0 A e n A e 2 2 A e 1 2 0 0 0 l-2v 4.27a (Air) = [D-] + [ D J (Ae) 4.27b where E = tangent Young's modulus v = tangent Poisson's ratio K = apparent tangent bulk modulus The constitutive relations can also be written in another equivalent form which i s adopted in INCOIL. or A o n ] A 0 2 2 ] ^ 3 3 ] Ao 1 2 + K B'+G' B'-G' B'-G' 0 B'-G' B'+G' B'-G1 0 B'-G' B'-G' B'+G' 0 1 1 1 1 1 0 1 1 1 1 0 1 1 1 1 0 0 0 0 0 0 A e l l I A e 2 2 j 3 3 I A e 1 2 4.28a [Ao] - ([D«] + [ D f ] ) [Ae] 4.28b where B' = 3 B 2(l+v) G , . _ J 2(l+v) B = tangent bulk modulus E = tangent Young's modulus v = tangent Poisson's ratio K a = apparent tangent bulk f l u i d modulus Equation 4.28 may be written in matrix notation as (Aa) = [D] (Ae) 4.29 where [D] is the constitutive matrix. 4.3.2 The Strain Displacement Matrix [B] The isoparametric elements, the geometry (r, z) and displacements (u, v) are both expressed by the same shape functions and are approximated as: c •> r - W ( 6 ) 4.32 = [ N ] ( « ' ) 4.33 [N] = N x 0 N 2 0 N 3 0 N 4 0 0 N , 0 N , 0 N q 0 N b ( 6 ) = ( > ! zl r 2 z 2 r 3 z 3 rk Zlf] ( 6 ' ) = [u! v x u 2 v 2 u 3 v 3 U l + v„] in which (fi) is the nodal coordinate vector, ( 6 ' ) is the incremental nodal displacement vector and the shape functions are the same as given in Section 4.2.2. and ri» Z i u., v. J J nodal coordinates in x and y directions respectively, incremental nodal displacements i n x and y directions respectively. The incremental strain vector can be expressed in terms of displacements as follows: Ae i A e 2 2 11 I Ae 33 A e 1 2 s 3u 3r 3v 3r u r 3v 3u ) 3r 3z 4.34 Substitution of u and v from Equation 4.33 into Equation 4.30 yields Ae 11 Ae 22 A e 3 3 A e 1 2 3NX 3r 0 Nl r 3Nj 3z~ 0 3Ni 3z~ 3Nj I F 3N2 3r~ 0 N 2 — 0 — r 3N2 3z~" 0 3N2 0 3N3 3r~ 3N3 TT~ 0 N3 r 3N3 37" 0 3N3 3z~ 0 3N3 w 52 W 0 r r 0 3z~ 0 31^ 3r~ vi j u 2 | v 2 u 3 v 3 UI+ 4.35 Equation 4.35 can also be written in matrix notation as (Ae) = [B] (6) 4.36 where [B] i s the strain displacement matrix However, the shape functions for isoparametric elements are defined with respect to the local coordinates s and t therefore cannot be differentiated directly with respect to the global x, y axes. In order to overcome this d i f f i c u l t y i t i s necessary to obtain a relationship between the derivatives of the two sets of coordinates and this can be achieved through the chain rule of partial differentiation. _8_ ds a_ at 3_ 3r 3_ 9z 4.37 where [J] -3r 3z 3s 3s 3r 3z_ 3t St is called the Jacobian matrix Hence the derivatives wvt x and y can be expressed as derivatives w»r«t s and t as follows 3_ 3r 3_ 3z [j] -1 3_ 3s 3_ 3t 4.38 A similar Gauss quadrature mentioned in Section 4.2.2. is employed to evaluate the above [B] matrix numerically. 4.3.3 The Stiffness Matrix [K] The stiffness matrix for the force displacement relatinship i obtained by the principle of v i r t u a l work. For a virtual nodal displacement vector (5)e> the external work done, W(ext), by the external force vector ( F ) caused by the virtual displacements is written as W(ext) = [s)Te (f) 4.39 The virtual strain vector caused by the virtual displacements vector i written as (Ai) = [B] (6) e 4.40 Hence the internal work done, W(int), caused by the virtual strain i s written as W(int) = / (Ae) T (Aa) dV 4.41a where V = volume of the element Substituting Equation 4.40 into 4.41a yields W(int) = J Y [&fe [ B ] T [D] [B] [6] dV 4.41b Applying the principle of virtual work Wext = Wint 4.42 Therefore (0 " / v [B] T [D] [B] [6] dV 4.43 (f) = [K] (6) 4.44 where / [ B ] [D] [B] dV is called the stiffness matrix for the element. The procedure of obtaining element stresses and strains i s the same as described in Section 6.2.3. 4.4 LOAD SHEDDING (PLAIN STRAIN) (Byrne 1983) Problems arise when any element within the solution domain violates the failure c r i t e r i a (Mohr Coulomb). That i s , for unloading of a shaft or tunnel, a plastic zone usually developes adjacent to the shaft. The extent of plastic zone and hence volume changes are only approximated since the stress redistribution has not been considered during the formation of the zone. The analysis predicts the stress path ABC instead of ABD on unloading. But the stress state at C (Figure 4.1) violates the failure criterion (Mohr Coulomb). If load shedding technique is used, the overstress can be distributed to the adjacent elements by applying an appropriate set of nodal forces described herein and brings the stress path BC back to the correct BD. The overstress, Ax, in the element can be removed by subtracting the computed stresses by Aa^i, Aa 2 2 and Aa i 2 amount as shown in Figure 4.21b. where Aa. . have the same notation as those in Section 4.2 and 4.3 (Ao) = Aa Aa Aa 11 22 = [T] Aa] Aa- 4.45 12 Ax 13 56 Fig.4.1 - Stresses associated with load shedding where (T) = A03 = minor principal stress = 0 Ao"i = major principal stress = 0 1 - 0 3 tan 2 (45 + -|-) AT13 = principal shear stress = 0 The derivation of these stresses changes and the transformation matrix i s in Appendix E, The removal of these overstresses can be achieved by applying a set of nodal forces which is obtained by the principle of virtual work. The incremental nodal force vector causes a virtual displacement vector. Hence the external work done, W(ext), can be written as W(ext) = ( s ) ^ ( A f ) e 4.46 The incremental virtual strain vector caused by the virtual displacement vector i s cos26 1 4 . cos29 2 2 sin26 L_ cos26 2 2 cos29 2 sin29 0 0 is called the transformation matrix Ui)e = [ B ] ( « ) e 4.47 58 Therefore, the internal work, done, W(int), i s W(int) = /. (Ae) T (Ac) t dA 4.48a A. 6 6 Substitution of Equation 4.43 to 4.48 yields W(int) - / ( 6 ) T [B] ( A C ) t dA 4.48b Applying the principle of virtual work W(ext) = W(int) 4.49 Hence , [Af ] g = t / A [ B ] T [Aa] e dA 4.50 where [B] is the strain displacement matrix in Section 4.2.2. and [Ao] e is the stress vector shown i n Equation 4.41. The failed element w i l l have a stress change of Ao"n> A622 a n d Ao"i2- However, the computed stresses may not l i e on the failure envelope due to the application of the nodal forces. Therefore iterations may be required to bring this to the assigned tolerance. The loading shedding technique presented herein gives the same results compared with INCOIL (a' = constant). However, the number of m iterations required to bring the computed stresses back to the failure envelope is less. With the incorporation of load shedding technique, the sand skeleton i s modelled as a non-linear elastic-plastic porous material. The s o i l skeleton is coupled with the pore fluids in the undrained model. For undrained conditions, this model allows stresses, pore pressure and deformations of o i l sand masses to be evaluated by f i n i t e elements. 60 CHAPTER 5 : COMPARISONS WITH EXISTING SOLUTIONS 5.1 INTRODUCTION It is important to check the validity of the developed model before any major application. Two types of comparisons are presented herein which are suitable for checking the stress-strain model of the s o i l skeleton and also the newly developed gas law model. For the validation of the analytical procedure, the computed results are compared with elastic and elastic-plastic closed form solutions. The gas law model is validated by comparing computed solutions with observed data, such as expansion of o i l sand cores and t r i a x i a l tests on gassy so i l s . 5.2 COMPARISONS WITH THEORETICAL RESULTS 5.2.1 Elastic Closed Form Solutions The theory for elastic closed form solution was f i r s t developed by Timoshenko (1941). The plane strain solutions of stresses and displacements in a thick wall cylinder are presented herein: Stresses a 2P. - b 2P (P - P ) a 2 b 2 i o . i o a — —r~> ? — + —TTV -T\—9 5.1a r b^ - a z (b^ - a z) r z a 2 P. - b 2P (P - P ) a 2 b 2 °Q b2 - a 2 (b 2 - a 2) r 2 5 ' l b Displacement /i n N N ( a 2 p- ~ b 2 P )r . _ (l-2v) (1+v) v 1 o' 0 - ^TZTTI a' (1+v) ( P i ~ P o ) a 2b 2 5.1c E (b z - a z) r where P i = P r e s s u r e o n t n e inner surface of cylinder P = Pressure on outer surface of cylinder o J E = Young's modulus v = Poisson's ratio and a, b, r are defined i n Figure 5.1 The response of unloading a thick wall cylinder is being investigated. The stresses and displacements predicted by the programme are in remarkably good agreement with the closed form solution as shown in Figure 5.2. Hence the analytical procedure i s validated 5.2.2 Elastic-Plastic Closed Form Solution For a tunnel or shaft problem, the i n i t i a l state of stress w i l l be the same throughout the domain. As the support pressure of the tunnel drops, yielding w i l l occur i f the strength of the s o i l is exceeded. Yielding develops on the inside face f i r s t as a plastic annular zone and extends radially outward i f the support pressure i s reduced further. Hence theie exists a plastic and an elastic zone in the domain concentrically. The solutions of stresses and displacements Fig.5.1 - Thick wall cylinder 63 O o o Radii (r/r0 ) E = 3000 MPa v = 1/3 i n i t i a l stress : ar = a, = 6000 kPa final stress : or = 2500 kPa inside radius : r„ = 1 m Fig.5.2 - Stresses and displacements around circular opening in an elastic material o o o 2 o 00 CO j_> o > 0) 64 Ed - Tl - J - J •Q— 0 0 0 p c l o s e d form o o o program 1 1 1 _l L 1 1 1 i 3000 kPa 1/3 v -E = V B 0 e 40° i n i t i a l s t r e s s f i n a l s t r e s s i n s i d e radius 4 6 Radii (r/r0) 8 10 or •= oe = 6000 kPa ar = 1500 kPa r 0 = 1m Fig.5.3 - S t r e s s e s and displacements around c i r c u l a r opening i n an e l a s t i c - p l a s t i c m a t e r i a l for the plastic and elastic zone w i l l be quite different. The stresses and displacements for the elastic zone are just an extension of Equation 5.1-by setting b to i n f i n i t y . They may be written as Stresses a 2 a = P + (P, - P ) -^ V 5.2a r o v i o/ r"1 a 2 aa = P - (P, - P ) —jr 5.2b 9 o v i o' r z Displacement E r v i o where P^ = Pressure on inner surface of cylinder P = Pressure on outer surface of cylinder o J E = Young's Modulus v = Poison's Ratio Different investigators, Gibson and Anderson (1951), Ladanyi (1963), Vesic (1972) and Hughes et a l developed closed form solutions of stresses and displacements for the plastic zone. Hughes et al presented the more elegant solutions which are presented herein: A l l the sand is assumed to f a i l with a constant ratio of principal stresses, so that N = tan 2 (45 + |) 5.3 The equilibrium equation that must be satisfied i s do a* - al - ^ + ^ i - 0 5.4 dr r Substituting for a' from Equation 5.3, integrating and using outer boundary conditions of a^ = at r = R. a' in = (1-N) Jin (|) 5.5 R where a' = radial stress at r within the plastic zone r r o R = radial stress at the outer boundary, R, of the plastic zone = P (1 - sin<j>) 5.6 p l-sinfr R - a [-2. ( l - sin*)] 2 s i n * 5.7 I Equation 5.5 governs the distribution of the radial effective stress within the plastic zone. Continuity of stresses and displacements between the plastic and elastic zone must be maintained. Hence Equation 5.6 is substituted into Equation 5.2c, the displacement, U at the elastic-plastic zone contact i s UR E" R Po S i n * 5.8 Hughes et al also show that the displacements, u, within the plastic zone ^ = ( | ) n + 1 ( ^ | ) 5.9 where j n = tan 2 (45 + u = dilation angle The response of unloading a tunnel is investigated. The comparisons i n Figure 5.3 show that the analysis of stresses and displacements in elastic-plastic materials predicted by the programme is i n good agreement with the closed form solutions. The minor discrepancies are due to the limit of Poisson's ratio and the coarseness of the mesh. An upper limit of 0.499 i s adopted i n the programme (MHANS) to maintain numerical sta b i l i t y , whereas the actual value should be 0.5. The agreement i n displacements and the extent of the plastic zone also confirm that the load shedding technique in the programme (MHANS) is working properly. This i s a satisfactory check o the programme in drained analysis. Unfortunately, the load shedding technique i n INCOIL cannot b successfully tested. The element type in this programme is QM-6. Non-equilibrium of stresses arise in QM-6 elements at high Poisson's ratio (v > 0.4) i f the geometry of the elements i s non-rectangular. Since the elements in a f i n i t e element mesh for modelling plane strain shaft problems are not rectangular, non-equilibrium of stresses arise before load shedding is required. 5.3 COMPARISONS WITH OBSERVED DATA 5.3.1 One Dimenional Unloading of Oilsand This i s an opportunity for the new gas law model to be checked with some f i e l d data. Since the unloading is 1-D, the validity of the schematic spring analogy model (Figure 1.1) can also be demonstrated. Unconfined oilsand core taken from dr i l l e d holes swells by 5 to 15% of the original volume (Dusseault 1980, Byrne et a l 1980). Maximum expansion potential generally cannot be reached because expansion stops when there is adequate intercommunication of gas voids to permit flow of gas out of the sample. This leads to disruption of the s o i l fabric. This i s the equilibrium saturation point i n which gas becomes mobile and the maximum gas saturation value is 15% (Amyx et a l 1960). However, the total amount of expansion i s impossible to predict and can only be measured for individual cores. The core liners are specifically designed oversized to prevent the jamming of the core within the barrel. Radial expansion is assumed to be completed when the o i l sand core i s brought up from the d r i l l e d hole. The liners and steel containers are assumed to be rigid and frictionless so that only axial expansion of the core i s allowed. Therefore, the oilsand core can be modelled as 1-D unloading after recovery. The i n i t i a l l y high stressed core specimen is modelled as shown in Figure 5.4a. Vertical stress i s reduced from 750 kPa to zero. The stresses, pore pressure and displacements are shown in Figure 5.4b. It can be seen that the change in total stress i s apportioned between the effective stess and the pore pressure as suggested by the spring analogy model. Case A: Gas saturation pressure = 100 kPa. I n i t i a l l y pore pressure i s above gas saturation pressure and saturation remains 100%. Therefore loads come off from the pore fluids while effective stress remains f a i r l y constant because pore f l u i d is the s t i f f e r phase as shown in Figure 5.4b. As pore pressure drops below the gas saturation pressure, gas starts to evolve which causes the pore fluid to become flexible. Stress change w i l l be taken up by the s o i l skeleton on further unloading u n t i l zero effective stress. Further reduction on boundary load at zero effective stress is entirely accommodated by the pore fluids. Case Bt Gas Saturation pressure = 300 kPa. Because of gas exsolution the pore fluids start as a flexible phase so effective stress drops to zero with no appreciable change of pore pressure on unloading. Reduction of boundary stress beyond this stage is entirely taken by the fluid phases as the so i l skeleton essentially has no stiffness at zero effective stress as shown i n Figure 5.4b. It may be seen from figure 5.4b that there are no appreciable displacements when the effective stress i s positive. The displacements essentially come from unloading at zero effective stress. The total displacements upon total removal of vertical stress l i e within the range of 5 to 15% of the original length. 70 unload to 450 kPa 300 kPa •Im = 1565 s -= 1000 • = - 65° R f -' 0.9 n = • 0.243 o n = = 0.045 w H = 0.2 o H = 0.02 w S = 100% n m A<f> = 0.5 = 22 0.25 o 5.4 - Model for one dimensional unloading of o i l sand w CO c '5 u o o _ \ - ubb = 100 kPa - u b b = 300 kPa j L i i i i 150 300 450 600 750 External Load (kPa) 900 g.5.4 - Response of o i l sand to one dimensional unloading 72 The results show that the new stress strain model has the excellent capability of predicting undrained response of o i l sand. 5.3.2 T r i a x i a l Tests on Gassy S o i l s (Sobkowicz 1982) Performance of gassy s o i l on laboratory t r i a x i a l tests are reported by Sobkowicz (1982). A review of his work shows that the gas law model is generally correct. An appraisal of the predictive capabilities of the gas law model i s made by comparing predicted and laboratory observed response of the immediate pore pressure, the immediate (short term) B value, the equilibrium pore presusre, the equilibrium (long term) B value, and value of saturation and displacements. Since only test 11 i n Sobkowicz's thesis i s documented in detail, a comparison between predicted and observed response for this test i s made to evaluate the val i d i t y of the undrained model. Analysis is performed using the same i n i t i a l condition as test 11: S = 99.75%, n = 32.28%, a = 1403.3 kPa, u = 652.3 kPa and 3 = 9E-6 kPa - 1. The s compressibilkity of solid i s comparable with that of water (3 = 5 x 10 - 7 kPa - 1) so that the B (short term) value w i l l not equal to one even for f u l l saturation. These components are converted into parameters to be read in by the programme. The conversion and parameters are shown in Appendix F. The unloading sequence i s shown in Table 5.1. During any phase of the isotropic unloading test, pore pressure responses are predicted from the knowledge of s o i l skeleton and fluid compressibilities which are a function or effective stress, pore pressure, saturation and porosity. For short term reponse, H i s set to zero in Equation 3.26b since there is no time for gas exsolution. H . , ^ = 0.02 and H „ , ^ = 0.86 are used for the air/water co2/water equilibrium response when gas exsolution i s complete. The comparisons are summarized in Table 5.1 and are presented graphically i n Figures 5.5 and 5.6. They include: 1) predicted undrained response by the present undrained model 2) measured undrained response (Test 11, Sobkowicz 1982) A careful examination of Table 5.1 and Figures 5.5 and 5.6 show that the predictive capability of the gas law model is remarkably good, especially for the long term undrained response. The minor discrepancies are due to the loss of gas from the sample as the result of gas diffusion and leakage through the membrane. The observed immediate pore pressure are higher than the predicted values because of the time elapse (15 to 30 seconds) between reducing the total stress and taking the f i r s t reading. Thus, the predicted short term B is always higher than the observed ones. It can be seen that the stress reduction i s apportioned between s o i l skeleton and pore fluids, depending on their compressibilities. The sample i n Test 11 was i n i t i a l l y saturated with respect to air in water and undersaturated with respect to carbon dioxide in water. On unloading, during the f i r s t few phases, as P 0 / . < P < P . , _ , a small amount of gas exsolves for H co2/water air/water' & (air/water) = 0.02. The change of effective stress and pore pressure TABLE 5.1a A Comparison of Computed and Measured Results (Test 11, Sobkowicz) Total Short Term B Long Term B Salw ra "t.» ori ( %,) Phase Stress (kPa) Predicted Measured Predicted Measured Predicted Measured A 1322.4 0.897 0.694 0.523 0.606 99.65 99.67 B 1220.5 9.843 0.69 0.477 0.433 99.50 99^54 C 1112.1 0.781 0.68 0.378 0.403 99.30 99.38 D 978.2 0.705 0.64 0.021 0.011 98.90 99.11 E 883.6 0.628 0.548 0.022 0.016 98.61 98.93 F 766.4 0.584 0.508 0.024 0.034 98.22 98.93 G 654.9 0.548 0.482 0.026 0.07 97.80 98.27 H 559.0 0.536 0.50 0.031 0.155 97.38 97.90 Jl 457.3 0.569 0.590 0.129 0.21 96.40 J2 1 1 1 1 95.15 J3 1 1 1 1 94.66 Predicted: Results predicted by Programme Measured: Results measured in Test 11 by Sobkowicz TABLE 5.1b A Comparison of Computed and Measured Results (Test 11, Sobkowicz) Phase Measured Strains (%) Horizontal Vertical Volumeteric Porosity (%) Predicted Measured A 0.112 E-l 0.112 E-l 0.336 E-l 32.30 32.30 B 0.275 E-l 0.275 E-l 0.825 E - l 32.33 32.33 C 0.490 E-l 0.490 E-l 0.147 32.38 32.36 D 0.929 E-l 0.919 E-l 0.2757 32.46 32.42 E 0.124 0.124 0.372 32.53 32.46 F 0.167 0.167 0.501 32.62 32.53 G 0.214 0.214 0.624 32.72 32.61 H 0.261 0.261 0.783 32.80 32.69 Jl 0.376 0.376 1.128 J2 J3 .5.5 - Comparisons of predicted and observed pore pressure 77 are roughly the same because the compressibilities of both s o i l skeleton and pore fluids are comparable. This characteristic i s similar to those of unsaturated s o i l s . On further unloading, as P < P „ , _ , a large amount of gas exsolves because the high solubility co2/water' ° ° a J (H „ = 0.86) of carbon dioxide in water. This causes a sudden coz increase i n f l e x i b i l i t y of the f l u i d phase and hence most of the load is transferred to the s o i l skeleton. When the effective stress in the skeleton approaches zero, the f l u i d once again becomes the s t i f f e r phase, hence the B value rises to one. This i s the typical behaviour of gassy s o i l on unloading. The predicted and measured displacements are in remarkably good agreement. This indicates that the input parameters (Appendix f) and the ratio of the parameters, = 0.6 K^, and n = 2m = 0.5 (Byrne and Cheung) are generally correct. CHAPTER 6 - STRESSES AROUND A WELLBORE OR SHAFT IN OIL SAND 6.1 INTRODUCTION The response of a wellbore in o i l sand upon unloading i s considered because i t i s an important problem in o i l recovery in o i l sand. In general, knoweldge of the stress solutions around a borehole i s of great importance in several situations: 1) borehole s t a b i l i t y 2) hydraulic fracturing 3) production or injection A theoretical solution for stresses around a wellbore was developed by Risnes et a l (1982) and the equations are presented herein. Validation of the programme (MHANS) for drained analysis i s made by comparing computed response with the closed form solutions developed by Risnes et a l . A linear elastic-plastic constitutive relationship i s used for the above validation. Upon validation of the programme, i t was used to study the behaviour of a wellbore in o i l sand upon unloading. Undrained and drained analyses were performed to obtain the short term and long term response respectively. In the undrained analysis, the gas exsolution is assumed to be very fast relative to the construction of wellbore. In the drained analysis, the pore pressure profile i s estimated by using Dupuit's theory (Section 6.3.3). 6.2 GENERAL MODEL DESCRIPTION The wellbore under consideration i s supported by f l u i d pressure. As a model, a vertical cylindrical hole through a horizontal layer of o i l sand i s considered. The geometry, f i n i t e element mesh, and i n i t i a l conditions of the problem are shown in Figure 6.1 and 6.2. Loading and geometry are assumed to be symmetrical around the well axis. Only radial displacement after the i n i t i a l overburden loading are considered. These correspond to the assumption of axisymmetric and plane strain conditions. The sand formation i s assumed permeable, isotropic, homogeneous and i n i t i a l l y f u l l y saturated. The material is assumed elastic-perfectly plastic and obeys Mohr Coulomb failure criterion. Only stress solutions for > at the elastic-plastic boundary w i l l be investigated. 6.3 THEORETICAL SOLUTIONS FOR STRESSES AROUND A BOREHOLE A closed form solution for streses around a well, using a linear elastic-plastic stress-strain relationship, can be obtained from Risnes et al (1982). The derivations of the stress solutions follow that of Risnes et a l (1982) with two additional assumptions. 1) Insitu state of stress is considered to be isotropic i n i t i a l l y , i.e. a = o. = a . r o z 2) The Mohr Coulomb failure criterion in a porous material i s f = a{ - 2S tana - 03 tan 2a. = 0 6.1 Fig.6.1 - Outline of the problem Fig. 6.2 76 elements , 78 nodes Finite element mesh for wellbore problem where S is the cohesion intercept (or apparent cohesion) a is the f a i l u r e angle, i.e. + -|— <f>' is the internal f r i c t i o n angle These symbols are also explained i n Figure 6.2b 6.3.1 Stresses 6.3.1.1 Stresses In Elastic Zone The stresses around a hole in an elastic thick wall cylinder, with porous material saturated with f l u i d , may be written as follows: R2, 'r = °ro + <0ro " 0 r i > F ^ R T [ l - ( ^ ) 2 ] o i o i ' 2(l-v) £n(R . R? » -°> " i ' *n(R o/R.)-R^* R °9 * °ro + (°ro " °ri> R ^ R T t 1 + o i ( P _ P ) 1-2* ^ o i ; 2(l-v) o i v o i ' [*n (^) - 1]) 6.3 Fig.6.2 b - Mohr Coulomb failure envelope 6.4 The procedure for obtaining these stress solutions is given in Appendix G. 6.3*1.2 Stresses in Plastic Zone As long as f = 0 (Equation 6.1), a plastic zone w i l l start to develop at the borehole wall, and then expanding i n size as support pressure is decreased. Equation 6.1 w i l l apply within the plastic plastic zone is considered, the elastic stress solutions at this boundary are given by Equation 6.2 to 6.3, with R = r = R., a = a o i rc r and P = P.. With the assumption of no fluid flow (i.e. P = P ) and c i r c o' , R Q » R^, the stress solutions from Equations 6.2, 6.3 and 6.4 may be written as: zone. If the stress state at the boundary between the elastic and a = a 6.4 rc rc = 2a - a 6.5 ro rc a = a 6.6 zc zo The elastic solutions (Figure 5.2) show clearly that the radial stress w i l l be the smallest at the boundary between elastic and plastic zones, and following the assumption (1) of i n i t i a l isotropic stress, the state of stress at the elastic-plastic boundary w i l l be a < a < an . rc zc 9c The stress solutions within the plastic zone may be derived by combining Coulomb failure criterion (6.1) and the equation of equilibrium. do a -a -* + JL-±=0 6.7 dr r The stress solutions for the plastic zone may be written as: For R. < r < R i c a = P. + ir^r; Zn |- + -^ (2S tana - -^r-) r i 2irhk R^ t 2Trhk [(I-)* - 1] 6.8 For R. < r < R i c °9 - p i + 2¥ht c1 + * n ir> + F <2S t a n a - iSk> [(t + 1) (I-)11 - l] 6.9 R i For R, < r < R. 87 °z " Pi + 2^hk <X + *N R7> + I (2S tan° - zSfc t ( t + L) < ! / 6.10 For R, < r < R b c a - (P + «a- £ n f-) + v l i l - + Cl-y)(l-2v) z i 2Trhk R ' 2irhk 1-v 6.11 <|)1 where t = tan 2a - 1, a = 45° + -j-u = fluid viscocity The procedure to derive Equations 6.7 to 6.11 is given in Appendix H. 6.3.1.3 Radius of the P l a s t i c Zone Radius of the inner Plastic Zone R, b At the boundary of inner and outer plastic zones, the tangential stress and vertical stress given by equations 6.9 and 6.11 are equal. Setting Equations 6.9 equal to 6.11 and r = R^ , yields. a x ( ^ V + a 2 = 0 i 6.12 88 where a i - £ (2S tana - [(t + 1) - v (t + 2)] a 2 = ( 1 . v ) J i J L - ( 1 ^ ) ( f f - P ) * 2irhk l - v v zo cr 7 ( 2 S t a n a ' 2$!k> ( 1 " 2 v ) Radius of the Entire Plastic Zone There are two requirements that must be satisfied at the boundary of elastic and plastic zones 1) Mohr Coulomb criterion must hold 2) Continuity of radial stress Inserting radial stress from Equation 6.8 and tangential stress from 6.3 into Mohr Coulomb failure criterion 6.1, the resulting equation for the radius of plastic zone R^ is 2 R R R R b x R c C £n ^ + b 2 £ n -2. + b 3 R 2 + b l + R 2 £ n |_ + b,. R 2 £ n _ | R R R R + b 6 in Y~ • in jp- + b 7 in ^ + b Q £ n _° + b g = 0 c i i c 6.13 89 where bx = (2S tan a - R ^ b 2 - - ^ r1 Ci R 2 • t 1 o B 3 * " 2(1=*) ( P o - V = 1 " 2 v 2 (1-v) 2uhk b 5 2irhk b - - b^ R 2 b8 = v + <PQ + V - 2 p i + 1 2 5 tan a t+2 uq i t 2TThkJ b 9 = - C 3 R 2 6.3.2 Stability It is noted that the radial stress component in Equation 6.8 consists of two r-dependent terms, one logarithmic and one to the power 90 of t. The last term w i l l become dominant when the exponent t has a value greater than about two. Setting C = 7- (2S tan a - ^rtjr) R ~ t 6.14 ° t 2irhk i If C is positive, radial stress in the plastic zone wi l l increase with r, and combined plastic-elastic solutions are possible. But when the flow rate q is large enough to cause C to become negative, radial stress w i l l decrease with increasing distance r, and combined solutions are not possible. Hence, there exists a stability criterion. C > 0 6.15 with the limit •X—T- = 2S tan a 6.16 2irhk This study concentrates on o i l sand which has S = 0. If the wellbore is supprted only by f l u i d pressure, equation 6.16 indicates that instability arises when flow into the wellbore occurs. 6.3.3 Pore Pressure Profile When steady-state conditions around the wellbore have been reached, the pore pressure i n the s o i l elements may be estimated i f the piezometric surface is known. Dupuit developed a theory which enables the quantity of steady-state seepage and the piezometric surface around a well to be evaluated. His theory is based on three assumptions: 91 1) the hydraulic gradient is equal to the slope of the free surface and i s constant with depth, 2) for small inclinations of the line of seepage the streamlines may be taken as horizontal. 3) The permeability of the s o i l Is constant With the terminology in Figure 6/2a. the flow when steady state conditions exist i s given by: Q - w k ( h 2 2 - h L 2) J^J^ 6.17 w and the location of the free surface is h 2 2 - h! 2 h 2 = h i 2 + . / T i ,—r- Jin (—) 6.18 1 £n(R/r ) v r ' w w There Is a controversy about the location of the piezometric surface predicted by Dupuit's theory, especially in the v i c i n i t y of the well. This is because the surface of seepage is omitted in Dupuit's prediction. But the problem i s modelled as a disk of sand below the seepage surface (Figure 6.3b). Only radial flow is assumed in the sand disk. Hence equation 6.16 i s accurate enough to estimate the pore pressure profile. h, ^^— f surface of seepage Flow 7"~ F.F. Mesh Fig.6. 3 - Idealised flow to a wellbore 93 6.4. Comparisons of Predicted Response and Closed Form Solution The response of unloading a borehole, with linear elastic-plastic porous material, i s investigated. A drained analysis was performed. The i n i t i a l and f i n a l conditions of the problem are shown i n Figure 6.4(a). When the f l u i d support pressure i s higher than the i n i t i a l pore pressure, no flow from the borehole into the sand formation i s assumed. Since the f i n a l f l u i d support (4100 kPa) i s higher than the i n i t i a l pore pressure, flow into the borehole i s not considered herein. The stress solutions computed by the programme are i n good agreement with the closed form solutions mentioned in Section 6.3, and shown i n Figure 6.4. The radius of the entire plastic zone i s small which shows that the borehole is stable at the fi n a l fluid support pressure of 4100 kPa. 6.5.1 Undrained Response The undrained non-linear elastic-plastic model wi l l now be used to study the response of a wellbore on unloading. The f i n i t e element mesh and i n i t i a l conditions are shown in Figure 6.1b. Only the long term undrained response w i l l be investigated i n this thesis because this condition is f e l t to be more r e a l i s t i c (t * 0) and i t was also shown that the long term undrained condition i s more c r i t i c a l than the immediate response (Sobkowicz, 1982) in terms of st a b i l i t y . The wellbore i s unloaded by decreasing the total stress at the wellbore wall, and the stress solutions and displacements are shown graphically in Figures 6.5. A careful examination of these figures indicates some interesting results: 94 W <=> 10 20 30 40 Radii (r/r0) 50 E = 60 MPa v = 0 . 4 5 + = 35° i n i t i a l s t r e s s f i n a l s t r e s s i n s i d e r a d i u s o u t s i d e r a d i u s Or = Or = l b = R = 0 8 = Oz 600 k P a 0 . 1 2 5 m 6 . 7 m = 4 0 0 0 k P a ( a ) F i g . 6 . 4 - - S t r e s s e s a r o u n d a w e l l b o r e i n a n e l a s t i c - p l a s t i c m a t e r i a l 95 O >< CL o o 00 o o C D CO CO CD u CO (D > o o o o CVJ o C d ° CO > o o CM O 0) V — I w ° o o closed form program -o-J L J I 10 20 30 40 Radii (r/r0) 50 (b) Fig.6.4 - Stresses around a wellbore in an elastic-plastic material Fig.6.5 - Undrained response of a wellbore in o i l sand on unloading 97 F i g . 6 . 5 - Undrained response of a w e l l b o r e i n o i l sand on unloading 98 Fig.6.5 - Undrained response of a wellbore in o i l sand on unloading 99 1) The support pressure can be reduced below the i n i t i a l pore f l u i d pressure and the wellbore Is s t i l l stable; instability is defined when large displacements start to occur at the wellbore wall. 2) Instability occurs at a support pressure of approximately 2500 kPa. 3) The size of the plastic zone remains small as long as the support pressure is higher than 2500 kPa (Figure 6.5 b). Once the support pressure drops below 2500 kPa, the size of the plastic zone increases rapidly (Figure 6.5 c). 4) Pore f l u i d pressure changes only occur in the plastic zone. The evaluation of the fluid pressure response e P depends on volumetric strain Ae = Ae + Ae . But in the v v v elastic zone, Ae e = - A e^ and Ae e = 0 so that Aee= 0, r 6 z v and hence no change in pore fluid pressure is predicted. 5) Once instability has been reached, a liquid zone with zero effective stress associated with a large plastic zone w i l l form adjacent to the wellbore. These zones wi l l extend into the sand formation rapidly upon further reduction of support pressure, leading to large displacements. 100 6.5.2 Drained Response For the drained condition, the pore fluid pressure is assumed to be known. Dupuit's theory described i n Section 6.4.3 i s adopted to estimate the pore pressure profile around the borehole in this study. Two typical pore pressure profiles are shown in Figure 6.7 with R = 100 and 150 m with the fl u i d support pressure fl u i d at 3200 kPa. There are some intermediate pore pressure profiles between support pressure of 3500 kPa to 3200 kPa, depending on the number of increments on unloading, but they are not shown here. When the f l u i d support pressure is above 3500 kPa, no flow from the wellbore into the sand formation is assumed. The wellbore i s unloaded i n the same manner as for the undrained analysis. The results of the stress solutons are shown i n Figure 6.7. A careful examination of these results indicate some interesting points: 1) To maintain borehole s t a b i l i t y , the support pressure cannot be reduced to less than the i n i t i a l pore pressure; instability is defined as large displacements start to occur at the wellbore wall. 2) The stress solutions only dif f e r by a few percent when the input pore pressure profiles are generated by using R = 100 and 150 m. Therefore, only one set of stress solutons is presented here in Figure 6.7. 3) Once st a b i l i t y has been reached, a liquid zone with zero effective stress associated with a large plastic zone w i l l extend into the sand formation rapidly upon further Fig. 6 . 6 - Pore pressure profile around a borehole 102 200 300 400 500 600 700 Support Pressure (kPa) (X101 ) (a) Fig.6.7 - Drained response of a wellbore in o i l sand on unloading 104 reduction of support pressure, leading to large displacements. 6.5.3 Implications of Undrained and Drained Analyses The analyses show that there are limits on the fluid support pressure reduction in order to maintain borehole s t a b i l i t y . Comparisons of both analysis are made at support pressure of 3500 kPa and 3200 kPa. 3500 kPa is the c r i t i c a l pressure below which instability occurs in drained condition. It i s noted that in Figure 6.9 the plastic zone in the undrained analysis i s much smaller than the one in drained analysis. At a support pressure of 3200 kPa, the borehole i s obviously unstable under drained conditions (Figure 6.9), showing a large plastic zone and consequently large displacements. But the borehole only exhibits a small plastic zone at this support pressure (3200 kPa) under undrained conditions (Figure 6.9). This i s because the pore pressure around the wellbore i s lower in the undrained case, which results in a higher effective stress. Consolidation i s the process which bridges the f u l l y undrained and drained conditions. An interesting point is that the pore pressure w i l l increase around the borehole during consolidation, leading to lower effective stress. Hence, the long term drained condition is less stable than the undrained condition. 6.6 Application to Oil Recovery For o i l production, the f i n a l support pressure must be reduced below the in-situ pore f l u i d pressure. Based on the undrained analyses with the wellbore supported only by flu i d pressure, the support g.6.8 - Comparisons of undrained and drained response of a wellbore in o i l sand at a support pressure of 3500kPa 106 Fig.6.9 - Comparisons of undrained and drained response of a wellbore in o i l sand at a support pressure of 3200kPa 107 pressure can be reduced below the in - s i t u pore f l u i d pressure and the wellbore is s t i l l stable. This allows the construction of the wellbore and i n i t i a l reduction of f l u i d support pressure below in-situ pore flui d pressure. However, for the drained condition, the wellbore becomes unstable which causes collapse of the well and hence no o i l production. Instability results i n the formation of large liquid and plastic zones (Figure 6.9) around the wellbore. Since the permeability i n the liquid and plastic zones are higher due to the expansion of sand skeleton, i t is desirable to have these zones around the o i l production well. This effectively increases the diameter of the well. To enhance o i l production and maintain wellbore s t a b i l i t y , a screen may be installed to provide effective support pressure after the liquid and plastic zones have formed. The three dimensional and viscous effects have not been considered in the analysis, however they may help in stabilizing the wellbore. 108 CHAPTER 7 : SUMMARY AND CONCLUSIONS A new stress-stress relationship for modelling the undrained response of o i l sand has been presented. An analysis which couples the s o i l skeleton and pore fluids i s used. The pore pressure changes are computed from the constraint of volume compatibility. Separate stress-strain models are required for both s o i l skeleton and pore fluids in this analysis. The conventional hyperbolic stress-strain model described by Duncan et a l i s adopted for the s o i l skeleton. The pore fluids stress-strain relationships are formulated on the basis of ideal gas laws. The developed model i s incorporated into a f i n i t e element programme for analysing the deformation behaviour of gassy soils (e.g. o i l sand). Upon validation i n Chapter 5, i t i s shown that the undrained model is capable of predicting the response of unsaturated to gassy s o i l s . Naylor has shown that this model can be used to predict the response of saturated so i l s . For non-rectangular QM-6 elements, equilibrium cannot be achieved for Poisson's ratio values greater than 0.4, but higher order element can remedy this. In the study of the response of wellbore i n o i l sand upon unloading, the fluid support pressure can be reduced below in-situ pore f l u i d pressure under undrained condition and the wellbore is s t i l l stable. However, for the drained conditin, the fluid support pressure cannot be reduced below the in-situ pore f l u i d pressure i n order to maintain wellbore s t a b i l i t y . For o i l production, the f l u i d support must be reduced below in-situ pore pressure which results in formation of large liquid and 109 plastic zones. It i s desirable to have these zones around the wellbore because the diameter of the well is effectively larger. To enhance o i l production and maintain wellbore s t a b i l i t y , a screen may be installed to provide effective support pressure after the liquid and plastic zones have formed. BIBLIOGRAPHY Atukorala, U.D. Finite Element Analysis of Fluid Induced Fracture Behaviour i n Oilsand. M.A.Sc. Thesis, U.B.C. 1983. Biot, M.A. General Theory of Three-Dimensional Consodilation. Journal of Applied Physics, Vol. 12, February 1941. Bishop, A.W. The Influence of an Undrained Change in Stress on the Pore Pressure in Porous media of Low Compressibility. Geotechnique, V.23, N.3. 1973. Bishop, A.W. The Influence of System Compressibility on the Observed Pore Pressure Response to an Undrained Change in Stress in Saturated Rock. Geotechnique, V.26, 1976. Brooker, E.W. Tar Sand Mechanics and Slopes Evaluation. 10th Canadian Rock Mechanics Symposium, Department of Mining Engineering, Queen's University, 1975. Byrne, P.M. and Cheung, H. Soil Parameters for Deformation Analysis of Sand Masses. Soil Mechanics Series No. 81, Department of C i v i l Engineering, U.B.C. 1984. Byrne, P.M. and Eldridge, T. A Three parameter Dilatant Elastic Stress-Strain Model for Sand. Soil Mechanics Series No. 57, Department of C i v i l Engineering, U.B.C. 1982. Byrne, P.M. and Janzen. W. Soilstress. Soil Mechanics Series No. 52, Department of C i v i l Engineering, U.B.C. 1981. Byrne, P.M. and Janzen, W. INCOIL. Soil Mechanics Series No. 80, Department of C i v i l Engineering, U.B.C. 1984. Chatterji, P.K., Smith, L.B., Insley, A.E. and Sharma, L. Construction of Saline Creek Tunnel in Athabasca O i l Sand. Canadian Geotechnical Journal, 1, 197 9. Christian, J.T. Undrained Stress Distribution by Numerical Methods. J.S.M.F.D., A.S.C.E., S.M.6, 1968. Cook, R.D. Concepts and Applications of Finite Element Analysis. 2nd Edition, John Wiley & Sons, 1981. Duncan, J.M. and Chang, C.Y. Nonlinear Analysis of Stress and Strain in Soils. J.S.M.F.D., A.S.C.E., S.M.5, 1970. Duncan, J.M., Byrne, P.M., Wong, K.S. and Mabry, P. Strength, Stress-Strain and Bulk Modulus Parameters for Finite Element Analysis of STresses and Movements in So i l Masses. Report No. UCB/GT/80-01, August 1980. Dusseault, M.B. Sample Disturbance i n Athabasca O i l Sand. Journal of Canadian Petroleum Technology, V.19, N.2, Dusseault, M.B. Undrained Volume and Stress Change Behaviour of Unsaturated Very Dense Sands. Canadian Geotechnical Journal, Volume 16, 1979. Dusseault, M.B. and Morgenstern, N.R. Shear Strength of Athabasca O i l Sands. Canadian Geotechnical Journal, V.15, N.2, 1978. Dusseault, M.B. and Morgenstern, N.R. Characteristics of Natural Slopes in the Athabasca O i l Sands, V.15, N.2, 1978. Florence, A.L. and Schwer, L.E. Axisymmetric Compression of a Mohr-Coulomb Medium Around or Circular Hole. International Journal for Numerical and Analytical Methods in Geomechanics, Vol. 2, 1978. Fredlund, D.G. Density and Compressibility of Air-water Mixture. Canadian Geotechnical Journal, V.3, 1976. Geertsma, J. Some Rock-Mechanical Aspects of O i l and Gas Well Completions. Paper EUR 38 presented at the 1978 SPE European Offshore Petroleum Conference and Exhibition, London, October 24-26. Hardy, R.M. and Hemstock, R.A. Shear Strength Characteristics of Athabasca Oil Sands. K.A. Clark Volume, Alberta Research Council, Information Series No. 45, 1963. Harr, Groundwater and Seepage. McGraw-Hill Harris, M.C., Poppen, S. and Morgenstern N.R. Tunnels in Oil Sand. Journal of Canadian Petroleum Technology, 1979. Harris, M.C. and Sobkowicz, J.C. Engineering Behaviour of O i l Sand. The O i l Sands of Canada-Venezuela, 1977. CIM Special Volume 17. Hughes, J.M.O., Wroth, CP. and Windle, D. Pressuremeter Tests in Sands. Geotechnique 27, 1977. Naylor, D.J. Discussion. Proceedings of the Symposium on the Role of Plasticity in Soil Mechanics, Cambridge, 13-15, p. 291-294. September 1973. Pasley, P.R. and Cheatham, J.B. Rock Stresses Induced by Flow of Fluids Into Boreholes. Soc. Pet. Eng., J., p. 85-94, March 1963. Risnes, R., Br a t l i , R.K. and Horsrud, P. Sand Stresses Around a Wellbore. Society of Petroleum Engineering Journal, Col. 22, No. 6, 1982. Smith, L.B. and Burn, P.M. Convergence-Confinement Method of Design for Shafts and Tunnels in Oilsands. Conference on Applied Oilsands Geoscience, Edmonton, Alberta, 1980. 113 Sobkowicz, J.C. The Mechanics of Gassy Sediments. Ph.D. Thesis, University of Alberta, C i v i l Engineering Department, 1982. Timoshenki, S. Strength of Materials. Vol. 2, 1941, Van Nostrand, New Yord. Todd, . Ground Water Hydrology. John Wiley and Sons, Inc. Vaziri, H. Forthcoming Ph.D. Thesis, Department of C i v i l Engineering, U.B.C. 1985. APPENDIX A The constitutive relationship may be written: (Aa') = [D'] (Ae) A.l in which (Ao 1) is the incremental stress vector, (Ae) is the incremental strain vector and [D*] is the incrmental effective stress strain matrix. The incremental strain vector i s related to the nodal displacements by: (Ae) = [B] (S) A.2 in which [B] is a matrix that depends on element geometry. By the principle of vi r t u a l work, the external work done by the virtual displacement is equal to the internal work done by the increment of vi r t u a l strains: (5) T ( f) - /A U e ) T (Ao') dA + J . (Ae")T (}) Au dA A.3 in which (f) = (Aa') = (Au) = the element force vector the element incremental effective stress vector the element incremental pore pressure vector Substituting for (Ae) and (Aa') from Equations A.l and A.2, (6) (f) = (6)T [B] T [D'] [B] (6) A o + (8)T [ B] T(*) Au A A.4 e in which A is the area of the element with unit thickness, e Rearrangement of Equation A.4 yields [K] (6) = (f) - ( K j (Au) A.5 in which [ K ] = [ B ] T [ D ] [ B ] A^ [ K ] = [ B ] T (1) A w L 1 ^o' e where [K] is the element stiffness matrix and (& w) is a load vector associated with the pore pressure. APPENDIX B Assumptions: 1) the volume^ of solids i s 1 unit, then the volume of voids is e units by the definition of void ratio, 2) the solids are incompressible. The undrained response of the element under a change on external pressure w i l l be Aa ' Skeleton (Ae ) - m v'SK B S K B.l Fluid (Ae )_ = B.2 v f in which (Ae ) , (Ae ), = volumetric strain: Skeleton, Fluid v SK v f ' Bg^, = bulk modulus: Skeleton, f l u i d Aa ' = mean effective stress change m Au = pore pressure change Since the s o i l skeleton and the pore fluids deform together when the conditions are undrained, i n additions to Equations B.l and 117 Compatibility e (A£ y) f - (1 + e) t^)^ B.3a or (Ae ) = - (Ae ) B.3b v f n v SK in which e is the void ratio and n is the porosity. Substituting B.3b into B.l, the pore pressure change, Au, may be written as: A u = K a < A ev> S K B ' 4 K f i n which K = — , the apparent bulk f l u i d modulus a n In f i n i t e element analysis Au = K (Ae ) a v = K a ( A e n + Ae 2 2 + Ae 3 3) 118 APENDIX C Poisson's ratio is equal to 0.5 in undrained analysis theoretically. This means that bulk modulus i s i n f i n i t e and numerical instability w i l l arise (this sentence does not sound right). Therefore the default of Poisson's ratio, v, i s always less than 0.5 (e.g. 0.495) to maiontain stability and accuracy. In the total stress model developed i n Section 3.2, i t i s the combined Poisson's ratio, for matrix [D] that controls the overall numerical st a b i l i t y , not just v. The elastic moduli i n matrix [D'] are related as G = E C.I 2(l+v) B = E C.2 3(l-2v) Assuming the elastic moduli in matrix [D] = [D'] + [ D F ] are related as G cb C.3 B cb C.4 119 in which suffic cb means combined, G = G , as pore f l u i d does not cb transmit shear. Just consider the direct stress terms i n the constitutive relationship A a n L + K a L*M+K a L*M+K a A e n Aa 1 2 s s L*M+K a L+K a L*M+K a Ae 2 2 Ao-33 L*M+K L a L*M+K a L +K a A E 3 3 C.5 4 M U. T E(l-V) i n w h i c h = (l+v)(l-2v) 1-v Substituting C.2 into C.5, yields, Aff H P+K a Q+K a Q+K a A e u Ao 2 2 = Q+K a P+K a Q+K a Ae 2 2 C.6 A 0 - 3 3 Q+K a Q+K a P+K a A E 3 3 120 v „ 3B(l-v) in which P = ( 1 + v / Q = 3 B V 1+v Adding the direct stress in Equation C.6 and rearranging T ( A q " + A g22 + A g33) = B(l-v) + 2 B V +K a (Ae^i + Ae 22 + ^£33) l-v 1+v = B C.7 cb Eliminating E from Equations C.3 and C.4, yields cb 3B c b - 2G Vcb = 6B°. + 2G 0 , 8 cb Substituting C.7 into C.8, rearranging 121 APPENDIX D In the incremental elastic method, two iterations are performed to obtain the tangent moduli of the s o i l skeleton. Hence two iterations are also employed to evaluate the tangent bulk fluid modulus in order to make the procedure compatible. The parameters at the end of previous increment w i l l be used in the analysis of the present increment, and then updated at the end of this increment. The tangent bulk flu i d modulus and the procedure to update the parameters are shown herein: Equation 3.24e and 3.26b are programmed as: e-n (1+e) + H n (1+e) + H „n (1+e) 0 l r w ww co2 w 3 f i - - L p + 3 N (1+e)] D.l w w J e-n (1+e) - n (1+e) + H,n, (1+e) + H n (1+e) „ _ 1 r o w b b ww 8f2 - e L + W 1 + e > ] D.2 + V w The parameters in Equation D.l and D.2 may be updated according to the following foaulae based on the assumptions: volume of solids i s 1 unit bitumen, water and s o i l solids are incompressible Ae = (1+e) Ae D.3 v e^ = e i + A e ^"^ S ce £ = S.e, D.5 f f 1 i n (1+e) + Ae -8 D.6 1+e + Ae n (1+e) < - rzzn D. 7 b 1+e+Ae n (1+e) w "w 1+e+Ae n' = r£—n D.8 P = u + Au D.9 123 APPENDIX E Refer to Figure 4.1b °1 = 9 + °12 cosec 29 E.l an + a22 03 = ^ ~ °12 cosec 28 E.2 x 1 3 = 0 E.3 The overstress in the element is removed by reducing AOj, A a 3 , A T 1 2 as follows AO! = Oi - a3 tan 2 (45° + |-) E.4 Ao-c, = 0 E.5 A x 1 2 = 0 E.6 These change in principal stresses can be expressed in terms of stresses i n x-y space: AO} + Aa 3 AO} - Aa 3 A o l l = 0 + * cos 29 E.7 124 Aoi + Aa-i Ao"i - Aero Aa 22 = cos 29 E.8 AOi - AO: A a 1 2 = cos 29 E.9 Equations E.7, E.8 and E.9 may be expressed in matrix form A c n Ao 22 Aa 1 2 r [ T ] Aoi Aa 3 Ati3 E.10 [ T ] is the transformation matrix [ T ] = 1_ cos 29 1_ cos 29 2 2 2 2 i 1_ cos 29 1_ cos 29 2 2 2 sin 20 2 sin 20 APPENDIX F Back Calculation of So i l Parameters Given the compressibility of s o i l 0 g = 9 E-b kPa - 1 and the effective stress o§ = 751 kPa, the bulk modulus K number can be B backcalculated by assuming a value of m. From Equation 3.12 B = *B P a ( p - ) m *' 1 a s a Substituting 3 , and P and assuming m = 0.25 (Byrne and Cheung) S cl F.2, = 665 Adopting the relationship K = 0.6 K (Byrne and Cheung) B E K E = 1108 Other parameters are depicted from Byrne and Eldridge on Byrne and Cheung reports. A complete set of s o i l parameters for Test 11 may b written as: \ = 1108 n = 0.5 h " 665 m = 0.25 Rf = 0.8 <t>' = 42° A<f>' = 8° n = 0.3228 (e = 0.4767) S = 0.9975 APPENDIX G Elastic Stress Solution If the f l u i d pressure i s included, the displacement u of an elastic material may be written as (X + 2G) ^ ( T i + ~) + 3 0 G.l dr dr r' dr where A = (1+v) (l-2v) G = E 2 (1+v) C 3 - 1 - -p-— which i s assumed equal to 1 Cb = sand matrix compressibility C, .= sand bulk compressibility b The pressure may be expressed by Darcy's law in radial form iP _ "q , dr 2irhKr The stresses are written as o = X e e + 2G e 6 + P G.3 r v r a. = X e e + 2G + P G.4 6 v 0 a = X e e + 2G e e + P G.5 Z V z e e e 6 6 Where e , e 0 and e are the elastic strain components and e = e + e r 9 z v v r + e e G.6 z Assuming the i n i t i a l loading cause a deformation only in the ve r t r i c a l direction, but no displacement i n the horizontal directions (e = en = 0), the i n i t i a l vertical strain e , is given by G.5 as r 9 zo a -P zo o :zo = X + 2G G , / Assuming only radial displacement after i n i t i a l loading, the strains are e = — G.8 r dr e u e9 = 7 G ' y e e = e z zo By solving Equation G.l with the boundary conditions a = a . when r = R. r r i i a = a when r = R, r ro i and combining the results with Equations G.8, G.9 and G.10, the stres solutions (Equations 6.2, 6.3, 6.4) can be found by inserting the result in Equations G.3 through G.5 APPENDIX H Plastic Stress Solutions (a < a < a„) r z 0 The conditions must be satisfied within the plastic zone 1) Equilibrium da a - o„ 2) Mohr Coulomb failure c r i t e r i a f = a. - a tan 2 a + (tan 2 a-1) P - 2S tan a = 0 H.2 0 r The flow rule associated with yield condition i s P , 9f ^ 2 e = X = - X tan z a r da r H.3a el - X ~ - - X H.3b 0 9a Q e P - \ | | - - 0 H.3c z da From Equations H.3a and H.3b, i t follows that e£ + tan 2 a = 0 H.4 The total strain components may be written as e r = e r + E r H ' 5 a e g = Eg + e^ H.5b e p e = e + e = e H.5c z z z zo = 0 because i t is assumed that there is only radial displacement after i n i t i a l loading. Combining Equations H.5a and H.5b, inserting into H.4, gives o e e o e + e. tan^ a = E + E q tan^ a H.6 r 0 r 0 Applying Hooke's law of elas t i c i t y for porous material, yields Ee® = a - v (a f l + a ) - (l-2v)P H.7a r r o z E E Q = a . - v ( o - + a ) - (l-2v)P H.7b 0 9 v r z Ee e = a - v (a + a Q) - (l-2v)P H.7c z z r 9 Substituting Equation H.5c into H.7c yields a = Ee + v (a + a.) + (l-2v)P H.8 z zo r 8 ' Combining Equations H.7 and H.8, inserting into the yield criterion Equation H.2 with the strain relation i n H.6, i t gives [tan1* a + 1 - v (tan 2 a + 1) 21 a = 2G £ + 2G e„ tan 2 a J r r 9 + [tan 2 a (tan 2 a - 1) - v (tan1* a - 1) + (tan 2 a + 1) (1 - 2 v)]p - [tan 2 a (l-v) - v] 2 Stan a + v (tan 2 a + 1) 2G E H.9 zo and [tan1* a + 1 - v (tan 2 a + l ) 2 ] o Q = 2Ge tan 2 a + 2Gen tan* a o r 6 + [v (tan1* a - i) - (tan 2 a - 1) + tan 2 a (tan 2 a + 1) (l-2v)]p - [v (tan 2 a + 1) - l ] 2Stan a + v tan 2 a (tan 2 a + 1) 2Ge H.10 zo Substituting Equations H.9 and H.10 into equilibrium Equation H.l, together with the strain-displacement relations G.8 and G.9, the displacement equation may be written as 2 r 2 — j + r - utan 2 a = 7^7 (- [tan 2 a (tan 2 a - 1) - v(tan l + - 1) dr + (tan 2 a + 1) (l-2v)] j^hk + ( t a f l 2 a + 1 ) < 1 _ 2 v) 2 S t a n a + v (tan 2 a - 1) 2G e 1 H . l l zoJ The displacement solution of Equation H.ll i s , tan 2a , . - tan 2a , _ „ .. 2Gu = A.^ + A 2r + Br H.12 and the corresponding strains are 2Ger = tan 2 .0 ^ r 1 ^ - 1 - tan 2 a A 2 r _ t a n 2 a + B H.13 0 0 A tan 2a . , -tan 2a . _ „ ,, 2GeQ = A ^ + A 2 r + B H.14 Where Aj and A 2 are constants of integrations which depend on i n i t i a l conditions. 134 B = [f2? L tan' z a ^ l-2v -i uq l-2v „ 0 v , - i + - — 2 — r J o u i - 7 — 2 — r 2 Stana •'a+l tan^a-l J 2irhk tan zct-l - v 2G e H.15 zo Substituting Equations H.13 and H.14 into Equations H.9 and H.10 together with Darcy's law for radial flow, gives P = P + -Hi- £n — * i 2TThk X n R i H.16 a = P, + -^rr In-I- - -\- (2Stana - ^ 3 - ) r i 2irhk R, t 2TThky o A l t + 2tan z a ^— r aa = p- + o^T * n I 7" (2Stana - tan 2 a -^tr") 9 I 2irhk t 2trhk H.17 u A l t + 2tan*a — r T H.18 a = P. + 0 .* Jin — 4Stan a - (tan z a + 1), z i ZTrhk R t L ^S_] + ( 1 7 ) ( 1 " 2 V ) (a -P ) - + v (tan 2a + 1), 2irhkJ 1-v zo o A l 2tan 2 a - r ' H.19 where t = tan 2 a - 1 T = tan 4 a + 1 - v (tan 2 a + l ) 2 The constant of Integration Aj can be found by inserting the boundary condition, o = P. when r = R., into H.17. r i i ' 136 APPENDIX J The quantity of steady-state seepage at any distance r from the centre of the well i s Q = k 2irrh ^ J«l x dr on setting the limit of integration, yields _ R , h 2 r hi w 1 and Q = IT k (h2 - h') A n I J.2 w Equating J . l and J.2, yields 2 r h § = (4 -h?) j^sbr-'C J.3 w integrating J.3 h 2 = C Jin r + D J.4 137 where D i s the constant of integration Substituting the boundary conditions r = and h = h^ in J.4, yields, 2 D = hi - C Jin r J.5 * XJ Substituting J.5 into J.4, yields 2 2 h l ~ h ? r h = hi + -—5-7— An — J.6 1 In R/r r w w
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- A stress-strain model for the undrained response of...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
A stress-strain model for the undrained response of oil sand Cheung, Ka Fai Henry 1985
pdf
Page Metadata
Item Metadata
Title | A stress-strain model for the undrained response of oil sand |
Creator |
Cheung, Ka Fai Henry |
Publisher | University of British Columbia |
Date Issued | 1985 |
Description | An efficient undrained model for the deformations analysis of oil sand masses upon undrained loading is presented in this thesis. An analysis which couples the soil skeleton and pore fluids is used. The soil skeleton is modeled as a non-linear elastic-plastic isotropic material. In undrained conditions, the constitutive relationships for the pore fluids are formulated based on the ideal gas laws. The coupling between the soil skeleton and the pore fluids is based upon volume compatibility. The undrained model was verified with the experimental results and one dimensional expansion of soil sand cores. Comparisons between computed and measured responses are in good agreement and suggest that this model may prove useful as a tool in evaluating undrained response of oil sand. The response of a wellbore in oil sand upon unloading was analyzed using the developed model. Such analysis are important in the rational design of oil recovery systems in oil sand. |
Subject |
Oil sands |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2010-05-27 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0062954 |
URI | http://hdl.handle.net/2429/25084 |
Degree |
Master of Applied Science - MASc |
Program |
Civil Engineering |
Affiliation |
Applied Science, Faculty of Civil Engineering, Department of |
Degree Grantor | University of British Columbia |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-UBC_1985_A7 C43.pdf [ 4.23MB ]
- Metadata
- JSON: 831-1.0062954.json
- JSON-LD: 831-1.0062954-ld.json
- RDF/XML (Pretty): 831-1.0062954-rdf.xml
- RDF/JSON: 831-1.0062954-rdf.json
- Turtle: 831-1.0062954-turtle.txt
- N-Triples: 831-1.0062954-rdf-ntriples.txt
- Original Record: 831-1.0062954-source.json
- Full Text
- 831-1.0062954-fulltext.txt
- Citation
- 831-1.0062954.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
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>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0062954/manifest