Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

A simple anisotropic bounding surface plasticity model for cyclic response of clays Seidalinov, Gaziz 2012

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata

Download

Media
24-ubc_2012_fall_seidalinov_gaziz.pdf [ 3.93MB ]
Metadata
JSON: 24-1.0073338.json
JSON-LD: 24-1.0073338-ld.json
RDF/XML (Pretty): 24-1.0073338-rdf.xml
RDF/JSON: 24-1.0073338-rdf.json
Turtle: 24-1.0073338-turtle.txt
N-Triples: 24-1.0073338-rdf-ntriples.txt
Original Record: 24-1.0073338-source.json
Full Text
24-1.0073338-fulltext.txt
Citation
24-1.0073338.ris

Full Text

A SIMPLE ANISOTROPIC BOUNDING SURFACE PLASTICITY MODEL FOR CYCLIC RESPONSE OF CLAYS by Gaziz Seidalinov A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF APPLIED SCIENCE in The Faculty of Graduate Studies (Civil Engineering) THE UNIVERSITY OF BRITISH COLUMBIA (Vancouver) November 2012 c
 Gaziz Seidalinov 2012 Abstract Natural clays are anisotropic in their in-situ state and have an undisturbed shear strength in excess of the remoulded strength. In addition, most of the structures founded on clay deposits must be designed to withstand cyclic loads such as seismic ground motions or ocean waves. When subjected to earthquake or wind induced cyclic loadings, clay exhibits complex response. A realistic modeling of clay response under irregular cyclic loading requires an ap- propriate description of the stress{strain relationship. This thesis extends the formulation of a Simple ANIsotropic CLAY plasticity (SANICLAY) model by incorporation of a bounding surface formulation for successful simulations of clay response under cyclic loading. The most important aspects of the proposed formulation are the position of a projection center and the ability to capture continuous stiness degradation. The proposed projection center is established in the instant of any stress reversal and it realistically re
ects the experimentally observed plastic strains. A damage parameter is also adopted to better simulate the con- tinuous stiness reduction during the course of applied cyclic loading. The proposed model is developed with the aim of maintaining the simplicity and yet including an adequate level of sophistication for successful modeling of the key features of clay response. The model formulation is presented in detail followed by its qualitative and quantitative comparison with experimentally observed clay response. The presented model validation demonstrates the capabilities of the model in capturing a number of important characteristic features of clay response in both monotonic and cyclic loadings. Limitations and recommendations for future work are discussed as well. ii Table of Contents Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ii Table of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iii List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iv List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . v List of Symbols . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vi Acknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . x 1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 2 Experimental and Theoretical Modeling of Clay Response . . . . . . . . 4 2.1 Experimental modeling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 2.1.1 Discussion on fast and slow undrained cyclic loadings . . . . . . . . . 4 2.1.2 Bacchus (1969) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2.1.3 Sangrey et al. (1969) . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 2.1.4 Takahashi et al. (1980) . . . . . . . . . . . . . . . . . . . . . . . . . 7 2.1.5 Zergoun (1991); Zergoun and Vaid (1994) . . . . . . . . . . . . . . . 8 2.2 Theoretical modeling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 2.2.1 Classical elastoplastic models . . . . . . . . . . . . . . . . . . . . . . 11 2.2.2 Plastic strains inside the yield surface . . . . . . . . . . . . . . . . . 13 2.2.3 Development of the bounding surface formulation . . . . . . . . . . . 14 2.2.4 Bounding surface models for monotonic response of clays . . . . . . . 22 2.2.5 Bounding surface models for cyclic response of clays . . . . . . . . . 26 iii Table of Contents 2.3 Missing features in reviewed models . . . . . . . . . . . . . . . . . . . . . . 32 3 Development of the Proposed Model . . . . . . . . . . . . . . . . . . . . . . 34 3.1 Reference Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 3.1.1 Formulation in the triaxial space . . . . . . . . . . . . . . . . . . . . 34 3.1.2 Formulation in the multiaxial space . . . . . . . . . . . . . . . . . . 38 3.2 Proposed Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42 3.2.1 Formulation in the triaxial space . . . . . . . . . . . . . . . . . . . . 42 3.2.2 Formulation in the multiaxial space . . . . . . . . . . . . . . . . . . 50 3.3 Sensitivity analysis on new model constants . . . . . . . . . . . . . . . . . . 55 3.3.1 Parameter h0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 3.3.2 Parameter ad . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58 4 Qualitative Model Performance . . . . . . . . . . . . . . . . . . . . . . . . . 61 4.1 The advantage of the proposed projection center . . . . . . . . . . . . . . . 61 4.2 Undrained cyclic triaxial loading . . . . . . . . . . . . . . . . . . . . . . . . 65 4.3 Circular stress path in the deviatoric plane . . . . . . . . . . . . . . . . . . 67 4.3.1 Response to cyclic circular stress path loading . . . . . . . . . . . . . 68 4.3.2 Eect of parameter C in circular stress path loading . . . . . . . . . 71 4.4 Generation of secant shear modulus and damping ratio variation with shear strain amplitude . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72 4.5 Generation of cyclic stress ratio versus number of cycles . . . . . . . . . . . 77 4.6 Generation of peak axial strain and residual pore water pressure versus number of cycles . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80 5 Quantitative Model Performance . . . . . . . . . . . . . . . . . . . . . . . . 85 5.1 Experiments required for model calibration . . . . . . . . . . . . . . . . . . 85 5.2 Reconstituted Georgia kaolin clay . . . . . . . . . . . . . . . . . . . . . . . . 87 5.2.1 Model calibration . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88 5.2.2 Model performance . . . . . . . . . . . . . . . . . . . . . . . . . . . 92 5.3 Structured Cloverdale clay . . . . . . . . . . . . . . . . . . . . . . . . . . . 98 iv Table of Contents 5.3.1 Model calibration . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98 5.3.2 Model performance . . . . . . . . . . . . . . . . . . . . . . . . . . . 102 5.4 Reconstituted Ariake clay (fast cyclic loading) . . . . . . . . . . . . . . . . . 111 5.4.1 Model calibration . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 5.4.2 Model performance . . . . . . . . . . . . . . . . . . . . . . . . . . . 112 6 Summary, Conclusions, and Recommendations for Future Research . . 115 6.1 Summary and conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115 6.2 Recommendations for further research . . . . . . . . . . . . . . . . . . . . . 117 Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118 v List of Tables 3.1 Model parameters for sensitivity analysis . . . . . . . . . . . . . . . . . . . . 56 3.2 Initial internal variables for sensitivity analysis . . . . . . . . . . . . . . . . . 56 4.1 Summary of model parameters and initial internal variables for simulations of undrained cyclic triaxial loading . . . . . . . . . . . . . . . . . . . . . . . . . 65 4.2 Summary of model parameters for simulations of circular stress path in the deviatoric plane . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 5.1 Model constants for the proposed SANICLAY-D-B model . . . . . . . . . . . 87 5.2 Initial internal variables for the proposed SANICLAY-D-B model . . . . . . 87 5.3 Model parameters for Georgia kaolin Clay . . . . . . . . . . . . . . . . . . . 91 5.4 Initial internal variables for Georgia kaolin clay . . . . . . . . . . . . . . . . 92 5.5 Model parameters for Cloverdale clay . . . . . . . . . . . . . . . . . . . . . . 103 5.6 Initial internal variables for Cloverdale clay . . . . . . . . . . . . . . . . . . . 103 5.7 Model parameters for Ariake clay . . . . . . . . . . . . . . . . . . . . . . . . 114 5.8 Initial internal variables for Ariake clay . . . . . . . . . . . . . . . . . . . . . 114 vi List of Figures 2.1 Comparison of accumulation of pore water pressure in fast (dash lines) and slow (solid lines) undrained cyclic triaxial loading tests with drainage allowed at the end of loading for a range of applied stress amplitudes, cy=2Suc, from Zergoun and Vaid (1994) c
Canadian Science Publishing or its licensors . . . 10 2.2 Radial mapping rule from Dafalias (1986b) . . . . . . . . . . . . . . . . . . . 18 2.3 Radial mapping rule from Dafalias and Herrmann (1986) . . . . . . . . . . . 19 2.4 Radial mapping rule from Anandarajah and Dafalias (1986) . . . . . . . . . 21 3.1 Reference model schematic diagram in the triaxial plane . . . . . . . . . . . 36 3.2 Proposed model schematic diagram in the triaxial plane . . . . . . . . . . . . 44 3.3 Proposed projection center . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47 3.4 Multiaxial representation of: (a{b) the bounding surface with N -cone, and (c{d) plastic potential with critical state failure envelope . . . . . . . . . . . 52 3.5 Parameter h0 sensitivity in undrained monotonic triaxial loading . . . . . . . 57 3.6 Parameter h0 sensitivity in undrained cyclic triaxial loading . . . . . . . . . 59 3.7 Parameter ad sensitivity in undrained cyclic triaxial loading . . . . . . . . . 60 4.1 Typical stress{strain path in undrained cyclic triaxial loading test after Sheu (1984) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 4.2 Projection center on rotational variable  in the triaxial plane . . . . . . . . 63 4.3 Stress{strain path simulations in undrained cyclic triaxial loading with the projection center on rotational variable : (a) Cp = 0, and (b) Cp = 0:6 . . . 64 4.4 Stress{strain path simulations in undrained cyclic triaxial loading with the projection center proposed in this thesis . . . . . . . . . . . . . . . . . . . . 64 vii List of Figures 4.5 SANICLAY family models performance in undrained cyclic triaxial loading: (a{b) SANICLAY, and (c{d) SANICLAY-D . . . . . . . . . . . . . . . . . . 66 4.6 Comparison of model performance with bounding surface eect . . . . . . . 67 4.7 Isotropic unloading, constant-p triaxial compression, and circular stress path: (a) deviatoric view, (b) shifted view . . . . . . . . . . . . . . . . . . . . . . . 69 4.8 Observed strains in the strain space (deviatoric view on the left, shifted view on the right) induced by circular stress path (clockwise) in the deviatoric plane 70 4.9 Eect of parameter C on induced strain path during circular stress path loading 71 4.10 Denition of equivalent shear modulus and the quantities used for the deter- mination of damping ratio . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73 4.11 Experimental data after Vucetic and Dobry (1988): (a) degradation of Gsec and (b) variation of damping ratio with 
 . . . . . . . . . . . . . . . . . . . 74 4.12 Simulations of Gsec=Gmax degradation with 
: (a) ad = 0 and (b) ad = 60 . . 75 4.13 Simulations of damping ratio variation with 
: (a) ad = 0 and (b) ad = 60 . 77 4.14 Experimental data after Zergoun and Vaid (1994) on cyclic stress ratio CSR required to reach "a = 3% in a number of cycles . . . . . . . . . . . . . . . . 78 4.15 Generation of CSR required to reach "a = 3% . . . . . . . . . . . . . . . . . 79 4.16 Development of: (a) peak axial strain "a, and (b) residual pore water pressure u versus number of cycles, experimental data from Zergoun and Vaid (1994) 81 4.17 Generation of: (a) peak axial strain "a, and (b) residual pore water pressure u versus number of cycles with dierent values of h0 and xed ad = 40 . . . . 83 4.18 Generation of: (a) peak axial strain "a, and (b) residual pore water pressure u versus number of cycles with xed h0 = 50 and varied values of ad . . . . . 84 5.1 Calibration of Cc and Cr in e{ log p plane based on isotropic consolidation test after Sheu (1984) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88 5.2 Calibration of Mc and Me in the triaxial plane based on undrained triaxial compression and extension tests starting from normally consolidated state, experimental results after Sheu (1984) . . . . . . . . . . . . . . . . . . . . . 89 viii List of Figures 5.3 Calibration of C based on undrained triaxial compression test on normally consolidated clay, experimental results after Sheu (1984) . . . . . . . . . . . 91 5.4 Initial hardening parameter calibrated to h0 = 50 based on undrained cyclic triaxial loading, experimental results after Sheu (1984) . . . . . . . . . . . . 92 5.5 Simulation of isotropic consolidation test, experimental results after Sheu (1984) 93 5.6 Simulation of CIUC and CIUE tests on normally consolidated clay, experi- mental results after Sheu (1984) . . . . . . . . . . . . . . . . . . . . . . . . . 94 5.7 Undrained cyclic triaxial loading stress paths with the digitized experimental results after Sheu (1984) (left) and model simulations (right) . . . . . . . . . 96 5.8 Undrained cyclic triaxial loading stress{strain paths with the digitized exper- imental results after Sheu (1984) (left) and model simulations (right) . . . . 97 5.9 Calibration of Cc and Cr in e{ log p plane based on oedometer test after Zer- goun (1991) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 5.10 Calibration of c in the triaxial plane based on undrained triaxial compression test on normally consolidated and overconsolidated specimens, experimental results after Zergoun (1991) . . . . . . . . . . . . . . . . . . . . . . . . . . . 100 5.11 Calibration of: (a) Si, and (b) ki based on oedometer test after Zergoun (1991)101 5.12 Calibration of C based on CIUC test on initially normally consolidated clay, experimental data after Zergoun (1991) . . . . . . . . . . . . . . . . . . . . . 102 5.13 Simulation of oedometer test on Cloverdale clay, experimental results after Zergoun (1991) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 104 5.14 Simulations of CIUC and CIUE tests on normally and overconsolidated clay, experimental results after Zergoun (1991) . . . . . . . . . . . . . . . . . . . . 105 5.15 Undrained cyclic triaxial loading stress path with the digitized experimen- tal results after Zergoun (1991); Zergoun and Vaid (1994) (left), and model simulations (right) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106 5.16 Undrained cyclic triaxial loading stress{strain path with the digitized exper- imental results after Zergoun (1991); Zergoun and Vaid (1994) (left), and model simulations (right) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107 ix List of Figures 5.17 Simulations of peak "a versus number of cycles in undrained cyclic triaxial loading for the range of available cyclic stress amplitudes . . . . . . . . . . . 110 5.18 Simulations of peak "a versus number of cycles in undrained cyclic triaxial loading for the range of available cyclic stress amplitudes with the original SANICLAY model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 5.19 Calibration of h0 and ad based on peak "a versus number of cycles in undrained cyclic triaxial loading, experimental results after Yasuhara et al. (1992) . . . 113 5.20 Simulation of cyclic stress ratio CSR versus number of cycles required to reach "DA = 5% . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114 x List of Symbols List of symbols is applicable for Chapters 3{6. ad rate of d evolution A destructuration distributing model parameter b similarity ratio of f = 0 to F = 0 C rate of  evolution Cc compression line slope in e{ log p plane Cp model parameter controlling position of projection center Cr unloading{reloading line slope in e{ log p plane Ce fourth-order elastic modulus tensor Cep fourth-order elastoplastic modulus tensor d damage parameter e void ratio ein initial e e deviatoric strain tensor ep plastic deviatoric strain tensor f = 0 loading surface F = 0 bounding surface G = 0 plastic potential G shear modulus Gs specic gravity of soil particles Gsec secant shear modulus Gmax maximum tangential shear modulus h bounding surface hardening function xi List of Symbols h0 initial value of h H Heaviside step function I identity tensor ki model parameter, controlling rate of destructuration K bulk modulus Kp plastic modulus Kp plastic modulus for stress states on bounding surface K0 earth pressure coecient at rest L loading index M critical state ratio in p{q plane Mc M in compression Me M in extension n unit tensor along r nc unit tensor along rc  N peak stress ratio on bounding surface p mean eective stress p image p on bounding surface pa atmospheric pressure pc p at projection center p0 size of bounding surface along p-axis p0;d destructured value of p0 p size of plastic potential along p-axis q deviatoric stress q image q on bounding surface qb closest qc projection onto bounding surface qc q at projection center qcy amplitude of uniform cyclic loading qn set of n internal variables q qc projection on r distance from projection center to image stress xii List of Symbols r back stress-ratio tensor s elastic nucleus factor s deviatoric stress tensor sc s at projection center s s on bounding surface sb projection of sc on bounding surface along nc Si initial structuration factor s projection of sc on Su undrained shear strength Suc undrained shear strength in compression u pore water pressure x model parameter required for  saturation  deviatoric back-stress ratio  deviatoric back-stress ratio tensor "a axial strain "pd destructuration plastic strain "DA double amplitude axial strain "pq plastic deviatoric strain "r radial strain "pv plastic volumetric strain " strain tensor "p plastic strain tensor  distance from current stress to image stress   generalized to multiaxial stress space  stress ratio K0  for K0-loading  Lode angle 
 shear strain  slope of unloading{reloading line in e{ ln p plane  slope of compression line in e{ ln p plane xiii List of Symbols  Poisson's ratio a axial stress e equivalent c r radial stress c consolidation stress  stress tensor   for stress states on bounding surface 1 major principal stress 2 intermediate principal stress 3 minor principal stress  shear stress cy shear stress amplitude of uniform cyclic loading  friction angle c friction angle in compression rF bounding surface gradient rf loading surface gradient hi Macauley brackets _ time derivative or rate  association with stress state on bounding surface : inner product of two tensors   outer product of two tensors tr trace operator xiv Acknowledgements I oer my sincerest gratitude to my supervisor, Dr. Mahdi Taiebat, who has supported me with his guidance throughout my research. His encouragement and high level of competence inspired me to explore the subject. I am truly indebted and thankful to him for the time he devoted to discussions on the subject matter and for his moral support. His sometimes critical comments helped me to become a better researcher. I attribute the level of this thesis to his eort. I would also like to thank Professor Yannis Dafalias for his help in formulating the motion of the projection center and for reviewing this thesis. Thanks to the Department of Civil Engineering faculty and sta members, the University of British Columbia, especially to the library sta for the resources provided. I would also like to thank all the people who have become my friends during my studies at UBC. I would like to thank my wife, Ainur, for her patience, love, and constant support. Thanks will not be enough for my parents, Ilyas and Karlygash, from whom I have been constantly away. Special thanks to my lovely brother, Arystan, who has been taking care of them and tried to ll my missing spot. Financial support provided by Bolashak International Scholarship of the President of the Republic of Kazakhstan and partial nancial support from the Natural Sciences and Engineering Research Council of Canada (NSERC) are gratefully acknowledged. xv 1. Introduction All stresses considered in this thesis are eective stresses. A constitutive model is a mathematical description of the stress{strain relationship for a material. The complexity of a constitutive model depends on the type of material to be simulated. Unlike man-made materials (e.g., concrete, steel), soil behavior is governed by the presence of water in addition to solid phase of soil particles. The complex nature of clays cannot be simulated satisfactory by simple models (e.g., von Mises, Mohr{Coulomb), as their behavior is known to be dependent on conning pressure and overconsolidation ratio (OCR), and they exhibit a signicant degree of anisotropy in the in-situ state. The dependence on conning pressure and to some extent on OCR can be simulated by the Modied Cam-Clay model (MCC) (Roscoe and Burland, 1966), which is based on Critical State Soil Mechanics (CSSM) theory developed at Cambridge University (Roscoe and Poorooshasb, 1963; Roscoe et al., 1963, 1958). The Simple ANIsotropic CLAY plasticity (SANICLAY) model developed by Dafalias (1986a), Dafalias et al. (2006) allows successive simulations of induced anisotropy by incorporation of rotational and distortional hardening of a so-called yield surface, which encloses elastic domain. Destructuration theory incorporated into SANICLAY is a further model extension (SANICLAY-D), proposed by Taiebat, Dafalias and Peek (2010), enabling one to simulate structure loss in sensitive clays. Experimental evidence on clays suggests that inelastic strains occur within the yield sur- face. This issue becomes fundamentally important for the loadings that involve stress paths inside the yield surface. Undrained cyclic loading on clays is accompanied by the reduction of mean stress, resulting in a stress path mainly contained within the yield surface. CSSM theory suggests that the reduction of mean stress is associated with plastic strains. As the classical elastoplastic models assume purely elastic strains for the stress states inside the yield surface, no change in mean eective stress can be simulated. A mechanism allowing 1 Chapter 1. Introduction reduction of mean stress for undrained cyclic loading is required. Iwan (1967) and Mroz (1967, 1969) independently, introduced similar models with several surfaces, enclosing each other, of progressively larger size. As the stress point reaches a larger surface, the stiness reduces, resulting in increased hysteretic work. A bounding surface framework, rst pre- sented for metals by Dafalias and Popov (1975) and then incorporated for soils by Dafalias (1986b), Dafalias and Herrmann (1986), and Anandarajah and Dafalias (1986), is another approach to simulate plastic deformations within the yield surface. Unlike in the models by Iwan and Mroz, stiness change within the yield surface is continuous in the bounding surface formulation. In addition, the formulation does not require storing in memory a nite number of surfaces for better simulation of stiness reduction for the stress states inside the yield surface. The bounding surface framework attracted the attention of researchers owing to its simplicity and has been incorporated into dierent models to simulate soil response. In the bounding surface formulation, the yield surface is replaced by a bounding surface, which plays a similar role as the yield surface in that it represents consolidation history. The bounding surface allows a certain degree of plastic deformation for the stress states inside it. The model proposed in this thesis (SANICLAY-D-B) incorporates the bounding sur- face formulation into the SANICLAY-D model, which is chosen to account only for isotropic destructuration sucient to simulate structure loss. The SANICLAY-D-B model is a contin- uation of the work presented by Taiebat, Dafalias and Kaynia (2010). The proposed model includes a modied projection center to avoid overdamping and a damage parameter for the reduction of stiness required to simulate the cyclic softening observed in experiments. The plastic modulus within the bounding surface depends on the relative distance between the current stress point and the image point, mapped from a projection center through the stress point on to the bounding surface. In addition to the distance dependence, the plastic modulus is also controlled by a hardening parameter, either constant or function. In the proposed formulation, the hardening parameter is a function of a newly introduced damage parameter. The projection center is located at the point of stress reversal at the instant of the reversal and further evolves to maintain its relative position within the evolving bounding surface. In general, the projection center can be anywhere within the bounding surface but never outside it. The number of constants that should be calibrated is 10 for reconstituted 2 Chapter 1. Introduction and 11 for structured clays, including two new model parameters related to the bounding surface. The thesis is presented in a manner that allows the readers, who are not familiar with the subject, to understand the proposed model formulation, its advantages and predictive capabilities. All stresses considered are eective stresses which govern the response of soil. In Chapter 2, a review of experimentally observed clay response in undrained cyclic and theoretical modeling of inelastic strains within the yield surface is covered. Theoretical mod- eling part also includes a brief introduction to classical elastoplastic models. The proposed SANICLAY-D-B formulation followed by a sensitivity analysis of new model parameters is presented to the reader in Chapter 3. A qualitative model performance is demonstrated in Chapter 4, consisting of:  a comparison of dierent projection center approaches reported in the literature and proposed for the current formulation,  a comparison of the proposed model with the previous SANICLAY family models in undrained cyclic triaxial loading and circular stress path in the deviatoric plane,  a model simulation of secant shear modulus and damping ratio variation with shear strain amplitude,  a model simulation of the cyclic stress ratio required to reach a certain strain level in a number of cycles,  a model simulation of attained axial strains and residual pore water pressure in a number of cycles. In Chapter 5, model calibration and quantitative performance are illustrated for three clays: Georgia kaolin clay (reconstituted), Cloverdale clay (structured) in undrained slow cyclic triaxial loading, and Ariake clay (reconstituted) in fast undrained cyclic triaxial loading. Finally, a summary, conclusions, and recommendations for future research are presented in Chapter 6. 3 2. Experimental and Theoretical Modeling of Clay Response In this chapter, a relatively detailed literature review is conducted on the experimentally observed response of clays in cyclic loading and a number of the proposed constitutive models for simulating that response. A bounding surface framework incorporated into the reviewed models appears to be an attractive approach for simulations of undrained cyclic loading. Finally, the features that need to be re
ected in cyclic loading simulations are discussed. 2.1 Experimental modeling 2.1.1 Discussion on fast and slow undrained cyclic loadings Imposed undrained cyclic loading on clays can be divided into two categories: slow and fast undrained cyclic loadings. Slow undrained cyclic loading tests are those in which the frequency of applied cyclic loading is slow enough to achieve equalization of pore water pressure in the whole specimen; otherwise, the tests are referred to as fast undrained cyclic loading ones in this study. It is important to understand advantages and disadvantages of both. The review of the experimental data is conducted on the tests performed in triaxial cell. During triaxial testing pore water pressure is usually measured on top and bottom platens of tested specimens. A measured pore water pressure is then deducted from applied total stresses to come up with eective stresses, governing the response of a tested material. In this way, positive pore water pressure is associated with reduction of eective stresses, while negative pore water pressure with increase of eective stresses. An element test requires the measured data on specimen boundaries to be representative of the entire specimen. The main advantage of slow undrained cyclic loading is that equalization of pore water pressure 4 2.1. Experimental modeling throughout the specimen is allowed, which ensures measurements of approximately average pressure in the entire specimen. Thus, interpreted eective stresses from slow undrained cyclic loading are more reliable. This is especially important for clays, where dissipation of pore water pressure requires a considerable amount of time. Fast cyclic loading tests miss this important feature and become less reliable in terms of representative eective stress measurements for the soil element. The rate eect of applied loading also plays an important role. The rate eect becomes more pronounced with increased plasticity of soil. Many investigators rely on the supposition that measured response is a result of the rate of applied loading only; however, in fast undrained cyclic loading, measured response is also a result of measuring un-equalized pore water pressure on the boundaries of the specimen, which is not representative of the material only. Of course, strain rate eect is present, but it is very dicult or even impossible to dierentiate soil response due to the rate eect from the measured response due to the pore water pressure inequalities within the specimen. Fast undrained cyclic loading tests invoke much higher measured strength than that in slow undrained cyclic loading. It is true that the increase of strength is due to the rate eect, but unreliable pore water pressure measurements are likely to erroneously increase measured strength as well. Thus, the strength in fast undrained cyclic loading is overesti- mated, rendering results that cannot be relied on. Although slow undrained cyclic loading does not account for the rate eect and therefore underestimates the strength of the soil, it gives conservative results. The emphasis is placed on slow undrained cyclic triaxial loading tests throughout the discussion of available experimental data. 2.1.2 Bacchus (1969) Undrained cyclic strain-controlled triaxial tests with constant strain amplitude were per- formed by Bacchus (1969) on commercially available Webb Mark 4 clay with plasticity index, PI = 26%. Tested specimens consist of 28% of clay. Bacchus (1969) was well aware of a need of pore water pressure equalization within cyclically loaded specimen; however, he adopted 0:2 Hz frequency as the appropriate one. From reported results on undrained cyclic triaxial loading on normally consolidated spec- 5 2.1. Experimental modeling imens the main ndings could be summarized as follows:  in all applied cyclic strain levels, during the rst initial cycles there was a decrease in mean stress p; the rate of reduction of p decreased with cyclic loading;  larger strain amplitude tests involved larger shear stresses to reach desired strains during rst few cycles; however, with an increasing number of cycles, lower shear stresses were attained for larger strain amplitudes in a given number of cycles; this was an indication of more rapid stiness loss for high strain amplitude tests;  small strain amplitudes stress{strain loops exhibited a subtle change in shape, and stress{strain loops became \S-shaped" with an increasing number of cycles; this be- havior might be attributed to the presence of a signicant amount of granular material. Undrained cyclic triaxial tests were also performed on samples isotropically consolidated to a range of preconsolidation stresses and with dierent overconsolidation ratios (OCRs). Specimens with larger OCRs exhibited an initial increase in p and reduction with further cycles. The extent of p increase, followed by its decrease, was larger for higher OCRs. In general, stress{strain loops for overconsolidated specimens were similar to those obtained for normally consolidated specimens with the dierence that less reduction in stresses and energy was observed during initial cycles. 2.1.3 Sangrey et al. (1969) Sangrey et al. (1969) paid a close attention to reliable pore water pressure measurements in undrained cyclic triaxial tests. They performed a series of one-way cyclic loading tests with constant stress amplitude on undisturbed Neweld clay. To ensure pore water pressure equalization within the specimens all tests were conducted with an axial strain rate of about 0:0002%=min. Their main nding was that for each OCR there existed a \critical cyclic stress level" (CCSL). The tests with cyclic stress amplitudes below CCSL attained equilibrium, and resultant strains were nite, whereas applied loading with stress amplitude above CCSL caused ultimate collapse with innite axial strains. While normally consolidated specimens 6 2.1. Experimental modeling subjected to undrained cyclic triaxial loading exhibited stress path migration towards the origin until the equilibrium was reached, overconsolidated samples showed an increase of p within the initial cycle. The increase of p was then followed by its reduction that was smaller for higher OCRs. 2.1.4 Takahashi et al. (1980) Takahashi et al. (1980) discussed observations from undrained cyclic triaxial tests with dier- ent frequencies on a low-plasticity sandy clay, Lower Cromer Till. The loading was applied in constant stress amplitude tests with sinusoidal waveform. In addition to conventional pore water pressure transducers, piezometer probes at mid-height of tested specimens were used. The period of applied cyclic loading was restricted to a minimum of 10 min (of about 0:0017 Hz maximum frequency), ensuring reliable measurements of pore water pressure at the mid-height of tested specimens. The following observations were made in examination of undrained cyclic triaxial tests with 30-min period and applied magnitude of (a  r) =2Suc = 0:75, where Suc is an undrained shear strength in compression:  a migration of the stress path towards the origin induced by the rise of pore water pressure within a constant volume specimen;  the stress path changed its shape as it approached the critical stress ratio; loops de- veloped in extension owing to a more dialative response than in compression;  as the stress path reached the critical stress ratio, axial strains development accelerated and larger magnitude of axial strains was observed in extension;  the stress path traveled beyond the critical stress ratio in extension as the rate of applied loading was increased; Takahashi et al. (1980) further aimed to investigate what in
uenced the rate of stress path migration towards the origin. From their ndings they concluded that the main factors were: cyclic stress ratio (CSR), consolidation history, and rate of applied loading. A larger CSR initially induced more rapid reduction of p. However, once the critical stress ratio was 7 2.1. Experimental modeling reached, p reduction occurred at a decreasing rate. A lower CSR, apparently, approached the critical stress ratio at lower p. As a result, overall reduction of p for lower CSRs was larger once the stress path was stabilized. Specimens sheared from high OCRs exhibited an increase in p, i.e., a dialative response during initial cycles, which was then followed by reduction of p with the magnitude of reduction inversely proportional to OCR. Takahashi et al. (1980) indicated that direct evidence on the eect of cyclic frequency was not available as the results were obtained by measuring pore water pressure at the boundaries; however, it is worth to elaborate on reported results of observed stress paths with dierent cyclic loading frequencies. In general, reduction of p during cyclic loading occurred to a lesser degree as loading frequency was increased. The rst series of tests were performed with cycle periods of 100 and 900 s, even though Takahashi et al. (1980) mentioned that the minimum cyclic period should be 10 min to ensure reliable pore water pressure measurements. Thus, the reported low reduction of p for smaller periods was partly due to not stabilized pore water pressures within the tested specimens. In the second series, soil response during cyclic tests with periods of 10, 30, 240, and 480 min were studied. The results clearly showed that as the loading period was increased the reduction of p was larger. Takahashi et al. (1980) ensured that pore water pressures at the peripheries of the specimens were measured accurately and radial non-uniformities were unlikely. They admitted that there might had been local irregularities of pore water pressure within the specimens. 2.1.5 Zergoun (1991); Zergoun and Vaid (1994) Zergoun (1991), Zergoun and Vaid (1994) justied the need of slow undrained cyclic triaxial loading for clays. They mentioned that non-uniformities within a specimen were dependent on the degree of end restraints and the rate of loading. As the center of specimen was believed to be the most representative of the whole specimen, the only way to obtain reliable pore water pressure measurements was to perform slow undrained cyclic triaxial loading tests with 0:05%/hr for strain-controlled and 60 kPa/hr for stress-controlled tests. Zergoun (1991), Zergoun and Vaid (1994) compared slow and fast undrained cyclic triaxial loading tests with constant stress amplitude on structured Cloverdale clay with PI = 24%. When comparing results for a particular stress amplitude in slow and fast cyclic loading the 8 2.1. Experimental modeling following observations were made:  the eect of stress amplitude was similar;  a much smaller number of cycles was required to induce a particular peak axial strain in slow cyclic loading; this was accompanied by more extensive p reduction of the stress path;  axial strains accumulation accelerated once the stress path nearly reached equilibrium for slow cyclic loading; in fast cyclic loading such acceleration occurred well before the equilibrium was attained;  larger values of pore water pressure were measured during slow cyclic loading; in Fig. 2.1 accumulation of pore water pressure reached a plateau in slow cyclic loading; referring to Fig. 2.1, fast cyclic loadings were accompanied by a pore water pressure increase during the rest periods at which applied stresses were held constant, while in slow cyclic loadings pore water pressure was unchanged; this indicated that equilibrium of pore water pressure was not achieved during fast cyclic loading and it did so during the rest periods involving further deformations; even though pore water pressure increased in rest periods after fast cyclic loading, residual pore water pressure was below that measured in slow cyclic loading for a given stress amplitude. From their ndings Zergoun (1991), Zergoun and Vaid (1994) conrmed that during fast undrained cyclic loading tests the eective stresses could not be adequately measured. In the tests with dierent stress amplitudes, Zergoun (1991), Zergoun and Vaid (1994) found a threshold stress amplitude. Axial strains decelerated if the applied stress amplitude was below the threshold and accelerate if it was above. Larger axial strains developed in extension for all stress amplitudes with a decreasing dierence between strains in compression and extension as the amplitude of applied loading was decreased. For stress amplitudes above the threshold, the following three stages of deformations were traced:  primary stage: axial strains develop at a decreasing rate;  secondary stage: strains develop at an approximately constant rate; 9 2.1. Experimental modeling Figure 2.1: Comparison of accumulation of pore water pressure in fast (dash lines) and slow (solid lines) undrained cyclic triaxial loading tests with drainage allowed at the end of loading for a range of applied stress amplitudes, cy=2Suc, from Zergoun and Vaid (1994) c
Canadian Science Publishing or its licensors  tertiary stage: strains develop at an increasing rate. For very high stress amplitudes the secondary stage was negligible, and the decreasing rate of strains during the primary stage was followed by an increasing rate in the tertiary stage. Two more important observations should be considered. The initial direction of applied cyclic loading played a signicant role. Two undrained tests with the same stress amplitude, one starting from compression and the other from extension, were compared. The cyclic test starting from extension exhibited more severe deformations than the test starting from compression. This implied the importance of an anisotropy for tested clay. The eect of a step increase in amplitude of cyclic loading was studied as well. In 10 cycles, the amplitude of loading was increased from cy=Suc = 0:63 to 0:67, where cy is an amplitude of applied shear stress. The induced strains, residual strains, and pore water pressure in a number of cycles simply corresponded to transformation of the response to the one that would had been obtained for cy=Suc = 0:67 with no change of amplitude of applied loading. 10 2.2. Theoretical modeling 2.2 Theoretical modeling Before going into detail of theoretical modeling of undrained cyclic loading, a review of general framework used in classical elastoplastic models is proposed to the reader. 2.2.1 Classical elastoplastic models Classical elastoplastic models assume the existence of a yield surface in a six dimensional stress space. Its general equation can be expressed as f (;qn) = 0 (2.1) where  is a stress tensor; and qn are n internal variables controlling the geometry of the yield surface. This yield surface encloses the region inside which any stress increment that does not cross the yield surface induces only elastic strains. Once a stress increment crosses the yield surface a plastic strain increment is induced in addition to the elastic one. A basic kinematic assumption is a decomposition of a total strain tensor " into elastic "e and plastic "p parts. The rate of " can be decomposed accordingly: _" = _"e + _"p (2.2) where a superposed dot stands for the material time derivative, i.e., rate. An input into a constitutive model could be either a strain increment or stress increment. General constitutive relations between stress increment and elastic strain increment are then given for the corresponding inputs, respectively, _ = Ee : _"e; (2.3) _"e = Ce : _ (2.4) where operator : is an inner product of two adjacent tensors; Ee and Ce are fourth order tensors of elastic moduli and compliances, respectively. Rates of plastic strain and internal variables are assumed to be given by Eqs. 2.5 and 11 2.2. Theoretical modeling 2.6. _"p = hLiNp (2.5) _qn = hLiqn (2.6) where L is a plastic scalar-valued multiplier, or loading index to be dened in the sequel of this section; hi are Macauley brackets such that hLi = L if L > 0, and hLi = 0 if L  0; Np and qn are tensorial and scalar/tensorial, respectively, functions of state variables  and qn. It is then common to deneNp = @g=@, where @g=@ is a gradient of a plastic potential given by: g (;qn) = 0 (2.7) An associative 
ow rule considers the yield surface and plastic potential coinciding with each other, while a non-associative 
ow rule assumes non-coincident yield surface and plastic potential. Then Eqs. 2.3 and 2.4 are re-written on the basis of Eqs. 2.2 and 2.5: _ = Ee :  _" hLi @g @  (2.8) _" = Ce : _ + hLi @g @ (2.9) The consistency condition requires that the stress point cannot lie outside the domain enclosed by the yield surface. Thus, once a stress increment crosses the yield surface, in- ternal variables are updated to maintain an updated stress point on the yield surface. The consistency condition is expressed in mathematical form as _f = @f @ : _ + @f @qn : _qn (2.10) Substituting Eq. 2.6 into Eq. 2.10 gives the loading index L = 1 Kp  @f @ : _  (2.11) 12 2.2. Theoretical modeling where Kp =  (@f=@qn) : qn is a plastic modulus. The constitutive relationships in Eqs. 2.8 and 2.9 can be re-written on the basis of Eq. 2.11 so that the the constitutive matrices are: _ = Eep _" =  Ee H (L) (E e : @g=@)  (@f=@) : Ee Kp + @f=@ : Ee : @g=@  : _" (2.12) _" = Cep : _ =  Ce + (@g=@)  (@g=@) Kp  : _ (2.13) where H (L) is the Heaviside step function, H (L) = 1 for L > 0 and H (L) = 0 for L  0; Eep and Cep represent an elastoplastic modulus and compliances matrices. Eqs. 2.12 and 2.13 are complete sets of constitutive relations for inputs of strain and stress increments, respectively, in a classical elastoplastic model. It is sucient to provide the yield surface f = 0, plastic potential g = 0 or any other appropriate function for Np, and evolution laws for plastic strain rate and rate of internal variables shown generally in Eqs. 2.5 and 2.6 2.2.2 Plastic strains inside the yield surface The observed experimental undrained cyclic response of clays in triaxial tests illustrated that the behavior of clay is not elastic within the yield surface. From what the question follows if it is necessary to account for inelastic strains within the yield surface in constitutive modeling. The answer depends on what type of loading is required to be simulated. For monotonic loading, concerned with peak and ultimate strength of clay or the stress states mainly on the yield surface, classical elastoplastic models with incorporated certain complexities are able to predict soil response with fair accuracy. If monotonic loading is concerned with simulated soil response for the stress states within the yield surface, model allowing inelastic strains within the yield surface should be considered. For earthquake loadings, the imposed stress path always involves loading and unloading processes. Even if the loading stress path is on the yield surface, unloading is directed towards the interior of the yield surface which is then followed by reloading. The necessity to simulate plastic deformations within the 13 2.2. Theoretical modeling yield surface becomes pronounced as the number of such loading-unloading-reloading cycles increases. Earthquake-induced deformations within the soil often involve a considerable number of these cycles, which does require simulation of inelastic soil response inside the yield surface for successive model predictions. As has been discussed in section 2.1.1, experimental testing requires that slow undrained cyclic loading should be performed for accurate pore water pressure measurements. In slow undrained cyclic loading tests rate eect is still present, but to a much lower degree than in fast undrained cyclic loading. Therefore, literature review is conducted below on rate- independent models with incorporated bounding surface formulation. The bounding surface framework is an attractive tool to model undrained cyclic loading imposed on clays used by many researchers. The bounding surface formulation can be incorporated into any classical elastoplastic model by adjusting a plastic modulus Kp and leaving all other relationships unchanged. First, development of the bounding surface formulations is discussed. Then, rate-independent elastoplastic models with incorporated bounding surface formulation for monotonic and cyclic loading simulations are reviewed. Not all of the reviewed models are referred to by the authors proposing them as bounding surface plasticity models; however, the concept adopted for plastic strains simulation within the yield surface is similar to that used in the bounding surface formulation. To keep a common denition of the yield surface, inside which only elastic strains are allowed, bounding surface term is used for the yield surface with inelastic strains allowed inside it. As in various publications, the bounding surface is referred to in dierent ways, a conventional yield surface term is to be used in such cases. 2.2.3 Development of the bounding surface formulation Attempts to simulate inelastic response within the conventional yield surface were made in the early 1950s; slip theory proposed by Batdford and Budiansky (1949) was adopted by Koiter (1953) for this purpose. Further development of models describing soil behavior in multiaxial stress space was conducted separately by Iwan (1967) and Mroz (1967, 1969). They adopted an idea of family of surfaces of progressively larger size, one enclosing the other, with constant work-hardening moduli for each surface. The inner-most surface enclosed a 14 2.2. Theoretical modeling truly elastic region, while the outer-most served as the conventional yield surface. The change of the work-hardening modulus was linear for the stress states on a single surface until a surface of larger size was reached; at this point the work-hardening modulus changed and corresponded to that of a larger surface at the current stress point until next surface of larger size was attained. In this way, piecewise constant work-hardening moduli allowed simulation of inelastic strains inside the outer-most surface. The main disadvantage of this approach is that a number of surfaces have to be stored in the model memory. The need of a simpler approach inspired development of the bounding surface formulation proposed by Dafalias and Popov (1975). Rules of plastic modulus variation were proposed based on uniaxial cyclic loading for metals and generalized to multiaxial stress space. Since its introduction, the bounding surface formulation has been used by various investigators. In a series of publications by Dafalias (1986a), Dafalias and Herrmann (1986), and Anandarajah and Dafalias (1986), the bounding surface formulation is reviewed for simulations of soil response. Because the current work is dedicated to cyclic response of clays, the publications mentioned will be discussed in detail. The bounding surface formulation introduces inelastic strains for the stress states inside the conventional yield surface through plastic modulus. All other relations assumed for the bounding surface are similar to those used in classical elastoplastic models. The conventional yield surface is replaced by the bounding surface. Within the bounding surface F = 0, plastic strains depend on the distance from the actual stress point ij to the corresponding image point ij on the bounding surface, where the latter is determined by a mapping rule. The degree of plasticity increases as this distance becomes smaller and reaches its limit once the stress point is on the bounding surface. A loading surface f = 0 passes through ij and is enclosed by the bounding surface so that both surfaces coincide once a current stress point is on the bounding surface. The bounding surface concept requires conformance of the following: 1. the mapping rule has to be such that a unique ij is dened for any stress state within the bounding surface, and @f=@ij = @F=@ij, where @f=@ij and @F=@ij are gradients of the loading surface and the bounding surface,respectively, for ij = ij to guarantee that the loading surface never intersects the bounding surface; the mapping 15 2.2. Theoretical modeling rule is not invertible in a sense that an innite number of the stress states within the bounding surface could correspond to a given ij. 2. plastic strains are allowed for the stress states within the bounding surface; thus, as ij must remain on the bounding surface, the plastic modulus for the stress states on the bounding surface is obtained from the consistency condition, _F = 0: Kp =  @F @qn : _qn (2.14) 3. The actual plastic modulus Kp is related to Kp through Euclidian distance:  = [(ij  ij) (ij  ij)]1=2 (2.15) such that Kp > Kp for  > 0, and Kp ! Kp for  ! 0. Dafalias (1986b) reviewed two dierent mapping rules; each of the two mapping rules deserves to be discussed in more detail. Original bounding surface mapping rule This mapping rule was originally proposed by Dafalias (1975). The loading and bounding surfaces, having kinematic hardening rules that can be independent, were expressed respec- tively as follows: f (ij  ij; qn) = 0 (2.16) F (ij  ij; qn) = 0 (2.17) where ij and ij are the corresponding kinematic hardening internal variables. Adopted kinematic and other hardening rules for the two surfaces should be such that surfaces never intersect. Then the mapping rule, with the assumption that bounding and loading surfaces had similar shapes, was given by: ij = b (ij  ij) + ij (2.18) 16 2.2. Theoretical modeling where b is a similarity ratio between the two surfaces. The proposed plastic modulus for any stress state was: Kp = Kp + h (jinj)  : nh(in  ) : ni (2.19) where  = (ij  ij); in is equal to  at the moment of initiation of new loading process; and n is a loading direction. The in was updated when (in : n  : n)  0. Notice that the rst term was always larger than the second for a loading process and vice versa for unloading as the terms changed signs. The use of tensorial quantities was advocated by the case of neutral loading simulations (see Dafalias, 1986b, for details). With the original bounding surface mapping rule, Dafalias (1986b) was aware that the problem of \overshooting" exists. This was the case for short unloading with new in after the loading process, followed by reloading with updated in. This would result in a very small denominator in Eq. 2.19, i.e., a large plastic modulus, as opposed to the case where no such small perturbation had occurred and stress{strain curve in reloading would overshoot the one in loading. However, the overshooting was not innite, as it was bounded by the bounding surface. Radial mapping rule Dafalias (1986b) referred to the work of Hashiguchi and Ueno (1977) as the pioneering work for the radial mapping rule, where the origin had been used as a projection center. The proposed radial mapping rule, shown in Fig. 2.2, required back-stress ratio , which was the projection center, to be always located within or on the bounding surface. The image point on the bounding surface was a projection from  through the current stress point. The projection center also served as the center of homology between the loading and bounding surfaces, such that @f=@ij points in the same direction as @F=@ij. This ensures that the two surfaces never intersect. The radial mapping rule was given by ij = b (ij  ij) + ij (2.20) where the value of b was determined by plugging Eq. 2.20 into F = 0, 1  b  1. Because of the two surfaces being homologous to each other no consistency condition was required 17 2.2. Theoretical modeling the 8in (symbolized by dmax in their paper) to the value it had when a minimum 8 was achieved in the past, if this 5 is reached again. A dis- continuous Kp may occur though, since the same 8 will be related to two 8in  Another attempt to avoid overshooting was made by Tseng and Lee (37) using the elastic unloading chord, but it appears to have been re- stricted only to the cases where elastic deformation takes place upon unloading/reloading. Radial Mapping Plasticity Formulation.—This particular formulation is characterized by the "radial mapping rule" as proposed by Dafalias (13), and deserves a detailed investigation since only its main features were described in Ref. 13. In addition, many existing models with dif- ferent names are in essence radial mapping plasticity formulations and will be discussed subsequently. Fig. 1 shows the concept. The bounding surface is schematically shown as a circle of center p,-, (kinematic hard- ening) and can be described analytically by Eq. 25b. The novelty here is the introduction of a second back-stress, a,y, always within or on F = 0, not as the center of an enclosed yield surface as in the two-surface approach, but as the projection center used to radially project the cur- rent stress cr,y onto <r,y on F = 0. Thus, a radial mapping from cr,y to cr,-,- is introduced, and Eq. 19 takes the form b(<r,y 7) + «« (29) The b is s i and can be determined by inserting the cr,-, from Eq. 29 into the analytical expression of F = 0. Observe that b is now variable, con- trary to the constant similarity ratio b entering Eq. 26. When er,y = <r,y the b = 1, i.e., the identity condition is satisfied, while when o^ = a,-, the b = 00 and &ij is indeterminate. Eq. 29 together with the assumption L,-, = Ljj = dF/d&ij, define indirectly a loading surface / = 0 shown by a dashed line in Fig. 1, which is homologous to F = 0 with a,, being the center of ELASTIC NUCLEUS LOADING SURFACE Jki -BOUNDING SURFACE F ( £ f M , , „ ) = !> FIG. 1.—Schematic Diagram of Bounding Surface and Related Concepts 977 Figure 2.2: Radial mapping rule from Dafalias (1986b) for the stress states within the bounding surface. The plastic modulus equation for any stress state was proposed as follows: Kp = Kp + Ĥ  hr  si = Kp + Ĥh b b 1  si 1 (2.21) where r is the distance from the projection center to the image point, Ĥ is a positive shape hardening function, s is an elastic factor, and r= = b= (b 1). When   r=s, i.e., for the stress states within the elastic nucleus, Kp =1, rendering purely elastic strains. For s = 1, the elastic domain degenerated into a single point, whereas for s ! 1, the elastic domain coincided with the bounding surface, rendering a conventional elastoplastic model. A particular radial mapping rule for the bounding surface formulations was presented by Dafalias and Herrmann (1986), as a continuation of previous work by Dafalias (1986b), for isotropic cohesive soils within the Critical State Soil Mechanics (CSSM) theory. In addition 18 2.2. Theoretical modeling to the radial mapping, the authors included Lode angle dependency in the proposed model, which enhanced predictive capabilities for more general stress states. For the formulation, I was a rst invariant of the stress tensor, and J was a second invariant of the deviatoric stress tensor. The proposed projection center in Fig. 2.3 was positioned on the hydrostatic line, represented by I-axis, and was related to the size of the bounding surface along this axis by: Ic = CpI0 (2.22) where Ic is the position of the projection center on the I-axis, Cp is a new model constant in the range of 0  Cp < 1, and I0 is the size of the bounding surface. Similar to Eq. 2.20, image stress for a given projection center was: I = b (I  CpI0) + CpI0 (2.23) The composite bounding surface, consisting of two ellipses and one hyperbola with con- tinuous tangents, as shown in Fig. 2.3, was used for better predictions at high OCRs. In addition, the use of R > 1 parameter, which determined the position of point I1 in Fig. 2.3 and related two major axes of the ellipse 1 as I1 = I0=R, was considered. Lode angle was denoted by . The plastic modulus for any stress state was the same as the one in Eq. 2.21. Dafalias and Herrmann (1986) elaborated on the shape hardening function Ĥ and considered the following equation: Ĥ = 1 + ein   g 2pa [zmh() + (1 zm)h0] (2.24) where ein is an initial void ratio;  and  are slopes of the normal compression and the rebound lines, respectively, in the e-ln p plane; g =  (@F=@p)2 + (@F=@q)2 1:2 was used to account for a jump of a bounding surface gradient in transition points between the three surfaces comprising the bounding surface; pa is an atmospheric pressure; z = J=J1; m = 0:02, rendering zm  1 to maintain dependency on h() rather than h0; h() is a shape hardening function of h, dependent on , and it required determination of hc and he for triaxial compression and extension; h0 = (hc + he) =2. 19 2.2. Theoretical modeling HYPERBOLA ELLIPSE 2 ELLIPSE 1 ELASTIC FIG. 1.—Schematic Illustration of Bounding Surface and Radial Mapping Rule in Stress Invariants Space it will be explained in the sequel. Eq. 30 has the following specific forms. For ellipse 1: F=(I-I0)[l+^I0)+(R-lfU) =0. For the hyperbola: 2 I- h N R\ N/ LJjL RJ \N R For ellipse 2: F = (7 - TI0)(I - (X + 2QI0] + p/2 = 0 T(Z + TF) T2 £ = - - = - — = — ; p = 0 . Z + 2TF' RA " N  (31)  (32) (33a) (33b) N Z(Z + 2TF') N i :; Z = - (1 + y - VT + 72) (33c) VT+f' R The dependence on ep occurs through loie7), and the dependence on a through the parameters N(a), R(<x), and A(a) according to Q(a) = g(a,c)Qc (34«) Qe c = ~ (34b) g(a,c) = — C >  , (34c) 1 + c - (1 - c) sin 3a 1270 Figure 2.3: Radial mapping rule from Dafalias and Herrmann (1986) Eects of hc and C parameters in triaxial compression were shown in detail. While the former controlled the degree of plastic deformation for the stress states within the bounding surface, the latter was successfully used to control contracting or dialative responses for overconsolidated stress states. Calibration of the parameters and model performance in undrained monotonic triaxial compression suggested that the projection center along the hydrostatic axis was an attractive feature of the model for monotonic loading simulations for high OCR states. Anandarajah and Dafalias (1986) further incorporated inherent and induced anisotropy into the bounding surface soil plasticity model. The orientation of the bounding surface was controlled by evolving anisotropy, represented by the second-order tensor aij in Fig. 2.4, during the loading process. The projection center was determined as CpI a 0 , where Cp was the same model parameter proposed by Dafalias and Herrmann (1986) and Ia0 was a size of the bounding surface; the radial mapping rule was also identical, given the projection center. Due to the bounding surface orientation being controlled by evolved anisotropy, the plastic modulus and the shape-hardening function were adjusted accordingly, with the latter 20 2.2. Theoretical modeling FIG. 1.—Schematic Illustration of Bounding Surface and Geometrical Interpreta- tion of Stress Variables was accounted for in a rather simple manner. Additionally, the concept of elastic nucleus was introduced in order to model the cyclic stabili- zation exhibited by clays. Some important modifications are incorporated in applying the con- cept of bounding surface for anisotropic clays. Referring to Fig. 1, the bounding surface is assumed to be oriented at any stage along a line (/"- axis) whose direction changes during the loading, representing the evo- lution of anisotropy. A symmetric second-order tensor 8°;, such that 8|8| = 3, defines the orientation of the bounding surface axis (/"-axis), in the stress space. According to the concept of bounding surface, every point in the stress space, such as point B in Fig. 1, is associated with a corresponding image point on the bounding surface via a suitable map- ping rule. The current theory, as in the isotropic formulation, employs a radial mapping rule (Fig. 1) where point C is obtained by radially pro- jecting point B from point A, which is assumed to be on the /"-axis. The projection center A is the counterpart of point Ic of the isotropic for- mulation, shown in figure 1 of Part II (11). The scalar distance between the origin and point.A is taken as CI"0 where C is a model parameter, referred to as the projection parameter. The gradient of the bounding surface at the image point C, defines the loading/unloading direction, as well as the direction of plastic strain rate for the actual stress point B. The plastic modulus Kp at the actual stress point a,y, is assumed to be a function of the plastic modulus Kp (bounding modulus) at the image point <r,y, and the scalar distance 8 between the image and the actual stress points. 1294 igure 2.4: R dial mapping rule from Anandar jah and Daf lias (1986) given by: Ĥ = 1 + ein   pa [z mh() + (1 zm)h0] h 9  @F=@ Ia 2 + 3  @F=@ Ja 2i (2.25) The model as tested in comparison with experimental data in undrained triaxial com- pression and extension tests on isotropically and K0-consolidated specimens. The use of the bounding surface with the projection center on evolving Ia-axis improved model perfor- mance during undrained shearing on K0-consolidated samples for the stress states within the bounding surface. In the radial mapping rule adopted by Dafalias and Herrmann (1986) and Anandarajah and Dafalias (1986) for isotropic and anisotropic cohesive soils, respectively, the projection center imposed a restriction on plastic str ins. Purely elastic strains occurred upon unloading until loading index L = ((@F=@ij) _ij) =Kp became positive as the denition of he plastic strains increment is _"ij = hLi @F=@ij, where Macauley brackets rendered null for negative L. This restriction might be appropriate for monotonic loading simulations, while unsatisfactory 21 2.2. Theoretical modeling results are expected for undrained cyclic loading simulations. 2.2.4 Bounding surface models for monotonic response of clays Although the study in this thesis is dedicated to simulations of cyclic clay behavior, it is also of interest to investigate incorporated bounding surface formulation proposed for monotonic loading simulations, as some important points and observations can be inferred. Whittle (1993) The bounding surface formulation was incorporated into the MIT-E3 model by Whittle (1993). The model was able to predict initial and induced anisotropy through a rotational hardening rule. A projection center was always positioned at the origin of the stress space for adopted radial mapping rule. An elastoplastic modulusHI was proposed based on K0-loading on normally consolidated clays to simulate parallel virgin consolidation line in compression plane for any radial consolidation path, strain-hardening, peak strength, strain-softening, and critical stress states. For the stress states within the bounding surface the actual elastoplastic modulus H was related to HI and value of H at rst yield (rst loading for stress states within the bounding surface) which is denoted by H0: H = HI +H0 g2 (2.26) where g2 is a mapping function describing the relative position of the current stress state. MIT-E3 model required calibration of 15 parameters, including initial void ratio. Model sim- ulations were then compared with experimental results in undrained compression (CK0UC) and undrained extension (CK0UE) tests on K0-consolidated clay. For OCR > 4, the model overestimated shear strength in compression and underestimated in extension. Jung and Yune (2011) Jung and Yune (2011) incorporated the bounding surface formulation into the ABSM model previously proposed by Banerjee and Yousif (1986). The main goal was to simulate small- strain behavior of a particular clay. This was achieved by the use of a nonlinear anisotropic 22 2.2. Theoretical modeling elastic model obtained from singular value decomposition of a compliance matrix. Further, they thoroughly investigated dierent projection center locations in adopted radial mapping rule. From stress probes in dierent directions, they concluded that the projection center at the last stress reversal rendered the most appropriate plastic strains. The plastic modulus for any stress state was modied by inclusion of an extra distance parameter to further enhance simulations at small strains. H = A + A0  (1 h)  =0 1 =0  + h  =1 1 =1  (2.27) where A is a plastic modulus at the image stress, A0 is a reference modulus at the apex of the bounding surface, h is a weighting factor, denition of  is the same as in previously described models, 0 is equal to  at initiation of a new loading process, 1 is a distance from the current stress point to the image stress point where the soil became stier prior to shearing. Computed stress{strain curves with 12 calibrated model parameters gave a good match for small strains observed in triaxial tests, while underestimated shear stiness at larger strains. Suebsuk et al. (2011) The bounding surface formulation was incorporated into Modied Cam-Clay model with de- structuration (MSCC) by Suebsuk et al. (2011). The origin of the stress space was adopted as the projection center in a radial mapping rule. Suebsuk et al. (2011) utilized ratio depen- dence of the plastic modulus via introduction of the image stress ratio, which was a function of current and image stress states:  = [(p=pj) + (q=qj)] (2.28) where (p; q) and (pj; qj) are current and image stress states, respectively. The plastic modulus was then given by the following equation: H = Hj + (hHj,i) (1 )2  (2.29) 23 2.2. Theoretical modeling where h is a model parameter, Hj is a plastic modulus at image point, Hj,i is equal to Hj at initiation of the loading process. The proposed model required calibration of 12 parameters: 6 parameters were the same as in the original Modied Cam Clay model, 5 parameters were required for destructuration simulations, and parameter h was needed for the bounding surface formulation. The model predictions were compared with undrained monotonic triaxial tests. Volumetric strains in the e-ln p plane induced by destructuration in consolidation process were closely matched by the model. Jiang et al. (2012); Ling et al. (2002) The model proposed by Ling et al. (2002) is an extension of the anisotropic model by Dafalias (1987), where rotational hardening was derived based on the critical state concept and deviatoric back-stress ratio  as a dissipating energy measurement. The associative 
ow rule and inherent anisotropy were used in the model; the bounding surface was a single ellipsoid in the multiaxial stress space. Isotropic, rotational, and distortional hardening rules were adopted for the bounding surface evolution. Distortional hardening was required to control the evolution of R. The projection center was always located on a back-stress ratio, representing stress- induced anisotropy, in the same way as proposed by Anandarajah and Dafalias (1986), while the plastic modulus was modied to reduce the number of model parameters as follows: Kp = Kp + 1 + ein   paRijRij   hr  si W (2.30) where Rij = @F=@ij and W is a model parameter replacing hc, he, z, and m, and denition of  is the same as in previously described models. In total, 12 model constants were required for simulations. Model calibration and perfor- mance were illustrated as compared with experimental data on isotropically and anisotrop- ically consolidated kaolin clay, San Fransisco Bay Mud, and Boston Blue clay to a range of OCRs. Overall model performance was enhanced for overconsolidated states; however, because of the associative 
ow rule, the model was not able to predict the softening observed in CK0UC tests on normally consolidated specimens. 24 2.2. Theoretical modeling Jiang et al. (2012) recently proposed an extension of the model discussed above. The non- associative 
ow rule proposed by Jiang and Ling (2010) was adopted in a new formulation. This 
ow rule considered the bounding surface and the plastic potential possessing the same back-stress ratio, controlling rotation of the two surfaces, and having dierent peak deviatoric stresses. For the plastic potential, such peak deviatoric stress represents the critical stress ratio at which innite shear strains occur without change in stresses, while peak deviatoric stress is used for simulations of the shear strength of soil. The projection center was not changed from the preceding model, while the plastic modulus was slightly modied in that the W power was omitted and the hardening parameter, h, was included in the second term of Eq. 2.30. With the newly proposed non-associative 
ow rule, softening during CK0UC stress path can be successfully simulated. Huang et al. (2011) Huang et al. (2011) proposed an anisotropic model with incorporated destructuration and the bounding surface formulation. The analytical description of the bounding surface was adopted from Ling et al. (2002), as was evolution of anisotropy fromWheeler et al. (2003). An exponential decay for a scalar valued variable, accounting for destructuration, was proposed with the rate of destructuration parameter and its initial value. The radial mapping rule was utilized in the bounding surface formulation with the projection center xed at the origin. The proposed plastic modulus considered degradation via introduced cumulative plastic deviatoric strain: Hp = Ĥp + pa " @F @p 2 +  @F @q 2#" 0 0     1 # (2.31)  =  0 exp ("ps ) (2.32) where "ps = R _"ps is a cumulative plastic deviatoric strain; ,  o, and  are model parame- ters; the  function was used to simulate degradation of stiness with accumulated devia- toric strain; and denitions of  and 0 are the same as in previously described models. A total of 13 model parameters were required for model simulations. Simulations were per- 25 2.2. Theoretical modeling formed for undrained triaxial compression on normally consolidated and overconsolidated clay. While shearing of normally consolidated specimens was simulated relatively well, lower shear strength was obtained for overconsolidated specimens than observed in experiments. The destructuration eect was captured to some extent in the compression plane. 2.2.5 Bounding surface models for cyclic response of clays The review of the bounding surface models that have been developed and evaluated for clay response in undrained cyclic loading is conducted below. Only some models are selected and presented in this discussion, as they illustrate some important features for cyclic loading simulations. Other models available in the literature are either somehow related to the ones discussed below or found to be irrelevant for elaboration in the current work (e.g., viscoplasticity models as the time eect is not being considered). Zienkiewicz et al. (1985) A simple model was proposed by Zienkiewicz et al. (1985) which incorporated the bounding surface formulation. The projection center xed at the origin was considered in the radial mapping rule. Zienkiewicz et al. (1985) distinguished the plastic modulus in unloading and reloading for the stress states within the bounding surface. Unloading was dened for the stress increments pointing towards the interior of the loading surface, passing through the current stress point and the origin, and reloading for the stress increments pointing outwards. Very simple plastic moduli for unloading and reloading were considered. The former was given by: HU = H cs L  0  
1 (2.33) and the latter by HL = H cs L  0   (2.34) where HcsL is a plastic modulus at the image stress, 
 and 
1 are model constants, and denitions of  and 0 are the same as in previously described models. In this way the plastic modulus was larger than HcsL for the stress states inside the bounding surface and was equal to HcsL for the stress states on the bounding surface. 26 2.2. Theoretical modeling The proposed model simulations were then compared with actual experimental data from undrained cyclic triaxial tests. It was clear that the model predicts considerable overdamp- ing as opposed to experimental data and a way to model plastic strains in unloading was required. Based on observations, Zienkiewicz et al. (1985) suggested the use of cumulative plastic deviatoric strain, as the proposed model was not able to simulate accumulation of plastic strains for clays once the stress path reached equilibrium in undrained cyclic triaxial loading. This was also done in the hope of achieving further reduction of stresses to simulate liquefaction for sands; however, it was still not accomplished with inclusion of cumulative plastic deviatoric strains. Zienkiewicz et al. (1985) expressed the need to simulate plastic strains in unloading to avoid overdamping. Al-Tabaa and Wood (1989); Stallebrass and Taylor (1997) Al-Tabaa and Wood (1989) extended the Modied Cam-Clay model by adding a kinemat- ically hardening yield locus within the Cam-Clay locus. The adopted translation of the yield locus guaranteed that the outward normal of the inner yield locus was in the same direction as the outward normal of the Cam-Clay locus. This avoided the possibility that the two loci would intersect. The response inside the yield locus was assumed to be elastic; once the stress point was on the yield locus inelastic strains took place. While the inner yield locus was allowed translating within the Cam-Clay locus, the latter was allowed only expanding or shrinking, as in the Modied Cam-Clay model. Expansion/shrinkage of both loci occurred any time plastic volumetric strains took place. The degree of plasticity for the stress states on the yield locus and inside the Cam-Clay locus was dependent on the distance from the current stress point to its conjugate point on the Cam-Clay locus. This enabled smooth transition of soil response from the initial stress point on the inner yield locus to the stress point on the Cam-Clay locus. This approach of smooth transition is somewhat similar to that proposed for soils by Dafalias (1986b); the dierence is that, instead of keeping in memory an additional memory parameter or having a projection center, the memory of the kinematic yield surface was required. The model was attractive in that it required only 2 more parameters in addition to those required for the Modied Cam Clay model. With the 6 model parameters close match of 27 2.2. Theoretical modeling volumetric and shear strains accumulation was achieved during one-way undrained cyclic triaxial loading consisting of 6 cycles. Stallebrass and Taylor (1997) further incorporated an additional kinematic hardening surface, history surface, which was always inside the outer surface and enclosed the inner yield locus. The history surface was required to account for an eect of recent stress history and yield at small strains. The translation laws, ensuring that the three surfaces never intersect, required that the centers of the surfaces moved in the direction of a vector from the current stress state to its conjugate point on the next surface. The surfaces translations occurred anytime plastic deformations took place. The proposed 3-SKH model required 8 model parameters. The model performance was illustrated for a single element simulations in monotonic and two cycles of undrained cyclic triaxial loading with overpredicted stiness reduction as opposed to observed experiments. Stallebrass and Taylor (1997) implemented the proposed model into nite element code and simulated the available centrifuge test of a stresses distribution underneath a foundation. The analysis showed relatively good model predictive capabilities under a low level of loads, while for increased loads the model underpredicted displacements underneath the foundation. Wathugala and Desai (1993) The Hierarchical Single Surface (HiSS) model approach by Desai et al. (1986) for isotropic material with associative 
ow rule was used by Wathugala and Desai (1993) as a basis model; the bounding surface formulation was incorporated for plastic strains generation inside the conventional yield surface. Wathugala and Desai (1993) suggested that other complexity features (e.g., non-associative 
ow rule and anisotropy) could be included in the proposed model if desired. The model considered virgin loading, unloading, and reloading separately. Virgin loading corresponded to the stress states on the conventional yield surface and increments of a stress directed outwards of the conventional yield surface. For virgin loading, the usual approach of nding the plastic modulus from the consistency condition was used. A reference surface, passing through the current stress point, was used to distinguish between unloading and reloading in the same way as loading surface is used in other models. 28 2.2. Theoretical modeling This surface was similar in shape to the conventional yield surface and also passed through the origin; both surfaces were homologous to each other with the center of homology lying at the origin. Based on laboratory tests, Wathugala and Desai (1993) proposed unloading to be essentially elastic, while for reloading the plastic modulus was a function of the plastic modulus for virgin loading such that the former approached the latter: HRL = HVLI1 +H VL I2 r1  1 ps r r2 (2.35) where HVLI1 and H VL I2 are plastic moduli at image stress and intersection of the conventional yield surface with the hydrostatic axis, respectively, r1 and r2 are model constants, ps and r are distances from the origin to the current stress state and to the image stress point, respectively. The ps=r was used to impose distance dependency of the plastic modulus for the reloading process. The model required calibration of 12 parameters based on isotropic consolidation and undrained triaxial compression/extension tests. Simulations of undrained cyclic triaxial test, consisting of 4 cycles of two-way loading, suggested that the model was able to capture gener- ation of plastic deviatoric strains within a considered number of cycles, while it overestimated damping as a unloading was simulated as elastic. The model was then implemented into nite elements code for pile installation, consolidation after installation, application of load tests, and nal cyclic loading simulation; the whole procedure of model implementation, calibration, and simulation was described by Desai et al. (1993). The results obtained gave a relatively good match to eld observations of pile segment displacement, shear transfer, and horizontal stress. Yu et al. (2007) A Clay And Sand Model (CASM) proposed by Yu (1998) is a unied model for simulations of soil behavior; CASM is based on the Modied Cam-Clay model with 2 extra parameters required to control the shape of the conventional yield surface. One of the parameters controled the position of the point where the conventional yield surface intersected with the CSL, while the other controled a \drop-like" shape. 29 2.2. Theoretical modeling Incorporation of the bounding surface formulation into CASM model required 5 extra parameters, rendering in total 12 parameters. The radial mapping rule assumed the projec- tion center at the origin of the stress space. The proposed model treated the plastic modulus for virgin loading, unloading, and reloading separately. The description of the virgin loading plastic modulus is omitted and interested readers are referred to the original paper. The plastic modulus for unloading process was given by: H = HU 1 1 
 (2.36) where HU is an unloading hardening parameter, representing the limit for H in unloading, 
 = p=pj = q=qj is a non-dimensional ratio between current and image stresses. Thus, for an unloading process with initial stress state on the bounding surface, HU  H  1. In reloading the plastic modulus was given by: H = Hj +HR  1   1 + "pq k (2.37) where Hj is equal to H at the image point, HR is a reloading hardening parameter, " p q is a cumulative plastic deviatoric strain, and k is a model parameter controlling the rate of a shakedown. In reloading, H = 1 for the stress states far away from the projected image point, and H ! Hj as 
 ! 1. The "pq allowed for stiness degradation as plastic deviatoric strains accumulated throughout loading. The model was further compared with bubble-type models proposed by Al-Tabaa and Wood (1989); Stallebrass and Taylor (1997) in simulations of one-way and two-way undrained cyclic triaxial loading with available experimental results. In terms of deviatoric strains, the proposed model illustrated better predictive capabilities and was satisfactory for volumetric strains. For two-way cyclic loading simulations, the model predictions of deviatoric strains were also satisfactory. Li and Meissner (2002) Li and Meissner (2002) proposed a two-surface plasticity model with an associative 
ow 30 2.2. Theoretical modeling rule and incorporated bounding surface formulation for simulation of clay cyclic undrained behavior. The bounding surface represented the consolidation history of clay, while the loading surface passing through current stress point was related to each loading process. Li and Meissner (2002) utilized the kinematic hardening rule of the bounding surface that allowed its translation in the stress space to simulate induced anisotropy. Expansion and contraction of the bounding surface were also permitted. The elastic domain was reduced to a single point at the last stress reversal. The plastic modulus within the bounding surface was dependent on the relative congurations of the loading and bounding surfaces and varied from its local to prescribed value on the bounding surface as the stress point approached the latter. The authors extended Masing's rule of stress reversal to a multiaxial case. They dened a memory center as the point where loading changed its direction so that it pointed towards the interior of the loading surface, i.e., the memory center was the point at the last stress reversal. The loading and bounding surfaces tangentially touched each other at the memory center and never intersected. The loading surface was homologous to the bounding surface and passed through the memory center and current stress point. Once the loading surface coincided with the bounding surface, the plastic modulus reached its prescribed value on the bounding surface, and the two surfaces evolved identically. If the stress path was reversed and directed towards the interior of the loading surface, the memory center was updated. A new loading surface evolved, starting from the point at the memory center, and the bounding surface was translated to the same point with its size unchanged so that the two surfaces were tangential. For the stress reversals occurring inside the bounding surface the latter indeed translated from its previous position to the memory center, whereas for the stress reversals on the bounding surface the latter stayed in the same position. To simulate cyclic degradation of clay stiness, Li and Meissner (2002) used a varying power parameter 
 enabling reduction of the plastic modulus: hm = Hm +  HM  Hm   0  (2.38) where Hm is a plastic modulus for the stress points on the bounding surface, HM is a local 31 2.3. Missing features in reviewed models value of the plastic modulus, and denitions of  and 0 are the same as in previously described models. The 
 was constant for monotonic loading, while for cyclic loading it increased hyperbolically with cumulative shear strain. It should be noted that the proposed HM required specication of the applied dynamic shear stress amplitude qd and monotonic shear stress at failure qf, which must be known prior to cyclic loading. This makes the model not applicable for simulation of irregular cyclic loading, where amplitude of shear stress cannot be specied in advance. The model required 12 parameters, including qd and qf, for cyclic loading simulations. Li and Meissner (2002) did not indicate if the clay tested in monotonic and undrained cyclic triaxial loadings was the same and prepared under similar conditions; however, a dierent set of parameters was used for monotonic and each cyclic loadings. Monotonic loading simulations closely matched experimentally observed clay responses. One-way and two-way undrained cyclic triaxial tests were well simulated in terms of reduction of p and accumulation of shear strains. 2.3 Missing features in reviewed models The presented literature review on experimental undrained cyclic response of clay suggests plastic strains should be re
ected in model simulations for the stress states inside the con- ventional yield surface. The stress path migration towards the origin was accompanied by plastic deviatoric/shear strains development. When the stress path almost stopped migrat- ing towards the origin, plastic deviatoric/shear strains kept developing with an increasing rate (Takahashi et al., 1980; Zergoun, 1991; Zergoun and Vaid, 1994). The bounding surface formulation is used by dierent researchers to simulate the cyclic response observed in the experiments. However, Zienkiewicz et al. (1985) noticed that model simulations with the plastic modulus possessing only the distance dependence lacked the feature of simulating in- creased rate of deviatoric/shear strain accumulation once the stress path stopped migrating towards the origin. They proposed cumulative plastic deviatoric strain for continuous reduc- tion of plastic modulus throughout applied cyclic loading. A similar approach for reducing the plastic modulus during simulations of undrained cyclic triaxial loading was used by Yu 32 2.3. Missing features in reviewed models (1998) for the isotropic model, and by Li and Meissner (2002) for the kinematic hardening two-surface model. Huang et al. (2011) used the same idea for monotonic loading simulations in an anisotropic model with destructuration. In addition to the observed continuous de- velopment of plastic deviatoric/shear strains during undrained cyclic triaxial tests, a certain cyclic stress amplitude existed below which the rate of deviatoric/shear strains development decreased, This trend was observed in all experiments (Bacchus, 1969; Sangrey et al., 1969; Takahashi et al., 1980; Zergoun and Vaid, 1994). Another important conclusion made by Zienkiewicz et al. (1985) was that with the pro- jection center xed at the origin simulated unloading exhibited purely elastic strains, while plasticity was in eect suggested by experiments. This resulted in overdamping produced in simulations of undrained cyclic triaxial loading. The projection center at the last stress reversal appears to solve a problem of overdamping in simulations of clay response. Jung and Yune (2011) concluded, based on stress probe experiments, that the projection center at the point of the last stress reversal predicted plastic strain increments closer to those observed in experiments; however, their focus was on application of the bounding surface formulation for simulations at small strains in monotonic loading rather than cyclic loading. Li and Meissner (2002), based on Masing's rule generalization, adopted the projection center xed at the last stress reversal and were able to obtain relatively good predictions of one-way and two-way undrained cyclic triaxial loading. The main disadvantage of their model was that it required the knowledge of the applied stress amplitude prior to cyclic loading. Some of the models reviewed treated the plastic modulus in loading, unloading, and reloading separately (e.g., Wathugala and Desai, 1993; Yu, 1998; Zienkiewicz et al., 1985). This typically required a larger number of parameters (e.g., Yu, 1998; Zienkiewicz et al., 1985) and treatment of unloading involving only elastic strains (e.g., Wathugala and Desai, 1993). The goal of this thesis is to incorporate the bounding surface formulation into an ex- isting model for better predictions of clay response in undrained cyclic loading imposed by earthquakes. In the development of the model, the position of the projection center and simulation of cyclic softening will be addressed based on the observed experimental evi- dence and the current trend in available constitutive models. The complexities incorporated 33 2.3. Missing features in reviewed models into the developed model should be adequate to maintain its simplicity for further model implementation. 34 3. Development of the Proposed Model A proposed model is developed by incorporation of the bounding surface formulation into an existing reference model. First, the reference model itself is described with some simplica- tions in the triaxial plane and then generalized to the multiaxial space. The proposed model is then presented in the same manner. The bounding surface constants in the proposed model are studied in a sensitivity analysis. 3.1 Reference Model A Simple ANIsotropic CLAY plasticity model with destructuration (SANICLAY-D) pro- posed by Taiebat, Dafalias and Peek (2010) is chosen as a reference model. The SANICLAY- D model incorporates induced anisotropy and destructuration mechanism. Some simpli- cations for SANICLAY-D, which are believed to be sucient for successive simulations in monotonic loading, are applied. 3.1.1 Formulation in the triaxial space The formulation of the reference model in the triaxial space is presented in terms of the following stress and strain quantities: p = 1 3 (a + 2r) ; q = (a  r) (3.1) "v = ("a + 2"r) ; "q = 2 3 ("a  "r) (3.2) 35 3.1. Reference Model where subscripts a and r denote the axial and radial directions, respectively, and subscripts v and q denote volumetric and deviatoric components, respectively. Elastic relations With the total strain rate additive decomposition, _" = _"e+ _"p, into elastic (subscript e) and plastic (subscript p) parts, the isotropic hypoelastic relations _"ev = _p K ; _"eq = _q 3G (3.3) are adopted. K and G are bulk and shear moduli expressed as: K = p (1 + e)  ; G = 3K (1 + e) 2 (1 + ) (3.4) where e is a void ratio and  is a slope of the rebound line in the e{ ln p plane. Flow rule and yield surface The volumetric and deviatoric plastic strain rates are given by the 
ow rule _"pv = hLi @g @p ; _"pq = hLi @g @q (3.5) where g is the plastic potential and L is a plastic multiplier or a loading index to be dened in the sequel of the formulation. hLi are Macauley brackets rendering hLi = L when L > 0 and hLi = 0 when L  0. The 
ow rule requires the plastic potential g = 0, which is proposed by Dafalias (1986a): g = (q  p)2  M2  2 p (p  p) (3.6) where p is the value of p at q = p, adjusted so that g = 0 passes through a given (p; q),  is a non-dimensional anisotropic variable or rotational hardening variable, M is a critical stress ratio, Mc and Me are critical stress ratios in compression and extension, respectively. M = Mc when the stress ratio  = q=p > , andM = Me when the stress ratio  = q=p < . 36 3.1. Reference Model One must have jj < M for real-valued p and q in Eq. 3.40. The plastic potential is of the shape of a rotated and distorted ellipse, as illustrated in Fig. 3.1. N N Mc Me α ∇g (p, q) p q Yield surface (f=0) Plastic potential (g=0) Figure 3.1: Reference model schematic diagram in the triaxial plane On the basis of Eqs. 3.5 and 3.40, the 
ow rule is given by _"pv = hLi @g @p = hLip M2  2 ; _"pq = hLi@g@q = hLi2p (  ) (3.7) The yield surface was assumed similar to the plastic potential by Dafalias (1986a), Dafalias et al. (2006), and Taiebat, Dafalias and Peek (2010) with M replaced by N and  by , where N represented peak q-stress, and  was an anisotropy variable for the yield surface. Jiang and Ling (2010) proposed the yield surface possessing the same anisotropy variable as the plastic potential. This approach is attractive because it eliminates storing the extra tensor. Following Jiang and Ling (2010), the expression for the yield surface is 37 3.1. Reference Model then given by f = (q  p)2  N2  2 p (p0  p) (3.8) where p0 is the size of the yield surface and represents the isotropic hardening variable, and  is the mentioned anisotropic variable. Clearly, one must have jj < N for real{valued p; q in Eq. 3.8. The yield surface is also shown in Fig. 3.1. Rate evolution equations for the internal variables Evolution laws for internal variables p0 and  are necessary in the sequel of the formulation. The p0 evolution is given by _p0 = hLip0 (3.9) p0 is dened to account for the destructuration mechanism as p0 = Si p0d + Si p0d (3.10) where Si  1 is an isotropic structuration factor, and p0d is a destructured value of p0 (or the value of p0 when Si = 1). p0d and Si are given, respectively, by p0d =  1 + e    p0d  @g @p  (3.11) Si = ki  1 + e    (Si  1) "pd (3.12) where  is a slope of the normal compression line in the e{ ln p plane, ki is a rate of de- structuration. Degradation of Si is considered through destructuration plastic strain rate _" p d dened as _"pd = hLi"pd = q (1 A) _"p2v + A _"p2q (3.13) A is a parameter distributing destructuration between plastic volumetric strain and plastic deviatoric strain increments; by default A = 0:5. 38 3.1. Reference Model The evolution of rotational hardening  is expressed as _ = hLi = hLi  1 + e    C  p p0 2 @g@p  j  xj b   (3.14) where the denition of  is self evident, C is a model constant controlling the rate of rotational hardening, and b = min (N;Me) is a bounding value of  to guarantee real values of p; q in Eqs. 3.40 and 3.8. Loading index and plastic modulus For the completion of the model, loading index L and plastic modulus Kp are obtained from the consistency condition _f = 0 applied to Eq. 3.8, which in combination with Eqs. 3.9 and 3.14 yields L = 1 Kp  @f @p _p+ @f @q _q  (3.15) Kp =   @f @p0 p0 + @f @   (3.16) 3.1.2 Formulation in the multiaxial space In multiaxial stress and strain spaces second-order tensors are denoted by boldface symbols. Triaxial stress quantities are generalized to multiaxial space, given a stress tensor  that can be decomposed into its hydrostatic pI and deviatoric s components, as follows: p = tr 3 ; s =   pI (3.17) where tr is a trace operator and I is an identity tensor; Triaxial strain quantities are similarly generalized to multiaxial space, given that a strain tensor " can be decomposed into its volumetric "vI=3 and deviatoric e parts: "v = tr"; e = " "vI 3 (3.18) A systematic generalization of the triaxial constitutive relations is based on the observa- 39 3.1. Reference Model tion that under triaxial conditions the following holds true: x = [(3=2) x : x]1=2 (3.19) where x is any deviatoric stress tensor and x is its triaxial counterpart; the symbol : denotes the trace of the product of adjacent tensors. On the basis of Eq. 3.19: q = [(3=2) s : s]1=2 (3.20) "q = [(2=3) e : e] 1=2 (3.21)  = [(3=2)  : ]1=2 (3.22) Elastic relations The strain tensor and its rate are decomposed into elastic and plastic counterparts, and the same is applied for "v and e. The multiaxial generalization of Eq. 3.3 is given by _"e = _s 2G + _p 3K I (3.23) Flow rule and yield surface The generalization of a 
ow rule is as follows: _"p = hLi @g @ ; (3.24) where g is the plastic potential and L is the loading index. The Eq. 3.40 for the plastic potential generalizes to g = 3 2 (s p) : (s p)  M2  3 2  :  p (p  p) = 0 (3.25) 40 3.1. Reference Model M is chosen to be Lode angle dependent: M = 2m (1 +m) (1m) cos 3; m = Me Mc (3.26a) cos 3 = p 6 trn3; n = r [(r) : (r)]1=2 (3.26b) where  is a Lode angle in the range from 0 to =3 corresponding to the stress-ratio denition of compression and extension, respectively, where eective stress-ratio (r) rather than r is used to dene . Given stress tensor , Eq. 3.25 yields the following solution for p: p = (3=2) (s p) : (s p) (M2  (3=2) : ) p + p (3.27) The gradient @g=@ in Eq. 3.24 is calculated based on Eqs. 3.25 and 3.26, as follows: @g @ = 3 (s p) + 1 3 p  M2  3 2 r : r  I+ @g @ @ @ (3.28a) @g @ = 6M2p (p0  p)  (1m) sin 3 (1 +m) (1m) cos 3  (3.28b) @ @ = 3 n3  (trn3)n 1 3 I (1 + tr(n2) trn3tr(n)) p [(3=2) (r) : (r)]1=2 sin 3 (3.28c) tr  @ @  = p  M2  3 2 r : r  +  @g @  3 [tr (n2) trn3tr (n)] [(3=2) (s p) : (s p)]1=2 (1 6tr2n3)1=2 (3.28d) Similar to the plastic potential, the yield surface is generalized to f = 3 2 (s p) : (s p)  N2  3 2  :  p (p0  p) = 0 (3.29) 41 3.1. Reference Model N is assumed to be Lode angle independent, which signicantly simplies the generalization of the yield surface gradient @f=@: @f @ = 3 (s p) + 1 3 p  N2  3 2 r : r  I (3.30) Rate evolution equations for the internal variables The p0 evolution in multiaxial space is given by _p0 = hLip0 (3.31) p0 is dened to account for destructuration mechanism as in the triaxial plane p0 = Si p0d + Si p0d (3.32) p0d and Si are generalized, respectively, to p0d =  1 + e    p0d tr  @g @  (3.33) Si = ki  1 + e    (Si  1) "pd (3.34) where Si is considered through _" p d and dened in multiaxial space as _"pd = hLi"pd = q (1 A) _"p2v + A (2=3) _ep : _ep (3.35) _ep = _"p  (1=3) _"pv I is a plastic deviatoric strain tensor. The evolution law for the plastic potential and the yield surface rotational hardening is generalized to  = hLi (3.36a)  =  1 + e    C  p p0 2 tr @g@  32 (r x) : (r x) 1=2  b  (3.36b) 42 3.2. Proposed Model b = r 3 2 min (N;Me)nx; nx = r x [(r x) : (r x)] 12 (3.36c) Loading index and plastic modulus For the completion of the model, loading index L and plastic modulus Kp are obtained from the consistency condition _f = 0 applied to Eq. 3.29, which in conjunction with Eqs. 3.31 and 3.36 yields the following: L = 1 Kp  @f @ : _  (3.37) Kp =   @f @p0 p0 + @f @ :   (3.38) 3.2 Proposed Model The bounding surface formulation is incorporated in this section into the reference model described in section 3.1. The resultant model (SANICLAY-D-B) is rst presented in the tri- axial plane and then generalized to multiaxial space. It is pertinent to say that the bounding surface approach determines only the plastic modulus dierent from the one in the refer- ence model, while other constitutive relations are unchanged. The plastic modulus then requires consideration of extra complexities associated with the bounding surface formula- tion. It might be necessary to incorporate additional complexities into the proposed for better model predictions. 3.2.1 Formulation in the triaxial space Similar to section 3.1, triaxial stress and strain quantities are used in formulation of the proposed model. 43 3.2. Proposed Model Bounding surface, plastic potential and loading surface In the SANICLAY-D-B model a bounding surface replaces the yield surface that served as a consolidation history in SANICLAY-D. The bounding surface is then given by F = (q  p)2  N2  2 p (p0  p) = 0 (3.39) where p0 is a size of the bounding surface, (p; q) are image stresses on the bounding surface (bar over symbols implies relevance to the bounding surface) to be explained in the sequel of the formulation, and  and N are the same as in the reference model. The plastic potential passes through the image stresses (p; q) and is expressed as G = (q  p)2  M2  2 p (p  p) (3.40) In all the relations involved in the reference model, f = 0 and g = 0 are replaced by F = 0 and G = 0, respectively, to render the proposed SANICLAY-D-B model. A loading surface f = 0 is a surface passing through current stresses (p; q) and is homol- ogous to the bounding surface with a projection center (pc; qc) being the center of homology (to be dened in the sequel of the formulation), as shown in Fig. 3.2. The analytical expres- sion of the loading surface is not required for SANICLAY-D-B, but it is helpful to visualize the loading surface in order to understand the loading, unloading and reloading processes explained later. Image stress The bounding surface formulation postulates that for any stress (p; q) a unique image stress (p; q) exists on the bounding surface (Dafalias, 1986b). The image stress (p; q) is obtained by a radial mapping of a stress point (p; q) from a projection center (pc; qc) onto the bounding surface as shown in Fig. 3.2. Given a projection center (pc; qc), image stresses are expressed as p = pc + b (p pc) (3.41) 44 3.2. Proposed Model N N Mc Me α ∇F ∇f (p, q) (p̄, q̄) (pc, qc) p q Plastic potential (G=0) Bounding surface (F=0) Loading surface (f=0) Figure 3.2: Proposed model schematic diagram in the triaxial plane q = qc + b (q  qc) (3.42) where b is a similarity ratio between the loading surface f = 0 and the bounding surface F = 0. The value b can be obtained by plugging Eqs. 3.41 and 3.42 into Eq. 3.39. This is solved term by term as follows: (q  p)2 = b2A1 + bB1 + C1 (3.43) 45 3.2. Proposed Model where A1 = [(q  p) (qc  pc)]2 (3.44) B1 = 2q [(q  p) (qc  pc)] (qc  pc) (3.45) C1 = (qc  pc)2 (3.46) Similarly, p(p0  p) = b2A2  bB2  C2 (3.47) where A2 = (p pc)2 (3.48) B2 = (2pc  p0) (p pc) (3.49) C2 = pc (pc  p0) (3.50) Then a quadratic equation is obtained from Eqs. 3.43 and 3.47 Ab2 +Bb+ C = 0 (3.51) where A = A1=  N2  2+ A2 (3.52) B = B1=  N2  2+B2 (3.53) C = C1=  N2  2+ C2 (3.54) The positive root corresponds to the image stress (p; q) projected through the stress point (p; q) from the projection center (pc; qc). Projection center Determination of the image stress on the bounding surface requires the projection center. The projection center at the last point of stress reversal can be used to simulate plastic 46 3.2. Proposed Model strains in unloading, required for better undrained cyclic loading simulations. The issue of simulations plastic strains in unloading during cyclic loading was discussed by Zienkiewicz et al. (1985). Jung and Yune (2011) showed from stress probing tests that simulations with the projection at the last stress reversal are favorable for monotonic loadings. The model proposed in this thesis, thus, considers projection center (pc; qc) established at the last point of stress reversal. Geometrically, stress reversal is considered in the proposed model any time stress increments ( _p; _q) are pointing into the interior of the loading surface (see Fig. 3.2). The projection center (pc; qc) is updated at the moment of the stress reversal. However, if the projection center (pc; qc) is simply chosen to be a xed point until the next stress reversal occurs, then there could be a case such that the projection center stays outside the bounding surface. This is true for the reference model which adopts isotropic and rotational hardening of the yield surface (the bounding surface in the proposed model) relatively to the origin. The projection center outside of the bounding surface violates the concept of the bounding surface suggesting the uniqueness of the image stress for a given projection center and stress state. Thus, the idea of the projection center at the last stress reversal should be modied for the proposed model to maintain the desired eect of plastic strains upon unloading and, at the same time, to keep the well-established bounding surface concept. The proposed idea is to maintain the relative position of the projection center while the bounding surface evolves. Consider (pc, qc) as an updated projection center at the moment of stress reversal. If one requires that the relative position of the projection center is kept unaltered owing to only a change of the size of the bounding surface, then the projection center (pc; qc) should evolve accordingly: _pc = _p0 p0 pc (3.55) _qc = _p0 p0 qc (3.56) The relative position of the projection center (pc; qc) should also be maintained upon rotation and distortion of the bounding surface. For this to be in eect, the following ratio 47 3.2. Proposed Model should be kept constant: qc  q qb  q = X (3.57) where, referring to Fig. 3.3, q = pc; qb is the closest to the projection center point on the bounding surface for a given pc and X is a constant ratio. N N α pc qc qα qb p q Bounding surface (F=0) Projection center Figure 3.3: Proposed projection center The denominator of Eq. 3.57 can be obtained from F = 0 in Eq. 3.39 as qb  q =  N2  2 pc (p0  pc)1=2 (3.58) Then the rate of qc due to rotation and distortion of the bounding surface is obtained based on Eq. 3.57, which yields _qc = " pc X pc (p0  pc) [(N2  2) pc (p0  pc)]1=2 # _ (3.59) Notice that the rst term of Eq. 3.59 considers the change of qc due to rotation of the bounding surface without distortion, while the second term adds the eect of distortion 48 3.2. Proposed Model induced by the bounding surface rotation. The two terms in conjunction with each other guarantee that the ratio X is maintained constant and the projection center is never outside the bounding surface. The rate of the projection center evolution due to isotropic, rotational and distortional hardening is expressed by combining Eqs. 3.55, 3.56, and 3.59: _pc = _p0 p0 pc (3.60) _qc = _p0 p0 qc + " pc X pc (p0  pc) [(N2  2) pc (p0  pc)]1=2 # _ (3.61) Loading index determining the stress reversal The loading index L is the same as that given in the reference model by Eq. 3.15 and is used in the proposed SANICLAY-D-B model to dene the point of the stress reversal. The stress reversal is assumed to take place whenever L = 0. The projection center is then updated to the point of the stress reversal. The use of L = 0 rather than L < 0 is required to update the projection center to the point of initiation of shearing preceded by unloading in isotropic consolidation. The neutral loading in the bounding surface formulation corresponds to the stress path sliding along the loading surface, which renders L = 0 and as a result zero plastic strains. With the proposed update of the projection center the possibility of the neutral loading is eliminated. However, the stress path with L ! 0+ would correspond to almost neutral loading which does not involve update of the projection center. Plastic modulus For the stress states on the bounding surface, the plastic modulus is obtained from the consistency condition _F = 0 as Kp =   @F @p0 p0 + @F @   (3.62) 49 3.2. Proposed Model Then, following Dafalias (1986b), the plastic modulusKp for any stress state (p; q) is assumed to be a function of plastic modulus Kp at the image point (p; q) and the Euclidian distances between (p; q) and (p; q) denoted by  and between (pc; qc) and (p; q) denoted by r: Kp = Kp + h p30  hr  si = Kp + h p30 hb=(b 1) si (3.63) where h is a positive hardening function, r= = b=(b1) with 1  b  1 dened previously, and s is an indirect measure of the radius of the elastic nucleus, which is homologous to the bounding surface with the projection center as the center of homology. The presence of p0 is required to maintain similar units for the second term as for Kp, i.e., stress to the third power. For s = 1, the elastic nucleus coincides with the bounding surface, implying that the stress state inside the bounding surface considers purely elastic deformations. The elastic nucleus degenerates to a point with s = 1. Splitting the total plastic modulus into its image part and an additional term is a key feature of the bounding surface formulation. The second term in Eq. 3.63 controls Kp in the following way. Within the bounding surface where b >> 1, the term takes large values relative to Kp, also resulting in high Kp; as b ! 1, the term assumes smaller values and Kp ! Kp. Thus, the plastic modulus for the stress states within the bounding surface decreases from innite for b >> 1 to nite values, approaching the plastic modulus at image stress as b ! 1. This allows a smooth transition from purely elastic strains to elastoplastic strains on the bounding surface. The default value of the elastic nucleus is chosen to be s = 1 for the current formulation. One could argue that it is necessary to have an elastic region at the initiation of a new loading process following reversal of stresses. However, this could be readily achieved by the model to a certain extent, since the projection center at the last stress reversal renders b ! 1 at initiation of a new loading process, resulting in Kp !1, which is equivalent to an initially purely elastic response. Hardening function Zienkiewicz et al. (1985) proposed a cumulative plastic deviatoric strain to induce reduction of plastic modulus Kp, as the use of only the distance dependent Kp was not sucient for 50 3.2. Proposed Model undrained cyclic triaxial loading simulations. Similarly, in the proposed SANICLAY-D-B model hardening function h in Eq. 3.63 is assumed to be a decaying function of its initial value h0 and new internal variable d, simulating damage, h = h0 1 + d (3.64) As d starts evolving, h decreases from initial value h0, leading to reduction of Kp. The proposed d evolution has a linear relationship with an increment of the plastic deviatoric strain _"pq, _d = ad _" p q (3.65) where ad is a new model parameter that controls the rate of d evolution. Internal variable d evolves any time plastic deviatoric strains _"pq occur. 3.2.2 Formulation in the multiaxial space Generalization of the proposed model is performed following the observation given by Eq. 3.19: q = [(3=2) s : s]1=2 (3.66) where, as for triaxial quantities, the bar over bold symbols implies the relevance to F = 0. Bounding surface and plastic potential The analytical expression for the bounding surface, given by Eq. 3.39, is generalized to F = 3 2 (s p) : (s p)  N2  3 2  :  p (p0  p) = 0 (3.67) where N is assumed to be independent of , as in the reference model. The plastic potential form in multiaxial space is given by G = 3 2 (s p) : (s p)  M2  3 2  :  p (p  p) = 0 (3.68) where M is Lode angle dependent and is expressed by Eqs 3.26a and 3.26b in the reference 51 3.2. Proposed Model model with r substituted for r. The bounding surface and plastic potential, corresponding to the triaxial stress state shown in Fig. 3.2, are plotted in the stress space in Figs. 3.4(a{b) and 3.4(c{d), respectively. The bounding surface is a rotated and distorted ellipsoid with the anisotropic tensor  being a central axis of the ellipsoid. Observe that N -cone is formed in the multiaxial space, as the parameter N is chosen to be independent of Lode angle in the proposed formulation. The N -cone is xed in the stress space with the space diagonal being a central axis of the cone. In the triaxial plane the gradients of the bounding surface at the intersection of the bounding surface with the stress ratio N in compression and extension are orthogonal to the p-axis. Similarly, the gradients of the bounding surface at the intersection of the bounding surface with N -cone are orthogonal to the stress space diagonal (11 = 22 = 33) in the stress space. The plastic potential passes through the image point on the bounding surface (see Fig. 3.2). Similar to the bounding surface, rotation of the plastic potential is controlled by the anisotropic tensor , which is also a central axis of the plastic potential. In the stress space, the shape of the plastic potential is not similar to the shape of the bounding surface as the critical state failure envelope, chosen to be Lode angle dependent, controls the shape of the plastic potential. Critical state failure envelope is xed in the stress space with the stress space diagonal being its central axis. The gradients of the plastic potential at the intersection of the plastic potential with the critical state failure envelope are orthogonal to the stress space diagonal. Image stress Given projection center tensor c can be decomposed into its hydrostatic pcI and deviatoric sc counterparts, multiaxial generalizations of Eqs. 3.41 and 3.42 are given by p = pc + b (p pc) (3.69) s = sc + b (s sc) (3.70) The value of b is determined by plugging Eqs. 3.69 and 3.70 into Eq. 3.67, which in 52 3.2. Proposed Model (a) (b) (c) (d) Figure 3.4: Multiaxial representation of: (a{b) the bounding surface with N -cone, and (c{d) plastic potential with critical state failure envelope 53 3.2. Proposed Model multiaxial space is solved term by term as follows: 3 2 (s p) : (s p) = b2A1 + bB1 + C1 (3.71) where A1 = 3 2 (s sc) : (s sc) 3 (p pc) (s sc) : + 3 2 (p pc)2 :  (3.72) B1 = 3 (s sc) : sc  3 (pcs+ (p 2pc) sc) : + 3pc (p pc) :  (3.73) C1 = 3 2 sc : sc  3pcsc : + 3 2 p2c :  (3.74) Eq. 3.47 is unchanged in multiaxial generalization. Then a generalized form of Eq. 3.51 is given by: Ab2 +Bb+ C = 0 (3.75) where A = A1=  N2  (3=2) : + A2 (3.76) B = B1=  N2  (3=2) : +B2 (3.77) C = C1=  N2  (3=2) : + C2 (3.78) Similar to the model formulation in triaxial plane, positive root is of interest for determina- tion of the image stress given by p and s. Projection center Projection center evolution due to only the change in size of the bounding surface is given in multiaxial space by _pc = _p0 p0 pc (3.79) _sc = _p0 p0 sc (3.80) 54 3.2. Proposed Model and pc evolution is unchanged in generalization. The multiaxial stress space equation equivalent to Eq. 3.57 is (sc  s) : (sc  s) (sb  s) : (sb  s) = X (3.81) where sb and s = pc are deviatoric tensors generalizing qb and q, respectively. sb is determined along unit tensor nc, which is expressed as nc = (rc ) = [(rc ) : (rc )]1=2 (3.82) where rc = sc=pc. Thus, it is necessary to nd a scalar-valued multiplier ab to solve for sb sb = s + ncab (3.83) ab can be determined by plugging Eq. 3.84 into Eq. 3.67 yielding ab =  (N2  (3=2) : ) pc (p0  pc) (3=2)nc : nc 1=2 (3.84) Given ab, Eq. 3.81 can be used to give the rate of sc due to rotational and distortional hardening of the bounding surface: _sc = " pc X pc (p0  pc) [(3=2) : ] 1=2 [(N2  (3=2) : ) pc (p0  pc)]1=2 # _ (3.85) The complete evolution of the projection center in multiaxial stress space is then given by a combination of Eqs 3.80 and 3.85 for sc: _pc = _p0 p0 pc (3.86) _sc = _p0 p0 sc + " pc X pc (p0  pc) [(3=2) : ] 1=2 [(N2  (3=2) : ) pc (p0  pc)]1=2 # _ (3.87) 55 3.3. Sensitivity analysis on new model constants Loading index and plastic modulus The loading index L in multiaxial space is the same as that given in the reference model by Eq. 3.37. A plastic modulus corresponding to the stress states on the bounding surface is general- ized to Kp =   @F @p0 p0 + @F @   (3.88) Given the multiaxial generalization of b and Kp, the plastic modulus Kp for any stress state is of the same form given in Eq. 3.63: Kp = Kp + h p30  hr  i = Kp + h p30 hb=(b 1) si (3.89) Hardening function For the generalization of hardening function h, only d evolution must be modied for mul- tiaxial space: h = h0 1 + d (3.90) _d = ad [(2=3) _e p : _ep]1=2 (3.91) where _ep is a plastic deviatoric strain increment tensor. SANICLAY-D-B model generalization is complete at this stage. 3.3 Sensitivity analysis on new model constants The proposed model has two constants corresponding to the bounding surface. Constant h0 is an initial value of h used to determine the plastic modulus within the bounding surface. Constant ad regulates the rate of evolution of the newly introduced damage parameter d. Sensitivity analyses on h0 and ad are performed in this section with the set of model constants and initial internal variables shown in Tables 3.1 and 3.2, respectively. Destructuration is not considered for illustration purposes. 56 3.3. Sensitivity analysis on new model constants Table 3.1: Model parameters for sensitivity analysis Model constant Designation Value Elasticity  0:03  0:2 Critical state Mc 1 Me 1 N 1  0:15 Rotational hardening C 5 x 1:7 Destructuration (no eect) ki 0 Bounding surface h0 * ad * * These constants are to be analyzed for sensitivity. Table 3.2: Initial internal variables for sensitivity analysis Model internal variable Designation Value Size of the bounding surface p0 200 kPa Initial damage parameter d 0 Structuration factor Si 1 Initial void ratio ein 0:7 Orientation of plastic potential and bounding surface  0 Initial stress state (p; q) (200; 0) kPa 3.3.1 Parameter h0 The mapping rule employed in the model formulation involves a plastic modulus dependent on a relative distance from the current stress state to the image stress. Eq. 3.89 shows that the plastic modulus is also controlled by hardening function h, which is determined by parameters h0 and d. In this section stress and stress{strain path variations with altered values of h0 are illustrated; ad = 0 is adopted so that h = h0 throughout imposed loading. Simulations of undrained triaxial compression starting from an isotropically consolidated state (CIUC) test on overconsolidated clay are performed. Undrained monotonic triaxial loading is simulated for overconsolidation ratio OCR = 1:5 and 3 with the results shown in Figs 3.5(a{b) and 3.5(c{d), respectively. The projection center is updated to the points of stress reversals. The corresponding image point during the simulations starting from OCR = 1:5 renders a gradient dictating contractive response. For lower h0 values, the 57 3.3. Sensitivity analysis on new model constants second term of Eq. 3.89 renders lower Kp, and more plasticity is involved in the process, i.e., a higher degree of contraction and reduction of p. With increased h0, Kp is increased and less contraction is simulated. The projection center for simulations starting from OCR = 3 dictates the gradient of the plastic potential imposing negative plastic volumetric strains, i.e., dialative response. The dialation results in increased p, as observed in simulations shown in Fig. 3.5(c{d). A lower h0 provides a stress path that is more favorable to increasing p than that with a higher h0. 0 0.2 0.4 0.6 0.8 1 1.2 0 0.2 0.4 0.6 0.8 p/p0,in q/ p 0 ,in OCR=1.5 240 120 60 30 h0=15 (a) 0 2 4 6 8 100 0.2 0.4 0.6 0.8 ε a  (%) q/ p 0 ,in 120 240 60 30 h0=15 (b) 0 0.2 0.4 0.6 0.8 1 1.20 0.2 0.4 0.6 0.8 p/p o,in q/ p o ,in OCR=3 240 30 h0=15 120 60 (c) 0 2 4 6 8 10 0 0.2 0.4 0.6 0.8 ε a  (%) q/ p 0 ,in 120 240 60 30 h0=15 (d) Figure 3.5: Parameter h0 sensitivity in undrained monotonic triaxial loading (ad = 0, p0;in: initial size of the bounding surface): (a{b) OCR = 1:5, and (c{d) OCR = 3 For illustration of h0 in cyclic loading, 6 cycles of undrained cyclic triaxial loading are simulated with four dierent h0 values starting from the normally consolidated state and are shown in Fig. 3.6. As in monotonic loading, a lower plastic modulus is involved for lower h0 values, resulting in a more contractive undrained stress path and larger axial strains. 58 3.3. Sensitivity analysis on new model constants The stress path migration towards the origin is not observed once OCR = 2 is reached since the projection center invokes image stresses at the critical state ratio where no plastic volumetric strains are involved. On the basis of the volumetric counterpart of Eq. 3.23, zero plastic volumetric strain correspond to _p = 0. The consequence of the imposed stress path and constant h0 is that the stress{strain curve traces an unchanged path. To have a systematic calibration approach, h0 can be determined with ad = 0 based on initial cycles. These cycles should not involve considerable accumulation of plastic deviatoric strains, since this would require consideration of damage parameter d as well. Model users should be aware of small changes in the simulated response during initial cycles once ad 6= 0 is utilized and consider this while calibrating h0. 3.3.2 Parameter ad To demonstrate the eect of the rate ad, loading involving considerable accumulation of plastic deviatoric strains is required. Therefore, only cyclic loading simulations are considered for sensitivity analysis on ad. 6 cycles of undrained cyclic triaxial loading starting from the normally consolidated state are imposed in simulations with dierent values of ad. h0 = 400 is considered for the simulations to achieve a wider range of h reduction during the course of applied loading. Within the initial cycles little dierence is observed in the results shown in Fig. 3.7, since accumulation of plastic deviatoric strain "pq is negligible. As the loading proceeds, damage parameter d evolves, with the rate determined by parameter ad. For ad = 0, the stress path reaches the state where it traces the same path in subsequent cycles, as also shown in h0 sensitivity simulations, which results in the stress{strain tracing an unchanged path as well. As ad is increased from zero, damage parameter d evolves during shearing, leading to continuous reduction of h. Observe that for ad > 0 the stress{strain curves keep evolving, with increasing axial strains in both sides even when the corresponding stress paths has reached equilibrium. The rate of axial strains development is controlled by the parameter ad. This illustrates the capability of the current model to simulate cyclic softening. Parameter h0 should be determined prior to calibration of ad. If only monotonic loading is of interest, then the damage parameter eect might not be required and one can set ad = 0. 59 3.3. Sensitivity analysis on new model constants 0 0.2 0.4 0.6 0.8 1 1.2−0.4 −0.2 0 0.2 0.4 p/p0,in q/ p 0 ,in h0=15 (a) −4 −3 −2 −1 0 1 2 3 4−0.4 −0.2 0 0.2 0.4 ε a  (%) q/ p 0 ,in (b) 0 0.2 0.4 0.6 0.8 1 1.2−0.4 −0.2 0 0.2 0.4 p/p0,in q/ p 0 ,in h0=30 (c) −4 −3 −2 −1 0 1 2 3 4−0.4 −0.2 0 0.2 0.4 ε a  (%) q/ p 0 ,in (d) 0 0.2 0.4 0.6 0.8 1 1.2−0.4 −0.2 0 0.2 0.4 p/p0,in q/ p 0 ,in h0=60 (e) −4 −3 −2 −1 0 1 2 3 4−0.4 −0.2 0 0.2 0.4 ε a  (%) q/ p 0 ,in (f) Figure 3.6: Parameter h0 sensitivity in undrained cyclic triaxial loading (ad = 0, p0;in: initial size of the bounding surface): (a{b) h0 = 15, (c{d) h0 = 30, (e{f) h0 = 60 60 3.3. Sensitivity analysis on new model constants 0 0.2 0.4 0.6 0.8 1 1.2−0.4 −0.2 0 0.2 0.4 p/p0,in q/ p 0 ,in ad=0 (a) −4 −3 −2 −1 0 1 2 3 4−0.4 −0.2 0 0.2 0.4 ε a  (%) q/ p 0 ,in (b) 0 0.2 0.4 0.6 0.8 1 1.2−0.4 −0.2 0 0.2 0.4 p/p o,in q/ p o ,in ad=15 (c) −4 −3 −2 −1 0 1 2 3 4−0.4 −0.2 0 0.2 0.4 ε a  (%) q/ p o ,in (d) 0 0.2 0.4 0.6 0.8 1 1.2−0.4 −0.2 0 0.2 0.4 p/p0,in q/ p 0 ,in ad=30 (e) −4 −3 −2 −1 0 1 2 3 4−0.4 −0.2 0 0.2 0.4 ε a  (%) q/ p 0 ,in (f) Figure 3.7: Parameter ad sensitivity in undrained cyclic triaxial loading (h0 = 400, p0;in: initial size of the bounding surface): (a{b) ad = 0, (c{d) ad = 15, (e{f) ad = 30 61 4. Qualitative Model Performance SANICLAY-D-B model formulated in the preceding Chapter is qualitatively compared with typical experimental response of clay. The advantage of the proposed projection center is demonstrated. Then the model is compared with the previous SANICLAY family models in undrained cyclic triaxial loading and circular stress path in the deviatoric plane. The proposed model performance is also illustrated in generation of the plots often used in practice for determination of dynamic characteristics of clay. 4.1 The advantage of the proposed projection center In this section the projection center on anisotropic variable  is qualitatively compared with the projection center proposed in this thesis. Undrained cyclic triaxial loading simulations with the stress amplitude qcy=p0;in = 0:39, where qcy is an amplitude of applied cyclic loading and p0;in is an initial size of the bounding surface, are performed for this purpose. For a reference to the reader, typical experimental results from undrained cyclic triaxial test performed by Sheu (1984) are shown in Fig. 4.1. The projection center on rotational variable  was initially proposed by Anandarajah and Dafalias (1986) and used by many researchers for monotonic and cyclic loading simulations. It requires the use of the following statement for the projection center in terms of triaxial quantities: pc = Cp p0 (4.1) qc = pc  (4.2) where Cp is a projection center parameter relating the position of the projection center on 62 4.1. The advantage of the proposed projection center −10 −8 −6 −4 −2 0 2 4 6 8 10 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 ε a  (%) q/ p 0 ,in   q cy/p0,in=0.39experiment (15 cycles) Figure 4.1: Typical stress{strain path in undrained cyclic triaxial loading test after Sheu (1984)  to the size of the bounding surface p0. The projection center is illustrated in the triaxial plane in Fig. 4.2. Generalization of Eqs. 4.1 and 4.2 to multiaxial space is given by c = Cp p0 (+ I) (4.3) The results of 15 cycles of simulations of undrained cyclic triaxial loading with the pro- jection xed on anisotropic variable  are shown in Fig. 4.3. The simulations start from an isotropically normally consolidated state with h0 = 50 and dierent values of Cp. The values of ad are adjusted to obtain approximately similar size of the stress{strain loops and are shown on each plot. Other model constants and internal variables are the same as in section 3.3. Parameter Cp = 0, used for simulations shown in Fig.4.3(a), renders the projection center xed at the origin. In this case the loading surface passes through the origin and the stress point. Observe that the simulated stress{strain curve involves elastic strains in unloading. Recall that unloading is dened for the stress increments pointing towards the interior of the loading surface. No plastic strains are expected in unloading according to Eqs. 3.24 and 3.37 since the projection center xed on anisotropic variable  invokes the stress increments pointing towards the interior of the loading surface. The plastic strains resume once the stress increments start pointing outwards the loading surface. In compar- ison with Fig. 4.1, simulated stress-strain path renders highly overdamped response. With 63 4.1. The advantage of the proposed projection center (p, q) (p̄, q̄) (pc, qc) N N Mc Me α p0 Cpp0 ∇F ∇f p q Bounding surface (F=0) Plastic potential (G=0) Loading surface (f=0) Figure 4.2: Projection center on rotational variable  in the triaxial plane Cp = 0:6 simulations shown in Fig.4.3(b) the projection center is shifted to p = 0:6 p0 on anisotropic variable . In this case, for a given stress state , the loading surface is of smaller size than that in simulations with Cp = 0 owing to the center of homology shifted from the origin. This results into reduced elastic part of the stress-strain path. However, overdamp- ing is still present while comparing the simulated response with typical experimental results shown in Fig. 4.1. Similar simulation of undrained cyclic triaxial loading is performed with the projection center proposed in this thesis and the results are shown in Fig. 4.4. The value of ad = 25 is chosen to obtain stress{strain loops of comparable size with the loops in simulations results obtained with the projection center on anisotropic variable . All other parameters are the same as before. Recall from the model formulation that the projection center is updated to the point of the stress reversal at the instant of stress increment pointing towards the interior of the loading surface. Then new loading surface starts evolving with the updated projection center from a point at the last stress reversal and stress increments point outwards the loading surface. Thus, the projection center at the last stress reversal guarantees that 64 4.1. The advantage of the proposed projection center −10 −8 −6 −4 −2 0 2 4 6 8 10 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 ε a  (%) q/ p o ,in Projection center on α q cy/p0,in=0.39 Cp=0 h0=50 ad=75 (a) −10 −8 −6 −4 −2 0 2 4 6 8 10 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 ε a  (%) q/ p o ,in Projection center on α q cy/p0,in=0.39 Cp=0.6 h0=50 ad=28 (b) Figure 4.3: Stress{strain path simulations in undrained cyclic triaxial loading with the pro- jection center on rotational variable : (a) Cp = 0, and (b) Cp = 0:6 plastic strains will occur almost immediately after the reversal, although with a relatively smaller magnitude that increases as the loading proceeds. The results obtained suggest that the problem of overdamping is eliminated with the projection center proposed in this thesis. −10 −8 −6 −4 −2 0 2 4 6 8 10 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 ε a  (%) q/ p o ,in Proposed projection center q cy/p0,in=0.39 h0=50 ad=25 Figure 4.4: Stress{strain path simulations in undrained cyclic triaxial loading with the pro- jection center proposed in this thesis 65 4.2. Undrained cyclic triaxial loading 4.2 Undrained cyclic triaxial loading SANICLAY developed by Dafalias (1986a); Dafalias et al. (2006) is a model of Cam-Clay class with incorporated stress-induced anisotropy. Its extension, SANICLAY-D, proposed by Taiebat, Dafalias and Peek (2010) can be successfully used to simulate clay destruc- turation. So far, SANICLAY family models have been used only for monotonic loading simulations. The proposed SANICLAY-D-B model is a further extension, enabling one to simulate undrained cyclic loading. Its qualitative performance is demonstrated in this section in undrained cyclic triaxial loading with the constant stress amplitude. Simulations are rst performed with the previous SANICLAY family models. A summary of model simulations and corresponding model abbreviations are given in Table 4.1. Using SANICLAY-D-B, one can adjust model parameters shown in the table to reproduce the results that would cor- respond to previous SANICLAY family models. The model parameters and initial internal variables not shown in the table are given previously in Tables 3.1 and 3.2, respectively. Table 4.1: Summary of model parameters and initial internal variables for simulations of undrained cyclic triaxial loading Abbreviation ki Si h0 ad SANICLAY 0 1 10100 0 SANICLAY-D 2 3 10100 0 SANICLAY-D-B 0 1 60 0 SANICLAY-D-B (with damage) 0 1 60 25 The rst two simulation results of initial 6 cycles of undrained cyclic triaxial loading are obtained with the original SANICLAY and SANICLAY-D models, shown in Figs. 4.5(a{b) and 4.5(c{d), respectively. During the initial quarter cycle the stress state is always on the bounding surface so that the plastic strains simulated by both models. Observe that SANICLAY-D model invokes larger plastic strains owing to the destructuration eect. Upon unloading, the stress path lies within the bounding surface rendering only elastic strains for both models until the stress path reaches the bounding surface again in extension. Note that if no rotational hardening had been considered, as in the Modied Cam-Clay model, the stress path would not have reached the bounding surface in extension. Once the stress path reaches the bounding surface, plastic strains resume in the simulations. As in the rst quarter 66 4.2. Undrained cyclic triaxial loading cycle, SANICLAY-D model exhibits more plastic strains induced by the destructuration mechanism. For the subsequent cycles the stress path lies entirely within the bounding surface, leading to purely elastic response for both SANICLAY and SANICLAY-D models. Neither anisotropy nor destructuration induce further reduction of p. 0 0.2 0.4 0.6 0.8 1 1.2 −0.4 −0.2 0 0.2 0.4 p/p0,in q/ p o ,in SANICLAY ki = 0 Si = 1 h0 = 10 100 ad = 0 (a) −2 −1 0 1 2 −0.4 −0.2 0 0.2 0.4 ε a  (%) q/ p o ,in (b) 0 0.2 0.4 0.6 0.8 1 1.2 −0.4 −0.2 0 0.2 0.4 p/p0,in q/ p o ,in SANICLAY-D ki = 2 Si = 3 h0 = 10 100 ad = 0 (c) −2 −1 0 1 2 −0.4 −0.2 0 0.2 0.4 ε a  (%) q/ p o ,in (d) Figure 4.5: SANICLAY family models performance in undrained cyclic triaxial loading: (a{b) SANICLAY, and (c{d) SANICLAY-D Simulation results obtained with the proposed SANICLAY-D-B model are shown in Fig. 4.6. Destructuration is not considered in order to illustrate only the eect of the bound- ing surface. The results shown in Fig. 4.6(a{b) are obtained without damage eect. Notice the reduction of p within the bounding surface owing to allowed plastic strains. However, as the stress path reaches OCR = 2, the image stress projected from the projection center through the stress point falls onto the critical state ratio. As a result, neither compres- sion nor dialation are involved, and the stress path traces similar path in each cycle. The h = h0 renders the plastic modulus purely dependent on the relative distance from the stress 67 4.3. Circular stress path in the deviatoric plane point to the image stress. As it was explained in section 3.3, this is accompanied by tracing unchanging stress{strain path. The damage eect in SANICLAY-D-B improves the model performance through reduction of hardening function h as can be seen in Fig. 4.6(c-d). Even though the stress path traces the same route at OCR = 2, the corresponding stress{strain path exhibits continuously reducing stiness. In this way, SANICLAY-D-B with damage eect can be successfully used to simulate continuously evolving plastic strains in undrained cyclic loading. 0 0.2 0.4 0.6 0.8 1 1.2 −0.4 −0.2 0 0.2 0.4 p/p0,in q/ p o ,in SANICLAY-D-B ki = 0 Si = 1 h0 = 60 ad = 0 (a) −10 −8 −6 −4 −2 0 2 4 6 8 10 −0.4 −0.2 0 0.2 0.4 ε a  (%) q/ p o ,in (b) 0 0.2 0.4 0.6 0.8 1 1.2 −0.4 −0.2 0 0.2 0.4 p/p0,in q/ p o ,in SANICLAY-D-B ki = 0 Si = 1 h0 = 60 ad = 25 (c) −10 −8 −6 −4 −2 0 2 4 6 8 10 −0.4 −0.2 0 0.2 0.4 ε a  (%) q/ p o ,in (d) Figure 4.6: Comparison of model performance with bounding surface eect: (a{b) SANICLAY-B and (c{d) SANICLAY-B with damage 4.3 Circular stress path in the deviatoric plane Undrained cyclic triaxial loading tests are commonly used for model calibration to subse- quently simulate soil response induced by earthquake loadings. However, in reality, earth- 68 4.3. Circular stress path in the deviatoric plane quakes involve a more complicated stress path than imposed in triaxial tests. A qualitative model simulation of strains induced by a circular stress path in the deviatoric plane is illus- trated in this section. 4.3.1 Response to cyclic circular stress path loading Model parameters and internal variables from section 3.3 are considered for the model sim- ulations. Dierent combinations of h0 and ad, shown in Table 4.2, are used to reproduce the original SANICLAY and the proposed SANICLAY-D-B models with and without damage eect. Each simulation consists of three stages starting from normally consolidated state. In the rst stage, isotropic unloading from p = 200 kPa to 150 kPa is induced, corresponding to the stress path A-B in Fig. 4.7. During the second stage, stress path B-C in Fig. 4.7, q is increased from 0 to 50 kPa in triaxial compression with constant p. The applied q is such that the stress point at the end of the second stage and subsequently applied cyclic stress path lie inside the bounding surface. For the last and main stage, the circular stress path in the deviatoric plane is applied clockwise, during which p and q are kept constant. To induce considerable plastic strains within the bounding surface, three cycles of the circular stress path are simulated. Table 4.2: Summary of model parameters for simulations of circular stress path in the deviatoric plane Abbreviation h0 ad SANICLAY 10100 0 SANICLAY-B 100 0 SANICLAY-B (with damage) 100 60 The results of the three simulations in the strain space are shown in Fig. 4.8 with the de- viatoric view on the left and the shifted view on the right for better illustration of volumetric strains. Observe that the strain path followed in simulations with the original SANICLAY model in Fig. 4.8(a{b) traces unchanging circular strain path within the three cycles as the simulated response is completely elastic. The results of simulations with the eect of the bounding surface with and without damage eect are shown in Figs. 4.8(c{d) and 4.8(e{f), respectively. The two sets of results 69 4.3. Circular stress path in the deviatoric plane σ33/σc  σ11/σc C A,B  σ22/σc Imposed stress path (a) A Space diagonal B σ33/σc C σ22/σc σ11/σc (b) Figure 4.7: Isotropic unloading, constant-p triaxial compression, and circular stress path: (a) deviatoric view, (b) shifted view dier in that more plastic strains are involved with the considered damage eect. Notice that this time plastic strains are induced throughout the whole process of the applied loading. The imposed image stress dictates positive plastic volumetric strains as can be seen from the strain path in the shifted view. In the deviatoric view the strain path is slightly elongated along "33-axis during the rst cycle. Such response is due to the distance dependence of the plastic modulus Kp. At the initiation of the circular stress path loading projection center is updated to the point C shown in Fig. 4.7. Then the relative distance from the current stress point to the image stress reduces during the rst-half cycle, rendering increasing plastic strains. Shortly after the rst half-cycle is complete, the projection center is updated and simulated plastic strains become initially negligible. Again the relative distance reduces as the loading proceeds and plastic strains increase until the projection center is updated shortly after the completion of the rst cycle. In this way larger degree of plastic strains approximately after the completion of rst half-cycle and rst cycle results into the strain path somewhat elongated along "33-axis. Notice that update of the projection center in subsequent cycles occurs at later stages and the elongation of the strain paths is rotated clockwise. One could argue that it might be unrealistic for clays to exhibit larger plastic strains at particular points during circular stress path. However, the circular stress path imposes the 70 4.3. Circular stress path in the deviatoric plane ε33 (%)  ε11 (%)  ε22 (%) SANICLAY: h0=1e10, ad=0 C (a) ε22 (%) ε33 (%) Space diagonal ε11 (%) C (b) ε33 (%)  ε11 (%)  ε22 (%) SANICLAY−B: h0=100, ad=0 C (c) ε22 (%) ε33 (%) Space diagonal ε11 (%) C (d) ε33 (%)  ε11 (%)  ε22 (%) SANICLAY−B: h0=100, ad=60 C (e) ε22 (%) ε33 (%) Space diagonal ε11 (%) C (f) Figure 4.8: Observed strains in the strain space (deviatoric view on the left, shifted view on the right) induced by the circular stress path (clockwise) in the deviatoric plane: (a{ b) original SANICLAY model, (c{d) proposed SANICLAY-D-B without damage, and (e{f) proposed SANICLAY-D-B with damage eect 71 4.3. Circular stress path in the deviatoric plane direction of the stress increments close to being orthogonal to the bounding surface gradient, which renders small values of L on the basis of Eq. 3.37. According to Eq. 3.24, the resulting plastic strains are also small. Therefore, even though plastic strains are abruptly reduced to zero upon the update of the projection center, the overall plastic strain magnitudes are not large enough to cause major deviations from expected clay response. It is more important that, unlike the original SANICLAY model, the bounding surface formulation enables one to simulate plastic strains for the stress states inside the bounding surface. 4.3.2 Eect of parameter C in circular stress path loading Apparently the rate of the rotational hardening of the bounding surface can control the induced strain path to some extent during the applied circular stress path loading. Simulation results of one cycle of circular stress path loading similar to that discussed above but with dierent values of C, and xed h0 = 50, and ad = 0 (i.e., no damage eect) are shown in Fig. 4.9. It should be noted that similar eect of C would apply for the stress states on the bounding surface. ε33 (%)  ε11 (%)  ε22 (%) SANICLAY−B: h0=50, ad=0 C=3 9 15 Figure 4.9: Eect of parameter C on induced strain path during circular stress path loading 72 4.4. Generation of secant shear modulus and damping ratio variation with shear strain amplitude Observe that, for lower values of C, the bounding surface does not rotate signicantly and it has to expand more to accommodate the imposed stress path. The less rotated and more expanded bounding surface gives more freedom for the stress path to be directed towards the interior of the loading surface and the projection center is updated at earlier stages. Then the updated projection center invokes the gradients of plastic potential mainly on the left side of the bounding surface in the deviatoric view. As a result, the induced increments of strain tensor, given by a Eq. 3.24, and strain path are pointing mainly along "33-axis. The resulting strains at the end of a cycle are above initial strains. As value of C increases the update of the projection center occurs at later stages and the resulting strain path becomes rounded with the strain path below initial strains at the end of a cycle. 4.4 Generation of secant shear modulus and damping ratio variation with shear strain amplitude Equivalent linear analysis is widely used in practice for site response analysis, where equiva- lent shear modulus and damping ratio for a given strain amplitude are required for represen- tation of soil properties to simulate non-linear soil response. Variation of the representative soil properties with applied shear strain amplitude in cyclic loading is obtained from lab- oratory tests. In this section, the eect of the model constants h0 and ad on simulated variation of representative soil properties with a shear strain amplitude is studied. Simula- tions of undrained cyclic loading in Direct Simple Shear (DSS) test are performed with the constraints corresponding to an element at the center of a DSS device. Model constants and initial internal variables are given in Tables 3.1 and 3.2, respectively. In order to come up with representative soil properties, a hysteresis loop at 5th cycle is considered following Seed et al. (1986). A schematic diagram of a typical hysteresis loop is shown in Fig. 4.10. The hysteresis loop is represented by a closed loop as explained by Hardin and Drnevich (1972). The shear modulus is measured as a secant shear modulus Gsec determined by the extreme points on the hysteresis loop, while the damping is proportional to the area inside the loop (Seed et al., 1986). It is typical to normalize the secant shear 73 4.4. Generation of secant shear modulus and damping ratio variation with shear strain amplitude modulus Gsec by a maximum tangential shear modulus Gmax for illustration of results. The damping ratio can be obtained from commonly used equation Damping ratio (%) = AL 4AT 100 (4.4) where AL is an area of a hysteresis loop; and AT is a stored equivalent elastic energy. Figure 4.10: Denition of equivalent shear modulus and the quantities used for the determi- nation of damping ratio Typical experimental data from Vucetic and Dobry (1991) on variation of secant shear modulus Gsec and damping ratio with shear strain 
 is shown in Fig. 4.11 to the reader for qualitative comparison with the simulated results. The curves are obtained from the experiments on soils with a range of plasticity index (PI). The simulated variation of secant shear modulus Gsec with shear strain amplitude 
 is shown in Fig. 4.12. Values of h0 = 50, 200, and 800 are adopted with ad = 0 and ad = 60 for the results shown in Figs. 4.12(a) and 4.12(b), respectively. Gsec starts reducing at 
 = 0:1% for the simulations obtained with the original SANICLAY model (dash line). With the proposed model adopting ad = 0, Gsec starts reducing already at 
 = 0:01%. Further reduction of Gsec is controlled by the parameter h0, giving a lower Gsec for a given 
 as h0 is decreased. Observe that as 
 ! 10% obtained values of Gsec tend to converge. This is true for any values of h0 adopted in the simulations as the stress path has to reach the bounding surface, where the plastic modulus and the degree of plastic strains reach their 74 4.4. Generation of secant shear modulus and damping ratio variation with shear strain amplitude 10−4 10−3 10−2 10−1 100 101 0 0.2 0.4 0.6 0.8 1 G se c/G m a x   PI=200 PI=100 PI=50 PI=30 PI=15 (a) 10−4 10−3 10−2 10−1 100 101 0 10 20 30 40 50 60 Shear strain (%) D am pi ng  ra tio  (% )   PI=200 PI=100 PI=50 PI=30 PI=15 (b) Figure 4.11: Experimental data after Vucetic and Dobry (1988): (a) degradation of Gsec and (b) variation of damping ratio with limiting values. For the results obtained with ad = 60, Gsec at lower values of 
 does not dier signicantly from those ones obtained with ad = 0 as the damage parameter d is still not evolved suciently. The eect of d becomes more pronounced as 
 is increased and the resulting Gsec variation curves shift down. The curves tend to converge for large values of 
 as for simulations with ad = 0; however, larger 
 needs to be achieved for the curves to converge. In comparison with the experimentally obtained Gsec degradation shown in Fig. 4.11(a), the results produced by SANICLAY-D-B model with damage eect are more plausible than those obtained with the original SANICLAY model. However, the proposed 75 4.4. Generation of secant shear modulus and damping ratio variation with shear strain amplitude 10−4 10−3 10−2 10−1 100 101 0 0.2 0.4 0.6 0.8 1 G se c/G m a x   ad=0 h0=800 h0=200 h0=50 SANICLAY (a) 10−4 10−3 10−2 10−1 100 101 0 0.2 0.4 0.6 0.8 1 Shear strain, γ, (%) G se c/G m a x   ad=60 h0=800 h0=200 h0=50 SANICLAY (b) Figure 4.12: Simulations of Gsec=Gmax degradation with 
: (a) ad = 0 and (b) ad = 60 SANICLAY-D-B model is not able to replicate Gsec degradation at strains 
 < 0:01% as mainly elastic strains are predicted by the model. It is interesting to note that the variation of Gsec could potentially serve as a basis for calibration of the bounding surface parameters h0 and ad. The former could be calibrated based on strain levels up to 
 = 0:03%, and the latter based on larger strains. The variation of damping ratio with shear strain 
 is illustrated in Fig. 4.13 with the similar values of h0 as before. The values of ad = 0 and ad = 60 correspond to the results shown in Figs. 4.13(a) and 4.13(b), respectively. SANICLAY model (dash line) predicts only elastic strains for 
 < 0:4% and, as a result, no damping variation since the stress path 76 4.4. Generation of secant shear modulus and damping ratio variation with shear strain amplitude apparently does not reach the bounding surface at the 5th cycle. The damping ratio starts increasing abruptly at approximately 
 = 0:4% and reaches about 60% of damping ratio at 
 = 10%. The proposed model simulates the damping ratio that starts increasing at 
 = (0:01{0:1)%, depending on the adopted value of h0. The damping ratio initially evolves faster for lower values of h0 in the simulated results obtained with both ad = 0 and ad = 60. The rate of damping ratio evolution then drops down slightly for lower values of h0, and the rate increases for higher values of h0. The reason is a reduction of damping, represented by AL in Fig. 4.10, with decreasing value of h0. This renders lower damping ratio on the basis of Eq. 4.4. The results obtained with ad = 60 show more signicant drop down of the rate of the increase of damping ratio at larger 
 as damage parameter d is considerably evolved at 5th cycles. Clearly, SANICLAY-D-B model with the damage parameter eect can better simulate experimentally obtained variation of damping ratio with shear strain 
, shown in Fig. 4.11(b), than the original SANICLAY model and SANICLAY-D-B model without the damage eect. Similar to Gsec degradation, the model cannot predict variation of damping ratio for 
 < 0:01% owing to simulated elastic strains only. Results obtained above are based on the model simulations with the xed set of model parameters, and their change may result in a dierent outcome. The parameters are typically calibrated based on experimental data, and each set corresponds to a particular soil. Some of the parameters can be related to particular soil properties. For example, parameter  can be correlated to PI (Wood, 1990):   0:6PI (4.5) This renders approximate PI = 25% for  = 0:15 used in the simulations. Referring to the above, the simulations can be approximately compared with the curves corresponding to PI = 30% in Fig. 4.11. It then follows that the model is able to successfully simulate the variation of the representative soil properties for 
 > 1%. Satisfactory model performance could be considered for 
 = (0:1{1)%, and much less satisfactory for 
 < 0:01%. For practical purpose, strains around well-designed geotechnical structures, such as retaining walls, foundations, and tunnels, are generally small (Clayton, 2011), and might be considered of about 
 = (0:01{1)%. However, for deformation analyses of dams and natural slopes, it 77 4.5. Generation of cyclic stress ratio versus number of cycles 10−4 10−3 10−2 10−1 100 101 0 10 20 30 40 50 60 D am pi ng  ra tio  (% )   ad=0 h0=800 h0=200 h0=50 SANICLAY (a) 10−4 10−3 10−2 10−1 100 101 0 10 20 30 40 50 60 Shear strain, γ, (%) D am pi ng  ra tio  (% )   ad=60 h0=800 h0=200 h0=50 SANICLAY (b) Figure 4.13: Simulations of damping ratio variation with 
: (a) ad = 0 and (b) ad = 60 is important to consider larger strains (
 > 1%). Thus, the proposed model is adequate for most geotechnical problems. 4.5 Generation of cyclic stress ratio versus number of cycles The eect of model parameters h0 and ad is considered in this section for generation of commonly used in practice plots of cyclic stress ratio (CSR) required to reach a particular 78 4.5. Generation of cyclic stress ratio versus number of cycles strain in a number of cycles. Simulations of undrained cyclic triaxial loading are performed for this purpose. CSR = cy=c = qcy=2c is typically used for sands, and CSR = cy=Su = qcy=2Su for clays, where c is a consolidation stress and Su is an undrained shear strength used for normalisation. The latter CSR is adopted for simulations presented in this section with Su = 64 kPa determined from undrained triaxial compression test simulation starting from isotropically normally consolidated state. Plots are generated for a given strain, representing adopted failure criterion, and results are obtained for a particular CSR and required number of cycles to reach this strain. In this section plots are generated for CSR required to reach peak "a = 3%. Experimental data after Zergoun and Vaid (1994) is shown in Fig. 4.14 to the reader for qualitative comparison with the simulated results. 0 20 40 60 80 100 0 0.2 0.4 0.6 0.8 1 1.2 Number of cycles Cy cli c st re ss  ra tio , q cy /2  S u   experiment Fitted curve Figure 4.14: Experimental data after Zergoun and Vaid (1994) on cyclic stress ratio CSR required to reach "a = 3% in a number of cycles Simulations are performed with the model parameters used in previous sections, except parameters h0 and ad. Each plot in Fig. 4.15 corresponds to simulations with values of h0 = 50; 200; and 800, and xed ad. The range of CSR = (0:4{1) is considered in simulations. Only results within 100 cycles are used for illustration of qualitative model performance. The dash line corresponds to the original SANICLAY model simulations in Fig. 4.15. Observe that "a = 3% can be reached in one cycle only with CSR = 1, corresponding to the stress level equivalent to undrained shear strength Su. The "a = 3% cannot be attained with CSR < 1, rendering number of cycles ! 1 as suggested by the horizontal dash line. 79 4.5. Generation of cyclic stress ratio versus number of cycles 0 20 40 60 80 100 0 0.2 0.4 0.6 0.8 1 1.2 Cy cli c st re ss  ra tio , q cy /2  S u   ad=0 h0=800 h0=200 h0=50 SANICLAY (a) 0 20 40 60 80 100 0 0.2 0.4 0.6 0.8 1 1.2 Number of cycles Cy cli c st re ss  ra tio , q cy /2  S u   ad=60 h0=800 h0=200 h0=50 SANICLAY (b) Figure 4.15: Generation of CSR required to reach "a = 3% SANICLAY-D-B simulation results shown in Fig. 4.15(a) are obtained with ad = 0, rendering h = h0 throughout cyclic loading. Notice that simulations with dierent values of h0 fall very close to each other for CSR = 1 and 0:9. Only simulations with h0 = 50 attain "a = 3% within 100 cycles for CSR < 0:9, and none of the simulations fall within 100 cycles for CSR  0:6. In Fig. 4.15(b), simulations are performed with ad = 60, rendering h decreasing from its initial values of h0 as damage parameter d evolves. This time fewer cycles are required to reach "a = 3% for all simulations with the same h0 values as in Fig. 4.15(a). It should be noted that the simulated damage eect becomes more pronounced with larger number of cycles 80 4.6. Generation of peak axial strain and residual pore water pressure versus number of cycles The results simulated with the original SANICLAY model are clearly far from the experi- mental data shown in Fig. 4.14, while simulations with the proposed SANICLAY-D-B model signicantly improve the results. One could argue that the use of the damage parameter d is not required as certain values of h0 and ad = 0 provide the results (Fig. 4.15(a)) that are close to those observed in experiments (Fig. 4.14). However, CSR versus number of cycles plot is generated only for a particular strain level, which neglects the strain history from before and after the strain level. Thus, dierent combination of h0 and ad can provide one with approximately similar CSR versus number of cycles plot, while at any other strain levels the simulated response would be completely dierent. In other words, a close qualitative match to experimental data obtained with, say, h0 = 800 and ad = 0 can be achieved with much higher h0 and ad 6= 0, which might be required for better simulations at strain levels other than "a = 3%. From the above simulations, the eect of h0 and ad on the generation of CSR versus number of cycles can be summarized:  as h0 is decreased, simulated soil response becomes softer and a lower plastic modulus Kp is involved for the stress states within the bounding surface; this results into smaller number of cycles required to bring the simulated response to a certain strain level;  with increased ad, the damage parameter d evolves more rapidly as plastic deviatoric strains accumulate throughout cyclic loading; increased damage parameter d reduces h from h0, and the plastic modulus Kp decreases, leading to larger plastic strains within the bounding surface. 4.6 Generation of peak axial strain and residual pore water pressure versus number of cycles In this section eect of model parameters h0 and ad are demonstrated for generation of peak axial strain "a and residual pore water pressure u versus number of cycles. Simulations of undrained cyclic triaxial loading with qcy=2Su = 0:5 of applied cyclic stress amplitude are performed with the model parameters used in previous sections, except parameters h0 and 81 4.6. Generation of peak axial strain and residual pore water pressure versus number of cycles ad. Typical experimental results from undrained cyclic triaxial loading tests with cy=Su = qcy=2Su = 0:69 reported by Zergoun and Vaid (1994) are shown in Fig. 4.16 as a reference for the reader. 100 101 102 −10 −8 −6 −4 −2 0 2 4 6 8 10 Pe ak   ε a (% )   q cy/2Su=0.69 experiment (a) 100 101 102 0 0.2 0.4 0.6 0.8 1 Number of cycles R es id ua l u/ σ e   q cy/2Su=0.69 experiment (b) Figure 4.16: Development of: (a) peak axial strain "a, and (b) residual pore water pressure u versus number of cycles, experimental data from Zergoun and Vaid (1994) The results obtained with the original SANICLAY model are shown with dash line in Fig.4.17. Observe that plastic strains are induced only during the initial loading as the stress state lies on the bounding surface and during the following reloading once the stress path hits the bounding surface in extension. In subsequent cycles, no plastic strains are simulated as suggested by the horizontal dash line. This is accompanied by residual pore water pressure 82 4.6. Generation of peak axial strain and residual pore water pressure versus number of cycles u accumulation only during simulated plastic strains, i.e., no further build-up of residual u. For the simulation results obtained with dierent values of h0 and xed ad = 40 shown in Fig. 4.17, parameter ad 6= 0 is chosen to induce accumulation of plastic strains for the stress path reaching equilibrium at OCR = 2. Observe that for lower values of h0 fewer number of cycles are required to induce a certain strain level, shown in Fig. 4.17(a). Accumulation of residual u in Fig. 4.17(b) can be interpreted to indicate when the stress path reaches equilibrium, where no further development of residual u is simulated. Clearly, the lower h0 values, the faster such equilibrium is reached since the plastic strains induce larger " p v resulting into faster p reduction. Note that no further peak "a would have been possible once OCR = 2 is reached, should ad = 0 had been utilized. In comparison with experimental data shown in Fig. 4.16, the simulations with the orig- inal SANICLAY model are not able to realistically replicate the observed response of clay. The proposed SANICLAY-D-B model with the evolving damage parameter d reproduces the observed response qualitatively well, although residual u is about twice smaller than observed in experiments owing to the equilibrium reached in simulations by the stress path at OCR = 2. In Fig. 4.18, the results are obtained with xed h0 = 50 and dierent values of ad. Simulations results obtained with ad = 0, render h = h0 throughout the applied cyclic loading. The resultant peak "a develops until the stress path reaches the equilibrium at OCR = 2. For ad 6= 0, the rate of damage parameter d evolution is larger with larger ad, which invokes faster reduction of plastic modulus Kp and fewer cycles are required to reach a certain strain level in Fig. 4.18(a). Observe in Fig. 4.18(b) that unlike simulations with dierent h0 values, variation of ad does not involve considerable dierence in accumulation of residual u, because eect of ad is more pronounced with increasing number of cycles. The proposed SANICLAY-D-B model simulation with no eect from damage parameter d cannot replicate experimental results shown in Fig. 4.16, while the simulations with evolving damage parameter d can reproduce accumulation of peak "a. However, the proposed model predicts about twice smaller residual u than that observed in the experiments, because the stress path migration is not allowed once the stress path reaches OCR = 2. The results of qualitative model performance illustrated in this section demonstrate the 83 4.6. Generation of peak axial strain and residual pore water pressure versus number of cycles ability of the model to simulate cyclic response of clays by means of new model parameters h0 and ad. Parameter h0 is necessary to simulate undrained cyclic response within initial few cycles and thus accumulation of residual u, while parameter ad enables one to simulate reduction of stiness with increasing number of cycles observed in experiments. 100 101 102 −10 −8 −6 −4 −2 0 2 4 6 8 10 Pe ak   ε a  (% ) h0=25 50 200100 400 ad=40 SANICLAY q cy/2Su=0.5 (a) 100 101 102 0 0.2 0.4 0.6 0.8 1 Number of cycles R es id ua l u/ σ c ad=40 qcy/2Su=0.5 200 400 h0=25 100 50 SANICLAY (b) Figure 4.17: Generation of: (a) peak axial strain "a, and (b) residual pore water pressure u versus number of cycles with dierent values of h0 and xed ad = 40 84 4.6. Generation of peak axial strain and residual pore water pressure versus number of cycles 100 101 102 −10 −8 −6 −4 −2 0 2 4 6 8 10 Pe ak   ε a (% ) q cy/2Su=0.5 h0=50 ad=80 40 20 10 5 2 0 (a) 100 101 102 0 0.2 0.4 0.6 0.8 1 Number of cycles R es id ua l u/ σ e h0=50 qcy/2Su=0.5 ad=80 0 (b) Figure 4.18: Generation of: (a) peak axial strain "a, and (b) residual pore water pressure u versus number of cycles with xed h0 = 50 and varied values of ad 85 5. Quantitative Model Performance This Chapter illustrates the proposed SANICLAY-D-B model performance in comparison with experimental data in monotonic and cyclic loading and is intended to serve as a guidance for model parameters calibration. Experiments required for calibration of model parameters are discussed. Then the model is tested quantitatively against experimental data on three dierent clays. 5.1 Experiments required for model calibration Abbreviations of the tests and their description used in this section are given below: CIUC undrained triaxial compression test on isotropically consolidated clay CIUE undrained triaxial extension test on isotropically consolidated clay CK0UC undrained triaxial compression test on K0-consolidated clay CK0UE undrained triaxial extension test on K0-consolidated clay The proposed SANICLAY-D-B model requires 11 model constants, whose roles are sum- marized in Table 5.3. Constants h0 and ad are novel in regards to SANICLAY-D model proposed by Taiebat, Dafalias and Peek (2010). Initial internal variables required for the proposed model with their denition are shown in Table 5.2. The reader is referred to the works by Dafalias et al. (2006) and Taiebat, Dafalias and Peek (2010) for comprehensive description of anisotropy and destructuration related parameters calibration, respectively. A brief description of the calibration procedure for the required constants is given as follows:  Calibration of  and  is based on, preferably isotropic, consolidation tests. The use of data from the K0-consolidation test is also acceptable. 86 5.1. Experiments required for model calibration  Calibration of  is based on initial K0-unloading stress path preceded by K0-consolidation test. Typical range of  = 0:2{0:3 for clays can be assumed, if the required stress path is not available.  Calibration of Mc and Me is based on triaxial compression and extension tests, respec- tively.  Calibration of N is based on CK0UC test on normally consolidated clay. If not avail- able, N can be calibrated along with C from CIUC or CIUE test on normally consoli- dated clay.  Calibration of C is based on the stress path from CK0UE test. If not available, C can be calibrated along with N from CIUC or CIUE test on normally consolidated clay. C = 3{20 is the typical range for most clays (Dafalias et al., 2006).  Calibration of x is based on the closed-form relation given by Eq. 27 in Dafalias et al. (2006).  Calibration of ki for structured clays is based on consolidation test to stresses, ensuring complete loss of structure. ki = 0 should be used for unstructured clays.  Calibration of h0 is based on stress and stress{strain paths from initial cycles of undrained cyclic loading, preferably on normally consolidated clay. If monotonic load- ing is of interest, the calibration can be performed in comparison with undrained monotonic loading test. At this stage ad = 0 should be used. Isotropic consolidation prior to cyclic loading guarantees that damage parameter d does not evolve. However, for K0-consolidated clay, calibrated h0 should be adjusted after ad is known to render h that is equal to the rst calibrated h0 prior to cyclic loading.  Calibration of ad is based on similar test data used for h0 calibration, but with much larger number of cycles. ad = 0 could be used for monotonic loading, if cyclic loading is not of interest. Below, model calibration and performance are shown for reconstituted Georgia kaolin clay (experimental data after Sheu (1984)), structured Cloverdale clay (experimental data 87 5.2. Reconstituted Georgia kaolin clay Table 5.1: Model constants for the proposed SANICLAY-D-B model Model constant Designation Description of its role Elasticity  Compressibility of overconsolidated clay  Poisson's ratio Critical state Mc Critical stress ratio in triaxial compression Me Critical stress ratio in triaxial extension N Shape of the bounding surface  Compressibility of normally consolidated clay Rotational hardening C Rate of evolution of anisotropy x Saturation limit of anisotropy, constant  = q=p path Destructuration ki Rate of destructuration Bounding surface h0 Initial hardening parameter ad Rate of damage evolution Table 5.2: Initial internal variables for the proposed SANICLAY-D-B model Designation Description of its role p0 Size of the bounding surface Si Structuration factor ein Initial void ratio  Orientation of plastic potential and bounding surface (p; q) Initial stress state in the triaxial plane after Zergoun (1991); Zergoun and Vaid (1994)) in slow cyclic loading, and reconstituted Ariake clay (experimental data after Yasuhara et al. (1992)) in fast cyclic loading. 5.2 Reconstituted Georgia kaolin clay Experimental data after Sheu (1984) on Georgia kaolin clay is used for model calibration. The soil tested was composed of 62% of clay and 38% of silt. The plasticity index for this clay was 20%. Kaolin clay was mixed with salty water for preparation of Georgia kaolin clay in order to reduce the eect of anisotropy through formation of a 
occulated structure. The mixture was poured into the consolidometer with two-way drainage and incremental vertical pressures were applied with the nal pressure of 205 kPa. Then block sample was trimmed into triaxial specimens. 88 5.2. Reconstituted Georgia kaolin clay 5.2.1 Model calibration Model calibration is based on available experimental results. Parameters for destructuration are set to the values corresponding to already destructured clay as the clay considered was reconstituted. Details of the calibration of model parameters for this clay are presented below. Parameters  and . In Fig. 5.1 the response of Georgia kaolin clay in isotropic con- solidation with three unloading-reloading sequences is presented in the e{ log p plane. The slopes of the normal compression and the rebound lines are Cc and Cr, respectively, in the e{ log p plane, whereas the required  and  are measured in the e{ ln p plane. The relation- ships between the slopes are as follows:  = Cc= ln 10 and  = Cr= ln 10. Then the values Cc = 0:28 and Cr = 0:085 measured in Fig. 5.1 render  = 0:121 and  = 0:037. 101 102 103 0.5 0.6 0.7 0.8 0.9 1 p (kPa) e   experiment C c =0.28 C r =0.085 Figure 5.1: Calibration of Cc and Cr in e{ log p plane based on isotropic consolidation test after Sheu (1984) Parameters Mc and Me. CIUC and CIUE stress paths on initially normally consolidated clay to a range of preconsolidation stresses c = 276, 345,414, and 552 kPa are shown in Fig. 5.2. From direct measurements at the ultimate stress states, critical state ratios 89 5.2. Reconstituted Georgia kaolin clay are estimated as Mc = 0:87 and Me = 0:86. For later calibration of x, friction angle in compression c might be required to come up with the K0 value; in this case, Mc is related to c by Mc = 6 sinc 3 sinc (5.1) 0 200 400 600 800 −400 −200 0 200 400 p (kPa) q (kP a)   M c =0.87 M e =0.86 414 345 276 Compression  σ c  (kPa): 552 414 345 276 Extension  σ c  (kPa): Figure 5.2: Calibration of Mc and Me in the triaxial plane based on undrained triaxial compression and extension tests starting from normally consolidated state, experimental results after Sheu (1984) Parameter N . Due to unavailability of CK0UC test, the value for N must be assumed. Looking ahead, N = Mc or N =Me would require parameter C to be very small relatively to its typical values for clay in the range of C = 3{20 mentioned previously. In addition, cyclic loading simulations would suggest that the value should be considerably smaller for better simulations, especially for high cyclic stress amplitudes. For these reasons, the decision is made to adopt N = 0:8. 90 5.2. Reconstituted Georgia kaolin clay Parameter . Due to unavailability of K0-unloading preceded by K0-consolidation test for Georgia kaolin clay, parameter  = 0:2 is assumed. This value has been calibrated for other clays by Taiebat, Dafalias and Peek (2010), and Dafalias et al. (2006). Parameter x. The parameter x for successful K0-loading simulation can be obtained using the following equation proposed by Dafalias et al. (2006): x = 2K0(1 =) B3K0 +  3 K0 + [2 (1 =)BM2c ] K0 M2c ; B =  2 (1 + ) 9 (1 2) (5.2) where  = 3=2 for K0-loading path; K0 is a measured value of the earth pressure coecient at rest; and K0 = 3 (1K0) = (1 + 2K0). Should K0 is not available, Jacky's equation, K0 = 1  sin, can be used. For Georgia kaolin clay,  = c = 22:3  from Eq. 5.1 yields K0 = 0:62, K0 = 0:51 and x = 1:69. Parameter C. Calibration of C requires all other constants calibrated in advance (Dafalias et al., 2006). In the proposed model novel bounding surface parameters h0 and ad are not required to be known for calibration of C. This is true since parameters h0 and ad are in eect for overconsolidated states, while C calibration is based on the stress path corresponding to normally consolidated state. From the available data, CIUC test is the most appropriate for C calibration, ensuring rotation of the bounding surface. The C = 3 is chosen as an optimum value from trial-and-error simulations shown in Fig. 5.3. Bounding surface parameters h0 and ad. Calibration of h0 and ad requires experimental data involving the stress path corresponding to overconsolidated states. As cyclic loading is of concern, the parameters should be calibrated based on cyclic loading test. To be consistent, calibration of the two parameters should be performed from the same test results. Parameter h0 is responsible for the initial part of the cyclic loading, when accumulation of damage parameter d is still not signicant to induce considerable reduction of h. Therefore, initial cycles from available undrained cyclic triaxial test are used for h0 calibration. h0 = 50 is calibrated from trial-and-error simulations with ad = 0 shown in Fig. 5.4. Once h0 is known, ad calibration is performed based on the same undrained cyclic triaxial test used for h0 calibration, but considering all available 15 cycles. For Georgia kaolin clay, ad = 7 is calibrated in a trial-and-error procedure. 91 5.2. Reconstituted Georgia kaolin clay 0 100 200 300 0 100 200 300 p (kPa) q (kP a)   experiment C=9 C=6 C=3 Figure 5.3: Calibration of C based on undrained triaxial compression test on normally consolidated clay, experimental results after Sheu (1984) Calibrated constants and initial internal variables for Georgia kaolin clay are shown in Tables 5.3 and 5.4, respectively. In Table 5.4, ein = 1:5 and p = 1 kPa are used for simulations of undrained cyclic triaxial loading. For simulations of consolidation, CIUC and CIUE tests, dierent values of ein and p are used, which are to be specied for the corresponding tests simulation. Table 5.3: Model parameters for Georgia kaolin Clay Model constant Designation Value Elasticity  0:037  0:2 Critical state Mc 0:87 Me 0:86 N 0:8  0:121 Rotational hardening C 3 x 1:69 Destructuration (no eect) ki 0 Bounding surface h0 50 ad 7 92 5.2. Reconstituted Georgia kaolin clay 0 100 200 300 400 500 −200 −100 0 100 200 p (kPa) q (kP a)   q cy=140.7 kPa experiment simulation (a) −3 −2 −1 0 1 2 3 −200 −100 0 100 200 ε a  (%) q (kP a)   h0=50 experiment simulation (b) Figure 5.4: Initial hardening parameter calibrated to h0 = 50 based on undrained cyclic triaxial loading, experimental results after Sheu (1984) Table 5.4: Initial internal variables for Georgia kaolin clay Model internal variable Designation Value Size of the bounding surface p0 1 kPa Structuration factor Si 1 Initial void ratio ein 1:5 Orientation of plastic potential and bounding surface  0 Initial stress state (p; q) (1; 0) kPa 5.2.2 Model performance Model constants calibrated previously (Table 5.3) are used in this section to simulate isotropic consolidation, CIUC, CIUE, and undrained cyclic triaxial tests on normally consolidated Georgia kaolin clay. Isotropic consolidation test The results of isotropic consolidation test simulation in comparison with experimental data are shown in Fig. 5.5. p0 = 75 kPa and ein = 0:88 are adopted for the simulation. Notice that plastic volumetric strains are involved upon reloading, even before the past maximum preconsolidation stress is reached. The magnitude of the plastic volumetric strains in reload- ing is in agreement with experimental results, suggesting that h0 = 50 calibrated based on 93 5.2. Reconstituted Georgia kaolin clay cyclic loading can be successfully used for monotonic loading simulations. 101 102 103 0.5 0.6 0.7 0.8 0.9 1 p (kPa) e   experiment simulation Figure 5.5: Simulation of isotropic consolidation test, experimental results after Sheu (1984) Undrained triaxial compression and extension tests Simulated stress and stress{strain paths of CIUC and CIUE loading on clay normally con- solidated to a range of preconsolidation stresses c = 276, 345, 414 and 552 kPa are shown in Fig. 5.6 against available experimental data. For the stress path in compression, simu- lations closely match experimental results, while in extension the \hook-type" stress path cannot be reproduced by the model because of the evolution of anisotropic variable  being a function of tr (@F=@) according to Eq. 3.36b, i.e., the stress path freezes as the critical stress ratio is reached (Dafalias and Taiebat, 2012; Taiebat and Dafalias, 2013). In terms of the stress{strain response, simulated path ultimately matches experimental data, although experimental results suggest larger plastic axial strains at earlier stages. Note that there is no eect from parameters h0 and ad for CIUC and CIUE on normally consolidated clay. 94 5.2. Reconstituted Georgia kaolin clay 0 200 400 600 800 −400 −200 0 200 400 p (kPa) q (kP a)   414 345 276 Compression  σ c  (kPa): 552 414 345 276 Extension  σ c  (kPa): (a) −20 −10 0 10 20 −400 −200 0 200 400 ε a  (%) q (kP a)   414 345 276 Compression  σ c  (kPa): 552 414 345 276 Extension  σ c  (kPa): (b) Figure 5.6: Simulation of CIUC and CIUE tests on normally consolidated clay, experimental results after Sheu (1984) 95 5.2. Reconstituted Georgia kaolin clay Undrained cyclic triaxial loading The main focus of the current study is on simulations of cyclic loading. For convenience, stress and stress{strain paths are shown separately in Figs. 5.7 and 5.8, respectively, with the digitized experimental data on the left and simulations on the right. The results of 22 cycles for qcy = 121:4 kPa, 15 cycles for qcy = 136 and 140:7 kPa, and 3 cycles for qcy = 165:5 kPa are available from the experiments. Simulations are not exhibiting the same amount of p reduction as suggested in experiments owing to the balance reached by the simulated stress path at OCR = 2 (see section 3.3). In the experiments the stress path falls well beyond the critical stress ratio. A similar response was observed by Takahashi et al. (1980). They found that if the drainage was allowed after the undrained loading, then the stress ratio nally reduced to the critical stress ratio. This suggests that unreliable pore water pressure measurements were the reason of the stress ratio exceeding the critical stress ratio. For the stress{strain path, the model is able to reproduce cyclic softening observed in the experiments for the whole range of cyclic stress amplitudes. It is no surprise that simulations for qcy = 140:7 kPa closely match experiments as h0 and ad calibration is based on this particular cyclic stress amplitude. For a slightly lower stress amplitude, qcy = 136 kPa, the simulated stress{strain path is also matching well the experiment. For stress amplitudes much dierent from 140:7 kPa, simulations slightly deviate from experiments, although the results can be considered satisfactory for simulations with the single set of model constants. Observe that the simulated strain levels are less than experimentally observed ones for 165:5 kPa cyclic stress amplitude and larger for qcy = 121:4 kPa. 96 5.2. Reconstituted Georgia kaolin clay 0 100 200 300 400 500 −200 −100 0 100 200 p (kPa) q (kP a)   q cy=121.4 kPa experiment (1−6,22 cycles) (a) 0 100 200 300 400 500 −200 −100 0 100 200 p (kPa) q (kP a)   q cy=121.4 kPa simulation (22 cycles) (b) 0 100 200 300 400 500 −200 −100 0 100 200 p (kPa) q (kP a)   q cy=136 kPa experiment (1−6,15 cycles) (c) 0 100 200 300 400 500 −200 −100 0 100 200 p (kPa) q (kP a)   q cy=136 kPa simulation (15 cycles) (d) 0 100 200 300 400 500 −200 −100 0 100 200 p (kPa) q (kP a)   q cy=140.7 kPa experiment (1−6,15 cycles) (e) 0 100 200 300 400 500 −200 −100 0 100 200 p (kPa) q (kP a)   q cy=140.7 kPa simulation (15 cycles) (f) 0 100 200 300 400 500 −200 −100 0 100 200 p (kPa) q (kP a)   q cy=165.5 kPa experiment (3 cycles) (g) 0 100 200 300 400 500 −200 −100 0 100 200 p (kPa) q (kP a)   q cy=165.5 kPa experiment (3 cycles) (h) Figure 5.7: Undrained cyclic triaxial loading stress paths with the digitized experimental results after Sheu (1984) (left) and model simulations (right) for cyclic stress amplitudes: (a{b) qcy = 121:4 kPa, (c{d) qcy = 136 kPa, (e{f) qcy = 140:7 kPa, and (g{h) qcy = 165:5 kPa 97 5.2. Reconstituted Georgia kaolin clay −15 −10 −5 0 5 10 −200 −100 0 100 200 ε a  (%) q (kP a)   q cy=121.4 kPaexperiment (1−6,19−22 cycles) (a) −15 −10 −5 0 5 10 −200 −100 0 100 200 ε a  (%) q (kP a)   q cy=121.4 kPasimulation (22 cycles) (b) −15 −10 −5 0 5 10 −200 −100 0 100 200 ε a  (%) q (kP a)   experiment (15 cycles) (c) −15 −10 −5 0 5 10 −200 −100 0 100 200 ε a  (%) q (kP a)   q cy=136 kPasimulation (15 cycles) (d) −15 −10 −5 0 5 10 −200 −100 0 100 200 ε a  (%) q (kP a)   experiment (1−7,13−15 cycles) (e) −15 −10 −5 0 5 10 −200 −100 0 100 200 ε a  (%) q (kP a)   q cy=140.7 kPasimulation (15 cycles) (f) −15 −10 −5 0 5 10 −200 −100 0 100 200 ε a  (%) q (kP a)   experiment (3 cycles) (g) −15 −10 −5 0 5 10 −200 −100 0 100 200 ε a  (%) q (kP a)   q cy=165.5 kPasimulation (3 cycles) (h) Figure 5.8: Undrained cyclic triaxial loading stress{strain paths with the digitized experi- mental results after Sheu (1984) (left) and model simulations (right) for cyclic stress am- plitudes: (a{b) qcy = 121:4 kPa, (c{d) qcy = 136 kPa, (e{f) qcy = 140:7 kPa, and (g{h) qcy = 165:5 kPa 98 5.3. Structured Cloverdale clay 5.3 Structured Cloverdale clay Experimental data after Zergoun (1991), Zergoun and Vaid (1994) on Cloverdale clay are used for model calibration. Soil tested was sensitive undisturbed marine clay consisting of 49% of clay and 45% of silt with PI = 24%. The sensitivity was a result of surface inltration following uplift above the sea level. The clay was block sampled from depth of (2:0{2:6) m below ground level. Then the triaxial specimens were trimmed from the same horizon and hydrostatically loaded under an eective stress of 200 kPa for 48 h. The consolidation stress beyond apparent precon- solidation stress destructured clay to some extend inducing reduced sensitivity. Prior to experiments, the specimens were left undrained under constant total stress for 24 h to allow pore water pressure generation due to arrest of secondary consolidation. 5.3.1 Model calibration Calibration of the model parameters given below is based on available experimental data. The clay considered requires calibration of the parameters ki and Si. Parameters  and . The response of Cloverdale clay in an oedometer test during one- dimensional compression and swelling is shown in Fig. 5.9. Lacking lateral stress mea- surements in this test, the data are plotted in the e{ log a plane, where a is an axial (vertical) stress. Cc and Cr measured in Fig. 5.9 are related to  and  as before, giving:  = Cc= ln 10 = 0:21 and  = Cr= ln 10 = 0:03. Note that while this indirect estimation of  is accurate, determination of  is only approximated, given the relative increase of the p=a ratio in swelling due to the increase of K0 for the overconsolidated state (Wood, 1990), i.e., rate of p and a change is no longer constant as in the loading process. Therefore,  calibration based on isotropic consolidation test, if available, is preferred (Papadimitriou et al., 2005). Parameters Mc and Me. The response of Cloverdale clay in CIUC tests on normally consolidated and overconsolidated specimens is shown in Fig. 5.10. The data is presented in the MIT stress plane, normalized by the value of the equivalent consolidation stress e, which is an equivalent stress on the normal consolidation line in the e{ log c plane corresponding to 99 5.3. Structured Cloverdale clay 100 101 102 103 0.6 0.8 1 1.2 1.4 1.6 σ a  (kPa) e   experiment C c =0.48 C r =0.07 Figure 5.9: Calibration of Cc and Cr in e{ log p plane based on oedometer test after Zergoun (1991) the same void ratio as the overconsolidated clay specimen (Zergoun, 1991). Given these tests, the friction angle in compression is estimated as c = 32 , rendering Mc = 1:29 according to Eq. 5.1. Me = 1:27 is chosen for better undrained cyclic triaxial loading simulations. Parameter N . Due to unavailability of a CK0UC test the value of N = 1 is assumed. Similar to the N preference for Georgia kaolin clay in section 5.2.1, this value is chosen owing to its appropriateness for calibration of the constant C and better undrained cyclic triaxial loading simulations. Parameter . Due to unavailability of K0-unloading stress path following K0-consolidation for Cloverdale clay the parameter  = 0:2 is assumed. Parameter x. Utilizing Eq. 5.2 along with Jacky's equation,  = 32 yields K0 = 0:47, K0 = 0:82 and x = 1:73. Parameters Si and ki. The best experiment for calibration of the initial internal variable Si and the value of destructuration parameter ki is an isotropic compression test on an intact sample until complete loss of structure occurs (Taiebat, Dafalias and Peek, 2010). The only 100 5.3. Structured Cloverdale clay 0 0.2 0.4 0.6 0.8 1 1.2 0 0.2 0.4 0.6 (σ a  + σ r ) /2 σ e (σ a  −  σ r) / 2 σ e   φ c =32° 1.0 1.3 2.0 3.8 Compression, OCR: Figure 5.10: Calibration of c in the triaxial plane based on undrained triaxial compression test on normally consolidated and overconsolidated specimens, experimental results after Zergoun (1991) available consolidation test is oedometric consolidation, which is used for the calibration of the parameters. Parameter Si can be estimated from the compression plane as the ratio of the axial stress (or mean stress, in the case of isotropic loading) at the end of the initial sti regime and the corresponding stress on the Intrinsic Compression Curve (ICC) at the same void ratio (Taiebat, Dafalias and Peek, 2010). The location of the ICC can be approximated by the dotted line in Fig. 5.11(a) and the ratio yields Si = 90=24 = 3:7 for Cloverdale clay. The parameter ki is calibrated to ki = 2 from trial-and-error simulations shown in Fig. 5.11(b). The initial value of p0 = 94 kPa has been used to successfully simulate the initial sti (elastic) part of the response. New model constants that still require calibration are set to h0 = 1e10 and ad = 0 to simulate sti response for overconsolidated state, since the calibration of ki involves normally consolidated state. Parameter C. CIUC is the most appropriate test from the available data for calibration of C, involving rotation of the bounding surface. From trial-and-error calibration procedure C = 3 is selected as an optimum value. Simulation results with dierent C values for the CIUC test on normally consolidated clay to p = 200 kPa are shown in Fig. 5.12. Bounding surface parameters h0 and ad. Model parameters h0 = 550 and ad = 68 are 101 5.3. Structured Cloverdale clay 100 101 102 103 0.6 0.8 1 1.2 1.4 1.6 σ a  (kPa) e   experiment ICC (approx.) 9024 Si=90/24=3.7 (a) 100 101 102 103 0.6 0.8 1 1.2 1.4 1.6 σ a  (kPa) e   experiment ki=2.0 ki=2.4 ki=1.6 (b) Figure 5.11: Calibration of: (a) Si, and (b) ki based on oedometer test after Zergoun (1991) 102 5.3. Structured Cloverdale clay 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 (σ a  + σ r ) /2 σ e (σ a  −  σ r) / 2 σ e   experiment C=3 C=9 C=6 Figure 5.12: Calibration of C based on CIUC test on initially normally consolidated clay, experimental data after Zergoun (1991) calibrated based on the available experimental cyclic loading with cy=Suc = 0:75 cyclic stress amplitude. The 11 calibrated model constants are shown in Table 5.5, and initial internal variables in Table 5.6. These parameters are to be used in simulations of the available experiments. 5.3.2 Model performance Model performance for Cloverdale clay is illustrated in this section. Simulations are con- ducted for one-dimensional consolidation, CIUC on normally consolidated and overconsoli- dated clay, and undrained cyclic triaxial tests on normally consolidated specimens. Oedometer test In Fig. 5.13 experimental results and model simulations are shown for an oedometer test on Cloverdale clay. Destructuration observed in laboratory testing can be well simulated by the model. In the reloading following initial unloading, plastic volumetric strains for the stress states before the preconsolidation stress is attained are predicted by the model in agreement 103 5.3. Structured Cloverdale clay Table 5.5: Model parameters for Cloverdale clay Model constant Designation Value Elasticity  0:03  0:2 Critical state Mc 1:29 Me 1:27 N 1  0:21 Rotational hardening C 3 x 1:73 Destructuration ki 2 Bounding surface h0 550 ad 68 Table 5.6: Initial internal variables for Cloverdale clay Model internal variable Designation Value Size of the bounding surface p0 94 kPa Structuration factor Si 3:7 Initial void ratio ein 1:48 Orientation of plastic potential and bounding surface  0 Initial stress state (p; q) (1; 0) kPa with the experiment. For unloading, simulated response closely matches experiment for a  20 kPa; however, the dialative nature of clay in K0-unloading cannot be predicted for a < 20 kPa, resulting in slightly underpredicted void ratio at the end of unloading. One might anticipate a bounding surface providing the desired dialation as the image point upon K0-unloading should be located so that the negative plastic volumetric strains are predicted. It is indeed so and the proposed model better simulates dialative response than the previous SANICLAY family models, but the degree of plastic strains with calibrated h0 and ad does not allow reaching dialation observed in the experiment. During parameters calibration, it was noticed that the parameter C is the factor that can control the degree of dialation through rotational hardening. For lower values of C, less rotation of the bound- ing surface renders more dialative response simulated for K0-unloading. Thus, adjusted rotational hardening in conjunction with the bounding surface formulation can potentially enhance model predictions in K0-loading. 104 5.3. Structured Cloverdale clay 100 101 102 103 0.6 0.8 1 1.2 1.4 1.6 σ a  (kPa) e   experiment simulation Figure 5.13: Simulation of oedometer test on Cloverdale clay, experimental results after Zergoun (1991) Undrained triaxial compression and extension tests Isotropic compression to p = 200 kPa is imposed prior to simulations of undrained triaxial compression and extension tests. Simulations of CIUC and CIUE tests on normally consol- idated and overconsolidated specimens are performed with the results shown in Fig. 5.14. Simulated response approximately follows the experimental results in terms of the stress path for OCR  2. As OCR increases the stress path deviates from experiment in that the dialative response observed in actual tests is not captured. Similarly, the stress{strain path is closely matching for OCR  2 throughout the whole process of undrained loading. For OCR = 3:8, initial stress{strain simulation exhibits larger plastic strains than desired and then falls close to the experiment. Undrained cyclic triaxial loading Experimental data, where slow undrained cyclic triaxial loading was performed by Zergoun (1991) to maintain the uniformity of pore water pressure, is used in this section. Prior to 105 5.3. Structured Cloverdale clay 0 0.2 0.4 0.6 0.8 1 1.2 0 0.2 0.4 0.6 (σ a  + σ r ) /2 σ e (σ a  −  σ r) / 2 σ e   1.0 1.3 2.0 3.8 Compression, OCR: (a) 0 2 4 6 8 10 0 0.2 0.4 0.6 ε a  (%) (σ a  −  σ r) / 2 σ e   1.0 1.3 2.0 3.8 Compression, OCR: (b) Figure 5.14: Simulations of CIUC and CIUE tests on normally and overconsolidated clay, experimental results after Zergoun (1991) 106 5.3. Structured Cloverdale clay simulations of undrained cyclic triaxial loading, isotropic consolidation to p = 200 kPa is simulated. The stress and stress{strain path in undrained cyclic triaxial loading are shown separately in Figs. 5.15 and 5.16, respectively, with the experimental results on the left and model simulations on the right. Unlike undrained cyclic triaxial loading tests on Georgia kaolin clay, the slow cyclic loading imposed on Cloverdale clay renders the stress path that does not travel beyond the critical stress ratio, as shown in Fig.5.15. The simulated stress path that freezes at OCR = 2 (see section 3.3) still deviates from the observed experiments, but to a smaller degree than the simulated stress path for Georgia kaolin clay. In terms of the stress{strain simulations, cy=Suc = 0:75 cyclic stress amplitude that is used for h0 and ad calibration is apparently well simulated. For cy=Suc = 0:70, the stress{strain response closely matches the experiment in terms of attained peak strains. Simulation of cy=Suc = 0:79 cyclic loading provides strains that are much smaller than the ones observed during the experiments. Finally, simulations in terms of peak "a versus number of cycles for the range of available cyclic stress amplitudes in undrained cyclic triaxial loading are shown in Fig. 5.17(b). Ex- perimental data is digitized and is shown in Fig. 5.17(a) for comparison . Considering that earthquake loading is of interest, a M7.5 earthquake corresponds to approximately 30 equiv- alent number of cycles imposed on ne-grained material according to Idriss and Boulanger (2008). Therefore, the results of simulations and experimental data are plotted within 30 cycles. Evolution of peak "a with number of cycles is closely simulated for cy=Suc = 0:75 and 0:69. Simulated cy=Suc = 0:79 cyclic loading requires about 5 extra cycles to reach the experimentally observed peak "a = 10%. For lower cyclic stress amplitudes, the picture becomes less satisfactory, simulating a softer clay response than desired. Cyclic stress am- plitude of cy=Suc = 0:62 reaches approximately "a = 10% in 30 cycles as simulated, while in experiments only "a = 2% is attained in 30 cycles. With a cy=Suc = 0:57, peak "a = 2% is simulated in 30 cycles, while less than 1% is reached in the experiments. A decreased rate of peak "a with a number of cycles is observed in the tests for cy=Suc = 0:54 approaching a plateau, whereas strains keep accumulating in simulation. Such plateau, however, is attained for lower stress amplitudes, although at a strain levels above that observed in the experi- ments. It should be noted that the data is plotted in logarithmic scale for the vertical axis 107 5.3. Structured Cloverdale clay 0 0.2 0.4 0.6 0.8 1 1.2 −0.4 −0.2 0 0.2 0.4 (σ a  + σ r ) /2 σ e (σ a  −  σ r) / 2 σ e   τ cy/Suc=0.70 experiment (5 cycles) (a) 0 0.2 0.4 0.6 0.8 1 1.2 −0.4 −0.2 0 0.2 0.4 (σ a  + σ r ) /2 σ e (σ a  −  σ r) / 2 σ e   τ cy/Suc=0.70 simulation (5 cycles) (b) 0 0.2 0.4 0.6 0.8 1 1.2 −0.4 −0.2 0 0.2 0.4 (σ a  + σ r ) /2 σ e (σ a  −  σ r) / 2 σ e   τ cy/Suc=0.75 experiment (13 cycles) (c) 0 0.2 0.4 0.6 0.8 1 1.2 −0.4 −0.2 0 0.2 0.4 (σ a  + σ r ) /2 σ e (σ a  −  σ r) / 2 σ e   τ cy/Suc=0.75 simulation (13 cycles) (d) 0 0.2 0.4 0.6 0.8 1 1.2 −0.4 −0.2 0 0.2 0.4 (σ a  + σ r ) /2 σ e (σ a  −  σ r) / 2 σ e   τ cy/Suc=0.79 experiment (4 cycles) (e) 0 0.2 0.4 0.6 0.8 1 1.2 −0.4 −0.2 0 0.2 0.4 (σ a  + σ r ) /2 σ e (σ a  −  σ r) / 2 σ e   τ cy/Suc=0.79 simulation (4 cycles) (f) Figure 5.15: Undrained cyclic triaxial loading stress path with the digitized experimental results after Zergoun (1991); Zergoun and Vaid (1994) (left), and model simulations (right) for cyclic stress amplitudes: (a-b) cy=Suc = 0:70, (c-d) cy=Suc = 0:75, and (e-f) cy=Suc = 0:79 108 5.3. Structured Cloverdale clay −3 −2 −1 0 1 2 3−0.4 −0.2 0 0.2 0.4 ε a  (%) (σ a  +  σ r) / 2 σ e   experiment (5 cycles) (a) −3 −2 −1 0 1 2 3 −0.4 −0.2 0 0.2 0.4 ε a  (%) (σ a  −  σ r) / 2 σ e   τ cy/Suc=0.70simulation (5 cycles) (b) −10 −8 −6 −4 −2 0 2 4 6 8 10 −0.4 −0.2 0 0.2 0.4 ε a  (%) (σ a  +  σ r) / 2 σ e   experiment (13 cycles) (c) −10 −8 −6 −4 −2 0 2 4 6 8 10 −0.4 −0.2 0 0.2 0.4 ε a  (%) (σ a  −  σ r) / 2 σ e   τ cy/Suc=0.75simulation (13 cycles) (d) −10 −8 −6 −4 −2 0 2 4 6 8 10 −0.4 −0.2 0 0.2 0.4 ε a  (%) (σ a  +  σ r) / 2 σ e   experiment (4 cycles) (e) −10 −8 −6 −4 −2 0 2 4 6 8 10 −0.4 −0.2 0 0.2 0.4 ε a  (%) (σ a  −  σ r) / 2 σ e   τ cy/Suc=0.79simulation (4 cycles) (f) Figure 5.16: Undrained cyclic triaxial loading stress path with the digitized experimental results after Zergoun (1991); Zergoun and Vaid (1994) (left), and model simulations (right) for cyclic stress amplitudes: (a-b) cy=Suc = 0:70, (c-d) cy=Suc = 0:75, and (e-f) cy=Suc = 0:79 109 5.3. Structured Cloverdale clay and therefore deviation from experimental data is somewhat exaggerated for lower cyclic stress amplitudes. Observe that at lower cyclic stress amplitudes the simulated curves (Fig.5.17(b)) lie above the ones experimentally obtained (Fig.5.17(a)) from initiation of cyclic loading. This suggests that either too large plastic strains are simulated at low strain levels or elasticity model constnats/relation used is/are not appropriate at low level of strains. Simulations results, shown in Fig. 5.18, obtained with the original SANICLAY model reveal that approximately the same peak "a are attained during initial cycles for the lowest cyclic stress amplitude. This eliminates the possibility of too large plastic strains induced by the model, implying that either adopted elasticity constants  and  or elasticity relation in the proposed model might be not appropriate. Changing parameter  would render inappropriate simulation of the consolidation process and would aect all other simulations, thus requiring re-calibration of some parameters. Parameter  can potentially enhance the simulated elastic response, but this would require its reduction to values far below the typical range for clays. Simulations with reduced  = 0:05 are shown in Fig. 5.17(c). Note that the overall simulated response does not change noticeably, except for the lowest cy=Suc, where elastic strains are governing strain during initial cycles. Thus, only simulations of very low cy=Suc can be controlled to some extent by elastic parameters. Undrained cyclic triaxial loading tests are also simulated with the original SANICLAY model for the whole range of available cyclic stress amplitudes with the results shown in Fig. 5.18. Model parameters and internal variables are the same as those used for the proposed SANICLAY-D-B model simulations. Parameters h0 and ad are not required for the original SANICLAY model. Clearly, no accumulation of peak "a can be simulated as the stress path falls inside the conventional yield surface, rendering purely elastic strains after the rst half of cycle for all cy=Suc. From the above simulated response of structured clay, the proposed SANICLAY-D-B model is considered to give adequate results within 30 cycles of applied cyclic loading, which is equivalent to a M7.5 earthquake (Idriss and Boulanger, 2008). Larger number of cycles might be considered with the reduced degree of satisfaction of the proposed model performance for low cy=Suc. 110 5.3. Structured Cloverdale clay 0 5 10 15 20 25 30 10−2 10−1 100 101 Number of cycles Pe ak  ε a  (% )   compression extension 0.54 0.570.620.69 0.79 0.75 0.47 τ cy/Suc=0.270.36 Experiment (a) 0 5 10 15 20 25 30 10−2 10−1 100 101  Number of cycles Pe ak  ε a  (% )   compession extension τ cy/Suc=0.270.360.470.540.570.620.69 0.750.79 SANICLAY−D−B (b) 0 5 10 15 20 25 30 10−2 10−1 100 101 Number of cycles Pe ak   ε  a  (% )   compression extension SANICLAY−D−B (ν=0.05) 0.36 τ cy/Suc=0.270.570.620.69 0.54 0.750.79 0.47 (c) Figure 5.17: Simulations of peak "a versus number of cycles in undrained cyclic triaxial loading for the range of available cyclic stress amplitudes: (a) experimental results from Zergoun and Vaid (1994), and model simulations with (b) calibrated model parameters, (c) calibrated model parameters and  = 0:05 111 5.4. Reconstituted Ariake clay (fast cyclic loading) 0 5 10 15 20 25 30 10−2 10−1 100 101 Number of cycles Pe ak   ε a  (% )   compression extension SANICLAY τ cy/Suc=0.270.47 0.36 0.79 0.75 0.69 0.62 0.57 0.54 Figure 5.18: Simulations of peak "a versus number of cycles in undrained cyclic triaxial loading for the range of available cyclic stress amplitudes with the original SANICLAY model 5.4 Reconstituted Ariake clay (fast cyclic loading) In order to illustrate if the model is able to reproduce the clay response in fast cyclic loading, although with less reliable experimental data (see section 2.1.1 for details), model parame- ters are calibrated approximately based on available data by Yasuhara et al. (1992). Two reconstituted marine clays, namely two Ariake clays, were used with the following index properties: PI = 69, Gs = 2:65; and PI = 72, Gs = 2:58. A clay slurry was consolidated under 55 kPa of vertical stress in a large consolidation vessel. Triaxial specimens were cut from the block samples formed in the consolidation vessel. The specimens were isotropically consolidated to c = 200 kPa with the drainage allowed through lter paper surrounding the specimen in the triaxial cell. Then isotropically consolidated specimens were subjected to undrained cyclic triaxial loading with constant stress amplitude. 5.4.1 Model calibration Because few experimental results are available, most of the model parameters are either approximated from known empirical relations or assumed as typical values for clays. A brief description of an approximate calibration procedure is given below. 112 5.4. Reconstituted Ariake clay (fast cyclic loading) Parameters  and . Given average PI  71% and Gs  2:62, Eq. 4.5 yields  = 0:41. The  = 0:05 is chosen simply as 1=8. Parameters Mc and Me. Yasuhara et al. (1992) obtained Mc = 1:68, corresponding to  = c = 41  from Eq. 5.1, based on monotonic undrained triaxial compression test with a strain rate of 0:05%=min. Such values for critical stress ratio M and friction angle  are questionable and perhaps are erroneous owing to the high rate of applied monotonic loading resulting in non-equalized pore water pressure distribution within the tested specimen (see section 2.1.1). However, all undrained cyclic triaxial tests were performed with high frequency loading of 1 Hz. Therefore, it was decided to adopt the reported Mc = 1:68 and assumed value of Me = 1:65 for simulations of high frequency cyclic loading. Parameter N . For Ariake clay, N = 1:68 is assumed. Parameter x. Based on Eq. 5.2 and the adopted value of Mc, x = 1:76 is determined. Parameter . For Ariake clay,  = 0:2 is chosen as this value is representative for most clays. Parameter C. C = 15 is calibrated from available CIUC test. Bounding surface parameters h0 and ad. For calibration of parameters h0 and ad, plots of peak "a versus number of cycles generated for qcy=2c = 0:3 and qcy=2c = 0:25 in undrained cyclic triaxial loading tests are used. This is an example of the parameters calibration based on the history of peak "a accumulation with number of cycles rather than already shown for Georgia kaolin and Cloverdale clay calibration of the parameters based on stress and stress{ strain paths. In trial-and-error simulations the most appropriate parameters are h0 = 1600 and ad = 80 with the nal simulation shown in Fig. 5.19. Calibrated model parameters are shown in Tables 5.7 and 5.8. 5.4.2 Model performance With the calibration of parameters completed, it is of interest to see how well the model can replicate the experimental data in cyclic stress ratio CSR versus number of cycles required to reach double amplitude axial strain "DA = 5%. The comparison of the experimental data with the simulated results is shown in Fig. 5.20. CSR = qcy=2c is adopted in accordance with the experiments. For CSR = 0:3, model simulates "DA = 5% reached in 10 cycles, while 113 5.4. Reconstituted Ariake clay (fast cyclic loading) 0 5 10 15 20 25 30 35 40 −5 −4 −3 −2 −1 0 1 2 3 4 5 Number of cycles Pe ak  ε a  (% )   experiment simulation 0.3 q cy/2σc=0.225 Figure 5.19: Calibration of h0 and ad based on peak "a versus number of cycles in undrained cyclic triaxial loading, experimental results after Yasuhara et al. (1992) in experiments the same strain level is reached in 6 cycles. As CSR decreases, simulated results fall closer to experiments. The two tted curves for simulation and experimental data points intersect at CSR = 0:24. For further reduced CSR, approximately the same number of cycles is required to reach "DA = 5%. In general, reasonable results are obtained with the calibrated model parameters. It is pertinent to say that the use of CSR versus number of cycles for calibration purposes is not suggested as it considers only one particular strain level that a given CSR reaches in a number of cycles and neglects the history of strains from before and after that particular strain level. As has been pointed out previously, an equivalent number of uniform cyclic loading representing a M7.5 earthquake is about 30 (Idriss and Boulanger, 2008). Considering this and the model simulations discussed previously, calibrated parameters can be used with satisfactory prediction of clay cyclic behavior for an earthquake magnitude of about M7.5 and possibly larger with reduced degree of satisfaction. Note that the model calibration is performed for 1 Hz frequency of cyclic loading, and the same model parameters cannot be used for loadings with dierent frequencies as the rate-independent model is not able to dierentiate between the applied frequencies. 114 5.4. Reconstituted Ariake clay (fast cyclic loading) Table 5.7: Model parameters for Ariake clay Model constant Designation Value Elasticity  0:05  0:2 Critical state Mc 1:68 Me 1:65 N 1:68  0:41 Rotational hardening C 15 x 1:76 Destructuration (no eect) ki 0 Bounding surface h0 1600 ad 80 Table 5.8: Initial internal variables for Ariake clay Model internal variable Designation Value Size of the bounding surface p0 37 kPa Structuration factor Si 1 Initial void ratio ein 2 Orientation of plastic potential and bounding surface  0 Initial stress state (p; q) (1; 0) kPa 0 10 20 30 40 50 0 0.1 0.2 0.3 0.4 0.5 Number of cycles Cy cli c st re ss  ra tio , q cy /2  σ c   experiment simulation Fitted curves Figure 5.20: [Simulation of cyclic stress ratio CSR versus number of cycles required to reach "DA = 5%, experimental results after Yasuhara et al. (1992) and tted power function y = aCSRb 115 6. Summary, Conclusions, and Recommendations for Future Research 6.1 Summary and conclusions The bounding surface formulation is incorporated into the SANICLAY-D model (Taiebat, Dafalias and Peek, 2010) to simulate plastic strains within the conventional yield surface during undrained cyclic loading. In addition to the eect of the bounding surface, damage parameter d is proposed to further reduce hardening function h required for plastic modulus Kp reduction. An evolving projection center is adopted in a way that allows it to keep its relative position within the bounding surface established at the moment of stress reversal. The eect of newly introduced initial hardening parameter h0 and rate of damage evo- lution ad is illustrated in the sensitivity analysis. The advantage of the proposed projection center is demonstrated in the qualitative model performance. The model is compared with previous SANICLAY family models in undrained cyclic triaxial loading, demonstrating the ability of the model to simulate plastic strains for the stress states inside the bounding sur- face and the advantage of damage parameter d. The bounding surface eect in generalized loading clearly illustrates the ability of the model to simulate plastic strains for circular stress path loading in the deviatoric plane with the stress path entirely lying inside the bounding surface. For better understanding of a role of newly introduced model parameters h0 and ad from practical point of view, the eect of the parameters on the generation of the following plots is illustrated:  secant shear modulus Gsec and damping ratio variation with shear strain amplitude; 116 6.1. Summary and conclusions  cyclic stress ratio CSR versus number of cycles required to reach a certain strain level;  peak axial strain "a and residual pore water pressure u versus number of cycles Calibration and model performance in comparison with experimental data is conducted for reconstituted and structured clays. The calibration process of the model parameters is shown and the model performance in monotonic and cyclic loadings is demonstrated. Destructuration of structured clay in monotonic loading can be well simulated by the de- structuration theory proposed by Taiebat, Dafalias and Peek (2010). It was reported by Idriss and Boulanger (2008) that 30 cycles of uniform cyclic loading are representative of a M7.5 earthquake for cyclic loading tests on clays. Thus, the number of cycles considered in simulations of cyclic loading is approximately 30. For reconstituted Georgia kaolin clay, sim- ulations of undrained cyclic triaxial loading closely match experimentally observed strains for higher cyclic stress amplitudes, whereas for lower cyclic stress amplitudes the model slightly overpredicts the observed strains. Reasonable simulation results are obtained for structured Cloverdale clay in undrained cyclic triaxial loading with wide range of available cyclic stress amplitudes, although for lower cyclic stress amplitudes the strains are also overpredicted. Finally, model simulation of Ariake clay in fast undrained cyclic triaxial loading illustrates that the new model parameters calibrated based on sucient experimental data can be suc- cessfully used to reproduce experimentally obtained cyclic stress ratio CSR versus number of cycles required to reach a certain strain level. The proposed model seems to have an adequate level of complexity for capturing many important characteristic features of clay response in both monotonic and cyclic loading. In the literature reviewed here (as well as those not discussed in the current work), none of the rate-independent bounding surface models demonstrated the ability to predict well a wide range of cyclic stress amplitudes. Some models were missing features of damage eect proposed in this work, others were not able to adequately simulate the damping, and yet others used model parameters that were not feasible to be calibrated for complex loading such as earthquake (e.g., cyclic stress amplitude parameter). 117 6.2. Recommendations for further research 6.2 Recommendations for further research The proposed SANICLAY-D-B model provides adequately detailed framework for predictions of stress-strain response of clays under irregular cyclic loading that may be induced by earthquake. The model has been validated in a number of element tests. The validation process can be continued for various experiments on dierent clays. Creation and collection of a more comprehensive experimental database on laboratory element tests on clays would provide the basis for further research. The model should also be properly implemented into a nite element or nite dierence framework. This would allow a detailed study on various problems of seismic wave prop- agation through clay deposits. In particular, the tool may be used for analysis of various foundations, dams, tunnels, and slopes subjected to earthquake loading. The model has been designed with particular attention on maintaining the simplicity of the formulation and the model calibration. However, if necessary, one may be able to address some of the shortcomings of the proposed model as discussed below. This would of course require a higher level of sophistication in the formulation and subsequently in the calibration process of the model. Although the model has signicantly improved from its previous versions, it still may not suciently reproduce variations of shear modulus and damping ratio observed in experiments at very small strains. This shortcoming may be due to the very simple elastic relation adopted in the model and may be solved by adopting a more sophisticated nonlinear elasticity formulation that provides hysteretic stress{strain relationships in small strain levels (e.g., Jung and Yune, 2011). Also the present model may overpredict the plastic strains at small cyclic stress amplitudes with increasing number of cycles. The hardening function hmay be the key issue for addressing some of the shortcoming of the model in wide range of cyclic stress amplitude and a more extensive study on this function may provide a further enhanced formulation, if needed. Finally the time dependency and creep have not been addressed in the presented formulation and may be considered in the future steps of the work. 118 Bibliography Al-Tabaa, A. and Wood, D. (1989), Experimentally based bubble model for clay., in `Pro- ceedings of the Third International Symposium NUMOG III', Elsevier Applied Science Publishers Limited, Niagara Falls, Canada, pp. 91{99. Anandarajah, A. and Dafalias, Y. (1986), `Bounding surface plasticity. III: Application to anisotropic cohesive soils.', Journal of Engineering Mechanics 112(12), 1292{1318. Bacchus, D. (1969), Cyclic deformation of a clay, Ph.D. dissertation, The University of Auckland, Auckland, New Zealand. Banerjee, P. and Yousif, N. (1986), `A plasticity model for the mechanical behaviour of anisotropically consolidated clay', International Journal for Numerical and Analytical Methods in Geomechanics 10(5), 521{541. Batdford, S. and Budiansky, B. (1949), A mathematical theory of plasticity based on the concept of slip, Technical Report TN 1871, NACA. Clayton, C. (2011), `Stiness at small strain: Research and practice', Geotechnique 61(1), 5{ 37. Dafalias, Y. (1975), On Cyclic and Anisotropic Plasticity: I) A General Model including Ma- terial Behavior under Stress Reversals. II) Anisotropic Hardening for Initially Orthotropic Materials, Ph.D. dissertation, University of California, Berkeley, United States. Dafalias, Y. (1986a), `An anisotropic critical state soil plasticity model', Mechanics Research Communications 13(6), 341{347. Dafalias, Y. (1986b), `Bounding surface plasticity. I: Mathematical foundation and hypoplas- ticity.', Journal of Engineering Mechanics 112(9), 966{987. 119 Bibliography Dafalias, Y. (1987), An anisotropic critical state clay palsticity model, in C. Desai, E. Krempl, P. Kiousis and T. Kundu, eds, `Constitutive laws for engineering materials: theory and application', pp. 513{521. Dafalias, Y. and Herrmann, L. (1986), `Bounding surface plasticity. II: Application to isotropic cohesive soils.', Journal of Engineering Mechanics 112(12), 1263{1291. Dafalias, Y., Manzari, M. and Papadimitriou, A. (2006), `SANICLAY: simple anisotropic clay plasticity model', International Journal for Numerical and Analytical Methods in Geomechanics 30(12), 1231{1257. Dafalias, Y. and Popov, E. (1975), `A model of nonlinearly hardening materials for complex loading.', Acta Mechanics 21, 173{192. Dafalias, Y. and Taiebat, M. (2012), Advantages and disadvantages of some rotational hard- ening rules in clay plasticity, in `Third International Workshop on Modern Trends in Geomechanics IW-MTG3, 4{5 September 2012'. Desai, C., Somasundaram, S. and Frantziskonis, G. (1986), `A hierarchical approach for constitutive modelling of geologic materials', International Journal for Numerical and Analytical Methods in Geomechanics 10(3), 225{257. Desai, C., Wathugala, G. and Matlock, H. (1993), `Constitutive model for cyclic behavior of clays. II: applications', Journal of Geotechnical Engineering 119(4), 730{748. Hardin, B. and Drnevich, V. (1972), `Shear modulus and damping in soils: Measurement and parameter eects.', American Society of Civil Engineers, Journal of the Soil Mechanics and Foundations Division 98(SM6), 603{624. Hashiguchi, K. and Ueno, M. (1977), `Elasto-plastic constitutive laws of granular materials', Proceedings of the International Conference on Soil Mechanics and Foundation Engineer- ing = Comptes Rendus du Congres International de Mecanique des Sols et des Travaux de Fondations (9), 73{82. 120 Bibliography Huang, M., Liu, Y. and Sheng, D. (2011), `Simulation of yielding and stress{stain behavior of shanghai soft clay', Computers and Geotechnics 38(3), 341{353. Idriss, I. and Boulanger, R. (2008), Soil Liquefaction during Earthquakes, Earthquake Engi- neering Research Institute (EERI), Oakland, California, USA. Iwan, W. (1967), `On a class of models for the yielding behavior of continuous and composite systems', Journal of Applied Mechanics 34(3), 612{617. Jiang, J. and Ling, H. (2010), `A framework of an anisotropic elastoplastic model for clays', Mechanics Research Communications 37(4), 394{398. Jiang, J., Ling, H. and Kaliakin, V. (2012), `An associative and non-associative anisotropic bounding surface nodel for clay', Journal of Applied Mechanics, Transactions ASME 79(3). Jung, Y. and Yune, C. (2011), `Reappraisal of the anisotropic bounding surface model of small-strain behavior for clays', KSCE Journal of Civil Engineering 15(3), 463{472. Koiter, W. (1953), `Stress-strain relations, uniqueness and variational theorems for elastic- plastic materials with a singular yield surface', Quarterly of Applied Mathematics 11, 350. Li, T. and Meissner, H. (2002), `Two-surface plasticity model for cyclic undrained behavior of clays', Journal of Geotechnical and Geoenvironmental Engineering 128(7), 613{626. Ling, H., Yue, D., Kaliakin, V. and Themelis, N. (2002), `Anisotropic elastoplastic bounding surface model for cohesive soils', Journal of Engineering Mechanics 128(7), 748{758. Mroz, Z. (1967), `On the description of anisotropic workhardening', Journal of the Mechanics and Physics of Solids 15(3), 163{175. Mroz, Z. (1969), `An attempt to describe the behavior of metals under cyclic loads using a more general workhardening model', Acta Mechanica 7, 199{212. Papadimitriou, A., Manzari, M. and Dafalias, Y. (2005), Calibration of a simple anisotropic plasticity model for soft clays, number 128, pp. 415{424. 121 Bibliography Roscoe, K. and Burland, J. (1966), On the generalised stress-strain behavior of \wet" clay, in Heymann and Leckie, eds, `Engineering Plasticity', Cambridge University Press, pp. 53{ 609. Roscoe, K. and Poorooshasb, H. (1963), `A theoretical and experimental study of strains in triaxial compression tests on normally consolidated clays', Geotechnique 13(1), 12{38. Roscoe, K., Schoeld, A. and Thurairajah, A. (1963), `Yielding of clays in states wetter than critical', Geotechnique 13(3), 211{240. Roscoe, K., Schoeld, A. and Wroth, C. (1958), `On the yielding of soils', Geotechnique 8(1), 22{53. Sangrey, D., Henkel, D. and Esrig, M. (1969), `The eective stress response of a saturated clay soil to repeated loading', Canadian Geotechnical Journal 6(3), 241{252. Seed, H., Wong, R., Idriss, I. and Tokimatsu, K. (1986), `Moduli and damping factors for dy- namic analyses of cohesionless soils', Journal of Geotechnical Engineering 112(11), 1016{ 1032. Sheu, W. (1984), Modeling of stress-strain-strength behavior of a clay under cyclic loading (Soil Dynamics), Ph.D. dissertation, Univeristy of Colorado, Boulder, Colorado, USA. Stallebrass, S. and Taylor, R. (1997), `The development and evaluation of a constitutive model for the prediction of ground movements in overconsolidated clay', Geotechnique 47(2), 235{253. Suebsuk, J., Horpibulsuk, S. and Liu, M. (2011), `A critical state model for overconsolidated structured clays', Computers and Geotechnics 38(5), 648{658. Taiebat, M. and Dafalias, Y. (2013), Rotational hardening and uniqueness of critical state line clay plasticity, in `Second International Symposium on Constitutive Modeling of Ge- omaterials: Advances and New Applications, 15{16 October 2012'. Taiebat, M., Dafalias, Y. and Kaynia, A. (2010), Bounding surface model for natural anisotropic clays, in P. Papanastasiou, E. Papamichos, A. Zervos and M. Stavropoulou, 122 Bibliography eds, `Ninth HSTAM International Congress on Mechanics, Proceedings of Vardoulakis mini-symposia', Limassol, Cyprus, pp. 43{47. Taiebat, M., Dafalias, Y. and Peek, R. (2010), `A destructuration theory and its application to SANICLAY model', International Journal for Numerical and Analytical Methods in Geomechanics 34(10), 1009{1040. Takahashi, M., Hight, D. and Vaughan, P. (1980), Eective stress changes observed during undrained cyclic triaxial tests on clay, in G. Pande and O. Zienkiewicz, eds, `Soils under Cyclic and Transient Loading', Vol. 1, 2, pp. 201{209. Vucetic, M. and Dobry, R. (1988), `Degradation of marine clays under cyclic loading.', Jour- nal of Geotechnical Engineering 114(2), 133{149. Vucetic, M. and Dobry, R. (1991), `Eect of soil plasticity on cyclic response', Journal of Geotechnical Engineering 117(1), 89{107. Wathugala, G. and Desai, C. (1993), `Constitutive model for cyclic behavior of clays. I: theory', Journal of Geotechnical Engineering 119(4), 714{729. Wheeler, S., Naatanen, A., Karstunen, M. and Lojander, M. (2003), `An anisotropic elasto- plastic model for soft clays', Canadian Geotechnical Journal 40(2), 403{418. Whittle, A. (1993), `Evaluation of a constitutive model for overconsolidated clays', Geotech- nique 43(2), 289{313. Wood, D. (1990), Soil behaviour and critical state soil mechanics, Cambridge University Press. Yasuhara, K., Hirao, K. and FL Hyde, A. (1992), `Eects of cyclic loading on undrained strength and compressibility of clay', Soils and Foundations 32(1), 100{116. Yu, H. (1998), `CASM: A unied state parameter model for clay and sand', International Journal for Numerical and Analytical Methods in Geomechanics 22(8), 621{653. Yu, H., Khong, C. and Wang, J. (2007), `A unied plasticity model for cyclic behaviour of clay and sand', Mechanics Research Communications 34(2), 97{114. 123 Bibliography Zergoun, M. (1991), Eective stress response of clay to undrained cyclic loading, Ph.D. dissertation, The University of British Columbia, Vancouver, British Columbia, Canada. Zergoun, M. and Vaid, Y. (1994), `Eective stress response of clay to undrained cyclic loading', Canadian Geotechnical Journal 31(5), 714{727. Zienkiewicz, O., Leung, K. and Pastor, M. (1985), `Simple model for transient soil loading in earthquake analysis. I: Basic model and its application.', International Journal for Numerical and Analytical Methods in Geomechanics 9(5), 453{476. 124

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

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>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.24.1-0073338/manifest

Comment

Related Items