xml version 1.0 encoding UTF8 standalone no
record xmlns http:www.loc.govMARC21slim xmlns:xsi http:www.w3.org2001XMLSchemainstance xsi:schemaLocation http:www.loc.govstandardsmarcxmlschemaMARC21slim.xsd
leader nam Ka
controlfield tag 001 001498284
003 fts
006 med
007 cr mnuuuuuu
008 041209s2004 flua sbm s0000 eng d
datafield ind1 8 ind2 024
subfield code a E14SFE0000574
035
(OCoLC)57724330
9
AJU6889
b SE
SFE0000574
040
FHM
c FHM
090
TA145 (ONLINE)
1 100
Sallam, Amr M.
0 245
Studies on modeling angular soil particles using the discrete element method
h [electronic resource] /
by Amr M. Sallam.
260
[Tampa, Fla.] :
University of South Florida,
2004.
502
Thesis (Ph.D.)University of South Florida, 2004.
504
Includes bibliographical references.
500
Includes vita.
516
Text (Electronic thesis) in PDF format.
538
System requirements: World Wide Web browser and PDF reader.
Mode of access: World Wide Web.
Title from PDF of title page.
Document formatted into pages; contains 228 pages.
520
ABSTRACT: The Discrete Element Method was first introduced by Cundall and Strack (1979) to model granular soils within the context of geotechnical engineering. The material is modeled as a random assembly of discrete elements. Each particle interacts with neighboring particles through contact forces that can be built up and broken at any time. The particles were modeled as discs in 2D or as spheres in 3D. Research studies have been conducted to improve the simulation of actual grain shapes. Ashmawy et al. (2003) developed the overlapping rigid clusters (ORC) method to accurately model irregular particle shapes. The idea relies on clumping a number of overlapping discs such that their coincides with that of the actual particle. In this dissertation, experimental verification program is presented. An experimental setup was built and modelgrains were manufactured in the laboratory. A numerical simulation for the experimental test was carried out.The numerical and experimental results were compared qualitatively and quantitatively. A good agreement was observed within small displacements ranges. However, results were considerably different at large displacements. Numerical results utilizing the ORC method were closer to the experimental results than those of discs. A sequential and operatorindependent procedure, which relies on the ORC concept, was developed. Identical inertial properties between the actual particle and the model were ensured. The new procedure was implemented for rounded and angular particles. The effect of particle shape and angularity on the strength and dilatancy characteristics of granular soils was investigated. A modified shape factor, which relies on the work introduced by Sukumaran and Ashmawy (2001), was developed. A series of pure shear testing simulations was performed on different shape and angularity particle groups.Angularity had a remarkable effect on strength and dilatancy properties compared to shape. The effect of interparticle friction on dilatancy was studied. An attempt was made to use an equivalent interparticle friction to model different particle shapes. It was concluded that there is no onetoone equivalency between interparticle friction and shape or angularity. Instead, the interparticle friction must be continuously altered as a function of confining pressure and void ratio to achieve the required effect.
590
Adviser: Ashmawy, Alaa.
653
dilation.
granular soils.
Discrete Element Method.
numerical modeling.
angularity.
690
Dissertations, Academic
z USF
x Civil Engineering
Doctoral.
773
t USF Electronic Theses and Dissertations.
4 856
u http://digital.lib.usf.edu/?e14.574
PAGE 1
Studies on Modeling Angula r Soil Particles Using the Discrete Element Method by Amr M. Sallam A dissertation submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy Department of Civil and Environmental Engineering College of Engineering University of South Florida Major Professor: Alaa Ashmawy, Ph.D. Manjriker Gunaratne. Ph.D. Sunil Saigal, Ph.D. Ashraf Ayoub, Ph.D. Carlos Santamarina, Ph.D. Autar Kaw, Ph.D. Date of Approval: November 12, 2004 Keywords: numerical modeling, discrete element method, granular soils, dilation, angularity Copyright 2004 Amr M. Sallam
PAGE 2
DEDICATION This work is dedicated to my wife, Sama r, my daughter, Nour, and my son, Seif. Samar has been a tremendous support. She sacrif iced her career to stay at home to take care of me and our kids. She was always ther e in the hard days throughout the research and writing processes. Without her love and support, this wo rk could not be completed.
PAGE 3
ACKNOWLEDGMENT I want to thank my supervisor, Dr Al aa Ashmawy, for his great ideas and support throughout the research. Dr Ashmawy was al ways there to suggest solutions for tough problems, to improve my writing quality, and to expand my academic knowledge in Geotechnical Engineering. I wa nt to thank my colleagues (Newel, Rory, Naim, Jessica, Jorge, Nivedita, Brian, and Whitney) for th eir encouragement. They are really good friends and I hope our friends hip will continue forever.
PAGE 4
i TABLE OF CONTENT LIST OF TABLES.............................................................................................................vi LIST OF FIGURES..........................................................................................................vii ABSTRACTÂ…. ..............................................................................................................xvii CHAPTER ONE: INTRODUCTION.................................................................................1 1.1 Modeling of angular particles................................................................................2 1.2 Effect of shape and angularity on dilation of granular soils..................................3 1.3 Organization of the study.......................................................................................4 1.3.1 Chapter one: introduction............................................................................4 1.3.2 Chapter two: literature review.....................................................................4 1.3.3 Chapter three: experimental ve rification of modeling angular sand particles using DEM...................................................................................5 1.3.4 Chapter four: the modified overlapping rigid clusters................................5 1.3.5 Chapter five: effect of par ticle shape and angularity on the dilation of granular soils.............................................................................6 1.3.6 Chapter six: conclusions and recommendations..........................................7 1.3.7 Appendix A: theory and background of DEM............................................7 CHAPTER TWO: LITERATURE REVIEW.....................................................................8 2.1 Discrete element method for modeling granular media........................................8 2.1.1 DEM versus FEM for modeling granular soils.........................................10 2.1.2 Earlier published work in DEM................................................................10 2.1.2.1 Fundamental investigations : monotonic and cyclic shear testing...................................................................................................11 2.1.2.2 Contact models is DEM..................................................................16 2.1.2.3 Coupled DEM and FEM.................................................................17
PAGE 5
ii 2.1.2.4 Experimental verification of DEM..................................................18 2.1.2.5 Other applications...........................................................................19 2.1.3 Modeling noncircular particles................................................................19 2.1.3.1 Modeling angular particle s as rigid clusters....................................24 2.1.3.2 Overlapping rigid cluste rs (ORC) technique...................................25 2.2 Dilatancy in granular soils...................................................................................29 2.2.1 Friction angle, dilation angle, and dilatancy.............................................31 2.2.1.1 Relation between fricti on and dilation angles.................................33 2.2.1.2 True and critical st ate friction angles..............................................37 2.2.2 Development of the theory of dilatancy....................................................38 2.2.2.1 RoweÂ’s theory.................................................................................38 2.2.2.2 BoltonÂ’s relative dilatancy index.....................................................40 2.2.3 Statedependent dilatancy..........................................................................43 2.2.3.1 Problems with unique relationship between d and .......................44 2.2.4 General expression for st ate dependent dilatancy.....................................46 2.3 Particle shape characterization............................................................................50 2.3.1 Particle shape charact erization techniques................................................51 2.3.2 The method introduced by Sukumaran and Ashmawy (2001)..................53 2.3.2.1 The distortion diagram....................................................................54 2.3.2.2 Shape factor.....................................................................................55 2.3.2.3 Angularity factor.............................................................................56 CHAPTER THREE: EXPERIMENTAL VERIFICATION OF MODELING ANGULAR PARTICLES USING DEM...............................................58 3.1 Introduction..........................................................................................................58 3.2 Background..........................................................................................................59 3.3 Overlapping rigid clusters (ORC)........................................................................60 3.4 Experimental program.........................................................................................62 3.4.1 Material.....................................................................................................63 3.4.2 Test setup...................................................................................................64 3.4.3 Image distortion effect...............................................................................68 3.4.4 Contact normal stiffness............................................................................70
PAGE 6
iii 3.4.5 Experimental testing procedures...............................................................72 3.4.6 Experimental results..................................................................................73 3.4.6.1 Fraser river sand with the rectangular piston..................................73 3.4.6.2 Michigan sand with the trapezoidal piston......................................74 3.5 Numerical simulation...........................................................................................75 3.6 Comparison and discussion.................................................................................83 3.6.1 Fraser river sand........................................................................................84 3.6.2 Michigan sand...........................................................................................88 3.7 Parametric study..................................................................................................92 3.7.1 Effect of interparticle friction....................................................................92 3.7.2 Effect of global damping...........................................................................98 3.7.3 Effect of contact normal stiffness............................................................104 3.8 Conclusion.........................................................................................................104 CHAPTER FOUR: THE MODIFIED OVERLAPPING RIGID CLUSTERS..............109 4.1 Introduction........................................................................................................109 4.2 The modified ORC method...............................................................................109 4.2.1 Operatorindependent partic le creation procedure..................................109 4.2.2 Kinetic compatibility equations...............................................................112 4.2.3 Procedure.................................................................................................113 4.3 Implementation of the modified ORC method..................................................114 4.3.1 Michigan sand particle # 1......................................................................115 4.3.2 Fraser river sand particle # 4...................................................................115 4.4 Conclusion.........................................................................................................119 CHAPTER FIVE: EFFECT OF PARTICLE SHAPE AND ANGULARITY ON THE DILATION OF GRANULAR SOILS.........................................122 5.1 Introduction........................................................................................................122 5.2 Background........................................................................................................122 5.2.1 Dilatancy in angular soils........................................................................122 5.2.2 Dilatancy and the discrete element method.............................................123 5.3 Modification to particle shape characte rization.................................................126
PAGE 7
iv 5.4 Particle groups with respec t to shape and angularity.........................................132 5.5 Numerical simulations.......................................................................................132 5.5.1 Sample preparation..................................................................................134 5.5.2 Sample consolidation to the desired confining pressure.........................134 5.5.3 Pure shear loading...................................................................................136 5.5.4 Simulation program.................................................................................138 5.6 Numerical results...............................................................................................140 5.6.1 Confining pressure (1.0MPa), voids ratio (0.162), and 2574 particles...................................................................................................140 5.6.1.1 Effect of angularity........................................................................140 5.6.1.2 Effect of particle shape..................................................................145 5.6.2 Confining pressure (500 KPa) voids ratio (0.25), and 2394 particles...................................................................................................152 5.6.2.1 Effect of angularity........................................................................152 5.6.2.2 Effect of particle shape..................................................................153 5.6.3 Confining pressure (250 KPa) voids ratio (0.35), and 2215 particles...................................................................................................158 5.6.3.1 Effect of angularity........................................................................158 5.6.3.2 Effect of particle shape..................................................................159 5.7 Critical state line and particle shape characteristics..........................................164 5.8 Effect of interparticle friction............................................................................166 5.8.1 Confining pressure (1.0 MPa) voids ratio (0.162), and 2574 particles...................................................................................................166 5.8.2 Confining pressure (500 KPa) voids ratio (0.25), and 2394 particles...................................................................................................167 5.8.3 Confining pressure (250 KPa) voids ratio (0.35), and 2215 particles...................................................................................................170 5.9 The equivalent interparti cle friction coefficient................................................175 5.10 Discussion........................................................................................................180 CHAPTER SIX: SUMMARY AND RECOMMENDATIONS.....................................182 6.1 Summary............................................................................................................182 6.2 Recommendations for future studies.................................................................184
PAGE 8
v REFERNCES ..............................................................................................................186 APPENDICES ..............................................................................................................196 Appendix A: theory and background of DEM........................................................197 A.1 Calculation cycle.......................................................................................197 A.2 Forcedisplacement law.............................................................................197 A.3 Motion.......................................................................................................200 A.4 Damping....................................................................................................201 A.5 Energy dissipation.....................................................................................203 A.6 Time step...................................................................................................203 A.7 Input parameters........................................................................................203 A.8 Contact constitutive models.......................................................................204 ABOUT THE AUTHOR.......................................................................................End page
PAGE 9
vi LIST OF TABLES Table 21. Literature survey of DEM up to 1992, Dobry and NG (1996)Â…Â…Â…Â…Â….12 Table 31. Maple wood propertiesÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….63 Table 32. Contact normal stiffness result s for Fraser River sand modelgrainsÂ…Â…..73 Table 41. Procedure of enforcing the kinetic compatibility equations for particle # 1 of Michigan sandÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….114 Table 42. Elements/Discs properties for particle # 1 of Michigan sandÂ…Â…Â…Â…Â…117 Table 43. Results of the ORC and the modified ORC techniques for Michigan sand particle # 1Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….117 Table 44. Elements/Discs properties fo r particle # 4 of Fraser River sandÂ…Â…Â…Â…120 Table 45. Results of the ORC and the modified ORC techniques for Fraser River sand particle # 4Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…120 Table 51. The modified shape factor fo r some standard geometric shapesÂ…Â…Â…...130 Table 52. Shape factors, global shape f actors, modified shape factors, and angularity factors for soil s used in the studyÂ…Â…Â…Â…...Â…Â…Â…Â…Â…Â…...131 Table 53. Shape and angularity fa ctors for different groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…..132 Table 54. Model dimensions and propertiesÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..135 Table 55. Simulation programÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…139
PAGE 10
vii LIST OF FIGURES Figure 21. Two circular particles in contactÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..21 Figure 22. Polygonshaped particlesÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..21 Figure 23. Ellipseshaped particles: (a) pa rticle assembly, (b) two ellipses in contact, (c) contact between an el lipseshaped particle and a wallÂ…Â…Â…..22 Figure 24. Ovalshaped particlesÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...24 Figure 25. Axisymmetrical particles: (a) sample particles, (b) contact between two multielement particlesÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..24 Figure 26. (a) Outline of sand particle (b) DEM disc element superimposed over sand particle, (c) a number of DEM disc particles are joined together is a rigid c onfiguration (cluster), (d) several possible combination of discs to form clustersÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..26 Figure 27. Disc elements inscribed within a particle outline to capture the shapeÂ…...26 Figure 28. Random assembly of eight circular particles (left) and the transformed equivalent particles (right)Â…Â…Â…Â…Â…Â…Â…...Â…Â…Â…Â…Â…Â…28 Figure 29. Dilation of dense sa nd during simple shear testÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…...32 Figure 210. Curved peak strength envelope on MohrCoulomb plotÂ…Â…Â…Â…Â…Â…Â…..32 Figure 211. Angles of friction and dilationÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….34 Figure 212. The sawtooth model for dilatancyÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...34 Figure 213. Definitions of Rowe Â’s geometry propertiesÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….39 Figure 214. Angle versus porosity for different materialsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…41 Figure 215. Deviation of triaxial dilatancy from biaxial stateÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â….41
PAGE 11
viii Figure 216. Variation of dilatancy with material state (data from Verdugo and Ishihara 1996). undrained re sponse of sand (a) with different densities (b) under di fferent confining pressuresÂ…Â…Â…Â…Â…Â…..46 Figure 217. Variation in the phase transformation stress ratio with material state (data from Verdugo and ishihara 1996)Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..46 Figure 218. The state paramete r and critical state lineÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…50 Figure 219. Relationship between peak fr iction angle and maximum dilation angleÂ….50 Figure 220. Illustration of WadellÂ’s procedure for eval uating particle roundnessÂ…..Â…54 Figure 221. Ideal geometric shape used to define the shape and angularity factorsÂ…...54 Figure 222. Distortion diagram for th e particle shown, Sukumaran (1996)Â…Â…Â…Â…Â…56 Figure 31. Fraser River sa nd modelgrains samplesÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..65 Figure 32. Michigan sand modelgrains samplesÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..65 Figure 33. The box alo ng with the pistonÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..66 Figure 34. The experimental setup for Fr aser River sand modelparticles with the rectangular pistonÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...66 Figure 35. The experimental setup for Michigan sand modelparticles with the trapezoidal pistonÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..Â…Â…...67 Figure 36. The two ruled Lsquare s used to locate modelgrainsÂ…Â…Â…Â…Â…Â…Â…Â…..67 Figure 37. Initial conditions for Fraser River sand modelgr ainsrectangular pistonÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...Â…Â…Â….68 Figure 38. Initial conditions for Mich igan sand modelgrainstrapezoidal pistonÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...Â…Â…Â…Â….69 Figure 39. Distortion calibration us ing graded tran sparent sheetÂ…Â…Â…Â…Â…Â…Â…Â…..70 Figure 310. Loading system used in dete rmining the contact normal stiffnessÂ…Â…Â…...71 Figure 311. Compression test on mode lgrain 5A of Fraser River sandÂ…Â…Â…Â…Â…Â….71
PAGE 12
ix Figure 312. Loaddisplacement curves for Fraser River Sand modelgrain # 8Z..Â…Â….72 Figure 313. Experimental in itial conditions (tests 15) for Fraser River sand modelgrainsÂ…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…76 Figure 314. Experimental results after imposing 50 mm displacement (tests 15) for Fraser River sand modelgrainsÂ…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….76 Figure 315. Experimental results afte r imposing100 mm displacement (tests 15) for Fraser River sand modelgrainsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….77 Figure 316. Experimental results after imposing 150 mm displacement (tests 15) for Fraser River sand modelgrainsÂ…Â…Â…Â…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…77 Figure 317. Experimental results after imposing 200 mm displacement (tests 15) for Fraser River sand modelgrainsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….78 Figure 318. Experimental results after imposing 250 mm displacement (tests 15) for Fraser River sand modelgrainsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….78 Figure 319. Experimental initial cond itions (tests 15) for Michigan sand modelgrainsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….79 Figure 320. Experimental results after imposing 50 mm displacement (tests 15) for Michigan River sand modelgrainsÂ…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…...79 Figure 321. Experimental results after imposing 100 mm displacement (tests 15) for Michigan River sand modelgrainsÂ…Â…Â…Â…Â…Â…Â…..Â…Â…Â…Â…Â…Â…...80 Figure 322. Experimental results after imposing 150 mm displacement (tests 15) for Michigan River sand modelgrainsÂ…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…...80 Figure 323. Experimental results after imposing 200 mm displacement (tests 15) for Michigan River sand modelgrainsÂ…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…...81 Figure 324. Experimental results after imposing 250 mm displacement (tests 15) for Michigan River sand modelgrainsÂ…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…...81 Figure 325. Experimental and numerical initial conditions fo r Fraser River sandÂ…Â…..83 Figure 326. Experimental and numerical initial conditions for Michigan sandÂ…Â…Â…..83 Figure 327. Experimental and numerical results for Michigan sand after imposing 50 mm displacementÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….85
PAGE 13
x Figure 328. Experimental and numeri cal results for Michigan sand after imposing 100 mm displacementÂ…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..86 Figure 329. Experimental and numeri cal results for Michigan sand after imposing 150 mm displacementÂ…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..86 Figure 330. Experimental and numeri cal results for Michigan sand after imposing 200 mm displacementÂ…Â…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..87 Figure 331. Experimental and numeri cal results for Michigan sand after imposing 250 mm displacementÂ…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..87 Figure 332. Experimental and numerical rotations for Fraser River sand after imposing 50 mm displacementÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…89 Figure 333. Experimental and numerical rotations for Fraser River sand after imposing 150 mm displacementÂ…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..89 Figure 334. Experimental and numerical rotations for Fraser River sand after imposing 250 mm displacementÂ…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..90 Figure 335. Experimental and numerical Vertical (Y) displacement for Fraser River sand after imposing 50 mm displacementÂ…Â…Â…Â…Â….Â…Â…Â…Â…Â….90 Figure 336. Experimental and numerical Vertical (Y) displacement for Fraser River sand after imposing 150 mm displacementÂ…Â…Â…Â…Â…Â….Â…Â…Â…...91 Figure 337. Experimental and numerical Vertical (Y) displacement for Fraser River sand after imposing 250 mm displacementÂ…Â…Â…Â…Â…Â…Â….Â…Â…...91 Figure 338. Experimental and numeri cal results for Michigan sand after imposing 50 mm displacementÂ…Â…Â…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…92 Figure 339. Experimental and numeri cal results for Michigan sand after imposing 150 mm displacementÂ…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..93 Figure 340. Experimental and numeri cal results for Michigan sand after imposing 250 mm displacementÂ…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..93 Figure 341. Experimental and numerical rotations for Michigan sand after imposing 50 mm displacementÂ…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…94
PAGE 14
xi Figure 342. Experimental and numeri cal rotations for Michigan sand after imposing 150 mm displacementÂ…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..94 Figure 343. Experimental and numerical rotations for Michigan sand after imposing 250 mm displacementÂ…Â…Â…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…..95 Figure 344. Experimental and numerical vertical (Y) displacement for Michigan sand after imposing 50 mm displacementÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â….Â…Â…..95 Figure 345. Experimental and numeri cal rotations for Michigan sand after imposing 150 mm displacementÂ…Â…Â…Â…Â…Â…..Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….96 Figure 346., Experimental and numer ical rotations for Michigan sand after imposing 250 mm displacementÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..96 Figure 347. Effect of interparticle fr iction (Fraser River sand) on particles rotation after imposing 50 mm displacementÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...97 Figure 348. Effect of in terparticle friction (Frase r River sand) on particles rotation after imposing 150 mm displacementÂ…Â…Â…Â…Â…Â…Â…...Â…Â…Â…..97 Figure 349. Effect of interparticle fr iction (Fraser River sand) on particles rotation after imposing 250 mm displacementÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….99 Figure 350. Effect of interparticle fr iction (Fraser River sand) on particles vertical (Y) displacement afte r imposing 50 mm displacementÂ…Â…Â…Â…...99 Figure 351. Effect of interparticle fr iction (Fraser River sand) on particles vertical (Y) displacement afte r imposing 150 mm displacementÂ…Â…Â…...100 Figure 352. Effect of interparticle fr iction (Fraser River sand) on particles vertical (Y) displacement afte r imposing 250 mm displacementÂ…Â…Â…...100 Figure 353. Effect of global damping (Fraser River sand) on particles rotation after imposing 50 mm displacementÂ…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….101 Figure 354. Effect of global damping (Fraser River sand) on particles rotation after imposing 150 mm displacementÂ…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…...Â…101 Figure 355. Effect of global damping (Fraser River sand) on particles rotation after imposing 250 mm displacementÂ…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…...102
PAGE 15
xii Figure 356. Effect of global damping (F raser River sand) on particles vertical (Y) displacement after im posing 50 mm displacementÂ…Â…Â…Â….Â…Â…Â….102 Figure 357. Effect of global damping (Fraser River sand) on particles vertical (Y) displacement after im posing 150 mm displacementÂ…Â…Â…Â….Â…Â…...103 Figure 358. Effect of global damping (F raser River sand) on particles vertical (Y) displacement after im posing 250 mm displacementÂ…Â…Â…Â…Â….Â…...103 Figure 359. Effect of contact normal s tiffness (Fraser River sand) on particles rotation after imposing 50 mm displacementÂ…Â…Â…Â…Â…Â…Â…...Â…Â…Â…..105 Figure 360. Effect of contact normal s tiffness (Fraser River sand) on particles rotation after imposing 150 mm displacementÂ…Â…Â…Â…Â…Â…...Â…Â…Â…Â…105 Figure 361. Effect of contact normal s tiffness (Fraser River sand) on particles rotation after imposing 250 mm displacementÂ…Â…Â…Â…Â…Â…...Â…Â…Â…Â…106 Figure 362. Effect of contact normal stiffness (Fraser River sand) on particles vertical (Y) displacement afte r imposing 50 mm displacementÂ…Â…Â…Â….106 Figure 363. Effect of contact normal s tiffness (Fraser River sand) on particles vertical (Y) displacement afte r imposing 150 mm displacementÂ…Â…...Â…107 Figure 364. Effect of contact normal s tiffness (Fraser River sand) on particles vertical (Y) displacement afte r imposing 250 mm displacementÂ…Â…Â…...107 Figure 41. Outline of particle # 1 of Michigan sand along with disc#1Â…Â…Â…Â…Â….111 Figure 42. Sequence of adding discs for particle # 1 of Michigan sand using the modified ORC methodÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….111 Figure 43. Outline of particle # 1 of Michigan sandÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…116 Figure 44. Generation of pa rticle # 1 of Michigan sa nd using the ORC technique....116 Figure 45. Outline of particle # 4 of Fraser River sandÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...118 Figure 46. Generation of particle # 4 of Fraser River sand using the ORC TechniqueÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...118 Figure 47. Sequence of adding discs for part icle # 4 of Fraser River sand using the modified ORC methodÂ…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…119
PAGE 16
xiii Figure 51. Schematic illust ration of bonded particleÂ…Â…Â…Â…Â…..................................124 Figure 52. Stress ratio and volumetri c dilation at normal stress of 100kPa (30,000 particle and aspect ratio 1: 3), after Ni et al. (2000)Â…Â…Â…Â…Â…..125 Figure 53. Effect of particle shape (normal stress = 100kPa) afterNi et al. (2000)Â…125 Figure 54. The shape factor versus the angularity factors for all particlesÂ…Â…Â…Â…..128 Figure55. The global shap e factor versus the angularity factors for all articlesÂ…Â…128 Figure 56. Sampling and rotation anglesÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…129 Figure 57. The modified shape factor versus the angularity factors for all particlesÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..129 Figure 58. Typical loadi ng conditions for strength and dilatancy testing: (a) triaxial; (b) direct shear; (c) si mple shear; (d) pure shearÂ…Â…Â…Â…Â….135 Figure 59. Stress path duri ng pure shear testing (within the context of critical state soil mechanics)Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…..137 Figure 510. Shear stressshear strain plot for circular particlesÂ…Â…Â…Â…Â…Â…Â…Â….Â…141 Figure 511. Volumetric strain shear strain plot for circ ular particlesÂ…Â…Â…Â…Â…Â…...141 Figure 512. Shear stresssh ear strain curves for differe nt angularity groupsÂ…Â…Â…Â…143 Figure 513. Volumetric st rainshear strain curves for different angularity groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…143 Figure 514. Maximum volumetric strain versus angularity factor for different angularity groupsÂ…Â…Â…Â…...Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….144 Figure 515. Maximum dila tion and residual shearing re sistance angles versus angularity factor for different angularity groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…144 Figure 516. Shear stressshear strain cu rves for different shape groupsÂ…Â…Â…Â…Â…...148 Figure 517. Volumetric strain shear strain curv es for different shape groupsÂ…Â…Â…..148 Figure 518. Maximum volumetric strain versus shape factor for different shape groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...149
PAGE 17
xiv Figure 519. Maximum dilation and resi dual shearing resistan ce angles versus shape factor for different shape groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….149 Figure 520. Shear stressshear strain cu rves for different modified shape factor groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….150 Figure 521. Volumetric stra inshear strain curves for different modified shape factor groupsÂ…Â…Â…Â…Â…Â…..Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….150 Figure 522. Maximum volumetric strain versus shape factor for different modified shape factor groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….151 Figure 523. Maximum dila tion and residual shearing re sistance angles versus shape factor for different modifi ed shape factor groupsÂ…Â…Â…Â…Â…Â…....151 Figure 524. Shear stresssh ear strain curves for differe nt angularity groupsÂ…Â…Â…Â…154 Figure 525. Volumetric strains hear strain curves for differe nt angularity groupsÂ… ..154 Figure 526. Maximum volumetric strain versus angularity factor for different angularity groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...Â…Â….155 Figure 527. Maximum dila tion and residual shearing re sistance angles versus angularity factor for different angularity groupsÂ…Â…Â…Â…Â…Â…Â…Â…...Â….155 Figure 528. Shear stressshear strain cu rves for different modified shape factor groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..Â…...156 Figure 529. Volumetric stra inshear strain curves for different modified shape factor groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….156 Figure 530. Maximum volumetric strain versus shape factor for different modified shape factor groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….157 Figure 531. Maximum dilation and resi dual shearing resistan ce angles versus shape factor for different modi fied shape factor groupsÂ…Â…Â…Â…Â…Â…...157 Figure 532. Shear stresssh ear strain curves for differe nt angularity groupsÂ…Â…Â…Â…160 Figure 533. Volumetric strain shear strain curves for diffe rent angularity groupsÂ…...160 Figure 534. Maximum volumetric strain versus angularity for different angularity groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…161
PAGE 18
xv Figure 535. Maximum dila tion and residual shearing re sistance angles versus angularity for different angularity groupsÂ…Â…Â…Â…Â…..Â…Â…Â…Â…Â…Â…Â…161 Figure 536. Shear stressshear strain curves for different modified shape groupsÂ…....162 Figure 537. Volumetric stra inshear strain curves fo r different modified shape groupsÂ…Â…Â….............................................................................................162 Figure 538. Maximum volumetric strain versus shape factor for different modified shape groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...163 Figure 539. Maximum dila tion and residual shearing re sistance angles versus shape factor for different m odified shape groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…...163 Figure 540. Critical state lines fo r different angularity groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â….165 Figure 541. Critical state lines fo r different shape groupsÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….165 Figure 542. Shear stressshear strain cu rves for different interparticle friction coefficients (circular particles)Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….Â…..168 Figure 543. Volumetric stra inshear strain curves for different interparticle friction coefficients (circu lar particles)Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..168 Figure 544. Maximum volumetric strain ve rsus interparticle friction coefficient (circular particles)Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..Â…Â…Â…Â…169 Figure 545. Maximum dilation and resi dual shearing resistan ce angles versus interparticle friction coefficient (circular particles)Â…Â…Â…Â…Â…...Â…Â…Â…169 Figure 546. Shear stressshear strain cu rves for different interparticle friction coefficients (circular particles)Â…Â…Â…Â…Â…Â…Â….Â…Â…Â…Â…Â…Â…Â…Â…Â…..171 Figure 547. Volumetric stra inshear strain curves for different interparticle friction coefficients (circu lar particles)Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..171 Figure 548. Maximum volumetric strain ve rsus interparticle friction coefficient (circular particles)Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..172 Figure 549. Maximum dilation and resi dual shearing resistan ce angles versus interparticle friction coefficient (circular particles)Â…Â…Â…Â…Â…Â…Â…Â…...172
PAGE 19
xvi Figure 550. Shear stressshear strain cu rves for different interparticle friction coefficients (circular particles)Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â….Â…Â…..173 Figure 551. Volumetric stra inshear strain curves for different interparticle friction coefficients (circu lar particles)Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..173 Figure 552. Maximum volumetric strain ve rsus interparticle friction coefficient (circular particles)Â…Â…Â…Â…Â…..Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…174 Figure 553. Maximum dilation and resi dual shearing resistan ce angles versus interparticle friction coefficient (circular particles)Â…Â…Â…Â…Â…Â…Â…Â…...174 Figure 554. Critical state lines for different interparticle friction coefficientÂ….Â…Â…..177 Figure 555. Critical state lines for circul ar particles (different interparticle friction coefficient) along with those for different angularity factor groups (f = 0.5)Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…178 Figure 556. Critical state lines for circul ar particles (different interparticle friction coefficient) along with those for different modified shape factor groups (f = 0.5)Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…..179 Figure A.1. Calculation cycle for the DEMÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…...198 Figure A.2. Forcedisplacement law in DEMÂ…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…Â…198
PAGE 20
xvii STUDIES ON MODELING ANGULAR SOIL PARTICLES USING DISCRETE ELEMENT METHOD Amr M. Sallam ABSTRACT The Discrete Element Method was first introduced by Cundall and Strack (1979) to model granular soils within the context of geot echnical engineering. The material is modeled as a random assembly of discrete elements. Each part icle interacts with neighboring particles through contac t forces that can be built up and broken at any time. The particles were modeled as discs in 2D or as spheres in 3D. Research studies have been conducted to improve the simulation of ac tual grain shapes. Ashmawy et al. (2003) developed the overlapping rigid clusters (O RC) method to accurately model irregular particle shapes. The idea relies on clumpi ng a number of overlapping discs such that their coincides with that of the actual particle. In this dissertation, experimental ve rification program is presented. An experimental setup was built and modelgrains were manufactured in the laboratory. A numerical simulation for the experimental test was carried out. The numerical and experimental results were compared quali tatively and quantitatively. A good agreement was observed within small displacements ranges. However, results were considerably
PAGE 21
xviii different at large displacement s. Numerical results utilizi ng the ORC method were closer to the experimental results than those of discs. A sequential a nd operatorindependent procedure, which relies on the ORC concept, wa s developed. Identical inertial properties between the actual particle and the model were ensured. The new procedure was implemented for rounded and angular particles. The effect of particle shape and a ngularity on the strength and dilatancy characteristics of granular soils was investig ated. A modified shape factor, which relies on the work introduced by Sukumaran and Ashmawy (2001), was developed. A series of pure shear testing simulations was performed on different shape and angularity particle groups. Angularity had a remarkable effect on strength and dilatancy properties compared to shape. The effect of interparticle friction on dilatancy was studied. An attempt was made to use an equivalent interparticle fricti on to model different particle shapes. It was concluded that there is no onetoone equi valency between interp article friction and shape or angularity. Instead, the interparticle friction must be conti nuously altered as a function of confining pressure and void ratio to ac hieve the required effect.
PAGE 22
1 CHAPTER ONE: INTRODUCTION Modeling granular materials has been a challenging topic for decades. The particulate nature of granular materials controls their engineering behavior. The continuum assumption has been used to idealize soils and rocks, and numerical solution techniques such as finite element, finite difference, and boundary element methods have been successfully used to model these materi als. The continuum mechanics models are phenomenological and are primarily concerned with the mathematical modeling of the observed phenomenon without giving detailed attention to the fundamental physical significance (Sitharam, 1999). Th e discrete nature makes th e constitutive relationship complex and needs an excessive number of para meters to be able to model the behavior accurately. Another alternative is to model the granular materials experimentally using photoelastic materials (De Josselin De Jong and Veruijit, 1950). This approach presented a physical basis for understanding the behavior of granular soils, but is experimental rather than numerical and is lim ited in terms of its ability to model real systems. The discrete element method (DEM) offers an alternative approach to simulate the behavior of granular materials. The part iculate nature is automatically simulated as particles are discrete and for ces are transferred through the c ontacts between particles. The method can handle a wide range of materi als constitutive behaviors, contact laws,
PAGE 23
2 and arbitrary geometries. The Discrete Element Method (DEM) for modeling granular soils within the context of civil engineeri ng was first introduced by Cundall and Strack (1979). The particles were modeled as a random assembly of discrete discs. Many research studies have since been conducted to improve the simulation of angular grain shapes. 1.1 Modeling of angular particles In the DEM scheme proposed by Cundall a nd Strack (1979), the granular particles were modeled as discs in 2D simulations a nd as spheres in 3D simulations. Modeling noncircular particles, in both 2D and 3D, has b een proposed in different research studies. Among those, are the u tilization of mathematical func tions to describe noncircular outlines, the approximation of the particle shape using polygons, and the combination several circular outlines into a cluster to form more complex shapes. Ashmawy et al (2003) proposed the use of overlapping rigid cl usters (ORC) to more accurately simulate angular particle shapes. The ORC method relies on clumping a number of overlapping discs such that the resulting outline coincides with and is almost identical to the actual particleÂ’s outline. Verification of the ORC method is essential in order to evaluate the accuracy and significance of the new procedure, which is the first part of the dissertation. An experimental validation program wa s carried out, which involves: 1. Designing and building an experiment al setup that allows tracking translations and rotations of modelgrains resulting from an external disturbance.
PAGE 24
3 2. Numerically simulating the experimental setup using PFC2D, with the same initial and boundary conditions and material properties. 3. Tracking the positions and rotations of modelgrains in both the experimental and numerical tests. 4. Qualitatively and quantitatively comparing the experimental and numerical results to verify the effectiveness of the ORC method. 5. Performing a parametric study to dete rmine the effect of interparticle friction, contact stiffness, and global damping on the numerical results. Although the ORC technique can accurately si mulate the granular particle shapes, the technique doesnÂ’t ensure that the center of gravity or the mass mo ment of inertia of both the actual and the created particles are iden tical. Modifications need to be made to ensure that the actual and the cr eated particles have identical inertial properties, which is the second part of the disserta tion. Another issue is that the elements within a particle were created manually without a certain ruled sequence. In order to generalize the ORC method, an operatorindependent and automa ted procedure for clumping the discs is required (this issue is also included in the seco nd part of the research). Having such an automated procedure will allow the creation of a large number of operatorindependent modelparticles in a short time compared to manual modelparticle creation, which will ensure more realistic and accurate numerical modeling. 1.2 Effect of shape and angularity on dilation of granular soils Dilation of granular soils plays an impor tant role in their stress and strain characteristics. The dependence of dila tion on both the density and the confining
PAGE 25
4 pressure has been studied and empirical and th eoretical correlations have been introduced (i.e. Rowe1962, Bolton 1986, Li and Dafalias 2000) The dependency of dilatancy on particle shape has been pointed out in many research studies Discrete element modeling has been used to explore the dependency of dilatancy on particle shape (i.e. Ting et al. 1995, and Ni et al. 2000). Ellipses with diffe rent aspect ratios, polygons, and clusters were used to simulate irregularly shaped pa rticles. Utilizing th e ORC method proposed by Ashmawy et al. (2003), angul ar shapes can be more accur ately simulated. The third part of the dissertation expl ores the effect of particle shape and angularity on the dilatancy of granular materials making us e of the ORC method. Particle shape is characterized using the shape factor and the angularity f actor introduced by Sukumaran (1996) and Sukumaran and Ashmawy (2001). 1.3 Organization of the study 1.3.1 Chapter one: introduction This chapter includes general introduction about modeling of granular soils using the discrete element method, modeling angular shapes, and soil dilatancy. The objectives of the research are then explained and the organizati on of the study is introduced. 1.3.2 Chapter two: literature review This chapter presents a literatur e review of the following topics: 1. The definition of discrete elemen t modeling, comparing DEM and FEM, and earlier published work on DEM a nd modeling noncircular particles.
PAGE 26
5 2. Dilatancy of granular materials, the development of the dilatancy theory, and state dependent dilatancy. 3. Conventional and modern methods used to characterize particle shapes and the method introduced by Suku maran (1996) and Sukumaran and Ashmawy (2003). 1.3.3 Chapter three: experimental verificati on of modeling angular sand particles using DEM This chapter includes the experimental and numerical work that has been conducted to verify the ability of the DEM to model angular par ticles using the ORC method. The chapter is organized as follows: 1. The experimental program (material, test setup, image distortion effect, procedure, and experimental results). 2. Numerical simulations for the experimental setup. 3. Qualitative and quantitative comparisons between the experimental and the numerical results. 4. Parametric study (the effect of in terparticle friction, contact normal stiffness, and global damping). 1.3.4 Chapter four: the modified overlapping rigid clusters This chapter includes the alteration that has been introduced to the ORC method to overcome the difference in center of grav ity and mass moment of inertia between the actual and the generated pa rticles. The chapter is organized as follows:
PAGE 27
6 1. The operatorindependent par ticle creation procedure. 2. The compatibility equations and their application to the particle creation procedure using the modified ORC method. 3. Verification of the proposed method usi ng an angular particle as well as a round particle. 1.3.5 Chapter five: effect of particle shape a nd angularity on the dilation of granular soils This chapter includes the numerical study performed to explore the dependency of the dilatancy of granular soil s on particle shapes and angu larity. The simulations were performed using PFC2D together with the ORC method to accurately model the angular shapes. The shape and angularity factors introduced by Sukumaran and Ashmawy (2001) are used to define the partic le shape and angularity groups. The chapter is organized as: 1. A background section that contains a br ief review of dilatancy for granular soils, dilatancy and DEM, and pa rticle shape characterization. 2. Modification of the particle shap e and angularity factors proposed by Sukumaran and Ashmawy (2001). 3. Separation of particle groups with resp ect to shape and angularity factors. 4. Numerical simulations for pure shea r testing in dr ained conditions (material properties, boundary conditions, and loading procedure). 5. Results for different groups of shape and angularity. 6. Control of interparticle friction angle to account for dilatancy.
PAGE 28
7 1.3.6 Chapter six: conclusions and recommendations In this chapter, conclusions are in troduced and recommendations for further future studies within the same context are stated. More ideas that can be considered an extension of the study and need to be investigated are also presented. 1.3.7 Appendix A: theory and background of DEM In this appendix, the detailed theory for the discrete element scheme developed by Cundall and Strack (1979) is introduced. The Calculation cycle, forcedisplacement law, motion, damping, energy dissipation, ti me step, input parameters, and contact constitutive models are explained briefly.
PAGE 29
8 CHAPTER TWO: LITERATURE REVIEW Numerical modeling, supported with an accu rate simulation of test conditions, soil parameters, constitutive models (FEM), and contact models (DEM), is a powerful tool to study physical phenom ena. Numerical modeling pr ovides the capability of performing extensive parametric studies because a large number of runs can be performed in a relatively short period of tim e. On the other hand, inaccuracies in determining soil parameters can affect the numerical modeling results dramatically. Experiments, however, reflect the actual beha vior of the material provided that all the field conditions have been accurately simula ted. The time consuming cost associated with building experimental setups are the main obstacles for the experimental work. Studying a certain phenomenon experimentally and confirming the results numerically may be the best approach and can reliably capture the response of the system with minimal cost. 2.1 Discrete element method for modeling granular media The particulate nature of the granular material usually governs the behavior of these materials. Pressuredependent shear strength and stiffness, dilatancy, pressure history, and continuously nonlinear stressstr ain response are essential properties when studying the behavior of granular media.
PAGE 30
9 The Discrete Element Method (DEM) for modeli ng granular soils was first introduced by Cundall and Strack (1979). The material is modeled as a random assembly of discrete elements, which is a better representation of the particulate nature of gr anular soils. Each particle interacts with nei ghboring particles through contact forces that can be built up and broken at any time. Diffe rent types of bonds can be applied between particles to simulate phenomena such as cementation. The particle assemb ly is governed by boundary and loading conditions. The DEM so lves the dynamic equilibrium equations for each element subjected to either body or boundary forces. The method is capable of analyzing multiple interacting deformable continuous, discontinuous or fracturing bodies undergoing large displacements and rotations. In the scheme developed by Cundall and Strack (1979), the interaction between particle s is monitored, contact by contact, and the instantaneous motion of the particles is co mputed accordingly. Deformation of the particles is not allowed. Inst ead, the particles are capable of "overlapping" at contact points as an alternative method to mode ling individual particle deformation. In DEM, the equilibrium contact forces and displacements are determined through a series of calculations tracing the moveme nt of the individual particles. These movements are the results of the propagati on, through the medium, of disturbances originating at the boundaries (Cundall and Strack, 1979). Th e discrete element method, in the version proposed be C undall and Strack (1979), is ba sed on the idea that a small enough time step should be chosen to ensure that, during a single tim e step, disturbances do not propagate from any disc further than its immediate neighbors. The detailed theory
PAGE 31
10 including motion, contact force calculations, kinematics, and damping is presented in detail in appendix A. 2.1.1 DEM versus FEM for modeling granular soils Cundall (2002) argued that Â“the case is made that continuum methods such as Finite Element Method (FEM) for rock and so il may be completely replaced by particle models such as DEM in 1020 yearsÂ”. While complex constitutive models and many parameters are required to model the complex behavior of a soil matrix in FEM, only few parameters are needed in the DEM because the behavior complexity arises as a manifestation of the discrete nature of th e system. Continuously nonlinear stress/strain behavior, dilation related to mean stress, transition from brittle to ductile behavior, hysteresis and memory, and non linear stress envelopes, appear automatically in a particle model (Cundall, 2002). The DEM mode l exhibits localization which is difficult to capture in a continuum mode l that uses a continuous mesh. 2.1.2 Earlier published work in DEM After the pioneer work of Cundall and Strack (1979) in developing the first discrete element scheme that can be used to model granular soils, researchers have started evaluating and improving the tec hnique. Some of the correspondi ng fields of study are as follows: 1. Fundamental investigation and applicat ion of the DEM in granular soils, cohesive soils, and powders. 2. Rock mechanics.
PAGE 32
11 3. Experimental and theoretical verification of the DEM. 4. Modeling different particle shapes. 5. Developing different contact models. 6. Coupled modeling methods. 7. Largescale and industr ial applications. Dobry and Ng (1992) presented an extensiv e survey of the rese arch studies that had been performed on the applications of the discrete element method from 1982 to 1992. The survey included 42 references in DEM applications. The results of the literature survey were retabulat ed and shown in tables (21) and (22). They developed a program called Â“CONBAL2Â” (CONTACT +T RUBAL in 2D) which is based on TRUBAL created by Strack an d Cundall (1978). Cyclic st raincontrolled loadings at constant volume on an isotropically conso lidated random assembly of grains were performed to simulate undrained loading, an d the results were comparable to these obtained by Seed and Idriss (1984). 2.1.2.1 Fundamental investigations: monotonic and cyclic shear testing Many research studies were performed to investigate the fundamental granular soil behavior using the DEM, especially in mo notonic and cyclic shear tests. Ni et al. (2000) studied the effects of the microproperties of granul ar material on its shear strength and shear stressstrain behavior. Three dimensional DEM simulations of the direct shear test were perfor med. Granular particles were simulated by means of pairs of bonded balls. The study revealed that: (i) the mobilized friction angle decreases with the increase of normal effective stress for both p eak and residual values, (ii) the particle
PAGE 33
Table 21. Literature survey of DEM up to 1992 [ Source: Dobry and Ng (1996)] Bathurst, R.J. and Rothenburg, L. (1988) Bathurst, R.J. and Rothenburg, L. (1989) Cundall, P.A. (1971) Cundall, P.A. (1974) Cundall, P.A. (1976) Cundall, P.A. and Strack, O.D.(1979) Â“GeotechniqueÂ” Cundall, P.A. and Strack, O.D.(1979) 3rd Num. Meth. Geomech Cundall, P.A. and Strack, O.D.(1979) Â“Report to NSFÂ” Cundall, P.A. and Strack, O.D.(1983) Cundall, P.A. (1988) Cundall, P.A. (1989) Huang, A.B. and Lee, J.S. (1989) IshibashiI., Agrarwal, T., and Ashraf, S.A. (1989) Iwashita, K., Tarumi, V., Casaverdi, L., Vermura, D., Meguro, K. and Hakuno, M. (1988) Iwashita, K. and Hakuno, M. (1989) Kishino, Y. (1988) Kiyama, H. anc Fujimura, H. (1983) Kuhn, M.R. and Mitchell, J.K. (1989) Ng, T. and Dobry, R. (1988) Ng, T. (1989) Oner, M. (1984) Petrakis, E., Debry, R., and Ng, T. (1989)23910111213141516172628293032333438394143 Periodic Boundar y 2D 3D Pol yg o n Rounded One size Two sizes Three sizes Four to ei g h t ?100 101 1000 1001 8500 Linear Lin. Pr. De p Sim p Mindlin Com. Mindlin S p ecial sch. Newton 2nd Not allowed Allowed D y namic Static C y clic Monotonic Soil Mech. Rock Mech. Grain flo w En g Proble m Type of problem Particle rotation Type of equilibriu m Type of loading Criterion Model size Contact law Relaxation scheme Boundaries Dimensions Particle shape Particle size distribution 12 <
PAGE 34
Table 21. Continued Petrakis, E. and Debry, R. (1989) Strack. O.D, and Cundall, P.A. (1978) Strack. O.D, and Cundall, P.A. (1984) Tarumi, Y. and Hakuno M. (1985) Thronton, C. and Randall, C.W. (1988) Ting, J.M. (1986) Ting, J.M. (1987) Ting, J.M. and Corkum, B.T.(1988) Â“1Â” Ting, J.M. and Corkum, B.T.(1988) Â”2Â” Ting, J.M. and Corkum, B.T.(1988) Â“3Â” Ti ng, J M ., C or k um, B T ., and Kauffman, C.R.(1989) Ushida, Y. and Hakuno, M. (1989) Uemura D. and Hakuno, M. (1987) Uemura D. and Hakuno, M. (1989) Walton, O.R. (1982) Walton, O.R. (1983) Walton, O.R. and Braun, R.L. (1988) Wu, W.Y. and Liu, Z.D.(1989) Yamamoto, T. and Hakuno, M. (1989) Zhang, Y. and Cundall, P.A. (1986)4451525354555657585960616263646566697071 Periodic Boundar y 2D 3D Pol yg on Rounded One size Two sizes Three sizes Four to ei g ht ?100 101 1000 1001 8500 Linear Lin. Pr. De p Sim p Mindlin Com. Mindlin S p ecial sch. Newton 2nd Not allowed Allowed D y namic Static C y clic Monotonic Soil Mech. Rock Mech. Grain flo w En g Proble m Particle rotation Type of equilibriu m Type of loading Type of problem Particle size distribution Model size Contact law Relaxation scheme Criterion Boundaries Dimensions Particle shape 13 <
PAGE 35
14 shape plays an important role in the particle movement, (iii) the particle size relative to the shear box had an effect on the residual bul k friction angle and the volumetric dilation, and (iv) the particle friction c ontributes significantly to the shear strength as the sample starts to dilate. Anandarajah (2000) modele d the behavior of kaolinite clay during onedimensional compression numerically using tw o dimensional DEM. Three interparticle forces are considered in the analysis: (i) the mechanical forces, (ii) the doublelayer repulsive forces, and (iii) Va n der Waals attractive forces. The numerical results were close to the behavior obse rved in the laboratory. Liu and Sun (2002) simulated a biaxial compression test on an assembly of circular rigid particles using the DEM by incorporating interparticle adhesive forces. The collapse behavior during isotropic compre ssion and biaxial shear were simulated by releasing the interparticle adhesive forces from the initial values to zero at constant mean compression pressure and constant stress rati o, respectively. Qualitatively similar results between the numerically simulated test and the triaxial test on unsaturated compacted clays were obtained. Mirghasemi et al (2002) reported the re sults of isotropic lo adingunloading and biaxial shear simulations on assemblies of polygonshaped particles with different particle angularity. The linear contact law was adopted. The average number of contacts in the assembly was found to be increasing wi th the increase of angularity. Also it was found that when the number of particle edges increases, the amount of rotation decreases. The assembly of angular particles was found to be more compressible. When more edgetoedge contacts are established, less voids space between particles remains and a more
PAGE 36
15 compacted assembly is then gained. The generated assemblies were subjected to a series of numerical biaxial loadingunloading simula tions and it was observed that in addition to initial density and coordination number, there are other factors, such as angularity, that influence the shear strength of angular mate rials. In the macroscopic level, they concluded that when the angul arity decreases, the maximum sh ear resistance occurs at lower axial strain. At large strains, the angle of intern al friction was found to be a function of angularity. Duri ng biaxial simulation, the round particles reached a steady state at lower values of strains than the a ngular particles. The maximum volume change of the angular particles was found to be greater than th at of the round particles. Sitharam et al (2002) presented numer ical simulation results of isotropic compression and triaxial static shear tests under drained and undrai ned stress paths on polydisperse assembly of loose and dens e spheres. The numerical simulation has periodic variation in space. They conclude d that the specific volumelogarithmic mean principal stress relationship is nonlinear a nd stress dependent. The compression paths were parallel at low stresses and tend to converge at high stresses. The initial arrangement of particles had an important influence on the behavior during shear tests. In drained tests, the initial state was found to have an effect on the presteadystate behavior of loose and dense assemblies at constant c onfining pressure, while at steady state their behavior was unique in terms of macroscopi c strength of the assembly. In undrained tests, the initial state influenced the behavior of the assembly.
PAGE 37
16 2.1.2.2 Contact models is DEM The contact model that governs the interp retable forces has received extensive interest in research studies. The two most common contact models utilized in the DEM are the linear elastic contact model and the He rtzMindlin contact model. Formulations for both models are explai ned in Appendix A. Mind lin and Deresiewicz (1953) introduced an investigation of the phenomena occurring at the contact surfaces of elastic spheres subjected to a variety of applied no rmal and tangential fo rces. Hertz (1951) considered only the forces normal to th e contact surfaces. After Hertz, some investigators studied the effect of forces tange ntial to the contact surface and revealed the necessity of taking slippage into account. Voyiad jis et al. (1992) proposed a formulation for an anisotropic distortional yield model fo r granular media. The model applied the numerical solution for the HertzMindlin prob lem of the contact of two elastic, rough spheres subjected to oblique compression. A se ries of hollow cylinder cyclic triaxial laboratory tests, in which both axial and to rsional loads were applied to glass bead specimens was performed and the results were confirmed using threedimensional numerical simulations performe d along the same stress path. Zhang and Whiten (1996) stud ied the different methods us ed to calculate contact forces between particles using spring and da mping models and unrea listic behavior was found for most of these methods. The nonline ar formula of Tsuji (1992) along with the particle separating when the force returns to zero rather than when the distance between the centers exceeds the sum of the radii was found to give realistic results. Zhang and Whiten (1999) presented a new calculation me thod for particle motion in the tangential
PAGE 38
17 direction in discrete element simulations. They stated that the dynamic friction model used for tangential calculations should incl ude sliding and rolling friction models. A simple coefficient of rolling friction was fo und to be insufficient for the purpose of determination of linear and a ngular motion in discrete el ement simulation. Zhang and Whiten (1999) presented a method for the ta ngential friction calc ulations. Detailed tangential calculations, sliding calculations, rolling calculations, updating particle motion, and the change of friction mode were intr oduced. The tangential calculation can be easily applied in conjunction with the e fficient normal calculation method proposed by Zhang and Whiten (1996). The simulation re sults appear to correspond well with the known properties of friction. 2.1.2.3 Coupled DEM and FEM Coupled FEM and DEM simulations have be en conducted such as that introduced by Barbosa and Ghaboussi (1992). They in troduced a method that can be used for modeling multiple interacting deformable bodies. The method was based on the assumption that deformation modes can be decoupled from rigid body motion. In this method, each deformable body was treated as an individual discrete unit, which was idealized by a finite element model. The discrete units may under go large displacements and rotations. The interactio ns between these units occurred through contact forces that were continually updated.
PAGE 39
18 2.1.2.4 Experimental veri fication of DEM Experimental verification of the ability of the DEM to capture the behavior of the granular soils has been of in terest for many years. Hryciw et al (1997) presented an automated video tracking and digital image anal ysis system that was developed to obtain soil particle displacement fields and velocities from smallsc ale laboratory experiments. A threedimensional simulation of a laborat ory plowing experiment was performed in which a onetoone correspondence was achieve d between the number of particles and their size distribution in simulation and physic al experiment. The video tracking data were used to evaluate the numerical simula tion in modeling the experiment kinematics. Good agreement was found between the experi mental and the simulated deformation patterns. Ng and Changming (2001) presented a co mparison between a granular material studied experimentally in simple shear a nd numerical modeling using DEM. Simple shear tests were conducted inside the magne tic core of magnetic resonance imaging (MRI) equipment. Spherical pha rmaceutical pills of were used as the granular material. The center of each pill was determined by the MRI. An experimental study was performed to determine the pharmaceutical pills properties required for the numerical modeling. Good agreement was found at both levels. Excellent agreement was found up to 10% shear strain. Beyond that, the peak shear points in the numerical results were slightly lower than those in the experimental results.
PAGE 40
19 2.1.2.5 Other applications Mcdowell and Harireche (2002) presented a DEM simulation for the behavior of sands in onedimensional compression. The yield stress was taken to correspond to the point of maximum curvature on th e plot of void ratio versus logarithm of stress. A single soil grain was modeled as an agglomerate of balls bounded together. Balls have been removed randomly from each agglomerate to si mulate flaws. The yield stress was found to be reducing with the increase of agglom erate size. The values of the yield stress predicted by the model were found to be lowe r than those of the experiments performed on silica sand. The authors attr ibuted this deviation to th e difference in shape between the silica sand particles and the agglomerates. 2.1.3 Modeling noncircular particles In the DEM scheme proposed by Cundall a nd Strack (1979), the granular particles were modeled as discs in 2D simulations and as spheres in 3D simulations (fig. 21). The resistance to rotation is much less for circul ar particles compared to that of the actual particles. Circular particles have an inherent tendency to roll. The macroscopic angle of shearing resistance for circular particles modeled using DEM is much less than that of the actual irregular particle assembly. Howeve r, the contact detection for the circular particles does not need complicated algorithms and also the direction of the contact normal force is always toward the center. As such, the computed normal contact forces never contribute to the moment acting on a pa rticle as each normal contact force always acts in the direction of the pa rticle centroid. Different particle shapes rather than
PAGE 41
20 discs/spheres were proposed, in many research studies, to be used in DEM simulations in order to improve the numerical simulation. Polygon particles were presented by Bar bosa and Ghaboussi (1992) and Matuttis et al (2000). Figure (22) shows a sample of polygonshaped particles. A more realistic representation of the particle behavior can be achieved using polygonshaped particles; however the method is comput ationally intensive. Ting et al (1993) realized th at ellipseshaped particles have fewer tendencies to rotate and also have a unique outward normal. They developed a numerical algorithm for the DEM using twodimensional ellipseshaped part icles (fig. 23). A robust algorithm for computing particletoparticle and particletowall contacts was derived. The results of the validation tests indicated that the ellipsebased DEM resulted in mechanical behavior that was similar to that of real so ils. Ellipseshaped particles, however, do not accurately represent the particle shape. Ng (1994) developed a new program called ELLIPSE2. The new program was used to explai n the stressstrain behavior of granular soils under monotonic drained compressive load ing and cyclic constant volume loading at various shear strains. A 2D random a ssembly of elastic rough quartz particles was used to simulate the granular soil particle s. Good agreement was observed between the numerical and the experimental results. The results confirmed that DEM simulations using an assembly of circular and elliptical particles give a better representation of the behavior of granular soils. Ting et al (1995) studied the influence of the contact of twodimensional ellipsebased particle on their interaction. The samples were isotropically compressed and then sheared in biaxial compression. In order to assess the relative
PAGE 42
21 Figure 21. Two circular particles in contact [Source: Cundall and Strack (1979)] Figure 22. Polygonshaped particles [Source: Matuttis et al. (2000)]
PAGE 43
22 importance of rolling and sliding, the contact deformation were separated into portions: (i) due to individual particle rotation, and (ii) due to particle translation. This decomposition demonstrated that particle rotation accounts for twice as much contact motion for round particles as particle transl ation. Lin and Ng (1997) developed a new 3Dellipsoidbased DEM code called ELLIPSE 3D. A numerical study on the mechanical behavior of monosized particle arrays using the developed ELLIPSE3D program was performed. Higher shear strengt h, larger initial modulus, more dilation, and less particle (a) (b) (c) Figure 23. Ellipseshaped particles: (a) particle assembly, (b) two ellipses in contact, (c) contact between an ellipseshaped particle and a wall [Source: Ting et al. (1993)]
PAGE 44
23 rotation were observed with the ellipsoid asse mbly during the triaxial test. The results demonstrated that using nonspherical particle s in discrete element modeling is essential to improve the simulations of granular materials. Potapov and Campbell (1998) proposed a comp utationally efficient model for the discrete element simulation of nonround particles. The idea is to approximate the ellipse by an oval shape whose boundary is determined by four circular arches of two different radii that are joined together in a continuous way (fig. 24). More complex shapes can be reproduced by changing the rad ii of the arches. The idea is computationally efficient because the determination of the contact and overlap between two circular arch segments is very similar to that of circles in the orig inal DEM code. The appl ication of the idea in 3D was not presented. A direct test of the methodÂ’s performance indi cated that it is very effective and less than two times slower than the corresponding circular models. Favier et al (1999) presented a method to model axisymmetrical particles as multisphere discrete elements. The particles are regenerated using overlapping spheres with fixed rigidity with resp ect to local coordinate system (fig. 25). The method can theoretically model any shape, however highly angular pa rticles cannot regenerated properly. Favier et al. ( 2001) used the proposed method to model the discharge of ellipseshaped particles through an orifice in a flat bottom h opper. The simulations were compared to the experimental results at th e same scale. A good agreement was noticed between the flow behavior pf the simulated and physical particle assemblies for all orifice sizes.
PAGE 45
24 Figure 24. Oval shaped particles [Source: Potapov and Campbell (2000)] Figure 25. Axisymmetrical particles: (a) sample particles, (b) contact between two multielement particles [Source: Favier et al. (1999)] (a) (b)
PAGE 46
25 2.1.3.1 Modeling angular particles as rigid clusters Jensen et al (1999) introduced the clus tering technique by comb ining a number of sphericalshaped particles in a semirigid configuration, whic h is a better representation of natural soil particles (fig. 26). Any number of particles to form a cluster as long as they rotate and translate as a rigid body. Th e interparticle contacts within the cluster are constrained to be linearelastic with a high stiffness. Thom as and Bray (1999) presented a different idea for clustering by imposing ki nematics restrictions on discs within a cluster to prevent relative translations and rotations. Although these techniques represent some improvement over earlier methods, the simulated particle outlines did not resemble those of actual particles. In addition, a substantial increase in computational time was incurred due to the high contact stif fness required within the cluster. 2.1.3.2 Overlapping rigid cl usters (ORC) technique Ashmawy et al (2003) proposed the use of Overlapping Rigid Clusters (ORC) technique to accurately simulate angular part icle shapes. A set of subroutines were introduced to PFC2D using ItascaÂ’s softwarespecifi c programming language, Fish. The builtin clump logic command was instrumented in formulating the new method. The clump logic allows temporary or permanent rigid bonding between several disk elements without detecting contacts or calculating contact forces be tween elements belonging to the same clump, which increases computati onal efficiency. The ORC method relies on clumping number of overlapping discs such th at the resulting outline coincides with and is almost identical to the actual particleÂ’s outline. First, twodi mensional outlines of a series of particles are obtained using a digital microscope or s canner. Overlapping circles
PAGE 47
26 Figure 26. (a) Outline of sand particle, (b) DEM disc element superimposed over sand particle, (c) a number of DEM disc par ticles are joined together is a rigid conifiguration (Cluster), (d) several possible combination of discs to form clusters [Source: Jensen et al. (1999)] Figure 27. Disc elements inscribed within a particle outline to capture the shape [Source: Ashmawy et al. (2003)]
PAGE 48
27 are then inscribed within the outline to captur e the shape, as shown in fig. (27). The number of overlapping circles that can accurate ly represent the actual particles depends on: 1. The degree of nonuniformity in the orig inal particle sh ape and angularity. 2. The desired level of geometric accuracy. 3. The required computation time limit. Typically, ten to fifteen disks will adequately capture the Â“trueÂ” shape of a particle. Due to overlapping, the density of each disk must be scaled to ensure that the mass of the particle remains pr oportional to the area, regardless of the nu mber of discs or level of overlap. An approximate method wa s proposed to scale the density as follows: p d p dA A (21) where d is the density of the disk elements, Ap is the area of the particle, Ad is the sum of the areas of the disk elements, and p is the density of the soil particle. Equation (21) ensures that the mass of the particle is pr oportional to the area but it neither guarantees that the moment of inertia nor the center of ma ss of the ORC particle is identical to that of the actual particle. In situations where dynamic moment s and rotational accelerations are significant, modifications n eed to be made by either adjusting the relative mass of the disk elements, or superimposing inertia bala ncing elements at specific points. The ORC method is implemented within PFC2D by means of a series of Fi sh functions that convert a particle assembly of discs into their corresponding angular particles as follows:
PAGE 49
28 1. Particles are automatically generated as discrete circular elements using any of the builtin particle s generation techniques. 2. The shape conversion algorithm is i nvoked to transform each circular particle into an angular equiva lent through repl acement with a corresponding set of disc elements, se lected randomly fr om a library of available particle shapes. 3. A random rotation between 0 and 360 is applied to each transformed particle to avoid preferre d orientations and subseque nt system anisotropy. Figure (28) shows a random assembly of circular particles that have been transformed to the corresponding angular shapes using the ORC method. In order for the system to remain computationally stable th roughout the solution, it is necessary that the computational time step, t, remains below a critical value of 2/ max, where max is the highest eigen frequency of the system. Because max changes every time the particle positions change, the eigen frequencies must be recalculated every time step by solving the global stiffness matrix, which is highly inef ficient. Instead, an approximate critical Figure 28. Random assembly of eight circ ular particles (left) and the transformed equivalent particles (right) [Source: Ashmawy et al. (2003)]
PAGE 50
29 time step, tc, is calculated from: i i ck / m min t (22) where mi and ki are the mass and equivalent stiffne ss of disk element i, respectively (Itasca, 1999). To ensure stability, the actual time step during computation is taken to be 0.8 tc. Because each angular particle consists of a clump of smaller disk elements with scaled (smaller) densities, the critical time st ep calculated from Eq. (22) is highly overconservative, and the corresponding computa tion time is unnecessarily long. To rectify this problem, the critical time step is adjust ed to correspond to the mass and stiffness of the clump, not the individual disk elements. The efficiency of the solution is thus preserved without compromising the nu merical stability of the system. To verify the ORC technique, Ashmawy et al. (2003) created particle assemblies of various degree of angularity. The assemb lies were subjected to a simulated undrained cyclic shear conditions to asses their lique faction susceptibility. The numerical model was able to capture the soil behavior closely. The findings indicated that at the maximum void ratio, the susceptibility to liquefaction is independent of voids ratio. The influence of particle morphology on liquefaction suscepti bility was found to be significant in the case of sands prepared at the same void ratio. 2.2 Dilatancy in granular soils The first referenced dilatancy was made by Professor Osborne Reynolds (1885). He stated, Â“Without attempting anything lik e a complete dynamical theory, which will require a large development of mathematics, I would point out the exis tence of a singular
PAGE 51
30 fundamental property of such granular medi a which is not possessed by known fluids or solids. I have called this unique property of granular masses Â“DilatancyÂ”, because the property consists in a definite change of bulk, consequent on a definite change of shape or distortional strain, any di sturbance whatever causing a change of volume and general dilationÂ”. Reynolds (1885) s howed that dense sands expa nd at failure whereas loose sands contract during shear to failure. This proves that part icle movements during deformation and failure are not necessarily in the direction of applie d shearing stress. MohrCoulomb failure criterion is the mo st commonly used failure criterion in geotechnical engineering. Closer investigatio ns revealed that soil behavior is more complex than that described by that model (Houlsby 1991). Terzaghi (1920) stated that Â“the fundamental assumption of the traditiona l earth pressure theo ries cannot, in fact, stand even superficia l examination. The fundamental error was introduced by Coulomb, who purposely ignored the fact that sand consis ts of individual grains and who dealt with the sand as if it were a homogeneous mass with certain mechanical properties. CoulombÂ’s idea proved very useful as a working hypothesis for the solution of one special problem of the earthpressure theor y, but it developed into an obstacle against further progress as soon as its hypothetical character came to be forgotten by CoulombÂ’s successorsÂ”. The stress stain diagram for dense sand exhibits a peak followed by a reduction in the stress at large strain, whereas for loose sa nds, no peak is usually observed (fig. 29). The volumetric strain correspond ing to shearing of the dense sand shows that dilation occurs as the test proceeds after a small initial compression. The magnitude of this
PAGE 52
31 dilation depends on the soil de nsity, confining pressure, gr ain size distribution, and particle shape, as described later in details. If tests of the same soil under the same density were carried under different levels of confining pressures, the corresponding dilations will be different. Th e result is that the failure en velope is curved in the MohrCoulomb plot, with the peak angle approaching the large strain fricti on angle at very high stress levels (fig. 210). Unde rstanding the rule of dilation is essential to explain the peak and large strain angles and the reduction of peak stre ngth with stress level. 2.2.1 Friction angle, dilation angle, and dilatancy The angle of friction ( \) can be expressed as the ratio of the shear stress to the stress normal to shearing plane as follows: \ 3 \ 1 \ 3 \ 1 \sin (23) On the other hand, the dilation angle, can be expressed as the ratio between a volumetric strain rate and a sh ear strain rate. For the case of plane strain, it can be defined as: 3 1 3 1sin (24)
PAGE 53
32 Figure 29. Dilation of dense sand during simple shear test [Source: Houlsby (1991)] Figure 210. Curved peak strength envelope on MohrCoulomb plot [Source: Houlsby (1991)]
PAGE 54
33 The angle of dilation is positive when the soil dilates and is negative when the soil contracts. The graphical representations of th e angle of friction and the angle of dilation are shown in fig. (211). The angle of dilation should strictly be defined in terms of the plastic components of the strain rates. For mo st soils, the elastic s tiffness is sufficiently high such that the elastic strains are much sm aller than the plastic strains and the total strain rate can be used as in the previous e quation. Dilatancy (d), however, is defined as the ratio of plastic volumetric strain increment to the plastic deviator strain increment as follows: sin d d d3 1 3 1 v (25) 2.2.1.1 Relation between friction and dilation angles If we consider a frictional block that slides on top of a flat plan e that has an angle of friction of \ cv, the ratio between the shear and the normal stress can be expressed as: \ cv \ ntan (26) where, \ cv is the critical state shearing resistance angle (constant volume). The simplest way to understand the relation between the fricti on and the dilation angles is to consider the sawtooth model shown in fig. (212). The teeth angle wi th respect to horizontal is whereas, the same angle of friction acts on the teeth of the saw. The ratio between the normal and the shear stresses and an expression fo r the friction, dilation, and critical state angles can easily be derived using simple statics as:
PAGE 55
34 Figure 211. Angles of friction and dilation [Source: Houlsby (1991)] Figure 212. The sawtoot h model for dilatancy [Source: Houlsby (1991)]
PAGE 56
35 \ cv \ \ ntan tan (27) \ cv \ (28) This type of relation is usually call ed Â“Flow RuleÂ” (Houlsby 1991). Another approach to correlate the friction and the di lation angles is that proposed by Taylor (1948), who suggested an ener gy correction to account fo r dilation. All friction relationships were expressed in terms of ener gy dissipation. In a simple shear setup, the dissipated work can be assumed to be proportional to both norm al and shear stresses. His results can be expressed as: tan tan tan\ cv \ (29) This relation is slightly different from the sawtooth model but still reveals the same trend. Another form is that pro posed by Taylor, followed by Skempton and Bishop (1950) separated the strength component due to friction from that due to expansion / dilation in the direct shear test as follows: Vr maxtan tan (210) where r is the residual shearing resistance angle, max is the angle of friction obtained by fitting CoulombÂ’s equation for direct shear test with c = 0, to the observed peak stress ratio, V is volume expansion per unit area, and is the equivalent shear displacement. Bishop (1954) developed a similar expression for the triaxial test:
PAGE 57
36 1 3 max \ 3 \ 1 r 2V V 2 1 45 tan (211) where 1 \ and 2 \ are the major and minor effective stresses, V / V is the rate of the unit volumetric change, and 1 is the rate of the major principal strain change. Newland and Allely (1957) considered the resultant di rection of movement during dilation and determined expressions for the reduced value of the angle of shearing resistance corrected for energy due to expansion, denoted by f, as follows: f max (212) 1V V tan (213) 1 max \ 3 \ 1 1 max \ 3 \ 1V V 1 V V tan (214) The value of the angle f is different from r. The difference is found to arise from the fact that the sample has also to expand against the compone nt of the reaction on the plane of contact resolved in the vert ical reaction, so it should be noted that f < r. A very similar approach was used in the deve lopment of the CamCl ay model by Schofield and Worth (1968). The final expression is very similar to TaylorÂ’s expression and can be summarized as: 2 sin 3 M sin 3 sin 6\ \ (215)
PAGE 58
37 2.2.1.2 True and critical state friction angles The true angle of friction between the mineral surfaces of the particles is a physical property that depends on the nature of the mineral, the properties of its surface roughness, and on the size of the load per particle (Bowden and Tabor 1954). The true friction angle is less than the constant vol ume shearing resistance angle (Rowe 1962). Caquot (1934) integrated force acting in two di rections at right angl es over the surface of a sphere on the assumption that in the ultim ate state, sliding o ccurs simultaneously on contact forces inclined in all the tangential di rections of a spherical surface. He derived the following expression: tan 2 1 tancv (216) Bishop (1954) also developed an approximate solution in the form: tan 3 10 tan 15 tancv (217) For a given material such as quartz, the angle of true fricti on depends on the size of the load per particle, which can vary by a factor of 104 between pebbles and siltsized particles, for a given average effect ive stress (Rowe 1962). As a result\ cv for quartz in water varies from 31 maximum for silt to 22 for pebbles. Proctor and Barton (1974) and Frossard (1979) reported that the interpar ticle fiction angle of saturated irregular quartz particles varies between 5 and 40. Bolton (1986) assumed typical value for the
PAGE 59
38 critical state shearing resistance angle for qua rtz and felspathic sands as 33 and 40, respectively. 2.2.2 Development of the theory of dilatancy 2.2.2.1 RoweÂ’s theory Rowe (1962) studied, both theoretically and experimentally, the properties of regular assemblies of spheres and rods and was able to relate the stress ratio to the strain rate ratio in terms of the geometry of packi ng. The materials used were steel, glass, and quartz. His findings can be summarized as follows: 1. Whatever the geometrical arrangement of the solids, th e stress ratio, at the peak strength and during the subsequent states of deformation, will follow the following rule: tan tan\ 3 \ 1 (218) Where; is the deviation angle at the cont act points from the major principal direction, can be seen as a geometrical pr operty of the packing (fig. 213), is the true angle of friction between the mineral surfaces, 1 \ and 3 \ are the major and minor effective principal stresses. 2. The energy ratio for a fixed orientati on of particle movement is given by the expression: tan tan E (219)
PAGE 60
39 3. Slippage occurs well past the peak when the stress ratio reaches unity. The slip plane observed in the laboratory and in the fiel d is not the cause of fa ilure at the peak, as assumed by Coulomb, but is a later result of failure. It is the plane at the instant of collapse and has nothing to do with the state of the particles at the peak strength. Rowe then applied the principle of least work to a random mass of irregular particles and established relations between the rate of dilatancy and the stress ratio for any given particle arrangement as follows: f 2 1 \ 3 \ 12 1 45 tan V V 1 (220) Figure 213. Definitions of RoweÂ’s geometry properties [Source: Rowe (1962)]
PAGE 61
40 where tan( f) may be written as k tan( ) where k is a coefficien t which increases with the degree of remolding and associated energy loss. Eliminating tan2(45f/2) from equations (218) and (220) and solving for the results can be expressed as: 1 \ 3 \ 1V V d 1 tan (221) Rowe introduced fig. (214), wh ich shows different values of satisfying equation (221). As proposed by Rowe (1962) there is a unique relation between the stress ratio and dilatancy (equa tions (218), (220), and (221)). This issue was discussed later by many researchers and will be explained later in details. 2.2.2.2 BoltonÂ’s relative dilatancy index Bolton (1986) collected data of the st rength and dilatancy of 17 sands in axisymmetric or plane strain tests at different densities and confining pressure combinations to experimentally correlate the angle of friction to the angle of dilation. He introduced a new relative dilatancy index, IR, which can be defined in terms of relative density, ID, and effective stress level, P\, as: 1 LnP 10 I I\ D R (222) Bolton rewrote the relative dilatancy i ndex, after a discussion by Tatsuoka (1987) regarding the effect of anis otropy and the behavior under lo w confining pressures, as follows:
PAGE 62
41 Figure 214. Angle versus porosity for different materials [Source: Rowe (1962)] Figure 215. Deviation of triaxi al dilatancy from biaxial state [Source: Schanz and Verneer (1996)]
PAGE 63
42 1 150 P Ln 5 I I\ D R (In case of P\ > 150 kPa) (223) 1 I 5 ID R (In case of P\ < 150 kPa) (224) Bolton (1986) implied that the relative d ilatancy can be correlated to dilatancy as follows: R max 1 VI 3 0 d d (225) In case of plane strain the e xpression should be written as: R max \ cv \ maxI 5 8 0 (226) Whereas, is case of triaxial st rain condition it was found that: R \ cv \ maxI 3 (227) Schanz and Verneer (1996) applied Rowe Â’s expression in two perpendicular planes in order to redefine the dilation angle in case of triaxial stress loading condition. As shown in fig. (215), if sliding on a certain plane is governed by the stress ratio 1/ 2 (mechanism A), sliding on a perpendicular plane should be governed by 1/ 3 (mechanism B). The following expression for th e angle of dilatancy in the triaxial test was deduced: 1 V 1 V2 sin (228)
PAGE 64
43 The authors supported the idea that the a ngle of dilatancy doesnÂ’t depend on of the testing method; triaxial or plane strain. RoweÂ’s expre ssion could be rewritten as: tr f tr tr f trsin sin 1 sin sin sin (229) The authors stated Â“RoweÂ’s stress dila tancy theory exhibits an appealing relationship between the friction angle and the di latancy angle for planar deformation. In factcv ps f However, this theory needs to be applied to triaxial conditions of stress and strain in order to obtain a relationship between the fr iction angle and the dilatancy angle. For this reason, the relationships th at given by Bolton (1986), Eqns (223) through (227), are now consideredÂ”. They showed that for triaxial testing conditions, BoltonÂ’s relative dilatancy coefficient, IR, can be correlated to the dilatancy angle as: R RI 7 6 I sin (230) 2.2.3 Statedependent dilatancy Rowe (1962) showed, based on the theory of least rate of internal work, that dilatancy (d) could be expressed as a function of the stress ratio, and the true angle of friction between the mineral surfaces of the particles as: C d d (231) Where C is a set of intrinsic material constants and is the stress ratio. Li and Dafalias (2000) argued that the previous equa tion worked well for cohesive soils. For granular soils, the applicability of the equa tion was found to be dependent on the density.
PAGE 65
44 It should be noted that Rowe (1962) obs erved the divergence between the proposed theory and the experimental results and poi nted out that a variable depending on the density and the stress history s hould be added to the stressd ilatancy relationship that he had derived earlier. Rowe attr ibuted the divergence to the r earrangement of the particle packing, a fact that was ignored in the origin al version of the stressdilatancy theory. Li and Dafalias (2000) attributed the lack of de nsity and stress level dependency in RoweÂ’s work to the minimization that has been done during imposing the minimum internal work hypothesis. This minimization made the stre ss ratio uniquely related to dilatancy and independent of the particles packing which c ontradicted the exact analytical conclusion deduced by Li and Dafalias (2000) for two different packings. 2.2.3.1 Problems with unique relationship between d and Li and Dafalias (2000) clar ified the problems associated with assuming a unique relationship between the dila tancy and the stress ratio. According to the concept of critical state soil mechanics, density and conf ining pressure, together, define whether the sand is loose or dense. To prove the defi ciency of the unique relation between d and they considered two specimens of the same sand; one is in the dense condition whereas the other is the loose condition. Subjected to shear loading and beginning from the same the loose sample contracted while the dense sample dilated (fig. 216). If the assumption that (d) in a unique function of is valid, the direction of plastic flow, and the undrained stress path, would be uniquely related to irrespective of the material state, which contradicts fig. (216).
PAGE 66
45 Another contradiction involves the transformation state. During the first stages of shearing, the sample passes through a transformation state at which = Md (the critical state ratio at transformation state) and d=0 (Is hihara et al. 1975) and then approaches the critical state line at which = M (the critical state ratio ) and d=0 (fig. 217). Applying the RoweÂ’s dilatancy equation will reveal that Md = M, which render the phase transformation intrinsic, and re sults in a unique phase transfor mation line, for a particular sand, at which the sand would change its phase from contractive to dilative irrespective of stress level and density. Li and Dafalias (2000) stated that there is another c ontradiction when analyzing undrained tests on dense sands. It is know n that dense sands tested in undrained conditions show a stress path that c onverges with a slope of approximately = M towards the ultimate state, whic h is the critical sate (dp\ = dq = d v = 0 and d q 0). During this stage, the rate of plastic volumet ric change tends towards a constant value, implying that dilatancy tends towards a zero valu e. If the assumption that (d) is a unique function of is valid, (d) would be a c onstant along this path. In other words, as shearing proceeds, the confining pressure will increase due to the restriction on volumetric change and the critical state would never be reached. The above observations lead to the conclu sion that the Rowe Â’s stressdilatancy relation works well only when the change in the material internal state is minor and can be neglected (Li and Dafalias 2000).
PAGE 67
46 Figure 216. Variation of dilatancy with material state (data from Verdugo and Ishihara (1996). Undrained response of sand (a) with different densities (b) under different confining pressures [Source: Li and Dafalias (2000)] Figure 217. Variation in the phase transformation stress ratio with material state (data from Verdugo and Ishihara (1996) [ Source: Li and Dafalias ( 2000 )]
PAGE 68
47 2.2.4 General expression for state dependent dilatancy Kabilamany and Ishihara (1990) provided experimental evidence showing that (d+ ) increases with the increase of shear de formation. Manzari and Dafalias (1997) presented a sand model in which a linear dependence of the phase transformation on was introduced. Li (1997) inve stigated the response of sand at the ultimate stress ratio and pointed out that the dilatancy is not re lated only to the stress ratio but is also a function of the plastic volumetric strain. Wan and Guo (1998) proposed a model with a dilatancy function modified from RoweÂ’s st ressdilatancy equati on. Cubrinovski and Ishihara (1998) also showed a dilatancy re lationship that depends on the material state represented by cumulative plastic shear strain. Wan and Guo (1999) proposed a modified st ressdilatancy mode l to describe the mechanical behavior of granular materials at different stress levels and densities. The developed model was based on the original well known RoweÂ’s dilatancy expression and was proposed as follows: f m cr f cr msin sin e e 1 sin e e sin sin (232) Where is a parameter, e and ecr are the current and the critical state void ratios. The factor cre / e is a measurement of the deviation of the current voids ratio from the mean stress dependent critic al voids ratio, which will provide a mean by which the dilation process can be corrected to include density, stress le vel, and deformation history
PAGE 69
48 dependence. Wan and Guo (1999) presented a procedure for determining the factor for different types of typical na tural and manufactured sands and introduced the constitutive law for the proposed model. Numerical simu lations were conducted and indicated that most of the basic behavioral features of gra nular materials could be reasonably captured. The proposed stressdilatancy model introduced a family of energy curves as a function of void ratio and stress level instead of the unique energy line proposed by RoweÂ’s (1962). Li and Dafalias (2000) introduced a genera l state dependent dilatancy expression in the following general form: C Q e d d (233) Where Q and C denote internal state variab les other than the voids ratio e. The formulation should satisfy the requirement that ; the dilatancy must be zero at both the critical state and the phase transformation point, that is: 0 C Q e e M d dcr and 0 C Q e e M d dcr (234) The wellknown state parameter was intr oduced by Been and Jefferies (1985) and can be defined as: cre e (235) Where e is the current voids ratio and ecr in the critical state voids ratio in the ep\ plane corresponding to the current p\ (fig. 218). The state parameter is a measure of how far and which side the material is from the critical state. Li and Dafalias (2000) made use of the state parameter after redefining it as follows:
PAGE 70
49 a c crP P e e e e (236) Where e, c and are the material constants determ ining the critical state line in the eP\ space (Li and Wang 1998) and Pa is the atmospheric pressure. The general form of their relation can be written as: M e d dm 0 (237) where d0 and m are two positive modeling parameters. It should be noted that the CamClay dilatancy expression (d = M Â– ) is a special case of equation (237) when m = 0 and d0 = M. Yang and Li (2004) suggested a unique linear relationship between the peak friction angle ( p), and the maximum dilation angle ( max) as follows (fig. 219): max \ cs \ p28 0 (238) Comparing this relation to that proposed by Bolton ( 1986) for the plane strain condition, that is: max \ cs \ p8 0 (239) it can be concluded that the excess friction angl e in triaxial condition is less than 40% of that in plane strain condition, which is c onsistent with the experimental findings.
PAGE 71
50 Figure 218. The state parameter and critical state line [Source: Yang and Li (2004)] Figure 219. Relationship between peak friction angle and maximum dilation angle [Source: Yang and Li (2004)]
PAGE 72
51 2.3 Particle shape characterization Particle shape characteristics affect the be havior of granular soil. Interlocking resistance, dilatancy, and the overall shear st rength behavior can vary dramatically with particle shapes. Many res earch studies have been condu cted to quantify the shape characteristics of granular soils and rocks. The shape characteristics have a noticeable impact on the dilatancy, which falls within the sc ope of this research. A critical review of conventional and modern techniques to charact erize particle shape is introduced in the following section. The two terms commonly used to charac terize particle shape are form roundness (Sukumaran and Ashmawy 2001). Roundness usua lly describes specific aspects of external morphology such as variations at the corners. On the other hand, form usually captures the overall geometry of a particle a nd reflects variations in the proportions of the particle. Surface texture is superimposed on the corners and reflects the properties of particle surfaces be tween corners. 2.3.1 Particle shape characterization techniques Both form and roundness are best defined in 3D. However, due to difficulty of 3D measurements, the majority of the studie s described in the literature rely on 2D projections of particle s. Most conventional methods of quantifying 3D form are based on the shortest, longest, and intermediate ort hogonal axes through formulas based on the axes ratios. Some other meas ures of 3D shape require the measurement of the volume of the particle. WadellÂ’s (1932) roundness is one of the conventional approaches widely
PAGE 73
52 used in 2D, Fig. (220). A corner is define d as the part of the outline with a radius of curvature equal to or less than the radius of curvature of the maximum inscribed circle. Each corner of the maximum projection outline is measured by finding the largest circle that will fit. The average roundness is defi ned as the ratio of the average of these diameters to the diameter of the maximum inscribed circle. Fourier analysis has been used as a mo re advanced technique to characterize particle shape. The idea is to unroll the tw o dimensional outline of the particle and fit a Fourier series. Anstey and Delmet (1973), Eh rlich (1970), and Ehrlich et al. (1980) used this approach to characterize the particles. The particle shape can be fit as closely as possible by including a suffici ent number of harmonic terms in a linear summation. Lower harmonics reflect largescale features of the outline, whereas higher harmonics reflect smallscale features. Although many re searchers have used this technique and have come up with particular combinations of the harmonic terms to express roundness and angularity, none of these solutions has gained universal acceptance (Sukumaran, 1996). A difficulty with such a re presentation is that if the radius intersects the outline twice, the representation is not possible (S ebestyn, 1969). A drawback of any of the Fourier series approaches arises when trun cating the Fourier series to obtain a more detailed representation. If the truncated coe fficients are used to reconstruct the outline, the new outline may not be continuous and may cross itself (Sukumaran,1996). Piper (1970) observed that a plot of the distribution of angles between chord lengths could help to estimate particle r oundness. A plot of the cumulative distribution of angles between adjacent chords gives a distribution, the varian ce of which increases
PAGE 74
53 with roughness. The drawback of distributiona l summaries is that it throws away all information on the relationship between ad jacent points on the perimeter (Sukumaran, 1996). Mandelbrot (1967) illustrated the use of fractals to determ ine the total length between two points along a free form. Orford and Whalley (1983) used fractal dimension to quantify the morphology of irregularshape pa rticles. A disadvantage of this technique is that it does not give any information about the overall particle shape. Also large changes in roughness are reflected as small ch anges in the fractal dimension (Sukumaran, 1996). 2.3.2 The method introduced by Sukumaran and Ashmawy (2001) Sukumaran (1996) presented a new techni que for quantifying particle shape and angularity. The true shape of the particle was approximated by an equivalent polygon. The true particle shape was characterized by co mparing it to an ideal shape (circle). The basic definitions required to describe the method are: 1. Shape Factor is related to the deviat ion of global particle outline from a standard outline or datum. 2. Angularity Factor is a measure of the number and sharpness of the corners. 3. Radial segment is a straight line connecting the particle centroid to a sampling point on the perimeter. 4. Sampling interval is the angle betw een two adjacent radial segments.
PAGE 75
54 Figure 220. Illustration of WadellÂ’s pr ocedure for evaluating particle roundness [Source: Sukumaran (1996)] Figure 221. Ideal geometric shape used to define the shape and angularity factors [Source: Sukumaran (1996)]
PAGE 76
55 2.3.2.1 The distortion diagram Consider a two dimensional projection of a particle that has been discretized by equal sampling intervals (fig. 221). A circle is discretized with the same sampling interval. Considering the circle chords as a datum and aligning the radial segments of both the circle and the particle the angle between correspondi ng chords can be measured. The angle sign depends on the direct ion of the particle chord vect or relative to that of the circle. These angles are cal led the distortion angles ( ). If the distortion angles are plotted against the cumulative sampling interval, the resulting diagram is the Distortion Diagram (fig. 222). The distortion diagra m is a mapping technique from which the particle shape can be fully reconstr ucted (Sukumaran and Ashmawy, 2001). 2.3.2.2 Shape factor Sukumaran and Ashmawy (2001) defined the particle shape in terms of the deviation of the particle from a circle by di viding the sum of the absolute values of ( ) by that corresponding to an elongated flat particle as follows: % 100 45 N SFN 1 i iParticle (240) Where N is the number of sampling points. For all practical particle shapes, values of shape factor will lie between zero and one, with the former value corresponding to a circle, and the latter value co rresponding to a flat plate.
PAGE 77
56 2.3.2.3 Angularity factor The angularity of a particle was defined in terms of the number and sharpness of the corners and can be obtained by dividing the sum of the difference between the internal angle ( ) of a spherical particle as appr oximated by an Nsided polygon and the corresponding internal angle of the particle (fig. 215). Instead of summing up the absolute differences in ( ) angles, the sum of the squares of the differences in internal Figure 222. Distortion diagram for the particle shown [Source: Sukumaran (1996)]
PAGE 78
57 angles were taken to amplify the influence of sharper corners. Th e angularity factor can be expressed as: % 100 N / 360 ) 180 ( 3 N / 360 180 AF2 2 N 1 i 2 2 iParticle (241) The previous expression will give a valu e of (AF = 0) for the angularity of a sphere.
PAGE 79
58 CHAPTER THREE: EXPERIMENTAL VERFICATION OF MODELING ANGULAR PARTICLES USING DEM 3.1 Introduction Discrete element modeling includes the entir e set of methods used to simulate the mechanical interaction between discrete bodi es. Different formulations have been introduced for a wide range of engineering applications (e.g., Cundall and Strack 1979; Shi and Goodman 1989; Mustoe 1992). Differ ent discrete element methods utilize different numeric solution schemes, particle deformability and rotation, and contact detection. The Distinct/Discrete Element Method was introduced by Cundall and Strack (1979) to model granular assemblies within th e context of geotechnical engineering. The method relies on a timestepping explicit scheme, where equa tions of motion are solved over a finite time interval to obtain incr emental displacements and rotations. Forcedisplacement laws are then invoked to calculate the resulting interparticle forces, which are fed back into the equations of motion ove r the next time increment (Itasca, 1999). The solution proceeds until a static or dynami c equilibrium threshold is reached. The numerical solution scheme is based on the idea that a small enough time step should be chosen to ensure that, during a single time step, disturbances do not propagate from any disc further than its immediate neighbors (Cundall and Strack, 1979).
PAGE 80
59 Until recently, limitations in computer speed and power represented a major obstacle to the widespread implementation of DEM in pr actice. In recent y ears, DEM has gained considerable momentum as a powerful tool in studying microstructural phenomena in granular assemblies. The first comp rehensive DEM software package, PFC2D, was introduced in 1995 by Cundall and his coworke rs, to model twodimensional circular elements (disks) with linear or HertzMindl in contact models (Itasca, 1999). A threedimensional version of the soft ware was later released. 3.2 Background In the DEM scheme proposed by Cundall a nd Strack (1979), the granular particles were modeled as discs in 2D simulations and as spheres in 3D simulations. Although computationally very efficient, circular partic les have an inherent tendency to roll, which reduces their interlocking resi stance. As a result, angularity induced dilation and particle interlocking are suppressed. Modeling noncirc ular particles, in both 2D and 3D, has been proposed using different approaches with varying degrees of success. One class of methods utilizes mathematical functions to describe noncircular outlines. Ting et al (1993) and Ng (1994) presented ellipseshape d particles for 2D simulations, while Lin and Ng (1997) introduced ellipsoidshaped particles for 3D Simulations. Ellipses and ellipsoids have a smaller tendency to rotate, but do not provide a clos e representation of the actual particle shapes. Potapov and Campbell (1998) pro posed ovalshaped particles, which are computationally more efficient than the ellipses. Fa vier et al (1999) presented an axisymmetricalshaped particle, which can model a variety of shapes. However, highly angular particles cannot be easily generated.
PAGE 81
60 Another approach relies on approximating the particle shape using polygons, as proposed by Barbosa and Ghaboussi (1992) and Ma tuttis et al (2000). A more realistic representation of the particle shape can be achieved using polygonshaped particles; however the method is computa tionally intensive. Anothe r approach is to combine several circular outlines into a cluster to form more complex shapes. Jensen et al (1999) introduced the clustering technique by combin ing a number of sphe ricalshaped particles in a semirigid configuration, which is a bett er representation of natural soil particles. Any number of particles can be linked together to form a cluster as long as they rotate and translate as a rigid body. The interparticle contacts with in the cluster are constrained to be linearelastic with a high stiffness. Thomas and Bray (1999) presented a different idea for clustering by imposing ki nematic restrictions on discs within a cluster to prevent relative translations and rotations. Al though these technique s represent some improvement over earlier methods, the simulated particle outlines did not resemble those of actual particles. In add ition, a substantial increase in computational time was incurred due to the high contact stiffne ss required within the cluster. 3.3 Overlapping rigid clusters (ORC) Ashmawy et al (2003) proposed the use of Overlapping Rigid Clusters (ORC) technique to accurately simulate angular pa rticle shapes. A set of subroutines was introduced to PFC2D using ItascaÂ’s softwa respecific programming language, FISH. The builtin clump logic command was instrument al in formulating the new method. The clump logic allows several disk elements to move as a rigid body without detecting contacts or calculating contact forces, which increases computational efficiency. Further
PAGE 82
61 description of the clump logic is given in Itasca (1999). The ORC method relies on Â“clumpingÂ” a number of overlapping disc elements, such that the resulting outline coincides with, and is almost identica l to the actual particleÂ’s outline. First, twodimensional outlines of a series of particles are obtained using a digital microscope or scanner. Overlapping circles are then inscribed within the outline to capture the shape, as shown in fig. (21). The number of overlapping circles that can accurately represent the actual particles depe nds on the degree of nonuniformity in the original particle shape and angularity, the desired level of geometric accuracy, and the required computation time limit. Typically, ten to fifteen disks will adequately capture the Â“trueÂ” shape of a particle. Due to overl apping, the density of each disk must be scaled to ensure that the ma ss of the particle remains propor tional to the ar ea, regardless of the number of discs or level of overlap (refer to Eqn. 21). The ORC method is implemented within PFC2D by means of a series of FISH fu nctions that convert a particle assembly of discs into their correspondi ng angular particles as follows: 1. Particles are automatically generated as discrete circular elements using any of the builtin particle s generation techniques 2. The shape conversion algorithm is i nvoked to transform each circular particle into an angular equiva lent through repl acement with a corresponding set of disc elements, se lected randomly fr om a library of available particle shapes 3. A random rotation between 0 and 360 is applied to each transformed particle to avoid preferred orientatio ns and subsequent system anisotropy.
PAGE 83
62 Figure (22) shows a random assembly of circular particles that have been transformed to the corresponding angu lar shapes using the ORC method. 3.4 Experimental program Verification of the ORC method is essentia l in order to evaluate the accuracy and significance of the new procedure. To veri fy the new technique, an experimental and numerical validation program was carried out by: 1. Designing and building an experiment al setup that allows tracking translations and rotations of modelgrains resulting from an external disturbance. 2. Numerically simulating the experimental setup using PFC2D, with the same initial and boundary conditions and material properties. 3. Tracking the positions and rotations of modelgrains in both the experimental and numerical tests. 4. Qualitatively and quantitatively comparing the experimental and numerical results to verify the effectiveness of the ORC method. 5. Performing a parametric study in order to determine the effect of interparticle friction, contact s tiffness, and global damping on the numerical solution. Fraser River sand and Michigan sand were chosen, as natural angular and rounded materials, respectively, for the experimental verification. Only two dimensional modelgrains were considered for this study in order to accurately track the particle motions.
PAGE 84
63 The scanned outlines of ten different Fraser River sand grains and five different Michigan sand grains were resized and pl otted on top of 25 mm thick maple wood. A thickness of 25 mm and an equivalent diamet er of 50 mm were chosen for the average discs dimensions in order to facilitate m odelparticles manufacturing and tracking. The actual particle area, Aap, was scaled as follows: ap 2A 25 Scale Area (31) 3.4.1 Material Modelgrains were manufactured from ma ple wood, which has reliable strength, density, and stiffness. Maple wood properties are summarized in Table (31). A scroll saw with hard carbon steel blade was used to cut the modelgrains by following the shape outlines plotted on top of the 25 mm maple w ood sheet. The upper and lower edges were smoothed using fine sand paper to prevent edge smearing or tearing. For each of the Fraser River sandÂ’s ten particle shapes, eight duplicate particles were manufactured for the testing program, whereas for each of th e Michigan sandÂ’s five particle shapes, fourteen particles were made. An alphanumer ic labeling scheme was used in order to Table 31. Maple wood propertiesDensity, Kg/m3406.6 Specific gravity0.659 Bending strength, MPa115 Stiffness, MPa14100 Hardness7290
PAGE 85
64 facilitate the recognition of different modelgrains during the digitized snapshots. For each particle, the center of mass and a direction marker were labeled in order to facilitate the tracking procedure during comparison, as shown in figs. (31) and (32) for some Fraser River and some Michigan sand particles, respectively. 3.4.2 Test setup The experimental boundary conditions were defined by a box and a piston, where the modelgrains were assembled. The i nner dimensions of the box were 500500 mm. The box inner faces were covered with lamina ted plastic sheets to reduce friction with modelgrains (fig. 33). A r ectangularshaped piston was used to disturb the Fraser River sand modelgrains, whereas a trapezoidals haped piston was used for Michigan sand modelparticles. Figures (34) and (35) s how the experimental se tup for Fraser River Sand and Michigan Sand, respectively. Both pistons were covered with laminated plastic sheets on all sides that are in contact with th e modelgrains. A stai nless steel ruled scale was attached to the top of the pistons to m easure their displacements as they penetrated into the box. A cast acrylic sheet was att ached to the top of the box to prevent the particles from Â“popping upÂ” duri ng the test. The cast acrylic sheet was thick enough to prevent any deflection that woul d have caused friction with th e top of the modelgrains. A digital camera with a 0.79 Mega pixel reso lution was used to capture the snapshots. The camera was centered and leveled at a pproximately 1 m above the box, and was supported by means of a wooden frame. Two ru led Lsquares were used simultaneously to locate the center point of each modelgrain for initial setup (fig. 36).
PAGE 86
65 Figure 31. Fraser River sand modelgrains samples Figure 32. Michigan sa nd modelgrains samples
PAGE 87
66 Figure 33. The box along with the piston Figure 34. The experimental setup for Fr aser River sand modelparticles with the rectangular piston
PAGE 88
67 Figure 35. The experimental setup for Michigan sand modelparticles with the trapezoidal piston Figure 36. The two ruled Lsquares used to locate modelgrains
PAGE 89
68 The modelgrains were distributed randomly within the box and the two Lsquares were used to determine the X and Y coordinates for each modelgrain. Th e initial direction of each modelgrain was then adjusted and veri fied using the two Lsquares. The initial condition of the modelgrains fo r both modelgrains types are shown in figs. (37) and (38). 3.4.3 Image distortion effect In order to accurately measur e the location of the particle s at various stages of the test, it is necessary to calibrate the capture d image. A 12.512.5 mm graded transparent sheet was, therefore, used to evaluate imag e distortion. The sheet was placed under the cast acrylic sheet on top of a few modelgrains as a support during calibration as shown Figure 37. Initial conditions for Fraser Rive r sand modelgrains rectangular piston
PAGE 90
69 in fig. (39). The image was analyzed using the following methods: 1. Counting the number of pixels within the center square, corner squares, and middle edge squares. 2. The squares sidesÂ’ lengths were extrac ted from the image and compared to the lengths of the center squaresÂ’ side lengths. 3. An image of the initial condition of th e modelgrains in the laboratory was digitized then the X and Y coordinate s for the modelgrains centers were extracted from the digitized image a nd compared to the experimental X and Y coordinates. The results confirmed that the distortion effect was minimal and can safely be neglected within the sc ope of this study. Figure 38. Initial conditions for Michigan sand modelgrains trapezoidal piston
PAGE 91
70 3.4.4 Contact normal stiffness The contact normal stiffness of the wooden modelgrains was determined experimentally using a compression loading machine (fig. 310). The compression test was performed on each of the 10 different shapes of the Fraser River sand modelgrains. The compression test performed on modelgrain 5A of Fraser River sand is shown in fig. (311). In order to study the effect of grai n angularity on the contact normal stiffness, the compression tests were performed in different di rections for the same modelgrain. Load displacement curves for modelgrain 8Z of Fr aser River sand is shown in fig. (312), which shows the effect of loading direction on the loaddisplacement results. The results for all the particles are summarized in table (32 ). Based on the experi mental results, an Figure 39. Distortion calibration using graded transparent sheet
PAGE 92
71 Figure 310. Loading system used in determining the contact normal stiffness Figure 311. Compression test on m odelgrain 5A of Fraser River sand
PAGE 93
72 average value of 106 N/m for the contact normal stiffness was used in the numerical model. The contact shear stiffness value was taken the same as the contact normal stiffness. 3.4.5 Experimental testing procedures The experimental test was set to the initia l conditions and the pi ston was fitted in place with zero displacement. The camera fr ame with the camera centered to the box was set in place. The piston was then moved su ccessively and snapshots were taken at 50 mm intervals. The test was terminated when th e piston could not be pushed further. The experimental test was repeated five tim es, ensuring identical initial and boundary conditions, in order to study the inherent the variability of the pro cedure due to slight changes in the initial conditions. Figure 312. LoadDisplacement curves fo r Fraser River Sand modelgrain # 8Z 0 5 10 15 20 25 30 051015202530354045Displacement, mLoad. N particle_8Z_1: Kn=1680000 N/m particle_8Z_2: Kn =850000 N/m particle_8Z_3: Kn = 634000 N/m Knav = 1.06X106 N/m
PAGE 94
73 3.4.6 Experimental results 3.4.6.1 Fraser river sand with the rectangular piston The resulted images for tests 1 through 5 for initial conditions along with images at displacements 50, 100, 150, 200, and 250 mm ar e shown in fig. (313) through fig. (318). The initial conditions, fig. (313), for all five tests were ensured to be identical. At 50 mm piston displacement an apparent Â“contact force chain,Â” which mostly connected the grains that mostly moved, developed and wa s almost identical for all five tests (fig. 314). At 100 mm displacement, slight variations in the force chain among the various tests Table 32. Contact normal stiffness results for Fraser River Sand modelgrains Particle Kn, N\m Type 1Z 9.11x 105Type 2Z 1.07x 106Type 3Z 8.38x 105Type 4Z 9.79x 105Type 5Z 1.00x 106Type 6Z 8.63x 105Type 7Z 1.36x 106Type 8Z 1.06x 106Type 9Z 9.16x 105Type 10Z 8.36x 105Maximum normal stiffness 1.36x 106Minimum normal stiffness 8.36x 105Average normal stiffness 9.93x 105
PAGE 95
74 began to appear (fig. 315). At 150 mm disp lacement, the force chain is very close in shape in all tests except test # 4 (fig. 316). After 200 mm displacement, the difference in force chains became more noticeable, with diffe rences in particles positions and rotations throughout the system (fig. 317). As th e displacement increased beyond this point, differences were noticed (fig. 318). The X and Y coordinates, along with the rotations of the modelgrains, for the five experiment al tests, were tracked by digitizing the experimental snapshots and were stored fo r performing the quantitative comparison with the numerical simulation. 3.4.6.2 Michigan sand with the trapezoidal piston For Michigan sand modelparticles, fig. (319) through fig. (324) show the resulted images for tests 15 for initial condition, displ acements 50, 100, 150, 200, and 250 mm, respectively. It should be noted that only during tests 2 and 3, a displacement of 250 mm could be achieved, while in tests 1, 4, and 5, the piston could not move further than 238240 mm. As shown in fig. (319), th e initial conditions for all five tests were ensured to be identical. The modelgrains of the first row were initially aligned around the trapezoidalshaped piston. More partic les displacements and rotations were noticed when the trapezoidalshaped piston was used, co mpared to the rectangular piston. At 50 mm piston displacement, particles positions and rotations were almost identical for all five tests (fig. 320). At 150 mm displacement, differences between the two groups began to appear but results fo r tests 1, 4, and 5 (group 1) can be considered identical ,whereas results for tests 2 and 3 (group 2) were different (fig. 322). For group 2, the piston was pushed to 250 mm, whereas for group 1, the piston could not be pushed
PAGE 96
75 further than 238240 mm. At this point, th e differences in particles positions and rotations throughout the system became mo re noticeable, and as the displacement increased further, large differences were not iced (figs. 323 and 324). The X and Y coordinates, for the Michigan sand modelgrai ns, along with the rotations, were tracked and recorded. It should be noted that although the same initial condition was ensured, reproducible results could not ac hieved. The results indicate that minor changes in the initial conditions may dramatically affect the overall behavior of the modelgrains assembly at large displacements. 3.5 Numerical simulation The Particle Flow Code, PFC2D, was used to numerically simulate the experimental setup for both Fraser River sand and Michigan sand. The box was simulated as four fixed rigid walls, whereas the piston was si mulated as rigid walls that moved vertically at constant velocity. The vertical moveme nt was slow enough to ensure a stable solution. The clumps were generate d at the same positions as the experimental setup using the ORC method desc ribed earlier. However, th e particles were prevented from rotating randomly and, instead, were a ligned in the same direction as the modelgrains. The initial condition for the nume rical models, along with the corresponding experimental initial conditions for Test 1, for both Fraser River sand and Michigan Sand, are shown in figs. (325) a nd (326), respectively.
PAGE 97
Figure 313. Experimental initial conditions (tests 15) for Fr aser River sand modelgrains Test 1 Test 2 Test 3 Test 4 Test 5 Figure 314. Experimental results after imposing 50 mm disp lacement (tests 15) for Fraser River sand modelgrains Test 1 Test 2 Test 3 Test 4 Test 5 76
PAGE 98
Figure 316. Experimental results after imposing 150 mm disp lacement (tests 15) for Fraser River sand modelgrains Test 1 Test 2 Test 3 Test 4 Test 5 Figure 315. Experimental results after imposing100 mm disp lacement (tests 15) for Fraser River sand modelgrains Test 1 Test 2 Test 3 Test 4 Test 5 77
PAGE 99
Figure 318. Experimental results after imposing 250 mm disp lacement (tests 15) for Fraser River sand modelgrains Test 1 Test 2 Test 3 Test 4 Test 5 Figure 317. Experimental results after imposing 200 mm disp lacement (tests 15) for Fraser River sand modelgrains Test 1 Test 2 Test 3 Test 4 Test 5 78
PAGE 100
Figure 320. Experimental results after imposing 50 mm di splacement (tests 15) for Michigan River sand modelgains Test 1 Test 2 Test 3 Test 4 Test 5 Figure 319. Experimental initial conditions (t ests 15) for Michigan sand modelgrains Test 1 Test 2 Test 3 Test 4 Test 5 79
PAGE 101
Figure 322. Experimental results after imposing 150 mm di splacement (tests 15) for Michigan River sand modelgrains Test 1 Test 2 Test 3 Test 4 Test 5 Figure 321. Experimental results after imposing 100 mm di splacement (tests 15) for Michigan River sand modelgrains Test 1 Test 2 Test 3 Test 4 Test 5 80
PAGE 102
Figure 324. Experimental results after imposing 250 mm di splacement (tests 15) for Michigan River sand modelgrains Test 1 Test 2 Test 3 Test 4 Test 5 Figure 323. Experimental results after imposing 200 mm di splacement (tests 15) for Michigan River sand modelgrains Test 1 Test 2 Test 3 Test 4 Test 5 81
PAGE 103
82 The density of the simulated grains wa s set to the maple wood density. The contact normal/shear stiffness was set to 106 N/m (the value obtained experimentally for the Fraser River sand modelgr ains). The simulated state, including grain positions and rotations, was saved at 50 mm piston movement intervals for further comparison with the experimental results. The resulted images for displacement increments, for both Fraser River sand and Michigan sand, will be shown later on. Because Fraser River sand particles are highly angular, and because the interparticle contact forces are lo w, the friction coefficient and contact stiffness do not have a significant effect on the displacement and ro tation of the particles. Instead, the displacement pattern is governed mostly by the shapes and interlocking mechanism between the grains. This is consistent with th e objective of this study, that is, to evaluate the ability of the ORC technique to model a ngular particles within a DEM framework. The effect of contact stiffness and interparti cle fiction coefficient was investigated and will be presented later in a parametric study. For further illustration of the particle sh ape, the experimental setup for Fraser river sand was modeled as simple discs in stead of the ORC method. The average diameter for the created discs was set 50 mm. The discs cente rs were taken the same as modelgrains centers. The same numerical te sting procedures were applied to the discs and positions and rotations were recorded and saved at 50 mm increments for further quantitative comparison with both the experi mental and the numerical modeling using ORC technique.
PAGE 104
83 Figure 325. Experimental and numerical initial conditions for Fraser River Experimental Test # 1 Numerical model Figure 326. Experimental and numerical initial conditions for Michigan sand Experimental Test # 1 Numerical model
PAGE 105
84 3.6 Comparison and discussion 3.6.1 Fraser river sand The numerical results, us ing both the ORC method and simple discs, along with the experimental results of Test 1 at displacements of 50, 100, 150, 200 and 250 mm are shown in fig. (327) through fig. (331). Fr om a qualitative standpoin t, it is clear that there is good agreement between the numeri cal results using th e ORC method and the experimental results within the small disp lacement range, whereas some deviations appeared when discs were used to simulate the particles. With a few exceptions, the grain chains are almost identical (figs. 327 and 328). The grain chain begins to deviate as the piston displacement increases, and the difference may be considered significant at large displacement, as seen in (figs. 330 and 331). It should be not ed that the numerical modeling results using the ORC method were closer to the experimental results than simple discs. However, it is important to note that large differe nces were noticed not only between numerical modeling and the experi mental tests, but also among the five Â“identicalÂ” experimental tests. It is also essential to note that the previous qualitative comparison was performed between the numerical simulation and experimental Test 1 only. It should also be noted th at the apparent Â‘force chai nÂ” that developed along the box may be a result of the experi mental setup. First, partic les with identic al shape were all aligned in the same directio n to facilitate the setup. Se cond, because all particles are of the same equivalent size. Particle arrangement takes an almost a row/column configuration. Because the box size is rela tively small compared to the size of the
PAGE 106
85 particles, boundaries may have a pronounced effect on the respons e. It is expected that a Â“realÂ” system will not behave similarly. For the purpose of this dissertation, the goal is to compare experimental and numerical results, so particle configuratio n is irrelevant as long as it matches in both the experime ntal and the numerical model. In order to evaluate the ability of the nume rical model to simulate the tests, it is necessary to compare the results in a quant itative manner. Specifically, if certain quantitative measures, such as particle rota tion and displacement, are compared between the experimental and numerical results, the numerical values shoul d ideally fall within the range of the experimental data. To th is end, a statistical comparison between the experimental results and the numerical resu lts was conducted. Histograms were drawn for cumulative vertical (Y) displacements a nd rotations of the modelgrains at 50 mm piston displacement intervals. The results of the five experimental tests and the numerical modeling using both ORC method and si mple discs were plotted. The results for particles rotations at piston displacements of 50, 150, and 250 mm are shown in fig. Figure 327. Experimental an d numerical results for Michigan sand after imposing 50 mm displacement Experimental Test # 1 Numerical model (ORC) Numerical model (Discs)
PAGE 107
86 Figure 329. Experimental and numerical re sults for Michigan sand after imposing 150 mm displacement Experimental Test # 1 Numerical model (ORC) Numerical model (Discs) Figure 328. Experimental and numerical resu lts for Michigan sand after imposing 100 mm displacement Experimental Test # 1 Numerical model (ORC) Nu merical model (Discs)
PAGE 108
87 Figure 330. Experimental and numerical resu lts for Michigan sand after imposing 200 mm displacement Experimental Test # 1 Numerical model (ORC) Numerical model (Discs) Figure 331. Experimental and numerical resu lts for Michigan sand after imposing 250 mm displacement Experimental Test # 1 Numerical model (ORC) Numerical model (Discs)
PAGE 109
88 (332) through fig. (334). The histogram s demonstrate that the numerical modeling results, using ORC method, lie within th e experimental resu lts except at large displacements, whereas that of using simple discs did not lie within the experimental results especially at moderate to large displacements. Similarly, the results for vertical displacement, shown in fig. (335) through fig. (337), indicate that the numerical modeling, using the ORC method, results can closely represent the experimental results, while simple discs resulted in inaccurate simu lation of the experiment al setup. Because the horizontal displacement s of particles are very small under the given boundary conditions, a similar comparis on would not be of much va lue and is therefore not reported within the context of this work. 3.6.2 Michigan sand The numerical results along with the experimental results of Test 1 at displacements of 50, 150, and 250 mm are shown in fig. (338) through fig. (340). Qualitatively, a good agreement between the numerical results using the ORC method and the experimental results within the small displacement range was noticed. Deviations began to appear between the experi mental and the numerical results at larger displacements. The differences between the numerical and experimental results became somewhat significant at large disp lacement, as seen in fig. (340). It should be noted that the deviations for Michigan sand simulation were less than that of Fraser River sand (as shown later in the quantitative comparison). The use of the trapezoidalshaped piston with Michigan sand involved more movements for the modelparticle s. The differences among the experimental tests are less th an those of Fraser River sand.
PAGE 110
89 Figure 332. Experimental and numerical rotations for Fraser River sand after imposing 50 mm displacement 0 10 20 30 40 50 60 70 70605040302010010203040506070 Rotation range, DegNo of Particles Test # 1 Test # 2 Test # 3 Test # 4 Test # 5 Numerical (ORC) Numerical (Discs) Fraser River Sand Rectangular Piston 50 mm DisplacementFigure 333. Experimental and numerical rotations for Fraser River sand after imposing 150 mm displacement 0 10 20 30 40 50 60 70 70605040302010010203040506070 Rotation range, DegNo of Particles Test # 1 Test # 2 Test # 3 Test # 4 Test # 5 Numerical (ORC) Numerical (Discs) Fraser River Sand Rectangular Piston 150 mm Displacement
PAGE 111
90 Figure 334. Experimental and numerical rotations for Fraser River sand after imposing 250 mm displacement 0 10 20 30 40 50 60 70 70605040302010010203040506070 Rotation range, DegNo of Particles Test # 1 Test # 2 Test # 3 Test # 4 Test # 5 Numerical (ORC) Numerical (Discs) Fraser River Sand Rectangular Piston 250 mm DisplacementFigure 335. Experimental and numerical Ve rtical (Y) displacement for Fraser River sand after imposing 50 mm displacement 0 10 20 30 40 50 60 70 3020100102030405060708090100110120130140150 Vertical Displacement, mmNo of Particles Test # 1 Test # 2 Test # 3 Test # 4 Test # 5 Numerical (ORC) Numerical (Discs) Fraser River Sand Rectangular Piston 50 mm Displacement
PAGE 112
91 Figure 336. Experimental and numerical Ve rtical (Y) displacement for Fraser River sand after imposing 150 mm displacement 0 10 20 30 40 50 60 70 3020100102030405060708090100110120130140150 Vertical Displacement, mmNo of Particles Test # 1 Test # 2 Test # 3 Test # 4 Test # 5 Numerical (ORC) Numerical (Discs) Fraser River Sand Rectangular Piston 150 mm DisplacementFigure 337. Experimental and numerical Ve rtical (Y) displacement for Fraser River sand after imposing 250 mm displacement 0 10 20 30 40 50 60 70 3020100102030405060708090100110120130140150 Vertical Displacement, mmNo of Particles Test # 1 Test # 2 Test # 3 Test # 4 Test # 5 Numerical (ORC) Numerical (Discs) Fraser River Sand Rectangular Piston 250 mm Displacement
PAGE 113
92 A statistical comparison between the expe rimental and the numerical results was conducted. Histograms were drawn for cumu lative vertical (Y) displacements and rotations of the modelgrains at 50 mm pist on displacement intervals. The results for rotations of the five experimental tests and the numerical mode ling were plotted and shown in fig. (341) through fig. (343). The histograms demonstrate that the numerical modeling results lie within the experimental results even at large displacements. The results for vertical displacement shown in fig. (344) through fig. (346) indicate that the numerical modeling results deviat ed from the experimental results and were not as close as rotations results especial ly for large displacements. 3.7 Parametric study 3.7.1 Effect of interp article friction The effects of interparticle friction a ngle has been studied through numerical Figure 338. Experimental an d numerical results for Michigan sand after imposing 50 mm displacement Experimental Test # 1 Numerical model
PAGE 114
93 Figure 340. Experimental an d numerical results for Michigan sand after imposing 250 mm displacement Experimental Test # 1 Numerical model Figure 339. Experimental an d numerical results for Michigan sand after imposing 150 mm displacement Experimental Test # 1 Numerical model
PAGE 115
94 Figure 341. Experimental and numerical rotations for Michigan sand after imposing 50 mm displacement Figure 342. Experimental and numerical rotations for Michigan sand after imposing 150 mm displacement 0 5 10 15 20 25 30 35 40 45 50 9080706050403020100102030405060708090 Rotation range, DegNo of Particles Test # 1 Test # 2 Test # 3 Test # 4 Test # 5 Numerical Michigan Sand Trapezoidal Piston 50 mm Displacement 0 5 10 15 20 25 30 35 40 45 50 9080706050403020100102030405060708090 Rotation range (bin), DegNo of Particles Test # 1 Test # 2 Test # 3 Test # 4 Test # 5 Numerical Michigan Sand Trapezoidal Piston 150 mm Displacement
PAGE 116
95 Figure 343. Experimental and numerical rotations for Michigan sand after imposing 250 mm displacement Figure 344. Experimental and numerical vertical (Y) displacement for Michigan sand after imposing 50 mm displacement 0 5 10 15 20 25 30 35 40 45 50 9080706050403020100102030405060708090 Rotation range (bin), DegNo of Particles Test # 1(240 mm) Test # 2 Test # 3 Test # 4(238 mm) Test # 5(240 mm) Numerical Michigan Sand Trapezoidal Piston 250 mm Displacemen t 0 5 10 15 20 25 30 35 40 45 50 3020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Test # 1 Test # 2 Test # 3 Test # 4 Test # 5 Numerical Michigan Sand Trapezoidal Piston 50 mm Displacement
PAGE 117
96 Figure 345. Experimental and numerical rotations for Michigan sand after imposing 150 mm displacement Figure 346. Experimental and numerical rotations for Michigan sand after imposing 250 mm displacement 0 5 10 15 20 25 30 35 40 45 50 3020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Test # 1 Test # 2 Test # 3 Test # 4 Test # 5 Numerical Michigan Sand Trapezoidal Piston 150 mm Displacemen t 0 5 10 15 20 25 30 35 40 45 50 3020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Test # 1(240 mm) Test # 2 Test # 3 Test # 4(238 mm) Test # 5(240 mm) Numerical Michigan Sand Trapezoidal Piston 250 mm Displacement
PAGE 118
97 Figure 347. Effect of interparticle friction (Fraser River sand) on particles rotation after imposing 50 mm displacement Figure 348. Effect of interparticle friction (Fraser River sand) on particles rotation after imposing 150 mm displacement 0 10 20 30 40 50 60 70 809080706050403020100102030405060708090 Rotation, degNo of Particles Average Experimental Numerical (f=1.00) Numerical (f=0.25) Numerical (f=0.05) Fraser River Sand Rectangular Piston Effect of interparticles Friction 50 mm Displacement 0 10 20 30 40 50 60 70 80 9080706050403020100102030405060708090 Rotation, degNo of Particles Average Experimental Numerical (f=1.00) Numerical (f=0.25) Numerical (f=0.05) Fraser River Sand Rectangular Piston Effect of interparticles Friction 150 mm Displacement
PAGE 119
98 modeling of Fraser River sand with the rectangular piston. The particles rotations for the numerical analysis using interparticle an internal friction coefficient of 0.05, 0.25, and 1.00 along with the average results for the five experimental tests ar e shown in figs. (347) through (349) for displacements 50, 150, and 250 mm, respectively. The results indicate that the interparticle friction has an effect on the particles rotations, which can be attributed to the small contact force that was encountered with the testing setup. Values between 0.251.00 were found to be appropriate for the numerical modeling, since they coincided with the experimental results in gene ral. Similarly, figs, (350) through (352) show the vertical displacements at the same displacements (50, 150, and 250mm), which confirmed that there is an effect for interparticle friction. 3.7.2 Effect of global damping In the discrete element model, energy dissipates through fr iction, contact, and global damping. As defined by Cundall and Strack (1979), global damping operates on the absolute velocities of the elements and it is intr oduced to the equation of motion calculations as shown in Appendix A. Global damping can be envisione d as the effect of dashpots that connect each element to the ground. The dashpots operate both on the velocity vector components and on the rotationa l velocity. In the experimental setup, the particles were allowed to move only in the XY plane. They were sliding on the boxÂ’s base which is covered with laminated plas tic sheets to minimize friction. The small friction with the base and the gravity effect were introduced by setting a global damping ratio to 1.4. The numerical test was repeated three times using global damping ratios of 0.70, 1.4, and 2.10 and the results, for both pa rticle rotations and vertical movement at
PAGE 120
99 Figure 349. Effect of interparticle friction (Fraser River sand) on particles rotation after imposing 250 mm displacement Figure 350. Effect of interparticle friction (Fraser River sand) on pa rticles vertical (Y) displacement after imposing 50 mm displacement 0 10 20 30 40 50 60 70 80 9080706050403020100102030405060708090 Rotation, degNo of Particles Average Experimental Numerical (f=1.00) Numerical (f=0.25) Numerical (f=0.05) Fraser River Sand Rectangular Piston Effect of interparticles Friction 250 mm Displacement 0 10 20 30 40 50 60 70 803020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Average Experimental Numerical (f=1.00) Numerical (f=0.25) Numerical (f=0.05) Fraser River Sand Rectangular Piston Effect of interparticles Friction 50 mm Displacement
PAGE 121
100 Figure 351. Effect of interparticle friction (Fraser River sand) on pa rticles vertical (Y) displacement after imposing 150 mm displacement Figure 352. Effect of interparticle friction (Fraser River sand) on pa rticles vertical (Y) displacement after imposing 250 mm displacement 0 10 20 30 40 50 60 70 80 3020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Average Experimental Numerical (f=1.00) Numerical (f=0.25) Numerical (f=0.05) Fraser River Sand Rectangular Piston Effect of interparticles Friction 150 mm Displacement 0 10 20 30 40 50 60 70 80 3020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Average Experimental Numerical (f=1.00) Numerical (f=0.25) Numerical (f=0.05) Fraser River Sand Rectangular Piston Effect of interparticles Friction 250 mm Displacement
PAGE 122
101 Figure 353. Effect of globa l damping (Fraser River sand) on particles rotation after imposing 50 mm displacement Figure 354. Effect of globa l damping (Fraser River sand ) on particles rotation after imposing 150 mm displacement 0 10 20 30 40 50 60 70 9080706050403020100102030405060708090 Rotation, degNo of Particles Average Experimental Numerical (D=2.10) Numerical (D=1.40) Numerical (D=0.70) Fraser River Sand Rectangular Piston Effect of Global Stiffness 50 mm Displacement 0 5 10 15 20 25 30 35 40 45 50 9080706050403020100102030405060708090 Rotation, degNo of Particles Average Experimental Numerical (D=2.10) Numerical (D=1.40) Numerical (D=0.70) Fraser River Sand Rectangular Piston Effect of Global Damping 150 mm Displacement
PAGE 123
102 Figure 355. Effect of globa l damping (Fraser River sand ) on particles rotation after imposing 250 mm displacement Figure 356. Effect of global damping (Fra ser River sand) on particles vertical (Y) displacement after imposing 50 mm displacement 0 5 10 15 20 25 30 35 40 45 50 9080706050403020100102030405060708090Rotation, degNo of Particles Average Experimental Numerical (D=2.10) Numerical (D=1.40) Numerical (D=0.70) Fraser River Sand Rectangular Piston Effect of Global Damping 250 mm Displacement 0 10 20 30 40 50 60 70 80 3020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Average Experimental Numerical (D=2.10) Numerical (D=1.40) Numerical (D=0.70) Fraser River Sand Rectangular Piston Effect of Global Stiffness 50 mm Displacement
PAGE 124
103 Figure 357. Effect of global damping (Fra ser River sand) on particles vertical (Y) displacement after imposing 150 mm displacement Figure 358. Effect of global damping (Fra ser River sand) on particles vertical (Y) displacement after imposing 250 mm displacement 0 10 20 30 40 50 60 3020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Average Experimental Numerical (D=2.10) Numerical (D=1.40) Numerical (D=0.70) Fraser River Sand Rectangular Piston Effect of Global Damping 150 mm Displacement 0 10 20 30 40 50 60 3020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Average Experimental Numerical (D=2.10) Numerical (D=1.40) Numerical (D=0.70) Fraser River Sand Rectangular Piston Effect of Global Damping 250 mm Displacement
PAGE 125
104 displacements of 50mm, 150mm and 250mm, are shown in figs. (353) through (358). The results for both global damp ing ratios of 2.1 and 1.4 were almost identical and very close to the experimental results, whereas mo re translations and rotations were observed when the global damping coefficient was set to 0.70. The results confirmed that using a damping ratio higher than the critical da mping ratio of 1.0 is satisfactory for the numerical modeling purposes. 3.7.3 Effect of contact normal stiffness The effect of contact stiffness has b een studied using numerical modeling of Fraser River sand with the rectangular pist on. Particles rotations for the numerical analysis using contact stiffness of 104, 105, and 106 kPa along with the average results for the five experimental tests ar e shown in figs. (359) through (361) for displacements 50, 150, and 250 mm. The results indicated that the contact st iffness had a minor effect on the particles rotations. A value of 106 N/m was found to be appropriate. It should be noted that the same value was determined fr om the compression test on the modelgrains and has been used throughout the numerical simulation of the experimental setup. Similarly, figs. (362) through (364) show the vertical displacements, which confirm the previous findings. 3.8 Conclusion An experimental validation of the ability of DEM to model angular particles created using overlapping rigi d clusters (ORC) was presen ted. The ORC technique creates a clump that contains circular particles in such a way that the outline of the clump
PAGE 126
105 Figure 359. Effect of contact normal stiffn ess (Fraser River sand) on particles rotation after imposing 50 mm displacement Figure 360. Effect of contact normal stiffn ess (Fraser River sand) on particles rotation after imposing 150 mm displacement 0 10 20 30 40 50 60 70 9080706050403020100102030405060708090 Rotation, degNo of Particles Average Experimental Numerical (k=1000000 N/m) Numerical (k=100000 N/m) Numerical (k=10000 N/m) Fraser River Sand Rectangular Piston Effect of Contact Stiffness 50 mm Displacement 0 10 20 30 40 50 60 70 9080706050403020100102030405060708090Rotation, degNo of Particles Average Experimental Numerical (k=1000000 N/m) Numerical (k=100000 N/m) Numerical (k=10000 N/m) Fraser River Sand Rectangular Piston Effect of Contact Stiffness 150 mm Displacement
PAGE 127
106 Figure 361. Effect of contact normal stiffn ess (Fraser River sand ) on particles rotation after imposing 250 mm displacement Figure 362. Effect of contact normal stiffn ess (Fraser River sand) on particles vertical (Y) displacement after imposing 50 mm displacement 0 10 20 30 40 50 60 70 9080706050403020100102030405060708090Rotation, degNo of Particles Average Experimental Numerical (k=1000000 N/m) Numerical (k=100000 N/m) Numerical (k=10000 N/m) Fraser River Sand Rectangular Piston Effect of Contact Stiffness 250 mm Displacement 0 10 20 30 40 50 60 70 80 3020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Average Experimental Numerical (k=1000000 N/m) Numerical (k=100000 N/m) Numerical (k=10000 N/m) Fraser River Sand Rectangular Piston Effect of Contact Stiffness 50 mm Displacement
PAGE 128
107 Figure 363. Effect of contact normal stiffn ess (Fraser River sand) on particles vertical (Y) displacement after imposing 150 mm displacement Figure 364. Effect of contact normal stiffn ess (Fraser River sand) on particles vertical (Y) displacement after imposing 250 mm displacement 0 10 20 30 40 50 60 3020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Average Experimental Numerical (k=1000000 N/m) Numerical (k=100000 N/m) Numerical (k=10000 N/m) Fraser River Sand Rectangular Piston Effect of Contact Stiffness 150 mm Displacement 0 10 20 30 40 50 60 3020100102030405060708090100110120130140150 Vertical Displacement range, mmNo of Particles Average Experimental Numerical (k=1000000 N/m) Numerical (k=100000 N/m) Numerical (k=10000 N/m) Fraser River Sand Rectangular Piston Effect of Contact Stiffness 250 mm Displacement
PAGE 129
108 coincides with that of the actual particle. An experimental set up was built and model particles were manufactured in the laboratory. The experi mental tests were modeled numerically using PFC2D. The experimental test result s of five identical tests were compared qualitatively and quantitatively w ith the corresponding numerical simulation. The ORC technique was found to be effectiv e in modeling the behavior of angular particle assemblies. While particle displacem ents and rotations of both the numerical and experimental systems were almost identical at small displacements, larger differences were observed at high displacements, even amo ng the five experiment al tests that have been performed ensuring identical initial and boundary conditions. This variability is attributed to minor changes in the initia l conditions, as well as other inherent uncertainties in the system, and cannot be avoided. However, the numerical model revealed that the ORC simulated system behavior is, in general within the range of the five experimental tests. The experimental validation was extended to two different types of model grains. A parametric study has been conducted to examine the effect of interparticle friction, contact stiffness, and global damping on the numerical results.
PAGE 130
109 CHAPTER FOUR: THE MODIFIED OV ERLAPPING RIGID CLUSTERS 4.1 Introduction Although the ORC method, described in chap ter three, ensures that the mass of the particle is proportional to the area, it neith er guarantees that the moment of inertia nor the center of mass of the ORC particle is id entical to that of th e actual particle. The differences in the mass properties will have minimal effect on the simulations of static problems. In order to accurately model dynamic problems, however, the mass properties should be identical for both the actual partic le and the numerical model. This chapter presents a modification for the Overlappi ng Rigid Clusters tec hnique that ensures identical matching for both the center of grav ity and the mass moment of inertia. It also introduces a sequenced particle creation pro cedure that ensures the minimization of the uncovered area of the original particle and is operator independent. 4.2 The modified ORC method 4.2.1 Operatorindependent particle creation procedure Twodimensional outlines of a series of particles are obtai ned using a digital microscope or scanner. The outline is then discretized into equal interval XY coordinates. The area, mass, center of gr avity, and mass moment of inertia about the normal axis passing through th e center of gravity for the particle are determined.
PAGE 131
110 The elements/discs were added to fill the par ticle area and coincide with the particle outline under the following sequence/conditions: The first element (Disc # 1) is added su ch that its center coincides with the particle center of gravity. The radius of disc # 1 is determined such that the Disc # 1 should be tangent to the particle outline at least at one point Figure (41) shows Disc # 1 for Michigan sand particle # 1. The next disc should be added in such a way that it minimizes the uncovered area of the actual particle. Overlapping is allo wed between the elements within the same particle. It is clear that th ere is a unique disc that satis fies the condition of minimizing the uncovered area, so the ope ratorindependency is ensured. Figure (42) shows the sequence of adding the discs to particle # 1 of Michigan sand satisfying the minimum uncovered area condition. The addition of the next discs should follow the same criterion of minimizing the uncovered area. The area ratio, defined in equation (41), should be checked after the addition of every disc to the particle. Once the area ratio reaches the desired value (an area ratio equal or more than 95% may be considered satisfactory for most particles), the element addition should be stopped. Defi ning the particle area as (Ap), the uncovered area of the particle as (Aunc), the area ratio (Ar) can be defined as: 100 1 p unc rA A A (41) Once the desired area ratio is satisfie d, discs addition should be stopped.
PAGE 132
111 Figure 41. Outline of particle # 1 of Michigan sand along with Disc#1 Particle Outline Element/Disc # 1 Uncovered Area R1(Xcg, Ycg) Particle Outline Element/Disc # 1 Uncovered Area R1(Xcg, Ycg) (1) (13) (12) (9) (8) (7) (6) (5) (4) (3) (2) (7) (10) (11) (1) (13) (12) (9) (8) (7) (6) (5) (4) (3) (2) (7) (10) (11) Figure 42. Sequence of adding discs fo r particle # 1 of Michigan sand using the modified ORC method
PAGE 133
112 4.2.2 Kinetic compatibility equations Four Â“compatibility discsÂ” are selected to satisfy inertial compatibility conditions. It is preferable to choose the discs that mo stly affect the mass and moment properties. For consistency, these will be the discs that have the maximum center coordinates in each of the four quadrants. The mass of the actual particle should be equal to that of the discs that constitute the clump. The mass balance equation can be expressed as: p i N 1 im m (42) Where N is the number of discs within the particle, mi is the mass of disc # i, and mp is the mass of the particle. The center of gravity of both the actual a nd the created clump s hould be identical. The eccentricity of the discs centers from th e actual particle center in both X and Y direction, should be vanished. This can be expressed as: 0 x mi i N 1 i (43) 0 y mi i N 1 i (44) Where xi is the distance in xdirection between the center of disc # i and the center of the particle and yi is distance in ydirection between the center of disc # i and the center of the particle. The mass moment of inertia for the cr eated and the actual particle, around a perpendicular axis passing through the center of gravity, should also be identical. This can be ensured as follows:
PAGE 134
113 p iz z N 1 iI I (45) where Izi is the mass moment of inertia for disc # i around the perpendicular axis passing through the particle center of gravity and Izp is the mass moment of inertia of the particle around its center of gravity. 4.2.3 Procedure The densities of the compatibility discs are treated as unknow ns of the kinetic compatibility equations. Solving the kinetic compatibility equations (Eqns 42 to 45) simultaneously, the densities of the compatibility discs can be determined. If the four densities that result from solving the equa tions are positive values, the procedure is terminated, and kinetic similitude is ensured. On the other hand, if one or more of the densities were found to be negative values, more iterations are required. In the next iteration the mass of the elements with negative densities should be set to small positiv e values (1 Â– 2 % of the total mass of the compatibility discs should be sufficient). The masses of the remaining compatibility discs, with positive density, should be decreased to ensure the overall mass compatibility. The reduction in the masses should be proportion to each mass. Four more compatibility discs should be chosen to invoke the compatibi lity equations for the next step. The four compatibility discs chosen for second iteration may be different from the four discs used in the first iteration. Also, some of the four compatibility discs used in the first iteration could be used in the next one. The procedure is terminated once positive values are
PAGE 135
114 obtained for the current four compatibility discs. Table (41) shows the successive iterations used for Michigan sand particle # 1. In most particle shapes, two to three itera tions are typically needed to satisfy the kinetic compatibility equations. An excel sp readsheet had been developed to solve the compatibility equations. Mass and inertia prop erties are calculated and the four kinetic compatibility equations are solved simultane ously. The spreadsheet allows performing iterations in case of negative unit densities are obtained until a stable solution is achieved. 4.3 Implementation of the modified ORC method Implementation of the modified ORC method was tested on two particles, one is rounded and the other is angular. For each part icle, the discs were inscribed within the particle outline using the crit erion described earlier for mi nimizing the uncovered area of Table 41. Procedure of enforcing the kinetic compatibility equations for particle # 1 of Michigan sand Disc # Initial unit densities Kinetic compatibility's unit densities Kinetic compatibility's relative masses Positive relative masses Positive unit densities 50.3053.1173712.6182324.5621.951 70.3050.4102724.627100.0000.015 80.3050.036257.677161.3380.023 90.3054.1103584.7052244.4732.573 20.3050.564671.418100.0000.018 60.3050.068453.489100.0000.015 110.3050.7895602.0674614.5910.976 120.3052.1951914.2791576.8490.591 51.951 2.342 92.573 2.808 110.976 0.434 120.591 1.153 1st iteration 2nd iteration 3rd iteration
PAGE 136
115 the actual particle. Michigan sand particle # 1 was chosen as an example for rounded particles, whereas Fraser River sand particle # 4 was chosen as an example for angular particle. 4.3.1 Michigan sand particle # 1 The outline of particle # 1 of Michigan sand is shown in fig. (43). The original ORC technique was used to clump discs/elemen ts in such a way that the outline of the clumped discs coincide with that of the part icle (Ashmawy et al 2003). The results are shown in fig.(44). The coordinates of cente r of gravity, for each element/disc relative to that of the actual particle, the radii, and the re lative densities are shown in the table (42). The modified ORC method was used to m odel the same particle. The suggested element addition sequence, by minimizing the uncovered area, was used. The results are shown in fig. (42). Disc numbering show s the discs adding sequence to the particle outline. Table (42) shows the elements properties. The mass, mass moment of inertia, and eccentricities in both X and Y directions are tabulated in table (43) for both the original ORC and modified ORC methods. The results s how that the modified ORC technique identically matched the original particle mass properties, whereas the original ORC method did not. 4.3.2 Fraser river sand particle # 4 The particle outline is shown in fig. (45). The results of using the original ORC method are shown in fig. (46 ), whereas those of using the modified ORC method are shown in fig. (47). The element geometrical data along with the re lative densities are
PAGE 137
116 Figure 44. Generation of particle # 1 of Michigan sand usin g the ORC technique Figure 43. Outline of particle # 1 of Michigan sand
PAGE 138
117 Table 42. Elements/Discs properties for particle # 1 of Michigan sand X, mmY, mmR, mm Unit Density X, mmY, mmR, mm Unit Density 1 2.9320.088.28 0.34 0.000.0020.480.31 2 17.3012.7416.06 0.34 13.592.8514.490.02 3 8.377.2927.60 0.34 14.7111.5910.620.31 4 0.0410.1722.97 0.34 1.3718.107.990.31 5 17.120.6418.24 0.34 12.4415.026.702.34 6 17.816.2815.71 0.34 0.597.3615.870.01 7 14.0215.709.98 0.34 8.237.1815.840.02 8 1.8316.5915.71 0.34 6.753.8516.360.02 9 8.932.3228.13 0.34 3.3218.845.732.81 4.907.8226.15 0.34 17.8912.938.200.31 11 0.340.9336.17 0.34 7.184.3813.350.43 12 18.190.3410.031.15 13 20.408.096.810.31 Disc No. Using The ORC techniqueUsing the modified ORC technique Table 43. Results of the ORC and the m odified ORC techniques for Michigan sand particle # 1 PropertyThe actual particle Using the modified ORC techniqueUsing the ORC technique Area, mm21963.501963.501963.50 Eccentricity in X direction, mm00 0.76 Eccentricity in Y direction, mm00 0.25 Mass moment of inertia, unit mass.mm26.41E+056.41E+051.51E+06
PAGE 139
118 Figure 45. Outline of particle # 4 of Fraser River sand Figure 46. Generation of pa rticle # 4 of Fraser River sand using the ORC technique
PAGE 140
119 shown in table (44). The mass properties for the created particle using both methods along with the original particle prope rties are shown in table (45). 4.4 Conclusion The overlapping rigid cluster (ORC) tec hnique was developed by Ashmawy et al. (2003) to accurately model the angular part icles within the cont ext of geotechnical engineering. Particles crea ted using the ORC technique we re found to be accurately matching the actual particles outlines. The method does not ensure identical center of (1) (3) (5) (6) (7) (7) (8) (2) (4) (12) (11) (9) (10) (12) (13) (14) (15) (16) (17) (18) (19) (20) (1) (3) (5) (6) (7) (7) (8) (2) (4) (12) (11) (9) (10) (12) (13) (14) (15) (16) (17) (18) (19) (20) Figure 47. Sequence of addi ng discs for particle # 4 of Fraser River sand using the modified ORC method
PAGE 141
120 Table 44. Elements/Discs properties for particle # 4 of Fraser River sand X, mmY, mmR, mm Unit Density X, mmY, mmR, mm Unit Density 1 5.3726.7410.890.360.000.0016.350.26 2 2.7913.6522.430.360.9710.9015.960.26 3 0.895.8527.240.363.4912.0914.850.03 4 23.091.118.470.3611.820.8412.470.26 5 6.4422.477.860.366.2322.919.030.26 6 4.8420.3717.700.366.3214.0412.860.26 7 11.752.668.550.360.453.6418.070.26 8 16.360.7615.710.366.3214.0412.860.26 9 2.135.4826.610.364.5421.546.330.03 10 7.210.9525.160.3613.3916.807.300.52 11 3.6012.9121.590.3614.7512.165.740.26 12 14.2311.038.300.365.7619.3711.520.03 13 13.5216.469.810.364.4621.297.120.26 14 3.538.3814.330.03 15 8.972.2312.930.25 16 6.3225.063.780.26 17 6.1226.925.922.41 18 14.891.243.610.26 19 17.894.936.281.47 20 16.7618.264.330.26 Disc No. Using The ORC techniqueUsing the modified ORC technique Table 45. Results of the ORC and the modi fied ORC techniques for Fraser River sand particle # 4 PropertyThe actual particle Using the modified ORC techniqueUsing the ORC technique Area, mm21963.501963.501963.50 Eccentricity in X direction, mm00 0.07 Eccentricity in Y direction, mm00 0.39 Mass moment of inertia, unit mass.mm26.82E+056.82E+054.84E+04
PAGE 142
121 gravity and mass moment of inertia for both th e original and the cr eated particles. A modification for the ORC technique was pr esented and verified for both rounded and angular particles. The modification ensures identical mass, center of gravity, and mass moment of inertia for both the actual and th e created particles. A conditional element addition sequence was proposed in order to facilitate the future automation of the procedure. The kinetic compatibility equati ons were imposed to ensure similitude with the original particle kinetic s. The method was implemented using Fraser river sand particle # 4 as a highly angular particle a nd Michigan sand particle # 1 as a rounded particle. The modified ORC can be consid ered an important step to improve the modeling of angular particle using DEM especially in dynamic applications.
PAGE 143
122 CHAPTER FIVE: EFFECT OF PARTICLE SHAPE AND ANGULARITY ON THE DILATION OF GRANULAR SOILS 5.1 Introduction The effect of particle shape characteristics on dilatancy of granular soil is studied in this chapter through 2D numerical simula tions of pure shear load ing on particles with different shapes and angularities. A brief in troduction about dilatanc y and particle shape characterization procedures is first introduced Model parameters, procedure, and testing program are included. The effect of interparti cle friction on the stre ngth and dilatancy is also explored. 5.2 Background 5.2.1 Dilatancy in angular soils The stress strain diagram in shear for dense sand exhibits a peak followed by a reduction in the stress at large strain, whereas for loose sands, no peak usually observed. The increase of volume during shearing of dens e sands shows that dilation occurs as the test proceeds, but after a small initial co mpression. The magnitude dilation depends on soil density and confining pressure (Houlsby, 1999). Many research studies have been performed to correlate dilation and peak shear ing resistance angles. The application of the well known sawtooth model (refer to chapte r two) reveals a di rect relation between
PAGE 144
123 the two angles (Eqn. 28). Other expre ssions include those by Taylor (1948), Skempton and Bishop (1950), Bishop (1954), Newland and Allely (1957), as Schofield and Worth (1968) are presented chapter tw o (refer to Eqns. 29 through 215). Rowe (1962) theory of dilatancy along with the state dependent dila tancy were discussed in details in chapter two. 5.2.2 Dilatancy and the discrete element method The discrete element method (DEM) was used, in many research studies, to explore fundamental granular soil properties such as stre ngth and dilatancy. Among those, Ting et al. (1995) presen ted the results of a DEM study on the influence of particle shape on strength and deformation of twodimensional ellipsebased particles. The samples were isotropically compressed to the required confining pressures and then loaded in biaxial compression. The aspect ratio for the ellipseshaped particles ranged from 1:1 to 1:3. The angle of sh earing resistance increased from 26 for an aspect ratio of 1:1 to 55 for an aspect ratio of 1:3. Volumetric change was found to be consistent with that of dense sands, with s lightly higher dilation observe d in the case of elongated particles. In order to assess the relative im portance of rolling and sliding, the contact deformation was separated into two portions; due to individual particle rotation and due to particle translation. This decomposition demonstrated that for round particles, particle rotation accounts for twice as much moti on as does particle translation. Ni et al. (2000) studied the effect of microproperties of granular materials on shear strength and dilation characteristics th rough 3D discrete numerical simulations of
PAGE 145
124 the direct shear test. Partic les were simulated by means of two permanently bonded balls. Particle shape was identified using a shape factor defined as: R r R SF (51) where R and r are the radii of the larger and smaller spheres, respectively (fig. 51). The shape factor was varied between 1 (sphere) a nd 2 (2 equal spheres). A sample of 30,000 particles was used throughout the study. The top platen was simulated by two layers of spheres whereas the bottom platen was simu lated by one layer. The top platen was subjected to a constant force in order to apply the de sired normal pressure. The lower half of the shear box was moved horizonta lly using a constant velocity until a displacement of 10mm was reached. The stress ratio and volumetric dilation were plotted versus shear displacement as shown in fig. (52). The experimental results showed a higher peak strength was achieve d at smaller deformation th an the numerical model. Also the dilation of the numerical model sample was more than three times that of the Figure 51. Schematic illustration of bonded particle [Source: Ni et al. (2000)]
PAGE 146
125 Figure 52. Stress ratio and volumetric dilation at normal stress of 100kPa (30,000 particle and aspect ratio 1:3) [Source: Ni et al. (2000)] Figure 53. Effect of particle shape (normal stress = 100kPa) [Source: Ni et al. (2000)]
PAGE 147
126 real sample. The authors attributed that to the relative size of the model sample compared to the actual sample size. Both the residual strength and volumetric dilation decreased with increasing number of partic les. The authors recommended a minimum number of particles of 30,000 for the simulatio n. The results shown in fig. (53) also indicate that the shear strength increases with increasing shape factor. Volumetric dilation increased by 80% by changing the shape factor from 1 to 2. The effect of interparticle friction on the overall dilatancy an d shear strength was studied by changing the interpar ticle friction angle from 25 to 35 The peak friction angle increased, while the effect was mi nimal on the residual friction angle. The volumetric dilation increased with increasin g interparticle fric tion angle but at a decreasing rate. The authors attributed the increase in dilation to the increase of the interface between particles and the subsequent increase in the shear zone thickness. 5.3 Modification to particle shape characterization Sukumaran (1996) and Sukumaran and Ashm awy (2001) used the particle's twodimensional projection to characterize its sh ape and angularity. The outline of a twodimensional projection of a particle can be quantified numerically by discretizing the perimeter. Thus the true shape of the partic le is approximated by an equivalent polygon (Fig. 215). The sampling angle is a critical factor in calculating the shape and angularity factor in the method proposed by Sukuma ran and Ashmawy (2001). Small sampling angles should be used if particle angularity is to be captured. However, if the overall shape is the main concern, the sampling angl e should be large enough to capture the first order morphology property (shape). In the study presented by Sukumaran and Ashmawy
PAGE 148
127 (2001), the sampling angle was assumed 9 for determining both the shape and the angularity factor. This may re sult in shape and angularity f actors that are interrelated, which indicates that they are not independent parameters. In order to check the uniqueness of shape and angularity factors, the shape factors of eight different soils are plotted against their correspond ing angularity factors. The results shown in fig. (54) indicate that the tw o factors are not independent of each other. A correlation factor of (R2 = 0.4008), by fitting a straight line to the data, was observed. Sukumaran (1996) defined the Global Shape Fact or using the same definition (Eqn. 240) of the shape factor but the sampling angle was 45 The 45 sampling axes were rotated at 9 to accurately capture the particle outline Figure (55) shows a plot of the global shape factors versus the angular ity factors for the same eight soils. The correlation factor decreased to (R2 = 0.1158). In order to be able to iden tify the effect of pa rticle shape and angularity separately on the dilatancy of granular soils, the shap e and angularity factor s should be totally independent. It is clear that increasing the sampling angle ha s a significant effect on the shape factor. To this extent, a sampling angle of 90 is proposed with a 9 rotation angle (fig. 56). In addition to separating th e effect of shape and angularity, the 90 intervals can be envisioned as two perpendicular axes, which represents a true global shape measurement. The modified shape factors were calculated according to the proposed sampling angle and plotted against the angularity factors. The result s are shown in fig. (57). It is clear that th e two factors are independent si nce almost no correlation is
PAGE 149
128 Figure 54. The shape factor versus the angularity f actors for all particles R2 = 0.40080 10 20 30 40 50 60 70 0102030405060708090 Shape FactorAngularity Factor R2 = 0.15880 10 20 30 40 50 60 70 01020304050607080 Global Shape FactorAngularity FactorFigure 55. The global shape factor versus the angularity factors for all particles
PAGE 150
129 R2 = 0.02060 10 20 30 40 50 60 70 0102030405060708090100 Modified Shape FactorAngularity FactorFigure 57. The modified shape factor versus the angularity factors for all particles 9 sampling angle 45 sampling angle 9 rotation angle 90 sampling angle 9 rotation angle Particle outline Circle outline 9 sampling angle 45 sampling angle 9 rotation angle 90 sampling angle 9 rotation angle Particle outline Circle outline Figure 56. Sampling and rotation angles
PAGE 151
130 observed (R2 = 0.0206). The modified shape factor was calculated for different standard geometric shapes and resulting values are sh own in table (51). Shape factors, global shape factors, modified shape factors, and angularity factors for the different soil types that have been used in the research are summarized in table (52). For the purpose of separating the effect of shape and angularity on dilatancy, the modified shape factor and the angularity factor were considered. It should be noted that using a sampling angle of 9 the smallest asperity that can be detected is an 18 asperity for both shape and angular ity factors. The use of a certain sampling angle may cause aliasing for some part icle shapes or even for some standard shapes. For example, the modified shape f actor for the sixpoint star was 33%, while a 100% modified shape factor was determined for ei ghtpointed and twelvepointed stars. ShapeModified SF, %ge Square & Circle0.00 Rectangle (side ratio 2:1) 64.11 Rectangle (side ratio 4:1)89.33 Rectangle (side ratio 10:1)98.73 Ellipse (aspect ratio 2:1)49.89 Ellipse (aspect ratio 4:1)81.33 Sixpointed Star33.33 Eightpointed Star0.00 Twelvepointed Star0.00 Table 51. The modified shape factor fo r some standard geometric shapes
PAGE 152
Table 52. Shape factors, global shape factor s, modified shape factors, and angularity factors for soils us ed in the study MaterialParticle #12345678910111213141516171819SF 43.5648.3745.7960.2121.0540.7850.6248.2443.8929.9936.8429.1838.7551.0835.7149.4440.0844.4538.38 Global SF34.6848.4046.8269.8111.3942.4833.3039.6845.7920.6623.3727.5628.0247.7422.3251.7027.3837.6432.48 Modified SF39.4366.7462.5993.2130.4445.0250.2937.8465.3729.2244.3536.9449.1254.7443.9472.4041.0547.5939.33 AF18.5215.2318.7718.257.989.7928.8431.9714.3612.2517.938.1421.6025.1519.8526.7431.1220.3924.87 SF40.7739.2551.2655.3252.8634.1143.8950.8736.5759.0347.24 Global SF43.3534.8247.8553.9848.6330.2137.6633.2822.4344.4048.36 Modified SF54.8754.9350.1373.8145.8438.6640.1536.9738.4748.6951.07 AF29.1413.6837.2822.8418.2614.9325.2439.9617.6650.7844.42 SF24.4926.0717.9339.4232.7236.4719.9031.1439.8627.29 Global SF21.5623.5116.8936.8434.6437.1514.6430.9736.7621.96 Modified SF53.4035.3846.8651.8132.7073.4256.0633.2641.0134.52 AF5.315.273.599.614.636.635.957.2010.6411.16 SF40.2141.0046.0955.7141.6848.9948.0938.5735.2345.7456.2238.5658.4135.38 Global SF32.4244.2144.7657.0134.3650.0940.1932.4126.2243.7452.0534.3846.4732.37 Modified SF45.8457.6859.8664.5439.0274.3240.0052.7443.6633.3757.7837.5554.1746.63 AF27.6816.9820.8022.2317.3721.8661.9426.3830.2226.8824.6329.5844.169.31 SF41.6541.1728.9338.5327.5240.5847.2742.9648.0142.6027.50 Global SF37.3738.7325.9334.8819.7133.2344.5343.7341.5945.5014.24 Modified SF40.7756.1437.7343.2139.3036.3864.9265.6362.4949.2132.02 AF16.0015.8310.9311.9911.0113.9016.1811.8024.479.4523.21 SF49.6562.4447.8176.6039.5939.7348.1554.3737.0750.3849.6148.3050.7948.6967.1154.1140.5053.2645.21 Global SF45.9562.5237.7565.3631.1427.5043.9351.5028.1545.1736.2741.0941.5430.2770.3948.7030.4454.5132.95 Modified SF44.0362.6450.9262.2040.9038.9664.9051.7042.1847.9942.9250.0354.8636.6970.6569.9444.1860.2451.34 AF35.6134.8627.5333.4121.3532.8625.1722.1817.3740.4726.6930.2333.0442.7719.6330.9815.3544.5818.10 SF42.3627.9841.2548.3827.2056.7546.7241.8433.7131.46 Global SF38.6919.4043.3142.2717.2156.9456.9429.6227.1120.75 Modified SF39.0234.3732.9845.3355.8956.2738.6536.5941.3644.38 AF8.2410.5213.0821.6612.2825.3413.9119.0912.9116.53 SF38.1826.3522.9633.9836.54 Global SF28.6820.4018.4129.1735.88 Modified SF39.3331.9529.9743.4350.43 AF12.438.687.0411.086.48Daytona Beach Sand Ottawa # 45 Ottawa # 20/70 Syncrude Tailings Sand Fraser River Sand Michigan Dune Sand Ottawa # 90 Ottawa 0.1 mm 131
PAGE 153
132 5.4 Particle groups with resp ect to shape and angularity Eight different soils that in clude 99 particles were used in the numerical model. The modified shape factors for all partic les varied between 29.22 and 93.21 with an average value of 48.32, whereas the angular ity factors varied between 3.59 and 61.94 with an average value of 20.73. The particles were divided into groups according to their shape and angularity factors. Particles belong to different soils could be classified in the same shape or angularity group. The propert ies of different groups are summarized in table (53). For each group, a subroutine was written to PFC2D using the program language Â“fishÂ” in order to create the corr esponding clumps. The clumps were created according to the overlapping rigid cluster method (Ashmawy et al., 2003). 5.5 Numerical simulations Typical loading conditions for evaluating soil strength and dilatancy properties Table 53. Shape and angularity f actors for different groups Particles Group SF range, %ge SF average, %ge AF range, %ge AF average, %ge AF130 7340.0005 107.50 AF232 7454.0020 2522.50 AF337 6050.0040 4542.50 SF130 3532.5005 2713.00 SF250 5552.5005 4425.00 SF370 7572.5007 2720.00 SF1(AF 1525)30 3532.5015 2520.00 SF2(AF 1525)50 5552.5015 2520.00 SF3(AF 1525)70 7572.5015 2520.00
PAGE 154
133 are shown in fig. (58). Only the stress increments are illustra ted in the figure. Among those, the triaxial test is the most commo n laboratory method used in both practice and research to determine soil strength parameters. The sample is initially subjected to a cell pressure and the veridical stress is increased until failure. As stated, dilatancy (d) can be defined, with respect to triaxial or biaxial testing conditions, as the ratio of plastic volumetric st rain increment to the plastic deviator strain increment. Conceptually, dilatancy is the ch ange of volume that co rresponds to shearing of granular soils (Reynolds, 1885). Howeve r, it is known that there is a volumetric change associated with the increase of confini ng pressure during the tria xial test. To this end, the volumetric strain measured during triaxi al testing of granular soils is attributed not only to changes in shear stre sses, but also to variations in mean total stress. In order to obtain volumetric strain that is solely due to shearing the soil (i.e. Dilatancy), either simple shear (fig. 57c) or pure shear (fig. 57d) testing conditions should be applied. On the other hand, results from direct shear test s (fig. 57b) are considerably biased due to boundary conditions and a stress path that are ha rd to interpret or analyze. Although the simple shear test (Fig. 57c) is difficult to perform, it has the a dvantage of applying pure shear loading conditions, so dilata ncy can be extracted directly from the test results. An identical stress path (at 45 angle) can be achieved in a triaxial configuration by simultaneously varying the vertical and hor izontal stresses in equal magnitudes and opposite directions. This is the pure shear loading sequence (Fig. 57d) that was used throughout the current simulations to study the effect of part icle shape and angularity on dilatancy.
PAGE 155
134 5.5.1 Sample preparation The particles were assumed to be of unit thickness (2D simulations). The particles were generated by radi us expansion in which a partic le size range, porosity, and sample size are specified. The number of the particles was determined to satisfy the sample dimensions, particle size range, a nd target porosity. The walls (boundary conditions) used to confine and load the samp le are generated longer than the actual size to account for large strains ex pected during the test After initial compaction, the lateral walls were given stiffnesses that are onetenth of the particle stiffness in order to simulate a soft boundary to ensure uniform stresses. The sample dimensions and properties are summarized in table (54). After the circular particles were created satisfying the target porosity, the angular substitution scheme was i nvoked to regenerate the angular particles according to the desired particle shape and angularity groups. 5.5.2 Sample consolidation to the desired confining pressure After the angular particles (clumps) were created, strains were determined by tracking the position of the walls. Stresses on each wall were computed by dividing the total force acting on the wall by the sample length/width. Stresses acting on the sample in each direction were obtained by taking the average of the stresses acting on each set of opposing walls. The sample was isotropically consolidated to th e desired confining pressure by moving the walls inward with a c ontrolled velocity. A servo Â“fishÂ” function was used to control the wall velocities in both horizontal and vertic al directions. The wall velocity is always proportional to th e difference between the current and target confining pressure. The conso lidation procedure usua lly involves volumetric change that
PAGE 156
135 Figure 58. Typical loading conditions for strength and dilatancy testing: (a) triaxial; (b) direct shear; (c) simple shear; (d) pure shear [Source: Ashmawy et al. (2003)] Table 54. Model dimens ions and properties Sample height, m12 Sample width, m6 Heighttowidth ratio2 Number of partic les2574 / 2394 / 2215 Minimum radius, m0.075 Maximum radius, m0.10 Confining pressure, MPa1.00 / 0.50 / 0.25 Density, kg/m31000 Normal stiffness, N/m5.00E+08 Shear stiffness, N/m5.00E+08 Interparticle friction coefficient0.50 Porosity0.14 / 0.20 / 0.26
PAGE 157
136 may be expansion or compression, depending on the desired confining pressure and porosity. The volume change corresponds to ap plying the confining pressure on circular particles was determined and saved. In order to ensure a unique porosity after the consolidation stage for all soils, samples that involve a volumetric change different from that of the discs were brought back to the same porosity by performing cyclic consolidation. After consolidation, the current conditions were stored as initial conditions in order to begin loading stage. 5.5.3 Pure shear loading The failure surfaces, within the context of critical state soil mechanics, in qp\space are shown in fig. (59). For sands, th e tension cutoff is a vertical wall at p/ = 0. Pure shear loading (confining pressure is c onstant throughout the test) was used in the numerical simulations. Initial conditions were expressed by (p0 \, 0, q0 = 0), where p0 \ is the confining pressure, 0 is the initial specific volume (1 +e), and q0 is the initial shear stress. The stress path for pure shear is an inclined path in the constant p0 \ plane, as shown in fig. (59). Pure sh ear loading is continued until the stress path hits the Hvorslev surface at (qpeak, peak) and then the strength decreases along HvorslevÂ’s surface within the same constant confining pressure plane, to the critical state line. It should be noted that, if the sample initial conditions is close enough to the CSL such that the stress path hits RoscoeÂ’s surface or the CSL, no peak will be defined, as shown later in some numerical simulations. In the numerical simulations, pure shear loading was applied by releasing the top and bottom platens from the ser vo control function and apply constant inward velocity
PAGE 158
137 Figure 59. Stress path duri ng pure shear testing (within the context of critical q v P\CSL NCLNo tension surface Hvorslevsurface Roscoe surface qpeakvpeakv0P0 \vCSL qCSL q v P\CSL NCLNo tension surface Hvorslevsurface Roscoe surface qpeakvpeakv0P0 \vCSL qCSL
PAGE 159
138 (compression loading). In order to ensure a stable solution, the velocity was increased gradually until the target velocity is achieve d. The top and bottom platens then continued to move at a constant velocity to the e nd of the test. Throughout the solution, the increase in the vertical stress was measured and target horizontal stress was set to decrease by the same amount, which ensures pure shear loading condition. The loading terminated at axial strain around 30% to ensure that critical state have been reached. In most tests, the critical state was reached befo re 30% axial strain but the loading continued to ensure that instability didnÂ’t occurred at very high strains. The confining pressure, shear stress, axial strain, and volumetric strain were tracked and saved using the Â“historyÂ” feature of the PFC2D for further analysis and comparison. 5.5.4 Simulation program The objective of the study is to explore th e dependency of dilatancy on particle shape and angularity. Pure shear testing was fi rst performed on the circ ular particles as a reference for comparison with different groups of angular particles. The pure shear test was then performed on groups AF1, AF2, AF3, SF1, SF2, SF3, SF1(AF 1525), SF2(AF 1525), and SF3(AF 1525) under the same conditi ons. Different sample states should be considered to evaluate the effect of particle shape and angularity on strength and dilatancy. The state paramete r introduced by Been and Jeffer ies (1986) is an established parameter to define the sample position with re spect to the critical state line. The state parameter determines how far the current c ondition is from the critical state. The numerical simulations were perfor med on the same shape and angularity groups with three different in itial state parameters. Conf ining pressures of 1.0MPa,
PAGE 160
139 500KPa, and 250KPa were used. The corresponding porosities were 0.14, 0.2, and 0.26, respectively. In all simulations, the interparti cle friction coefficient was set to 0.5, which corresponds to an angle of in terparticle friction of 25.56 (within the recommended range of 5 Â– 40 ). The performed numerical simulations were summarized in table (55). The effect of friction angle was also studi ed in order to explore whether changing the friction coefficient, for a ci rcular particles, can be utili ze as an alterna tive to modeling angular particles. For that purpose, the test was performed on the circular particle using friction coefficients in the range of 0.25 Â– 2. 50, which correspond to interparticle friction angles of 14 68 GroupP/= 1.0 MPa & n = 0.14 P/= 500 KPa & n = 0.20 P/= 250 KPa & n = 0.26 AF1 AF2 AF3 SF1N/AN/A SF2N/AN/A SF3N/AN/A SF1 (AF 1525) SF2 (AF 1525) SF2 (AF 1525) Discs (f = 0.25) Discs (f = 0.50) Discs (f = 0.75) Discs (f = 1.00) Discs (f = 1.25) Discs (f = 1.50) Discs (f = 1.75) Discs (f = 2.00) Discs (f = 2.25) Discs (f = 2.50) Table 55. Simulation Program
PAGE 161
140 5.6 Numerical results 5.6.1 Confining pressure (1.0MPa), voids ratio (0.162), and 2574 particles Figure (510) shows the shear stressshear strain results for pure shear test performed on circular particles. The result s define a peak shear strength of 500 kPa which occurred at a shear strain of 4.0% foll owed by a decrease in th e shear strength to a residual value of 400 kPa. The volumetric stra in was plotted against the shear strain and shown in fig. (511). The plot clarifies th at dilation occurs thr oughout the test with the maximum rate of dilation being at low shear strains (3.0% 8.0%). Beyond the peak strength, the rate of dilation (dilatancy) begins to decrease until the critical/steady state is reached, at which the rate of dilation is zer o. Figure (511) shows that zero rate of dilation occurred at shear strain of 13% a nd remain almost constant. The maximum volumetric strain was in the range of 3.2% to 3.72%. The maximum rate of dilation was determined by performing regression analysis to the portion of vol umetric strainshear strain curve that corresponds to shear strain of 3.0% to 8. 0% and a value of 0.329 was determined. Moving average was used in presenting the numerical results. 5.6.1.1 Effect of angularity The shear stressshear strain and volumetric strainshear strain curves for groups AF1, AF2, and AF3 are shown in figs. (512) and (513). The results show that the peak strength increases with incr easing angularity factor. The peak strength in groups AF2 and AF3 occurred at smaller strains than in th e circular particles. The residual/critical state shear strength increased compared to circular particles. Groups AF2 and AF3
PAGE 162
141 Figure 510. Shear stressshear strain plot for circular particles 0 100 200 300 400 500 600 0102030405060 Shearl Strain, %Shear Stress, KPa Confining Pressure = 1.0 MPa Voids Ratio = 0.162Figure 511. Volumetric strainshear stra in plot for circular particles 0 1 1 2 2 3 3 4 4 0102030405060 Shear Strain, %Volumetric Strain, % Confining Pressure = 1.0 MPa Voids Ratio = 0.162
PAGE 163
142 showed approximately the same critical state shear strength, which implies that the effect of angularity on critical state sh ear strength decreases at high angularity factors. The plot of volumetric strainshear stra in (fig. 513) shows that th e volumetric strain increased dramatically for group AF3 (alm ost twofold), while about 25 % increase was noticed for group AF2, which clarifies the ef fect of angularity on dilation. The large increase in the dilation for group AF3 compared to that of group AF2 can be explained by recalling that both of the groups has almost the same crit ical state shear stre ngth but group AF3 has higher peak strength (fig. 512). The excess strength of AF3 over AF2 (the difference between the two peaks) is due to the differen ce in dilation indicated by the increase of volumetric strain for group AF3. The volumet ric strain of group AF1 was approximately the same as that of circular particles, which can be attributed to the very low angularity of group AF1. The maximum dilation angle Â“ Â” (the slope of the vol umetric strainshear strain curve is equal to sin ) in case of AF1was lower than that of circular particles, which can attributed to the stre ss strain behavior of both of them. The circular particles showed a clearly defined peak followed by a residual strength, which is 20% lower than the peak strength, so a reasonable maximum d ilation angle is observed. However, in case of AF1, the increase of the peak over the re sidual strength was 12% so a lower maximum dilation angles was observed. An almost id entical maximum dilation angle was observed in case of AF2, while, a noticeable in crease was noticed for group AF3. The angularity factor was plotted against the maximum volumetric strain. The results are shown in fig. (514), which shows a remark able increase of maximum volumetric strain with angularity.
PAGE 164
143 Figure 512. Shear stressshear strain curves for different angularity groups 0 100 200 300 400 500 600 700 800 0102030405060 Shear Strain, % Shear Stress, KPa AF1 AF2 AF3 Discs Confining Pressure = 1.0 MPa Voids Ratio = 0.162Figure 513. Volumetric strainshear strain curves for different angularity groups 0 1 2 3 4 5 6 7 8 9 0102030405060 Shear Strain, %Volumetric Strain, %ge AF1 AF2 AF3 Discs Confining Pressure = 1.0 MPa Voids Ratio = 0.162
PAGE 165
144 Figure 514. Maximum volumetric strain versus angularity factor for different angularity groups 3 4 5 6 7 8 9 051015202530354045Angularity Factor (AF), %Maximum Volumetric Strain, % Confining Pressure = 1.0 MPa Voids Ratio = 0.162Figure 515. Maximum dilation and residua l shearing resistance angles versus angularity factor for different angularity groups 0 10 20 30 40 50 60 70 051015202530354045 Angularity Factor (AF), %Angles, degree Maximum dilation angle Residual shearing resistance angle Confining Pressure = 1.0 MPa Voids Ratio = 0.162
PAGE 166
145 The correlation between maximum volumetric st rain and angularity factor (at SF 40 54) can be written as: ) AF ( 02 0 max ve 52 3 ) ( [R2 = 0.9787] (51) A plot of the angularity factor versus bot h the maximum dilation angle and the residual shearing resistance angle is shown fig. (515). The maximum dilation angle decreased for group AF1 and increased for groups AF2 and AF3 compared with that of the circular particles. The residual shearing resistance increased at a decreasing rate and then remained constant ( res)AF2 = ( res)AF3. The results show that dilatancy is strongly affected by the angularity. Specifically, the angularity fact or has a solid correlation w ith strength and dilatancy properties of granular soils. 5.6.1.2 Effect of particle shape The shear stressshear strain and volumetric strainshear strain relationships for groups SF1, SF2, and SF3 along with the circular particles are plotted in figs. (516) and (517). The results show an increase of both the peak and residual strengths compared with circular particles. The volumetric st rain was almost identical for group SF1 and circular particles, while it was higher fo r groups SF2 and SF3. The maximum volumetric strain of group SF3 was slightly more than that of group SF2. According to the difference in the shape factor between groups SF2 and SF3, it was expected that group SF3 dilates much more than group SF2, which didnÂ’t happen. This may be attribute to the difference in angularity between the pa rticle that constitute the two groups (as
PAGE 167
146 discussed later). The maximum dilation angles or the three groups were almost the same, which implies that the effect of shape is not clea rly defined. It is cl ear that the results are not consistent and there is no clear trend can be drawn, unlike the observations for different angularity groups. This may be attri buted to the fact that the shape factors for groups AF1, AF2, and AF3 are in the same rang e. In contrast, the angularity factor for different shape groups is not necessarily within the same range. Therefore, there may be a combined influence of shape and angularity. The shape factor was plotted against the ma ximum volumetric strain for the three groups along with the circular particles. Figure (518) shows that the maximum volumetric strain increases with the incr ease of shape factor, but a sh arp decrease in the increasing rate was noticed at high shape factors. The effect of particle shape on both the maximum dilation and the residual shearing resistance angles is shown in fig. (519). The results show that the eff ect of shape is not as clear as that of the angularity. As stated before, this may be attributed to interference of the influence of shape and angularity. To overcome this in terference and clearly identify the effect of particle shape on dilatancy, three more shape groups were created such that the range of angularity factors was identical for all of them. These groups are SF1(AF 1525), SF2(AF 1525), and SF3(AF 1525). The results for these grou ps are shown in figs (520) through (523). The peak strengths were almost identical for the three shape groups with an increase of 25% of that of the circular particles. However, the residual strengths were higher for higher shape factors. The results imply that the effect of shape factor on the peak strength is minimal compared with that of angularity factor.
PAGE 168
147 The maximum volumetric strain increased w ith increasing the shape factor (fig. 522). The maximum volumetric strain can be expressed in terms of shape factor (AF 1525) as follows: ) SF ( 0083 0 max ve 57 3 ) ( [R2 = 0.9279] (52) Figure (223) shows that the maximum dilation angle increased slightly for SF1(AF1525) and decreased for both SF2(A F1525) and SF3(AF1525). It can be correlated to the shape factor (AF 1525) as follows: 258 19 SF 2022 0 SF 0038 02 max [R2 = 0.998] (53) Residual shearing resistance angles were plotted versus the shape factor and the results are shown in fig. (523). It can be correlated to the shap e factor (AF 1525) as follows: ) SF ( 005 0 residuale 482 21 [R2 = 0.9531] (54) The shape factor groups with constant av erage angular shape (t he modified shape factor groups) had more clear effect on dila tancy and strength properties. The steady state for particles with very high shape factor was reached after large strains, while for lower shape factors steady state was es tablished at lower strain values.
PAGE 169
148 Figure 516. Shear stressshear strain cu rves for different shape groups 0 100 200 300 400 500 600 700 010203040506070 Shear Strain, % Shear Stress, KPa SF1 SF2 SF3 Discs Confining Pressure = 1.0 MPa Voids Ratio = 0.162Figure 517. Volumetric strainshear strain curves for different shape groups 0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 01020304050607080 Shear Strain, %geVolumetric Strain, % SF1 SF2 SF3 Discs Confining Pressure = 1.0 MPa Voids Ratio = 0.162
PAGE 170
149 Figure 518. Maximum volumetric strain versus shape factor for different shape groups 3.0 3.5 4.0 4.5 5.0 5.5 6.0 01020304050607080Shape Factor (AF), %Maximum Volumetric Strain, % Confining Pressure = 1.0 MPa Voids Ratio = 0.162Figure 519. Maximum dilation and residua l shearing resistance angles versus angularity factor for different shape groups 0 5 10 15 20 25 30 35 40 01020304050607080 Shape Factor (SF), %Angle, degree Maximum dilation angle Residual shearing resistance angle Confining Pressure = 1.0 MPa Voids Ratio = 0.162
PAGE 171
150 Figure 520. Shear stressshearstrain curves fo r different modified sh ape factor groups 0 100 200 300 400 500 600 700 0102030405060708090 Shearl Strain, %Shear Stress, KPa SF1 (AF 1525) SF2 (AF 1525) SF3 (AF 1525) Discs Confining Pressure = 1.0 MPa Voids Ratio = 0.162Figure 521. Volumetric strainshear strain curves for different modified shape factor groups 0 1 2 3 4 5 6 7 8 0102030405060708090 Shear Strain, %Volumetric Strain, % SF1 (AF 1525) SF2 (AF 1525) SF3 (AF 1525) Discs Confining Pressure = 1.0 MPa Voids Ratio = 0.162
PAGE 172
151 Figure 522. Maximum volumetric strain ve rsus shape factor for different modified shape factor groups Figure 523. Maximum dilation and residua l shearing resistance angles versus angularity factor for different m odified shape factor groups 3.0 3.5 4.0 4.5 5.0 5.5 6.0 6.5 7.0 7.5 01020304050607080Shape Factor (SF), %Maximum Volumetric Strain, % Confining Pressure = 1.0 MPa Voids Ratio = 0.162 0 5 10 15 20 25 30 35 40 01020304050607080 Shape Factor (SF), %Angle, degree Maximum dilation angle Residual shearing resistance angle Confining Pressure = 1.0 MPa Voids Ratio = 0.162
PAGE 173
152 5.6.2 Confining pressure (500 KPa), voids ratio (0.25), and 2394 particles 5.6.2.1 Effect of angularity The shear stressshear strain and volumetric strainshear strain curves for groups AF1, AF2, and AF3 are shown in figs. (524) and (525). The peak shear strength was observed neither for the circular particle nor for the angularity groups (the sample initial condition is very close to the CSL, as expl ained later). The resi dual strength was higher for group AF1 than for circular particles, wh ile the increase for AF2 and AF3 was very close and higher than group AF1. The volumet ric strain increased gradually with the increase in angularity. Group AF2 needed 50% shear strain to reach a constant volumetric strain. Figures (526) and (5 37) show the depende ncy of the maximum volumetric strain, the maximum dilation angle, and the residual shearing resistance angle on particle angularity. Th e correlation between maxi mum volumetric strain and angularity factor (SF 40 54) can be written as: 186 1 AF 1525 0 AF 0012 0 ) (2 max v [R2 = 0.9824] (55) The maximum dilation angle can be correlate d to the angularity factor (SF 4054) as follows: ) AF ( 0386 0 maxe 2157 3 [R2 = 0.8142] (56) The residual shearing resistance angle increased at a decreasing rate with angularity (fig. 527). The co rrelation of the residual shear ing resistance angle to the angularity factor (SF 4054) may be written as follows:
PAGE 174
153 049 25 AF 6063 0 AF 0097 02 residual [R2 = 0.9871] (57) 5.6.2.2 Effect of particle shape The shear stressshear strain and volumetric strainshear strain relationships for groups SF1(AF 1525), SF2(AF 1525), and SF3(AF 1525) along with the circular particles are plotted in figs. (528) and (529). Again, the peak shear strength could not be defined for the circular particle nor for the modified shape groups. The behavior of groups SF2(AF 1525) and SF3(AF 1525) was almost identical and the steady state shear strength was not reached even after al most 60% shear strain. No clear trend was observed for the volumetric strain of the th ree groups. All of them exhibited more volumetric strain than that of circular partic les. The shape factor was plotted against the maximum volumetric strain for the three groups along with the circular particles. Figure (530) shows that the maximum volumetric strain increases with the shape factor according to the following correlations: ) SF ( 0164 0 max ve 07 1 ) ( [R2 = 0.9606] (58) or 0195 1 SF 0309 0 ) (max v [R2 = 0.9865] (59) The maximum dilation angle increased fo r SF1(AF 1525) and then decreased for both SF2(AF 1525) and SF3(AF 1525) as shown in fig. (531), which confirm the same trend observed under [ p/ = 1.0 MPa and n = 0.14 ]. Figur e (531) shows th at the residual angle increased for SF1(AF 1525) and almo st identical increase was observed for SF2(AF 1525) and SF3(AF 1525). A correlation can be assumed between the residual shearing resitance angle and the shape factor (AF 1525) as follows:
PAGE 175
154 Figure 524. Shear stressshear strain curves for different angularity groups 0 50 100 150 200 250 300 350 400 010203040506070 Shear Strain, %ge Shear Stress, KPa AF1 AF2 AF3 Discs Confining Pressure = 500 KPa Voids Ratio = 0.25Figure 525. Volumetric strainshear strain curves for different angularity groups 0 1 2 3 4 5 6 7 8 010203040506070 Shear Strain, %Volumetric Strain, % AF1 AF2 AF3 Discs Confining Pressure = 500 KPa Voids Ratio = 0.25
PAGE 176
155 Figure 526. Maximum volumetric strain versus angularity factor for different angularity groups 0 1 2 3 4 5 6 051015202530354045Angularity Factor (AF), %Maximum Volumetric Strain, % Confining Pressure = 500 KPa Voids Ratio = 0.25Figure 527. Maximum dilation and residua l shearing resistance angles versus angularity factor for different angularity groups 0 5 10 15 20 25 30 35 40 45 50 051015202530354045Angularity Factor (AF), %Angle, degree Maximum dilation angle Residual shearing resistance angle Confining Pressure = 500 KPa Voids Ratio = 0.25
PAGE 177
156 Figure 528. Shear stressshear strain curves fo r different modified sh ape factor groups 0 50 100 150 200 250 300 350 400 0102030405060 Shear Strain, %Shear Stress, KPa SF1 (AF 1525) SF2 (AF 1525) SF3 (AF 1525) Discs Confining Pressure = 500 KPa Voids Ratio = 0.25Figure 529. Volumetric strainshear strain curves for different modified shape factor groups 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 0102030405060 Shear Strain, %Volumetric Strain, % SF1 (AF 1525) SF2 (AF 1525) SF3 (AF 1525) Discs Confining Pressure = 500 KPa Voids Ratio = 0.25
PAGE 178
157 Figure 530. Maximum volumetric strain ve rsus shape factor for different modified shape factor groups 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 01020304050607080Shape Factor (SF), %Maximum Volumetric Strain, % Confining Pressure = 500 KPa Voids Ratio = 0.25Figure 531. Maximum dilation and residua l shearing resistance angles versus angularity factor for different mo dified shape fact or groups 0 5 10 15 20 25 30 35 40 45 50 01020304050607080 Shape Factor (SF), %Angle, degree Maximum dilation angle Residual shearing resistance angle Confining Pressure = 500 KPa Voids Ratio = 0.25
PAGE 179
158 ) SF ( 0046 0 residuale 165 25 [R2 = 0.9038] (510) 5.6.3 Confining pressure (250 KPa), voids ratio (0.35), and 2215 particles 5.6.3.1 Effect of angularity The shear stressaxial strain and volumetric strainshear strain curves for groups AF1, AF2, and AF3 are shown in figs. (532) and (533). The peak shear strength was not observed for the circular particle nor for the modified shape groups except for group AF3, which shows a slightly defined peak st rength. The residual strength increased for all group with almost identical increase for groups AF2 and AF3. The volumetric strain increased dramatically for AF1 and AF2 and increased even further for group AF3. Group AF2 needed 60% shear strain to reach a constant volumetric strain. Figures (534) and (535) show the de pendency of the maximum volumetric strain, the maximum dilation angle, and the residual shearing resi stance angle on particle angularity. The maximum volumetric strain increased at a d ecreasing rate as the modified shape factor increased. The maximum dilati on angle may be correlated to the angularity factor (SF 4054) as follows: 852 1 ) AF ( 2185 0max [R2 = 0.806] (511) The residual angle increased with almost identical values observed for AF2 and AF3. A correlation of the residual angle to the angularity factor (SF 4054) may be written as: 755 32 ) AF ( 4159 0 ) AF ( 0065 02 residual [R2 = 0.9959] (512)
PAGE 180
159 5.6.3.2 Effect of particle shape The shear stressshear strain and volumetric strainshear strain relationships for groups SF1(AF 1525), SF2(AF 1525), and SF3(AF 1525) along with the circular particles are plotted in figs. (536) and (537). No peak shear strength was defined except for group SF1(AF 1525), where a slightly de fined peak strength is observed. Group SF3(AF 1525) underwent continuous increase in shear strength througho ut the test. The volumetric strain behavior for groups SF1(AF 1525) and SF2(AF 1525) was almost identical and exhibited about thre e times the final volumetric stra in of circular particles. Group SF3(AF 1525) reached the steady state with respect of volumetric strain after a shear strain of 65%. The modified shape factor was plotted ag ainst the maximum volumetric strain for the three groups along with the circular particles. Figure (5 38) shows that the maximum volumetric strain increases with the shape factor (AF 1525) according to the following correlations: 8678 0 ) SF ( 046 0 ) (max v [R2 = 0.9583] (513) Figure (539) shows that the maximum dilation angle increased for SF1(AF 1525) and then decreased for both SF2(AF 1525) and SF3(AF 1525) according to the following correlation: 1649 3 ) SF ( 1203 0 ) SF ( 0019 02 max [R2 = 0.842] (514) Figure (539) shows a clear trend fo r the increase of the residual shearing resistance angle with the increase of sh ape factor, which can be expressed as:
PAGE 181
160 Figure 532. Shear stressshear strain curves for different angularity groups 0 50 100 150 200 250 01020304050607080 Shear Strain, % Shear Stress, KPa AF1 AF2 AF3 Discs Confining Pressure = 250 KPa Voids Ratio = 0.35Figure 533. Volumetric strainshear strain curves for different an gularity groups 1 0 1 2 3 4 5 6 7 01020304050607080 Shear Strain, %Volumetric Strain, % AF1 AF2 AF3 Discs Confining Pressure = 250 KPa Voids Ratio = 0.35
PAGE 182
161 Figure 534. Maximum volumetric strain versus angularity for different angularity groups 0.0 1.0 2.0 3.0 4.0 5.0 6.0 051015202530354045Angularity Factor (AF), %Maximum Volumetric Strain, % Confining Pressure = 250 KPa Voids Ratio = 0.35Figure 535. Maximum dilation and residua l shearing resistance angles versus angularity factor for different angularity groups 0 10 20 30 40 50 60 051015202530354045 Angularity Factor (AF), %Angle, degree Maximum dilation angle Residual shearing resistance angle Confining Pressure = 250 KPa Voids Ratio = 0.35
PAGE 183
162 Figure 536. Shear stressshear strain curves for different modified shape groups 0 50 100 150 200 250 051015202530354045 Shear Strain, %Shear Stress, KPa SF1 (AF 1525) SF2 (AF 1525) SF3 (AF 1525) Discs Confining Pressure = 250 KPa Voids Ratio = 0.35Figure 537. Volumetric strainshear strain cu rves for different modified shape groups 1 0 1 2 3 4 5 01020304050607080 Shear Strain, %Volumetric Strain, % SF1 (AF 1525) SF2 (AF 1525) SF3 (AF 1525) Discs Confining Pressure = 250 KPa Voids Ratio = 0.35
PAGE 184
163 Figure 538. Maximum volumetric strain ve rsus shape factor for different modified shape groups 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 01020304050607080Shape Factor (SF), %Maximum Volumetric Strain, % Confining Pressure = 250 KPa Voids Ratio = 0.35Figure 539. Maximum dilation and residua l shearing resistance angles versus angularity factor for different modified shape groups 0 10 20 30 40 50 60 01020304050607080 Shape Factor (SF), %Angle, degree Maximum dilation angle Residual shearing resistance angle Confining Pressure = 250 KPa Voids Ratio = 0.35
PAGE 185
164 761 30 ) SF ( 1581 0residual [R2 = 0.9945] (515) 5.7 Critical state line and pa rticle shape characteristics The previous results confirmed that the strength and dilatancy properties of granular soils depend on particle shape and angularity. The changes observed on dilatancy characteristics of different groups are not attribut ed solely to particle angularity or shape. Confining pressure and the corre sponding initial void ratio contributed to the behavior as well. To isolate the effect of particle shape an d angularity, the critical state lines for different groups were identified in the vp/ space. The specific volume was plotted against the confining pressure for ci rcular particles and fo r different angularity groups, and the results are shown in fig. (540). The results show that the cr itical state lines for differe nt angularity groups were shifted above that of circular particles, whic h confirms that the dilatancy increases with angularity. The results also explained the reas oning behind the absence of a peak in some stressstrain curves when c onfining pressures of 250KPa or 500 KPa were applied. The initial conditions were very cl ose to the critical st ate line such that th e stress path during pure shear intersected with Hvorslev surface ve ry close to the critical state line. The critical state lines for di fferent modified shape factor groups are shown in fig. (541). The critical state line for the first and second groups were al most identical, while that of the third group was shifted considerably. Again, the result confirmed the dependency of dilatancy on particle shape. It can be concluded that the angularity factor had a more remarkable effect on dilatancy than the shape factor.
PAGE 186
165 Figure 540. Critical state lines for different angularity groups 1.15 1.20 1.25 1.30 1.35 1.40 1.45 1.50 20030040050060070080090010001100 Confining pressure, KPaSpecific volum e Initial conditions Circular particles Group AF1 Group AF2 Group AF3Figure 541. Critical state lines for different shape groups 1.15 1.20 1.25 1.30 1.35 1.40 1.45 1.50 20030040050060070080090010001100 Confining pressure, KPaSpecific volum e Initial conditions Circular particles Group SF1(AF 1525) Group SF2(AF 1525) Group SF3(AF 1525)
PAGE 187
166 5.8 Effect of interparticle friction Interparticle friction is r ecognized as an important pr operty that affects strength and dilatancy properties of the granular so ils. The interparticle friction angle was correlated to the critical state shearing resistance angle in different studies (Caquot, 1934; Bishop, 1954). Since the critical state shearing resistance along with the dilation constitute the overall shearing resistance of th e granular soils, the interparticle friction angle is envisioned as a critic al property to the global shear strength. A numerical study of the effect of interparticle friction angl e on the dilatancy and strength properties of angular soils was performed. The input parame ters are those shown in table (54). The interparticle friction coefficient was changed in the range 0.25 2.50 (step = 0.25). The study was performed under the same confining pressures used earlier in this chapter (250KPa, 500KPa, and 1.0 MPa). 5.8.1 Confining pressure (1.0 MPa), voids ratio (0.162), and 2574 particles The shear stressshear strain and volumetric strainshear strain curves are shown in fig. (542) and (543). The results show that the peak strength increases with increasing interparticle friction, with the peak occurring at small st rains. The rate of increase decreases for higher interparticle frictions. The critical state shear strength slightly increased with increasing interpar ticle friction and was almost identical for interparticle friction coefficients higher than 1.75. The volumetric strainshear strain plot (fig. 543) shows that the volumetric strain increases with the increase of the interparticle friction but again, with decreasing rate for higher values. For interparticle frictions
PAGE 188
167 higher than 2.0, almost iden tical maximum volumetric stains were observed. The maximum dilation angle dependent dramatically on interparticl e friction, but was almost identical for values higher than 1.50. Th e effect of interparticle friction on maximum volumetric strain, maximum dilation angle, and residual shearing resistance angle are summarized in figs. (544) and (545). Th e behavior can be expressed by fitting trend lines to the data as follows: 47 5 ) f ln( 289 2 ) (max v [R2 = 0.9842] (516) 213 36 ) f ln( 411 21max [R2 = 0.9734] (517) 7743 0 max) f ( 4 31 [R2 = 0.9622] (518) 117 24 ) f ln( 337 4residual [R2 = 0.9819] (519) 1888 0 residual) f ( 883 23 [R2 = 0.9755] (520) Where f is the interparticle friction coefficient (f = tan ) and is the interparticle friction angle. 5.8.2 Confining pressure (500 KPa), voids ratio (0.25), and 2394 particles The shear stressshear strain and volumetric strainshear strain curves are shown in fig. (546) and (547). The results show that the peak was not defined for friction coefficient 0.25, 0.50, 0.75, and 1.00, wher eas for higher friction coefficients, peak strengths are observed and they increased with increasing the interparticle friction. The rate of increase decreases for higher interp article frictions. Th e critical state shear strength slightly increased with increasing interparticle friction and was almost identical
PAGE 189
168 Figure 542. Shear stressshear strain curves for different interparticle friction coefficients (circular particles) Figure 543. Volumetric strainshear strain curves for different interparticle friction coefficients (circular particles) 0 100 200 300 400 500 600 700 800 0102030405060 Shear Strain, %Shear Stress, KPa f = 0.25 f = 0.50 f = 0.75 f = 1.00 f = 1.25 f = 1.50 f = 1.75 f = 2.00 f = 2.25 f = 2.50 Confining Pressure = 1.0 MPa Voids Ratio = 0.162 0 1 2 3 4 5 6 7 8 9 10 0102030405060 Shear Strain, %Volumetric Strain, % f = 0.25 f = 0.50 f = 0.75 f = 1.00 f = 1.25 f = 1.50 f = 1.75 f = 2.00 f = 2.25 f = 2.50 Confining Pressure = 1.0 MPa Voids Ratio = 0.162
PAGE 190
169 Figure 544. Maximum volumetric strain versus interparticle friction coefficient (circular particles) 2 3 4 5 6 7 8 0.00.51.01.52.02.53.0Interparticle friction coefficientMaximum Volumetric Strain, % Confining Pressure = 1.0 MPa Voids Ratio = 0.162Figure 545. Maximum dilation and residua l shearing resistance angles versus interparticle friction coefficient (circular particles) 0 10 20 30 40 50 60 0.00.51.01.52.02.53.0 Interparticle friction coefficientAngle, degree Maximum dilation angle Residual shearing resistance angle Confining Pressure = 1.0 MPa Voids Ratio = 0.162
PAGE 191
170 for interparticle friction coefficients higher than 1.50. The volumetric strainshear strain plot (fig. 547) shows that the volumetric st rain increases with increasing interparticle friction but again, the rate of increase decreases for higher values. For friction coefficient of 0.25, contraction was observed. Th e maximum rate of dilatancy increased dramatically for lower interparticle frictions and was almost identical for values higher than 1.50. The effect of interparticle fr iction on maximum volum etric strain, maximum dilation angle, and residual shearing resist ance angle are summarized in figs. (548) and (549), which can be expressed as follows: 1442 3 ) f ln( 4807 2 ) (max v [R2 = 0.9675] (521) 239 0 ) f ln( 2242 0max [R2 = 0.9883] (522) 2273 0 max) f ( 083 30 [R2 = 0.9384] (523) 552 30 ) f ln( 3864 6residual [R2 = 0.9539] (524) 5.8.3 Confining pressure (250 KPa), voids ratio (0.35), and 2215 particles The shear stressshear strain a nd volumetric strainshear strain relationships are shown in fig. (550) and (551). The results show th at the peak was not defined for all friction coefficients. The critical state shear stre ngth slightly increased with increasing interparticle friction and was almost identical for interparticle friction coefficients higher than 1.25. The volumetric stra inshear strain plot (fig. 551) shows that the volumetric strain increases with the increase of the inte rparticle friction but at a decreasing rate for higher values. At f = 0.25, contraction is observed, while a slight decrease and then
PAGE 192
171 Figure 546. Shear stressshear strain curves for different interparticle friction coefficients (circular particles) 0 50 100 150 200 250 300 350 400 450 0102030405060 Shear Strain, %Shear Stress, KPa f = 0.25 f = 0.50 f = 0.75 f = 1.00 f = 1.25 f = 1.50 f = 1.75 f = 2.00 f = 2.25 f = 2.50 Confining Pressure = 500 KPa Voids Ratio = 0.25Figure 547. Volumetric strainshear strain curves for different interparticle friction coefficients (circular particles) 1 0 1 2 3 4 5 6 7 0102030405060 Shear Strain, %Volumetric Strain,% f = 0.25 f = 0.50 f = 0.75 f = 1.00 f = 1.25 f = 1.50 f = 1.75 f = 2.00 f = 2.25 f = 2.50 Confining Pressure = 500 MPa Voids Ratio = 0.25
PAGE 193
172 Figure 548. Maximum volumetric strain versus interparticle friction coefficient (circular particles) 1 0 1 2 3 4 5 6 0.00.51.01.52.02.53.0 Interparticle friction coefficientMaximum Volumetric Strain, % Confining Pressure = 500 KPa Voids Ratio = 0.25Figure 549. Maximum dilation and residua l shearing resistance angles versus interparticle friction coefficient (circular particles) 10 5 0 5 10 15 20 25 30 35 40 0.00.51.01.52.02.53.0 Interparticle friction coefficientAngle, degree Maximum dilation angle Residual shearing resistance angle Confining Pressure = 500 KPa Voids Ratio = 0.25
PAGE 194
173 Figure 550. Shear stressshear strain curves for different interparticle friction coefficients (circular particles) 0 50 100 150 200 250 300 010203040506070 Shear Strain, %Shear Stress, KPa f = 0.25 f = 0.50 f = 0.75 f = 1.00 f = 1.25 f = 1.50 f = 1.75 f = 2.00 f = 2.25 f = 2.50 Confining Pressure = 250 KPa Voids Ratio = 0.35Figure 551 Volumetric strainshear strain curves for different interparticle friction coefficients (circular particles) 2 0 2 4 6 8 10 010203040506070 Shear Strain, %Volumetric Strain,% f = 0.25 f = 0.50 f = 0.75 f = 1.00 f = 1.25 f = 1.50 f = 1.75 f = 2.00 f = 2.25 f = 2.50 Confining Pressure = 250 MPa Voids Ratio = 0.35
PAGE 195
174 Figure 552. Maximum volumetric strain versus interparticle friction coefficient (circular particles) 1 0 1 2 3 4 5 6 7 8 0.00.51.01.52.02.53.0 Interparticle friction coefficientMaximum Volumetric Strain, % Confining Pressure = 250 KPa Voids Ratio = 0.35Figure 553. Maximum dilation and residua l shearing resistance angles versus interparticle friction coefficient (circular particles) 10 0 10 20 30 40 50 60 0.00.51.01.52.02.53.0Interparticle friction coefficientAngle, degree Maximum dilation angle Residual shearing resistance angle Confining Pressure = 250 KPa Voids Ratio = 0.35
PAGE 196
175 increase in volume was observed at f = 0.5. The maximum rate of dilatancy increased noticeably for lower values of interparticle fr ictions. The effect of interparticle friction on maximum volumetric strain, maximum dilati on angle, and residual shearing resistance angle are summarized in figs. (552) and (5 53), which can be expressed as follows: 2592 3 ) f ln( 1077 3 ) (max v [R2 = 0.9540] (525) 1908 0 ) f ln( 2695 0max [R2 = 0.9340] (526) 757 38 ) f ln( 504 7residual [R2 = 0.9473] (527) 5.9 The equivalent interparti cle friction coefficient The critical state lines corresponding to circular particles using different interparticle friction coefficients are shown in fig. (554). The criti cal state lines shifted up and became more steep with the increase of the interparticle friction coefficient. Contraction was observed for (f = 0.25) under confining pressures of 250KPa and 500KPa because the initial conditions define a co ntractive soil. The results show that the numerical simulations are ve ry sensitive to the selected interparticle friction coefficient/angle. As stated in many research studies, the main obstacle to the DEM is that simulations are time consuming especially wh en actual particle shapes are used along with large number of particles. Unless the DEM is able to simulate millions of particles with their actual shapes, no actual scale proble ms can be accurately simulated. Another alternative is to change some properties of circ ular particles to control their behavior such that they reproduce the behavior of angular particles. The ci rcular particles, although not
PAGE 197
176 a true representation of the particles, are the most efficient shape to model. Computation time is very small compared to the most eff ective alternative schemes that utilize other shapes. Interparticle friction coefficient is one of the most significant properties in DEM simulations. As shown in fig. (554), the in terparticle friction angl e affects the position of the critical state line. It was observed also that the cr itical state lines for different angularity and shape groups were sh ifted up with respect to that of circular pa rticles. To this extent, it may be possible to adjust th e interparticle friction angle to account for the particle shape characteristics. The target in terparticle friction coeffi cient can be extracted by replotting the critical st ate line for the required angularity or shape group on top of the critical state lines shown in fig (563) and choose the interparticle friction that matches the critical stat e line for the group. The procedure was applied for the three angularity groups (AF1, AF2, and AF3) and the results are shown in fig. (555). The results show that the critical state line for each angularity group varies as a function of friction values, so it would not be accurate to choose an equivalent value of friction to account for the angularity. One option is to continuously alter the friction coefficient as a function of ep/ state in order to match the response of the angular particles. The procedure was repeated for the th ree shape groups (SF1(AF 1525), SF2(AF 1525), and SF3(AF 1525)) and the result s are shown in fig. (556).
PAGE 198
177 Figure 554. Critical state lines for different interparticle friction coefficient 1.15 1.20 1.25 1.30 1.35 1.40 1.45 1.50 020040060080010001200 Confining pressure, KPaSpecific volum e Initial conditions f = 0.25 f = 0.50 f = 0.75 f = 1.00 f = 1.25 f = 1.50 f = 1.75 f = 2.00 f = 2.25 f = 2.50
PAGE 199
178 Figure 555. Critical state lines for circular particles (different interparticle friction coefficient) along with those for different angularity factor groups (f = 0.5) 1.15 1.20 1.25 1.30 1.35 1.40 1.45 1.50 020040060080010001200 Confining pressure, KPaSpecific volum e Initial conditions f = 0.25 f = 0.50 f = 0.75 f = 1.00 f = 1.25 f = 1.50 f = 1.75 f = 2.00 f = 2.25 f = 2.50 AF1 AF2 AF3
PAGE 200
179 Figure 556. Critical state lines for circular particles (different interparticle friction coefficient) along with those fo r different modified shape fact or groups (f = 0.5) 1.15 1.20 1.25 1.30 1.35 1.40 1.45 1.50 020040060080010001200 Confining pressure, KPaSpecific volum e Initial conditions f = 0.25 f = 0.50 f = 0.75 f = 1.00 f = 1.25 f = 1.50 f = 1.75 f = 2.00 f = 2.25 f = 2.50 SF1(Af 1525) SF2(Af 1525) SF3(Af 1525)
PAGE 201
180 5.10 Discussion The effect of particle shape and angul arity characteristic s on strength and dilatancy properties was studie d by performing discrete elem ent numerical simulations of the triaxial test under pure shear conditions. A modificat ion of the particle shape characterization procedure proposed by Sukumaran and Ashmawy (2001) was introduced. A total of 99 particles were di vided into different particle shape and angularity groups. The modeled pure shear test was performed on circ ular particles as a reference and on different shape and angul arity groups under three different initial conditions. The shear stressshear strain and volumetric strainshear strain curves were developed for circular and different par ticle group. The dependency of maximum volumetric strain, maximum dilatio n angle, and the resulted re sidual shear st rength angle on shape and angularity factors was explored The critical state lines for shape and angularity groups along with the co rresponding lines for circular particle were drawn. It is concluded that the critical state lines, for different shape and angularity groups, were shifted up compared with that of the circular particles, which indicates an increase of volumetric strain and hence dilation. The effect of interparticle friction coefficient on the behavior was studied. Again, the critical state lines showed significant dependency on interparticle frictions. An attempt was made to use an equivalent interparticle friction to model different particle shapes. It was conc luded that there is no onetoone equivalency between interparticle friction and shape or angularity. Instead, the interparticle friction must be continuously altered as a function of confining pressure and void ratio to achieve the required effect.
PAGE 202
181 It is noted that the numerical simulation resulted in relatively high dilation rates and volumetric strains compared to typical experi mental values. This is attributed to two causes. First, because the simulations are twodimensional, the tendency of particles to interlock is higher than in 3D due to the restriction enforced on the third degree of freedom. Second, because the sample is re latively small in size, the deformation is uniform across the sample and no shear banding or strain localization occurs. In addition, the pure shear loading conditions further suppr esses strain localiza tion. If a larger assembly is simulated, it is likely that strain localization wi ll occur and smaller volumetric strains will be observed.
PAGE 203
182 CHAPTER SIX: SUMMARY AND RECOMMENDATIONS A thorough study of modeling angular partic les using DEM has been performed. The Particle Flow Code (PFC2D) was used to perform the numerical simulations. The study can be divided into three main parts. Th e first part is the e xperimental verification of the ability of DEM to model angular part icles using the ORC method. The second part is the modification of the ORC method to ensure an identical center of gravity and mass moment of inertia for both actual and created particles and to define a userindependent particle creation procedure. Th e third part of the study is st udying the effect of particle shape and angularity on the di latancy of angular soils. 6.1 Summary An experimental setup was built and m odel particles were manufactured and modeled numerically. The experimental test re sults of five identical tests were compared qualitatively and quantitativel y with the corresponding numer ical simulation. The ORC technique was found to be e ffective in modeling the beha vior of angular particle assemblies. While particle displacements and rotations of bot h the numerical and experimental systems were almost identical at small displacements, larger differences were observed at high displacements, even amo ng the five experiment al tests that have been performed ensuring identical initial and boundary conditions. This variability is
PAGE 204
183 attributed to minor changes in the ini tial conditions, as well as other inherent uncertainties in the system, and cannot be a voided. The ORC simulated system behavior is within the range of the five experiment al tests. The experi mental validation was extended to a different type of sands. The experimental setup was altered to induce a different mode of disturbance to the m odel grains. A parame tric study has been conducted to study the effect of interparticl e friction, contact st iffness, and global damping on the rotations and displacements of the model particles. The results indicated that the global damping is a critical para meter in DEM modeling and must be carefully selected. The effect of interparticle fric tion and normal contact stiffness on the response was also evaluated. A modification for the ORC technique was presented and verified for both rounded and angular particles. The modificat ion ensures an identi cal mass, center of gravity, and mass moment of inertia for both the actual a nd the created particles. A conditional element addition sequence was propos ed in order to faci litate the future automation of the procedure. The compa tibility equations were imposed to the compatibility elements/discs to determine their modified density. The method was verified manually using Fraser river sand par ticle # 4 as a highly angular particle and Michigan sand particle # 1 as a rounded pa rticle. The modified ORC method was found to be applicable for both sand particles and sh ould be applicable fo r almost all particle shapes, which may need more iteration to reac h the final densities. The modified ORC can be considered a step forward in impr oving the modeling of angular particle using DEM especially in dynamic applications.
PAGE 205
184 Discrete element numerical simulations of pure shear conditions were performed for different particle shape and angularity groups. The modified shape factor and the angularity factor were used to separate pa rticles into shape and angularity groups. The stressstrainvolume curves were presented an d analyzed. The critical state lines for shape and angularity groups along with the corresponding lines fo r circular particles were determined. The critical state lines, for different shape and angularity groups, were shifted up compared with that of the circular particles, which indicates an increase of volumetric strain and hence dilation. The effect of interparticle friction coefficient on the behavior was studied. Again, the critical state lines showed significant dependency on interparticle frictions. An attempt was made to use an equivalent interparticle friction to model different particle shapes. It was conc luded that there is no onetoone equivalency between interparticle friction and shape or angularity. Instead, the interparticle friction must be continuously altered as a function of confining pressure and void ratio to achieve the required effect. 6.2 Recommendations for future studies Further research is required to evalua te and compare the interparticle contact forces obtained from the expe rimental tests and the numeric al simulation. Different experimental setups that involve high interp article contact forces should be constructed and verified numerically. The photostress t echnology may present a viable approach to perform the contact force comparison. The automation of the modified ORC method is essential to facilitate and expand the number of created particle. The modi fied ORC method should be used to run
PAGE 206
185 dynamic problems and the results should be co mpared to the original ORC method to address the increase of simulations accuracy The extension of the ORC method to the three dimension actual particle simulations can be considered a cutting edge future research topic. The work should be extended to larger num ber of particle shape and angularity groups under various initial conf ining pressures and porosities. There is a need for a fundamental theory that explains par ticle movement during shear and how the interlocking resistance develops The equivalent interparticle friction coefficient is an interesting research topic for the fu ture of modeling angular particles.
PAGE 207
186 REFERENCES Anandarajah, A., 2000, Â“Numerical simulation of onedimensional behavior of kaoliniteÂ”, Geotechnique, Vol. 50, No. 5, pp. 505519. Anstey, R. L. and Delmet, D.A., 1973, Â“Fourie r analysis of zooeci al shapes and fossil bryozoansÂ”, Bulletin of the Geological Society of America, No 84, pp. 17631764. Ashmawy, A.K., Sukumaran, B., and Hoang, A. V., 2003, Â“Evaluating the influence of particle shape on liquefaction behavi or using Discrete Element MethodÂ”, Proceedings of the thirteenth intern ational offshore and polar engineering conference (ISOPE 2003) H onolulu, Hawii, May 2003. Barbosa, R., and Ghaboussi, J., 1992, Â“Discr ete finite element methodÂ”, Engineering Computations, Vol. 9, pp 253266. Barrett, P.J., 1980, Â“The shape of rock particles, a critical reviewÂ”, Sedimentology, Vol. 27, pp. 291303. Bathurst, R.J. and Ruthenburg, L. 1989, Â“Inve stigation of micromechanical features of idealized granular assemblie s using DEMÂ”, Proc. of 1st US Conference on Discrete Element Methods, Golden, CO, Engineering Computations (1992),Vol. 9, pp 199210. Been, K. and Jefferies, M.G. 1985. Â“A state parameter for sandsÂ”. Geotechnique, Vol. 35, No. 2, pp. 99112. Bolton, M.D. 1986. Â“The strength and dilatanc y of sandsÂ”. Geotechnique, Vol. 36, No. 2, pp. 6578. Campbell, C. S., and Brennen, C. E., 1985, Â“Computer simulation of granular shear flowsÂ”, Journal of Fluid Mech anics, Vol. 151, pp. 167188.
PAGE 208
187 Caquot, A. 1934. Â“Equilibre des massifs a frottement interneÂ”. Stabilit des Terres Pulverentes et Cohrentes, Paris: Gauthier Villars. Cundall, P. A., 1978, Â“Ball: A computer prog ram to model granular media using the distinct element methodÂ”, Technical note TNLN13, Advance Technology group. London: Dames and Moore. Cundall, P. A., 1971, Â“ A computer model for simulating progressive large scale movements in block rock systemsÂ”, Pr oceedings of the Symposium of the International Society of Rock Mechanics, Nancy, II, Article 8. Cundall, P. A., 1974, Â“A computer model fo r rockmass behavior using interactive graphics for the input and output of ge omechanical dataÂ”, Report for Contract number DACW4574C006, for the Missouri River Divisi on, US Army Corps of Engineers, University of Minnesota. Cundall, P. A., 1976, Â“Explicit finitediffere nce methods in geomechanicsÂ”, Proceedings of the 6th Conference of Numerical Methods in Geomechanics, pp. 132150. Cundall, P. A., and Strack O. D. L., 1979, Â“A discrete numerical model for granular assembliesÂ”, Geotechnique, Vol. 29, pp. 4765. Cundall, P. A., and Strack O. D. L., 1979, Â“The development of constitutive laws for soil using distinct element method Â”, Proc. 3th Conference of Nume rical Methods in Geomechanics, Aachen, pp. 289298. Cundall, P. A., and Strack O. D. L., 1979, Â“T he distinct element method as a tool for research in granular media Â”, Part II, Report to NSF, Dept of Civil and Mineral engineering, Univ. of Minnesota. Cundall, P.A., and Strack O.D.L., 1983, Â“M odeling of microscopic mechanisms in granular materialÂ”, In Mechanics of Granular Materials: New Method and Constitutive Relations, (Eds. J.T. Jenkins and M. Satake), Elsevier, Amsterdam, pp. 137149. Cundall, P.A., 1988, Â“Computer simulation of de nse sphere assembliesÂ” In Mechanics of Granular Materials: New Method and Cons titutive Relations, (Eds. J.T. Jenkins and M. Satake), Elsevier, Amsterdam, pp. 113123. Cundall, P.A., and Strack, O.D.L. 1989, Â“N umerical experiment on localization in frictional materialsÂ”, Ing. Archiv, Vol. 59, pp. 184159.
PAGE 209
188 Ehrlich, R., Brown, B.J., Yarus, J.M., and Prygocki, R.S., 1980, Â“The origin of shape frequency distributions and the relationship between size and shapeÂ”, Journal of Sedimentary Petrology, Vol. 50, pp. 475484. Favier, J., AbbaspourFard, M, Kremmer, M ., and Raji, A., 1999, Â“Shape representation of axisymmetrical, nonsphe rical particles in discrete element simulation using multielement model particlesÂ”, Engineer ing Computations, Vol. 16, pp. 467480. Favier, J. F., AbbaspourFar d, M. H., and Kremmer, M ., 1999, Â“Modeling nonspherical particles using multisphere discrete elementÂ”, Journal of Engineering Mechanics, Vol. 127, No. 10, pp. 971976. Houlsby, G.T. 1991. Â“How the dilatancy of soils affects their be haviorÂ”, European Conference on Soil Mechanics and Founda tion Engineering, Vol. 4, pp. 11891202. Hryciw, R., Raschke, S., Ghalib, A., Horner, D., and Peters, J., 1997, Â“Video tracking for experimental validation of discrete el ement simulations of large discontinuous deformationsÂ”, Computers and Geotec hnics, vol. 21, No. 3, pp. 235253. Huang, A.B. and Loe, J.S., 1989, Â“A DEM pa rticulate calibration chamberÂ”, Proceedings of the 1st US Conference on Discrete Element Methods. Golden, CO. Ishibashi, I., Agarwal, T., and Ashraf, S.A ., 1989, Â“Anisotropic behavior of glass spheres by a discrete element method and laborator y experimentsÂ”, Proceedings of the 1st US Conference on Discrete Element Methods. Golden, CO. Iwashita, K., Taruni, V., Casaverde, L., Ve rmura, D, Meguro,K. and Hakuno, M., 1988, Â“Granular assembly simulation for ground collapseÂ”, In Mechanics of Granular Materials, (Eds. J.T. Jenkins and M. Sata ke), Elsevier, Amsterdam, pp. 125132. Ishibashi, I.and Hakuno, M., 1989, Â“A DEM simulation for cliff collapseÂ”, Proceedings of the 1st US Conference on Discrete Element Methods. Golden, CO. Jensen, R., Bosscher, P., Plesha, M., and Ed il, T., 1999, Â“DEM simulation of granular mediastructure interface: effect of surface roughness and particle shapeÂ”, International Journal for Numerical and Analytical Methods in Geomechanics, Vol. 23, pp. 531547.
PAGE 210
189 Kato, S., Yamamoto, S., and Nomani, S., 1999, Â“DEM simulation of biaxial compression test with results of triaxial compression test for unsaturated soilÂ”, Journal of Applied Mechanics, JSCE, Vol. 2, pp. 419426. Kishino, Y., 1988, Â“Disc model analysis of gr anular mediaÂ”, In Mechanics of Granular Materials, (Eds. J.T. Jenkins and M. Sata ke), Elsevier, Amsterdam, pp 143152. Kiyama, H. and Fujimura, H., 1983, Â“Applicatio n of CundallÂ’s discrete block method to gravity flow analysis of rock like granular materialÂ”, Proceedings of the Japanese Society of Civil Engineering (JSCE), Vol. 333, pp. 137146. Krumbein, W.C., 1941, Â“Measurement and geological significance of shape and roundness of sedimentary particlesÂ”, Jour nal of Sedimentary Petrology, No. 11, pp. 6472. Kuhn, M. and Mitchell, J., 1989, Â“The modeli ng of soil creep with the discrete element methodÂ”, Proceedings of the 1st US Conference on Discrete Element Methods. Golden, CO. Kuhn, M. and Mitchell, J., 1993, Â“New pe rspectives on soil creepÂ”, Journal of Geotechnical Engineering, vol. 119, No. 3, pp. 507522. Li, X.S., and Dafalias, Y.F. 2000 Â“Dilata ncy for Cohesionless soilsÂ”. Geotechnique, Vol. 50, No. 4, pp. 449460. Lin, X., and NG, T.T., 1997, Â“A threedimens ional discrete element model using arrays of ellipsoidsÂ”, Geotechnique, Vol. 47, No. 2, pp. 319329. Liu, S. H., and Sun, D. A., 2002, Â“Simulati ng the collapse of unsaturated soil be DEMÂ”, International Journal for Numerical and Analytical Methods in Geomechanics, Vol. 26, pp. 633646. Mandelport, B.B. 1967. Â“How long is the cost of Great Britain? Statis tical self similarity and fraction dimensionÂ”. Science, Vol. 166, pp. 636638. Manzari, M.T., and Nour, M.A. 2000. Â“Significa nce of soil dilatancy in slope stability analysisÂ”. Journal of Geotechnical and Geoenvironmental Engineering, Vol. 126, No. 1, PP 7580.
PAGE 211
190 Matuttis, H.G., Luding, S., and Herrmann, H. J., 2000, Â“Discrete element simulations of dense packings and heaps made of spheri cal and nonspherical particlesÂ”, Powder Technology, Vol. 109, No 1, pp. 278292. McDowell, G. R. and Harireche, O., 2002, Â“D iscrete element modeling of yielding and normal compression of sandÂ”, Geotechnique Vol. 52, No 4, pp. 299304. Mindlin, R.D., and Deresiewicz, H., 1953, Â“E lastic spheres in contact under varying oblique forcesÂ”, National Conference of the Applied M echanics Division, Minnesota., June 1820, ASME, pp. 327344. Mirghasemi, A. A., Rothenburg, L., and Matyas E. L., 2002, Â“Influence of particle shape on engineering properties of twodi mensional polygonshaped particlesÂ”, Geotechnique, Vol. 52, No 3, pp. 209217. Mirghasemi, A. A., Rothenburg, L., and Maty as, E. L., 1997, Â“Numerical simulations of assemblies of twodimensional polygonshap ed particles and e ffects of confining pressure on shear strengthÂ”, Soils & F oundations, Vol. 37, No 3, pp. 4352. Mustoe, G.G.W., 1992, Â“A generalization of the discrete element methodÂ”, Engineering Computations, Vol. 9, pp. 181190. Newland, P.L. and Allely, B.H. 1957. Â“Vol ume change in draine d triaxial tests on granular materialsÂ”. Geotechnique, Vol. 7, No. 1, pp. 1728. Nakai, T., Hoshikawa, T., and Hinokio, M. 1997. Â“Drained and undrained behavior of sand under cyclic loadingÂ”. Proceedings of Seventh Interna tional Offshore and Polar Engineering Conference, Vol. 1, pp. 901906, Honolulu, USA. Ng, T.T. and Dobry, R., 1988, Â“Computer modeling of granular soils under cyclic loading: a particulate mechanics appro achÂ”, Report CE8803, Department of Civil Engineering, Rensselaer Po lytechnic Institute, Troy, NY. Ng, T.T., 1989, Â“Numerical simulation of gr anular soil under monotonic and cyclic loading: a particulate mechanics appro achÂ”, PhD Dissertation, Department of Civil Engineering, Rensselaer Po lytechnic Institute, Troy, NY. Ng, T.T., 1994, Â“Numerical simulation of gra nular soil using elliptical particlesÂ”, Computers and Geotechnics, Vo l. 16, No. 2, pp. 153169.
PAGE 212
191 Ni, Q., Powrie, W., Zhang, X., and Harkness, R., 2000, Â“Effect of particle properties on soil behavior:3D numerical modeling of shear box testsÂ”, ASCE Geotechnical Special Publications, Numerical Methods In Geotechnical Engineering. 96, pp. 5870. Oner, M., 1984, Â“Analysis of fabric change during cyclic loading of granular soilsÂ”, Proceedings of the 8th World Conference on Earthquake Engineering, pp. 5562. Orford, J.D. and Whalley, B.W. 1983, Â“The use of fractal dimension to quantify the morphology of irregular shaped particle sÂ”. Sidementology, Vol. 30, pp. 655668. Petrakis, E., Dobry, R., and Ng, T., 1989, Â“S mall strain response of random arrays of elastic spheres using a nonlinear dis tinct element procedure, in wave propagationÂ”, Granular Media, (Eds Ka ramanlidis and R.B. stout), AMD, Vol 101, ASME, New York, pp 1727. Petrakis, E., and Dobry, R., 1989, Â“Micromech anical behavior and modeling of granular soilÂ”, Report to Air force Office of Scientific research, Department of Civil Engineering, Rensselaer Polytechnic Institute, Troy, NY. Piper, D.J.W., 1971, Â“The use of DMac penc il follower in routine determination of sedimentary parametersÂ”, Data proce ssing in biology and geology (ed. J. T. Cubill). pp. 97103. Potapov, A., and Campbell, C., 1998, Â“A fast model for the simulation of nonround particlesÂ”, Granular Matter Vol. 1, pp 914. Powers, M.C., 1953, Â“A new roundness scale fo r sedimentary particlesÂ”, Journal of Sedimentary Petrology, No. 2, pp 117119. Reynolds, O. 1885. Â“On the dilati on of media composed of rigi d particles in contact, with experimental illustrationsÂ” Philosophical Magazine, Series 5, Vol. 20, pp. 469481. Rothenburg, L., and Bathurst, R. J., 1989, Â“A nalytical study of induced anisotropy in idealized granular materials Â”, Geotechnique, Vol. 39, No 4, pp. 601614. Rowe, P.W. 1962. Â“The stress dilatancy relation for static equilibrium of an assembly of particles in contactÂ”. Proceed ings of the Royal Society, Series A, Vol. 269, pp. 500527.
PAGE 213
192 Rothenburg, L., and Bathurst, R. J., 1992, Â“Micromechanical features of granular assemblies with planner elliptical particle sÂ”, Geotechnique, Vol. 42, No 1, pp. 7995. Schanz, T. and Vermeer, P.A. 1996. Â“Angl e of friction and dilatancy of sandÂ”. Geotechnique, Vol. 46, No. 1, pp. 145151. Schofield A. N. and Worth, C.P. 1968. Â“Cr itical state soil mech anicsÂ”. McGraw Hill, London. Sebstyan, G.S. 1969. Â“On pattern recognition with application to silhouettesÂ”. PhD Thesis, Massachusetts Institute of Technology, MA. Seed, H. B., Wong, R. T., Idriss, I. M., and Tokimatsu, K., 1984, Â“Moduli and damping factors for dynamic analysis of cohe sionless soilsÂ”, Report UCB/EERC84/I4, University of California, Berkeley. Sitharam, T., Dinesh, S., and Shimizu, N., 2002, Â“Micromechanical modeling of monotonic drained and undrained shear beha vior of granular media using threedimensional DEMÂ”, International Journal for Numerical and Analytical Methods in Geomechanics, Vol. 26, pp. 11671189. Skempton, A.W. and Bishop A.W. 1950. Â“Meas urement of shear strength of soilsÂ”. Discussion by A.W. Bishop. Geotec hnique, Vol. 2, No. 2, pp. 113. Skempton, A.W. 1954. Â“The pore pressure co efficient A and BÂ”. Geotechnique, Vol. 2, No. 2, pp. 113. Strack O. D. L. and Cundall, P. A, 1978, Â“T he distinct element method as a tool for research in granular mediaÂ”, Part I, Report to NSF, Department of Civil and Mineral engineering, University of Minnesota. Minneapolis, MN. Strack O. D. L. and Cundall, P. A, 1984, Â“F undamental studies of fabric in granular materialsÂ”, Interim report to NSF, CEE8310729, Dept of Civil and Mineral engineering, University of Minnesota. Minneapolis, MN. Stroud, M.A. 1971. Â“The behavior of sand at low stress levels in the simple shear apparatusÂ”. Ph.D. Thesis, University of Cambridge, Cambridge, United Kingdom. Sukumaran, B., and Ashmawy, A. K., 2001, Â“Quantitative characterization of the geometry of discrete part iclesÂ”, Geotechnique, Vol. 51, No. 7, pp. 619627.
PAGE 214
193 Sukumaran, B. 1996. Â“Study of the effect of pa rticle characteristics on the flow behavior and strength properties of pa rticulate materialsÂ”. Ph.D. Thesis, Purdue University, West Lafayette, IN. Sukumaran, B. and Ashmawy, A.K. 2001. Â“Qua ntitative characterizat ion of the geometry of discrete particlesÂ”. Geotechni que, Vol. 51, No. 7, pp. 619627. Taylor, D.W. 1948. Â“Fundamentals of soil mechanicsÂ”. Wiley, New York. Tarumi, Y. and Hakuno, M.A., 1989, Â“A DEM for sand liquefactionÂ”, Proceedings of the 1st US Conference on Discrete Element Methods. Golden, CO. Terzaghi, K. 1920. Â“Old earthp ressure theories and new test resultsÂ”. Engineering News Record, Vol. 85, No. 14. (Reprinted in Fo rm Theory to Practice in Soil Mechanics 1960, New York: J. Wiley and Sons). Thornton, C. and Randall, C.W., 1988, Â“Applicat ions of theoretical contact mechanics to solid particle system simulationsÂ”, In Mechanics of Granular Materials: New Method and Constitutive Relations, (Eds. J.T. Jenkins and M. Satake), Elsevier, Amsterdam, pp. 133142. Ting, J., Khawaja, M., Meachum, L., and Ro well, J., 1993, Â“An ellipsebased discrete element model for granular materialsÂ”, International Journal for Numerical and Analytical Methods in Geomechanics, Vol. 17, pp. 603623. Ting, J.M., Corkum, B.T., Kauffman, C.R., and Greco, C.A., 1986, Â“Application of the distinct element method in geotechnical engineeringÂ”, Proceedings of the 2nd International Symposium on Numerical Models in Geomechanics, Ghent. Ting, J.M., Corkum, B.T., Kauffman, C.R., and Greco, C.A., 1987, Â“Discrete numerical modeling of soil: validation and applicati onÂ”, Publication 8703, Department of Civil Engineering, University of Toronto, Toronto, Canada. Ting, J.M., and Corkum, B.T., 1988, Â“Strengt h behavior of granul ar materials using discrete numerical modelingÂ”, Proceedings of the 6th International Symposium on Numerical Models in Geomechanics, Innsbruck, pp. 305310. Ting, J.M., and Corkum, B.T., 1988, Â“Soilstructure interaction by discrete element modelingÂ”, Proceedings of the Canadian Society of Civil Engineering (CSCE), Calgary, Canada, pp. 196215.
PAGE 215
194 Ting, J.M., and Corkum, B.T., 1988, Â“Applica tion of DEM in geotechnical engineeringÂ”, Proceeding of the 3rd International Conference on Computational Civil Engineering, Vancouver, Canada, pp. 578594. Ting, J.M., Corkum, B.T., Kauffman, C.R., and Greco, C.A., 1987, Â“Discrete numerical model for soil mechanicsÂ”, Journal of Ge otechnical Engineering, ASCE, Vol. 3, pp 379398. Ting, J., Meachum, L., and Rowell, J., 1995, Â“Effect of particle shape on the strength and deformation mechanisms of ellipseshaped granular assemblagesÂ”, Engineering Computations, Vol. 12, pp. 99108. Uckida, Y. and Hakuno, M.A., 1989, Â“A DEM si mulation for debris flowÂ”, Proceedings of the 1st US Conference on Discrete Element Methods. Golden, CO. Uemuro, D. and Hakuno, M., 1987, Â“Granular assembly simulation with CundallÂ’s model for dynamic collapse of structural f oundationÂ”, Proceedings of the Japanese Society of Civil Engineering (JSCE), Vol. 4, pp. 181190. Uemuro, D. and Hakuno, M., 1989, Â“A DEM simulation for pile penetrating into groundÂ”, Proceedings of the 1st US Conference on Discrete Element Methods. Golden, CO. Voyiadjis, G., Thiagarajan, G., and Petrak is, E., 1992, Â“Anisotropi c distortional yield model for granular mediaÂ”, Microstruc tural Characterization In Constitutive Modeling of Metals and Granular Media, ASME, Vol. 32, pp. 119134. Wadell, H., 1932, Â“Volume shape and roundness of rock particlesÂ”, Journal of Geology, Vol. 40, pp. 443451. Walton, O.R., 1982, Â“Explicit particle dyn amics model for granular materialsÂ”, Proceedings of the 4th Int. Conference on Numerica l methods in Geomechanics, Edmonton, pp.12621269. Walton, O.R., 1983, Â“Particledynamic calcula tions of shear flowÂ”, In Mechanics of Granular Materials: New Method and Cons titutive Relations, (Eds. J.T. Jenkins and M. Satake), Elsevier, Amsterdam, pp 327338. Walton, O.R. and Braun, R.L., 1988, Â“Vis cosity and temperature calculations for assemblies of inelastic frictional discsÂ” Journal of Rheology, Vol. 30, pp. 949980.
PAGE 216
195 Wan, R.G., and Guo, P.J. 1999. Â“A pressure and density dependent dilatancy model for granular materialsÂ”. Soils and Founda tions, Vol. 39, No. 6, pp. 110. Worth, C.P. 1958. Â“The behavior of soils a nd other granular media when subjected to shearÂ”. Ph.D. Thesis, University of Cambridge, United Kingdom. Worth, C.P. and Bassett, R.H. 1965. Â“A st ressstrain relationship for the shearing behavior of sandÂ”. Geotechnique, Vol. 15, No. 1, pp. 3256. Wu, w.Y. and Liu, Z.D., 1989, Â“The simulati on of fabrics for granular material with DEMÂ”, Proceedings of the 1st US Conference on Discrete Element Methods. Golden, CO. Yamamoto, T. and Hakuno, M.A., 1989, Â“A Dem simulation for tunnel excavationÂ”, Proceedings of the 1st US Conference on Discrete Element Methods. Golden, CO. Yang, J. and Li, X.S. 2004. Â“Statedependent strength of sands from the perspective of unified modelingÂ”. Journal of Geotechni cal and Geoenvironmen tal Engineering, Vol. 130, No. 2, pp. 186198. Yudhbir and Abedinzadeh, R., 1991, Â“quantification of pa rticle shape and angularity using image analyzerÂ”, Geotechnical Tes ting Journal, Vol. 14, No 3, pp 296308. Zhang, Y., and Cundall, P.A., 1986, Â“A nume rical simulation of slow deformationsÂ”, Proc. Symp. On Mechanics of Particulat e Media, Tenth US National Congress on Applied Mechanics, University of Texas, Austin. TX. Zhang, D., and Whiten, W.J., 1996, Â“The calcu lation of contact forces between particles using spring and damping modelsÂ”, Po wder Technology, Vol. 88, pp. 5964. Zhang, D., and Whiten, W.J., 1999, Â“A new calculation method for particle motion in tangential direction in discrete element simulationsÂ”, Powder Technology, Vol. 102, pp. 235243.
PAGE 217
196 APPENDICES
PAGE 218
197 Appendix A: theory a nd background of DEM A.1 Calculation cycle The calculations performed in the discrete element me thod alternate between the application of NewtonÂ’s second law to the discs and a force displacement law at the contacts. Using NewtonÂ’s law, the motion (d isplacement, velocity, a nd acceleration) of a particle can be determined, wh ereas the forcedisplacement law is used to determine the contact forces from displacements/overlaps, Cundall and Strack (1979). The calculation cycle is shown in fig. (A1). A.2 Forcedisplacement law Two discs are taken to be in contact if the distance between their centers, D, is less than the sum of their radii, Rx, Ry, fig. (A2). According to fig. (A2), the normal and tangential components of the rela tive velocities can be give n, Cundall and Strack (1979), as: i i i i ie y x e X n (A1) y y x x i i i i iR R t y x t X s (A2)
PAGE 219
198 Appendix A (continued) Figure A.1. Calcula tion cycle for the DEM Figure A.2. Forcedisplacement law in DEM
PAGE 220
199 Appendix A (continued) The time integration of the velocity gives the relative displacement components in the normal and tangential direction as: t e y x t n ni i i (A3) t R R t y x t s sy y x x i i i (A4) Applying the forcedisplacement law using the contact stiffness properties and the derived normal and tangential displacements, the normal and shear forces increments can be calculated as follows: t e y x K n K Fi i i n n n (A5) t R R t y x K s K Fy y x x i i i s s s (A6) The force increments are then added to the sum of all forces in order to determine the final forces after this particular time step: n 1 N n N nF F F (A7) s 1 N s N sF F F (A8)
PAGE 221
200 Appendix A (continued) A coulombtype friction law is then inco rporated by comparing the magnitude of the shear/tangential force, Fs, with the maximum possible shear strength: C tan F FN max S (A9) Where, is the smaller of the interparticle friction angles of the two discs in contact and C is the smalle r of their cohesions. If Fs is larger that Fs,max, Fs is set equal to Fs,max preserving the sign of Fs. The resultant forces and moments can then be determined. Using NewtonÂ’s second law, th e linear and angular acc elerations can be determined and the cycle continues. A.3 Motion In order to determine the linear and ro tational velocity used in the forcedisplacement law, the current resu ltant force and moment at time tN is assumed to be constant during the period from tN1/2 to tN+1/2. NewtonÂ’s law can then be applied as: xi i xF x m (A10) x x xM I (A11) Knowing the moment of inertia, Ix and the mass, mx, and assuming that the linear and rotational accelerations are constants duri ng the time step, the linear and rotational velocity can be determined as follows:
PAGE 222
201 Appendix A (continued) t m F x xN x xi 2 1 N i 2 1 N i (A12) t I MN x x 2 1 N x 2 1 N x (A13) The new values for the velocities are used to update the positions and rotation of the discs as: t x x x2 1 N i N i 1 N i (A14) t2 1 N x N x 1 N x (A15) A.4 Damping Two forms of damping are introduced, fr iction damping and viscous damping. Friction damping occurs during sliding when th e absolute value of the shear force at any point is Fsmax. Viscous damping can be expressed in two types of damping, contact and global damping. Contact damping, which can be envisioned as resulting from normal and shear contact dashpots, opera tes on the relative velocities at the contacts. A damping force should be included in the force summa tion. The damping force will affect the moment summation as well:
PAGE 223
202 Appendix A (continued) t m D F x xN x xi xi 2 1 N i 2 1 N i (A16) t I MN x x 2 1 N x 2 1 N x (A17) The contact damping force in normal and tangential directions can be expressed as: i 2 1 N i i n n N ne y x c n c D (A18) 2 1 N y y x x i 2 1 N i i s s N sR R t y x c s c D (A19) n nK c (A20) s sK c (A21) The global damping, which may be e nvisioned as the effect of dashpots connecting each particle to ground, operates on the absolute velocity of the discs and is included in the calculation of motion: i xi xi i xx C D F x m (A22)
PAGE 224
203 Appendix A (continued) x x x xC M I (A23) xm C (A24) xI C (A25) A.5 Energy dissipation Energy is dissipated in the discrete el ement model through friction, contact, and global damping in order that the assemblies reach a state of equi librium for all conditions. A.6 Time step To keep the numerical scheme stable, the time step should be less than the critical time step. The critical time step can be estimated on the basis of single degree of freedom, SDOF, system of mass (m) connected to ground by spring of stiffness (k) as: k m 2 t (A26) A.7 Input parameters The input parameters for the discrete elem ent codes can be divided into geometric and physical properties. The ge ometrical data may include pos itions and orientations of straight rigid boundaries (usually strainc ontrolled boundary) and positions and radii of discs. The physical Properti es such as disc density, c ohesion, interparticle friction
PAGE 225
204 Appendix A (continued) coefficient, shear contact stiffness, and norma l contact stiffness should be included. Also contact damping, global damping, and time step are very important inputs. A.8 Contact constitutive models Simple constitutive models at the contact are used to simulate the overall behavior of an array of particles. Three main categ ories of constitutive modeling are applied, a stiffness model, a slip model, and a bonding model. Among these models, the contactstiffness models will considered here in more details. Two contactstiffness models are available, the linear elastic and HertzMindlin contact model. It is also possible for the user to input a userdefined contact model that is different from the builtin contact stiffness models. The linear elastic contact model is defined by the normal and shear stiffness, kn and ks, for balltoball and balltowall. The contact normal stiffness can be expressed as: B n A n B n A n nk k k k k (A27) The contact shear stiffness can be expressed as: B s A s B s A s sk k k k k (A28)
PAGE 226
205 Appendix A (continued) The HertzMindlin contact model is a nonlinear contact formulation based on the approximation of the theory of Mindlin and De resiewicz (1953). It is applicable only to the case of spheres in contact. It doesnÂ’t reproduce the contin uous nonlinearity in shear. Also it cannot be used between two balls that are joined by contact bond Â“not defined for tensionÂ”. The contact no rmal stiffness can be expressed as: n \ nU 1 3 R 2 G 2 k (A29) The contact shear stiffness can be expressed as: 3 / 1 n i 3 / 1 \ 2 sF 2 R 1 3 G 2 k (A30) Where, Un is the sphere overlap, and Fn is the magnitude of the normal contact force. For balltoball contact the multipliers are given by: B A B A \R R R R 2 R (A31) B AG G 2 1 G (A32) B A2 1 (A33)
PAGE 227
206 Appendix A (continued) For balltowall contact, the multipliers are given by: ball \ R R (A34) ballG G (A35) ball (A36
PAGE 228
ABOUT THE AUTHOR Amr Sallam received a B.Sc. and a M.S. in Civil Engineering from University of Alexandria, Egypt in 1994. He starte d teaching Soil Mechanics and Foundation Engineering in Alexandria University an d the Egyptian Air Force Academy in 1994. Besides teaching, he joined FACB Consulting Group in 1995 as a geotechnical engineer. He continued teaching and cons ulting until joining the Ph.D. pr ogram in the University of South Florida in 2001. While in the Ph.D. program at the University of South Florida, Mr. Sallam was active in the Geotechnical Society and the Amer ican Society of Civil Engineering. He worked as a research assistant on various pr ojects. He was the laboratory instructor for two full semesters. He also worked as a teac hing assistant for the Soil Mechanics class. He coauthored six publications in time doma in reflectometry (TDR) and discrete element method (DEM) and gave several presentations at professional meetings and conferences.
