SIMULATION OF PROGRESSIVE DAMAGE DEVELOPMENT IN BRAIDED COMPOSITE TUBES UNDER AXIAL COMPRESSION by Carla Jane McGregor B.Sc. (Civil Engineering), Lakehead University, 2002 A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF APPLIED SCIENCE in THE FACULTY OF GRADUATE STUDIES r. Civil Engineering THE UNIVERSITY OF BRITISH COLUMBIA April 2005 © Carla Jane McGregor, 2005 Abstract ABSTRACT Composite materials have recently been investigated as an attractive alternative to traditional materials used as energy absorbers in the automotive industry. One of the factors inhibiting the implementation of these lightweight energy absorbers is the lack of accurate, robust, and physically meaningful material models capable of simulating damage growth in composite materials under a variety of applied loads. A plane-stress continuum damage mechanics based model for composite materials, CODAM (COmposite DAMage), recently developed at the University of British Columbia (Williams, 1998; Williams et al., 2003; Floyd, 2004) has proved to be capable of predicting the complete tensile behaviour of composites from initiation of damage to complete failure. The current study focuses on extending the model to capture the complete compressive response of composite materials as well as their behaviour under load reversals. Refinements to the model are based on the experimentally observed compressive failure mechanisms presented in the literature. In particular, the mechanical consequences of kinking and kink band broadening, in conjunction with matrix cracking and delamination, are represented in the model. The model parameters defining the compressive response are related to experimentally observed behaviour, thus maintaining the physical basis of CODAM. These improvements contribute to the development of CODAM as a predictive tool in the analysis of gross damage growth in composite materials. To demonstrate the capability of the model, the predicted response of braided composite panels under tensile loads, compressive loads and load reversals are compared to experimentally observed behaviours. The model is also evaluated on its ability to predict the damage propagation and energy absorption in axial crushing of braided composite tubes. Model predictions for both panel and tube specimens correlate well with the experimental results. The successful tube crushing predictions suggests that CODAM could offer an attractive design aid for the future incorporation of lightweight composite energy absorbers into crashworthiness structures. ii Table of Contents T A B L E OF CONTENTS Abstract ii List of Tables v List of figures vi Nomenclature xv Acknowledgements . xviii Chapter 1 Introduction 1 1.1 Motivation 1 1.2 Goals and Objectives 3 1.3 Outline of the Thesis 3 Chapter 2 Background and Literature Review 5 2.1 Introduction 5 2.2 Composite Tube Crushing - Experimental 6 2.2.1 Test Methodologies 6 2.2.2 Crushing Modes and Mechanisms 7 2.2.3 Energy Absorption Calculations 14 2.2.4 The Effect of the Reinforcing Fibre 14 2.2.5 The Effect of the Matrix Material 18 2.2.6 The Effect of the Fibre Architecture 19 2.2.7 The Effect of the Geometry of the Specimen 23 2.2.8 The Effect of the Fibre Volume Fraction 26 2.2.9 The Effect of the Testing Speed. 26 2.3 Composite Tube Crushing - Numerical 28 2.3.1 Finite Element Codes 30 2.3.2 Material Models 35 2.3.3 Delamination Modeling 42 2.4 Summary 44 Chapter 3 Numerical Implementation of a Compressive Damage Model 45 3.1 Introduction 45 3.2 Compressive Failure in Composite - Background 45 3.2.1 General Compressive Failure Mechanisms 45 3.2.2 Experimental Investigation of Failure Mechanisms 50 3.2.3 Predicting the Compressive Response of Fibre-Reinforced Composites 58 Table of Contents 3.3 Compressive damage modelling using CODAM 63 3.3.1 Overview of the CODAM material model 63 3.3.2 Overview of the Analog Model 66 3.3.3 Incorporation of the Analog Model into CODAM. 75 3.4 Verification of CODAM 82 3.4.1 Analytical Solution 84 3.4.2 Numerical Solution 91 3.5 Summary 97 Chapter 4 Calibration and Validation of Compressive Model 99 4.1 Introduction 99 4.2 Panel Material Properties and Experimental Details 99 4.3 CODAM Parameters 103 4.4 Tensile Simulation Results 112 4.5 Compressive Simulation Results 115 4.6 Strain Reversal Simulation Results 117 4.7 Summary 119 Chapter 5 Tube Crushing Simulations 120 5.1 Introduction 120 5.2 Experimental Details and Energy Balance of Tubes 120 5.3 General Simulation Details 123 5.4 Simulations with Plug Initiator 128 5.4.1 Single Ply '. 128 5.5 Simulations without Plug Initiator 144 5.5.1 Single Ply 144 5.5.2 Double Ply 150 5.5.3 Four Ply 153 5.6 Summary 156 Chapter 6 Conclusions and Future Work 159 6.1 Summary and Conclusions 159 6.2 Future Work 161 References 163 Appendix A Implementation of Compressive Damage Model 171 Appendix B LPT and Ply Discount Calculations For Braided Material 173 i v List of Tables LIST OF TABLES Table 2.1 List of various materials studied in composite tube crushing investigations. 17 Table 2.2 List of researches that have studied circular, square and/or rectangular cross-sections 25 Table 2.3 List of various testing rates investigated by researchers 28 Table 2.4 List of experimental and numerical studies on the axial crushing of composite tubes 30 Table 3.1 Guidelines for remaining lamina stiffness values (at matrix damage saturation) for use in LPT 82 Table 3.2 COD AM parameters for the single element simulations 83 Table 4.1 Average mechanical properties of the Fortafil #566/Hetron 922 composite constituents (DeTeresa et al, 2001) 99 Table 4.2 Experimentally measured initial modulus and peak stresses (DeTeresa et al., 2001) 103 Table 4.3 COD AM parameters in both the local x and _y-directions showing upper and lower bound assumptions 104 Table 5.1 Summary of experimentally measured tube crushing properties 123 Table 5.2 Consistent System of Units 124 Table 5.3 Average damage parameters of side and corner elements 136 Table 5.4 Predicted variation of energy absorbed by the complete failure of corner and side elements for different fracture energies 139 Table 5.5 Average damage parameters of side and corner elements 148 Table 5.6 Average damage parameters of side and corner elements 153 Table 5.7 Average damage parameters of side and corner elements 156 v List of Figures LIST OF FIGURES Figure 2.1 Principal scheme of the energy absorption devices in an automobile (Kerth andMaier, 1994) : 5 Figure 2.2 Schematic showing the general load-displacement features characteristic of progressive axial tube crushing 8 Figure 2.3 Two trigger mechanisms investigated by Thornton (1990) and double chamfered tubes studied by Song et al. (2002) 9 Figure 2.4 Schematic showing the two main failure modes, (a) progressive folding and (b) progressive crushing (Hull, 1991) 10 Figure 2.5 Details of the two extreme modes of progressive crushing, (a) splaying and (b) fragmentation proposed by Hull (1991) 10 Figure 2.6 Transverse shearing proposed by Farley and Jones (1992a). This failure mode is equivalent to Hull's fragmentation mode 11 Figure 2.7 Lamina bending proposed by Farley and Jones (1992a). This failure mode is equivalent to Hull's splaying mode 12 Figure 2.8 Local buckling proposed by Farley and Jones (1992a). This failure mode is equivalent to Hull's progressive folding mode 12 Figure 2.9 Brittle fracturing proposed by Farley and Jones (1992a). This failure mode is a mix of transverse shearing and lamina bending 13 Figure 2.10 Schematic of failure mechanisms observed in the crush zone of 3-D braided tubes failing in the splaying mode, showing (a) compressive fibre breakage, (b) kink banding, (c) bending failure, (d) shear slippage due to bending, (e) shear deformation and transverse cracking in braider yarns, (f) resin fracture and debonding (g) crack propagation (Chiu et al., 1997). 14 Figure 2.11 Typical crush load-displacement curves of braided composite tubes studied by Chiu et al. (1999). (K = Kevlar, C = carbon) 16 Figure 2.12 (a) Biaxial and (b) triaxial braid architectures (Karbhari et al., 1997) 21 Figure 2.13 Comparison of specific energy absorbance traces for biaxial and triaxial braided circular tubes (Karbhari et al., 1997) 21 Figure 2.14 Sketch of a triaxially braided perform showing braiding angle and axial/braiding yarns (Chiu et al., 1998) 22 v i List of Figures Figure 2.15 Variation of crush stress and SEA of tubes with (a) varying braiding angle and (b) varying axial yarn content (Chiu et al., 1998) 23 Figure 2.16 Details of the finite model to predict the crushing response of circular tubes used by Farley and Jones (1992b) 32 Figure 2.17 Comparison of predicted and measured energy absorption for (a) Kevlar/epoxy and (b) graphite/epoxy Farley and Jones (1992b) 32 Figure 2.18 (a) Typical failure behaviour in 20-ply carbon/PEEK composite tubes and response of (b) 3D model and (c) axisymmetric model (Tay et al. 1998)... 33 Figure 2.19 Schematic showing (a) the finite element model and (b) comparison of predicted load-displacement profile to the experimental for circular tubes (Beard and Chang, 2002) 34 Figure 2.20 Schematic showing (a) the finite element model and (b) comparison of predicted load-displacement profile to the experimental for square tubes (Beard and Chang, 2002) 34 Figure 2.21 Predicted failure behaviour of filament wound graphite fibre reinforced tube (Kerth and Maier, 1994) 37 Figure 2.22 Comparison of numerical (a) force-displacement profile and (b) internal energy profile to experimentally measured trends (Kerth and Maier, 1994). 38 Figure 2.23 Representation of damage progression in a laminate shell element by successive failure of individual plies defining the laminate (Caliskan, 2002). 38 Figure 2.24 Examples of stress-strain profiles that result from various Ej and SLIM, values (Xiao et al., 2003) 39 Figure 2.25 Mixed model used by Morthorst and Horst (2004) to capture the in-plane and delamination properties of a laminate 40 Figure 2.26 Schematic of general strain softening material model showing (A) damage initiation point and (B,C) unloading along a path defined by a reduced elastic modulus 41 Figure 2.27 Illustration of mesh sensitivity observed when a strain-softening material model is defined, showing (a) rectangular panel with a zone of decreased strength subjected to tensile load, localization of failure into a band of elements in (b) mesh A density and (c) mesh B density, and (d) dependency of the peak loads and post-peak response on the mesh density (Bazant and Planas, 1998) 41 vn List of Figures Figure 2.28 Typical stress-strain response defined for springs in cohesive failure delamination model (Fleming, 2001) 43 Figure 3.1 Global Euler buckling of a thin plate under a remote applied stress 46 Figure 3.2 Transverse tensile failure of a unidirectional lamina 47 Figure 3.3 Delamination of outer plies in an idealized laminate comprised of four plies. 47 Figure 3.4 Two modes of elastic microbuckling proposed by Rosen (1965) 48 Figure 3.5 Details of kink banding with fibres aligned with remotely applied load 49 Figure 3.6 Graphical representation of differences between out-of-plane and in-plane kinking 49 Figure 3.7 Specimen geometry, showing the location of the strain gauges and the distribution of damage along the width of the specimen after compressive loading (Sivashanker et al., 1996) 52 Figure 3.8 Progression of stress (as estimated by measured product of gauge strain and elastic modulus) at location A as the damage tip approaches and passes.... 53 Figure 3.9 Summary of the various damage mechanisms observed in multidirectional laminate (a) general through thickness out-of-plane failure mechanisms showing kinking, delaminations and off-axis matrix damage (b) region CI, showing out-of-plane kinking (c) close up of kink band (d) region C2, showing delamination and kink band formation in 45° layer (e) region C3, showing kinking in 0° layer extending into the 45° layer (Sivashanker, 2001) 55 Figure 3.10 Experimental study done by Quek et al. showing (a) details of the specimen dimension and (b) loading directions (Quek 2004a) 57 Figure 3.11 Typical failure mechanism observed in braided composites due to compression showing global out-of-plane tow buckling and inter-tow matrix damage (Quek, 2004a) 57 Figure 3.12 Degree of crimp observed in non-tensioned (a) and tensioned (b) braided panels (West and Adams, 1999) 58 Figure 3.13 Uniaxial numerical compressive stress-strain response compared to the experimental data for varying degrees of imperfection. Initial modulus is well predicted, but the peak stress is over-predicted (Quek et al. 2004b)... 62 Figure 3.14 Generalized normalized modulus curve. 64 Figure 3.15 Generalized damage growth curve 65 viii List of Figures Figure 3.16 General one-dimensional stress-strain response, assuming previous damage growth and normalized modulus curves and that F equals the strain in the direction of interest 65 Figure 3.17 Schematic of overall analog model showing the two main boxes, the slider and the lock-spring element 69 Figure 3.18 (a) Schematic of Laminate Box, showing the arrangement of the T/C springs and fuses in a parallel infinite series. Repeating units (RU) have the same stiffness values, but each successive fuse has a higher failure load (b) Force-displacement response of a single RU. (c) Force-displacement response of laminate box consisting of an infinite number of RU 70 Figure 3.19 (a) Schematic of Rubble Box, showing the arrangement of gap and spring elements (parallel infinite series) and the slider element below. Repeating units (RU) have the same stiffness values, but each successive RU has a higher gap size (b) Force-displacement response of a single repeating unit showing inactivation before gap closure, (c) Force-displacement response of rubble box 72 Figure 3.20 Response of individual boxes and analog model in force-displacement domain showing elastic and damaging regions. Possible unloading paths (or reloading paths with prior damage) are also shown as dashed lines 73 Figure 3.21 Idealized tensile damage growth in a representative volume of a [ 9 0 / 0 / 9 0 J . laminate and force-displacement response predicted by the analog model. 74 Figure 3.22 Idealized compressive damage growth in a [ 9 0 / 0 / 9 0 J s RVE and resulting force-displacement response predicted by the analog model 74 Figure 3.23 Schematic showing the erroneous bi-linear unloading path and corrected linear unloading path for a one-dimensional case assuming modified Bazant scaling 76 Figure 3.24 Schematic comparing the calculation of elemental post-peak strain in (a) the previous version of COD AM and (b) after the correction 77 Figure 3.25 Schematic showing the progression of failure mechanisms in compression. Matrix damage is assumed to initiate first, followed by fibre damage in the kink band 80 Figure 3.26 Analytical Tensile Response (x=y): (a) Normalized modulus curve and (b) Stress-strain curve for the single element simulations 85 Figure 3.27 Graphical representation of various energy densities for a given element subjected to a given pre-failure tensile strain (x=y) 86 ix List of Figures Figure 3.28 Breakdown of various energies (from stress-strain profile) for an element, showing three components of the total energy density 87 Figure 3.29 Analytical Compressive Response (x=y): (a) Normalized modulus curve and (b) Stress-strain curve for the single element simulations 88 Figure 3.30 Graphical representation of energy development in compression (x=y) 89 Figure 3.31 Graphical representation of energy development in compression (x=y) 90 Figure 3.32 Breakdown of various energies (from stress-strain profile). 91 Figure 3.33 Loading applied to a single (a) shell element and (b) brick element in the local and global co-ordinate systems 92 Figure 3.34 Numerical stress-strain profile for 5 mm shell and solid element (Tension). Results agree with the analytical solution 93 Figure 3.35 Numerical output of total energy density for shell and solid element (Tension). Agrees with the analytical solution 94 Figure 3.36 Numerical stress-strain profile for 5 mm shell and solid element (Compression). Agrees with the analytical solution 94 Figure 3.37 Numerical output of total energy density (Compression). Agrees with the analytical solution 95 Figure 3.38 Stress-strain profile for strain reversal simulation. Loaded in compression, then unloaded and loaded to failure in tension 96 Figure 3.39 Numerical output of total energy density (Compression - Tension) 97 Figure 4.1 Schematic of basic 2-dimensional braid architecture, showing the four loading directions in the panel tests 100 Figure 4.2 Schematic showing the testing configuration for tension and compression. 101 Figure 4.3 Typical failure characteristics in tension (a) and compression (b). Local x-direction (DeTeresa et al., 2001) 102 Figure 4.4 Schematic showing the orientation of the tows in the four loading directions (shown for one axial-braider-braider tow group). The off-axis tows (tows not aligned with the loading direction) intersect at least one edge 102 Figure 4.5 Example of the plateau stress observed in compression (local x-direction). Apparent strain calculated using machine displacement (DeTeresa et al., 2001) 103 x List of Figures Figure 4.6 Fibre bridging observed when a crack propagates parallel to a fibre tow, forming a large damage (or process) zone (Salvi et al., 2003) 106 Figure 4.7 Examples of compressive damage heights observed in (a) the local y-direction of [0/±30]2 carbon/vinyl ester panels(DeTeresa et al., 2001) and (b) the local x-direction, (c) the local y-direction and (d) 45° from the local x-direction of [0/±60]2 carbon/vinyl ester panels(Xiao, 2003b) 110 Figure 4.8 Tensile simulation details. Overall dimensions are 38 x 162 mm. Ramped displacement is applied to top nodes, while bottom nodes are fixed 113 Figure 4.9 Example of the failure behaviour (a) and nominal stress-average strain response (b) of the panel in tensile simulations, shown here for the upper bound local x-direction. Predicted initial modulus and peak stress are similar to upper bound experimentally measured values 114 Figure 4.10 Comparison of upper and lower bound predictions of tensile peak stress with experimentally measured values 115 Figure 4.11 Compressive simulation details. Overall dimensions are 38 x 22 mm. Ramped displacement is applied to top nodes, while bottom nodes are fixed. 116 Figure 4.12 Example of the failure behaviour (a) and nominal stress-average strain response (b) of the panel in compressive simulations, shown here for the upper bound local x-direction. Note that compressive stresses and strains are shown as positive values. Predicted values are similar to the experimental values 116 Figure 4.13 Comparison of upper and lower bounded predictions of compressive peak stress with experimentally measured values 117 Figure 4.14 Stress-displacement profile of specimen first compressed before being pulled in tension to failure. Panel still has significant strength in tension. Note that the reported displacement is measured from the machine. (DeTeresa et al. 2001) 118 Figure 4.15 Numerical prediction of the nominal stress-average strain response during a strain reversal (compression-unload-tension). The peak stresses, unloading behaviour, and energy dissipation are reasonably predicted 119 Figure 5.1 Typical failure characteristics of braided tubes without a plug initiator from the front (a and b) and from the bottom (c). A significant amount of loose material gathers within the tube (DeTeresa et al., 2001) 121 Figure 5.2 Typical failure characteristics of braided tubes crushed with a plug initiator. Tube tears at the corner and the four fronds bend away from the plug (Xiao, 2003c) .' 122 xi List of Figures Figure 5.3 Schematic of the basic components in the tube simulation. The tube contacts either the plug or the rigid plate 124 Figure 5.4 (a) Displacement applied to the orthotropic tube, (b) Distribution of longitudinal and (c) circumferential strains along the length of the tube due to the applied displacement 128 Figure 5.5 Variation of (a) longitudinal and (b) circumferential strains at the bottom of the tube. Terms in brackets refer to strains normalized with respect to the applied global strain of -0.02 129 Figure 5.6 Comparison of predicted upper and lower bound force-displacement profiles to the experimentally measured trace for single ply tube with plug. Experimental results provided by (Xiao, 2003a) 130 Figure 5.7 General failure behaviour predicted for single ply braided tubes crushed with a plug initiator 131 Figure 5,8 Comparison of predicted (a) peak forces and (b) Specific Energy Absorption (SEA) to. experimental values. Error bars refer to upper and lower bound values 131 Figure 5.9 Fringe plot of longitudinal strains in single ply tube crushed with plug. Absolute minimum = -0.001, absolute maximum = -0.70 132 Figure 5.10 Fringe plot of circumferential strains in single ply tube crushed with plug. Absolute minimum = 0.0, absolute maximum = 0.74 132 Figure 5.11 Fringe plot of in-plane shear strains in single ply tube crushed with plug. Absolute minimum = -0.15, absolute maximum = 0.026 133 Figure 5.12 Variation of (a) longitudinal strains and (b) circumferential strain at a cross section 165 mm from the bottom of the tube 133 Figure 5.13 Growth and decay of various (a) global and (b) dissipative energy measures predicted by LS-DYNA as the tube crushes with lower bound parameters. 134 Figure 5.14 Fringe plot of the (a) local x-direction and (b) local y-direction damage parameters during initial stages of tube crushing with the plug initiator... 135 Figure 5.15 Variation of the SEA with tensile fracture energies (local y-direction) for single ply with plug tube simulations 137 Figure 5.16 Variation of the SEA with compressive fracture energies (local x-direction) for single ply with plug tube simulations 138 xu List of Figures Figure 5.17 Assumed distribution of the energy absorbed in the tube simulation with upper bound parameters, showing the components of energy absorbed by complete failure of the elements 140 Figure 5.18 Variation of SEA with changes in the longitudinal plateau stress for single ply with plug tube simulations 141 Figure 5.19 . Variation of SEA with changes in the coefficient of friction for single ply with plug tube simulations 142 Figure 5.20 Comparison of the overall failure behaviour predicted by different mesh densities for single ply with plug tube simulations 143 Figure 5.21 Comparison of the force-displacement profiles predicted for different mesh densities 144 Figure 5.22 Comparison of (a) peak forces and (b) SEA values predicted for different mesh densities for single ply with plug tube simulations 144 Figure 5.23 Comparison of two experimental force-displacement profiles for single ply tube without plug. Experimental results provided by (Xiao, 2003a) 146 Figure 5.24 Comparison of predicted upper bound force-displacement profiles to the experimentally measured trace for single ply tube without plug. Experimental results provided by (Xiao, 2003a) 146 Figure 5.25 Comparison of predicted lower bound force-displacement profiles to the experimentally measured trace for single ply tube without plug. Experimental results provided by (Xiao, 2003a) 147 Figure 5.26 General failure behaviour as predicted for single ply braided tubes crushed without a plug initiation showing morphology from the (a, c) side and (b, d) bottom for 15 mm and 130 mm crush distances 147 Figure 5.27 Comparison of predicted (a) peak forces and (b) Specific Energy Absorption (SEA) to experimental values. Error bars refer to upper and lower bound values 148 Figure 5.28 Variation of SEA with changes in the longitudinal plateau stress for single ply without plug tube simulations 149 Figure 5.29 Variation of SEA with changes in the longitudinal plateau stress for single ply without plug tube simulations 150 Figure 5.30 Comparison of two experimental force-displacement profiles for double ply tube without plug. Experimental results provided by (Xiao, 2003a) 151 xm List of Figures Figure 5.31 Comparison of predicted upper bound force-displacement profiles to the experimentally measured trace for double ply tube without plug. Experimental results provided by (Xiao, 2003a) 151 Figure 5.32 Comparison of predicted lower bound force-displacement profiles to the experimentally measured trace for double ply tube without plug. Experimental results provided by (Xiao, 2003a) 152 Figure 5.33 General failure behaviour as predicted for double ply braided tubes crushed without a plug initiator showing (a) side and (b) bottom morphology 152 Figure 5.34 Comparison of predicted (a) peak forces and (b) Specific Energy Absorption (SEA) to experimental values. Error bars refer to upper and lower bound values 153 Figure 5.35 Comparison of predicted upper and lower bound force-displacement profiles to the experimentally measured trace for four ply tube without plug. Experimental results provided by (Xiao, 2003a) 154 Figure 5.36 General failure behaviour as predicted for four ply braided tubes crushed without a plug initiation showing (a) side and (b) bottom morphology.... 155 Figure 5.37 Comparison of predicted (a) peak forces and (b) Specific Energy Absorption (SEA) to experimental values. Error bars refer to upper and lower bound values 155 Figure 5.38 Comparison of measured and predicted peak force values of single ply, double ply and four ply tubes. Error bars refer to upper and lower bound values 157 Figure 5.39 Comparison of measured and predicted SEA values of single ply, double ply and four ply tubes. Error bars refer to upper and lower bound values 157 Figure A.l Variation of plateau stress with effective strain (local x-direction) 172 xiv Nomenclature N O M E N C L A T U R E GENERAL NOMENCLATURE 1, 2, 3 Lamina principal directions x, y, z Laminate principal directions (Local co-ordinate system) X, Y, Z Global co-ordinate system a Acceleration D Diameter E Modulus (Ek)(, (Ek )f Kinetic energy, initial and final {Epo,) Potential energy, initial (Eabs)f Absorbed/dissipated energy, final (Eheat)/ Energy lost due to heat generation, final Ej Strain at peak stress (MAT58), i=x-tension, x-compression, y-tension, y-compress, shear Ef Fibre elastic modulus Em Matrix elastic modulus G Composite shear modulus G/ Fibre shear modulus Gm Matrix shear modulus k Interlaminar shear strength m Mass Q Lamina stiffness tensor Q Transformed lamina stiffness tensor SLIMj Ratio of the minimum stress to the maximum stress (MAT58), i=x-tension, x-compression, y-tension, y-compress, shear S Side length / Thickness v Velocity Vf Fibre volume fraction Vm Matrix volume fraction Vv Void volume fraction w Kink band width /3 Kink band inclination angle y Shear strain yy Composite yield strain in shear Ah Crush distance s Strain xv Nomenclature V Poisson's ratio vf Fibre Poisson's ratio vm Matrix Poisson's ratio a Normal Stress Remote stress Ultimate compressive strength T Shear stress Ty Composite yield stress in shear <t> Fibre rotation </> " Initial fibre misalignment ANALOG MODEL SPECIFIC NOMENCLATURE F Force F'f ' Fc Fibre fuse strength of T/C fuses in the Laminate box F'L ' Fc ' 1 m Matrix fuse strength of T/C fuses in the Laminate box K Yielding force of slider element Stiffness value of springs in the rubble box Stiffness values of fibre T/C springs in the Laminate box K >Km Stiffness values of matrix T/C springs in the Laminate box L Laminate box n Number of repeating units in the rubble box or laminate box R Rubble box S Displacement 8* .,Sgn Gap distances in the rubble box COD AM SPECIFIC NOMENCLATURE Superscripts C Compression T Tension Subscripts c Relating to characteristic size e Relating to an element General E Instantaneous secant modulus En Initial undamaged modulus xvi Nomenclature E N o r m a l i z e d instantaneous secant modulus Ems N o r m a l i z e d secant modulus at matr ix damage saturation F Damage potential function (effective strain function) Fmt' Fms Effec t ive strain parameters, in i t ia t ion and saturation o f damage i n matr ix Ffl, Ffs Effec t ive strain parameters, in i t ia t ion and saturation o f damage i n fibre Gf Fracture energy per unit area Ga Ini t ia l undamaged modulus hc Character is t ic height he E lement height k S c a l i n g factors required to scale fracture energy n Ra t io o f characteristic height to the element height K , L , M Strain interaction terms S, T , U Shear strain interaction terms ye Recoverable elastic strain energy density (per unit vo lume) yf Fracture energy density (per unit vo lume) Yjf, B a n d broadening energy density (per unit vo lume) ym To ta l energy density (per unit vo lume) yu Pre-peak energy diss ipat ion density (per unit vo lume) °'pirneau Compress ive plateau stress co Damage variable cofi Magn i tude o f the damage parameter attributed to fibre damage coms Magn i tude o f the damage parameter attributed to matr ix damage Acknowledgements A C K N O W L E D G E M E N T S I would like to acknowledge the help of many people during my study. Firstly, I would like to thank my advisors, Dr. Reza Vaziri and Dr. Anoush Poursartip for their expertise and encouragement during the course of this project. Their guidance and advice was very helpful. I would also like to thank the members of .the UBC Composites Group for their stimulating discussions and insightful suggestions. I am especially grateful to Navid Zobiery, whom has been abundantly helpful. I thank General Motors Corporation for the opportunity to work on such an interesting project and extend special thanks to Dr. Xinran Xiao for her guidance and support throughout the project. The financial support from the Natural Sciences and Engineering Research Council (NSERC), University of British Columbia, Dr. Reza Vaziri and Dr. Anoush Poursartip was greatly appreciated as I worked on this project. I would like to thank my friends who have provided a stimulating and very enjoyable environment in which to learn and grow. Lastly, and most importantly, I wish to thank my family. Their endless support, encouragement, and love is invaluable. To them I dedicate this thesis. xvin Chapter 1 Introduction C H A P T E R 1 INTRODUCTION The prediction of the response of fibre reinforced composites to severe loads remains a difficult task because of the complex nature of the failure mechanisms and because of the lack of accurate and robust predictive tools. Consequently, expensive experimental testing programs are often required before a composite part can safely be incorporated into a structure. Accurate predictive tools in this area would reduce the time and money spent on testing these composite parts and potentially act as a design aid for future applications using composite materials. 1.1 Motivation The use of composites in industrial applications is increasing. Composites, in the context of this thesis, refer to a heterogeneous material composed of continuous fibres embedded in a polymer matrix. Areas of increased composite use include aerospace (wings, fuselage, propellers), automotive (body parts, bumpers), marine (boats, off-shore structures), and sporting equipment (skis, snowboards, golf clubs, tennis rackets, bicycles). Composites are attractive alternatives to conventional materials in large part because of their high specific stiffness and strength (specific tensile strength of typical carbon/epoxy composites are 4-6 times larger than steel and aluminium). In addition to their superior material properties, designing with composites is extremely flexible. The material can be specifically tailored and optimized to meet the design requirement of the structure. With all of these advantages, there still exists a general reluctance to incorporate composite materials into common engineering applications. This reluctance is attributed to many factors, including uncertainties about performance, cost, and a lack of simple design guidelines. Recently, a considerable amount of attention has centred on using composite materials in the automotive industry because of their inherent weight saving attribute, which translates into better fuel efficiency. One particular application is the front rail structure, responsible for absorbing front-end crash energy. Typically, metallic tubular structures are employed, which 1 Chapter 1 Introduction absorb the energy via progressive buckling involving significant plastic deformations. Composite tubular structures, on the other hand, absorb energy through a complex combination of fracture mechanisms, including matrix cracking, fibre breakage, and delamination. The actual sequence and interaction of these failure processes is highly dependent on material and geometric properties, as well as boundary conditions. Given the compressive nature of the applied crushing loads, damage is dominated by compressive failure mechanisms. Although tensile failure in composites is generally better understood, the area of compressive failure in composites is still developing. Matrix cracking (or yielding), fibre kinking and delamination are generally accepted as being the main failure mechanisms in composites under compressive loads. The effects of imperfections, such as fibre waviness or discontinuities, which are not detrimental to the tensile behaviour, are magnified in compression. Fibre waviness has been shown to lead to premature damage growth in composites due to progressive microbuckling of fibres (West and Adams, 1999). This premature damage growth generally results in a compressive strength that is lower than the tensile strength. In order to efficiently integrate composite structures into applications dominated by compressive loads, it is crucial to be able to accurately predict their compressive behaviour. This involves accounting for the various failure mechanisms and possible interactions between these mechanisms. Two distinct types of constitutive material models, micromechanical and macromechanical, are typically used to predict the response of composites. Micromechanical models use detailed schemes that discretize the individual fibres and matrix, and describe the overall composite behaviour based on the strain and stress states of the constituents. These types of models can accurately capture the initial stiffness, peak stress and post peak response, but are inevitably computationally expensive when used to analyze real sized structures. Macromechanical models are typically defined at the level of a representative volume of element (RVE) and attempt to represent the effect of damage growth in the individual constituents by smearing the response over the RVE. These types of models tend to be much more computationally efficient, but often require parameters that are difficult or impossible to obtain through experiments, rendering them non-physical and limiting their predictive capabilities. To the best of the author's knowledge, a purely 2 Chapter 1 Introduction predictive material model capable of efficiently simulating damage progression in composites under a variety of loading conditions does not exist. 1.2 Goals and Objectives The definite lack of numerical models capable of capturing, the crushing response of composite tubular structures prompted an investigation into the effectiveness of CODAM, a material model previously developed at the University of British Columbia, in predicting the behaviour of these structures. CODAM is a continuum damage mechanics based model for composite materials incorporated into the non-linear finite element code LS-DYNA (Williams, 1998; Williams et al., 2003; Floyd, 2004). CODAM has previously been shown to be capable of predicting the tensile behaviour of composites (Floyd, 2004), but until recently its ability to simulate compressive damage growth and failure was not evaluated. Initial investigation indicated that an accurate compressive model needed to be developed and incorporated into CODAM. The objective of this thesis is to enhance the capabilities of the CODAM constitutive model by incorporating a compressive damage model based on experimentally observed compressive failure mechanisms in composites. The model is evaluated based on the accuracy of the predicted strengths of flat braided panels subjected to tensile loads, compressive loads and load reversals. The model is further validated by comparing the predicted damage propagation and energy absorption in the axial crushing of braided composite tubes to the experimentally measured values. A supplementary goal is to study the sensitivity of the tube crushing analysis to variations in the model parameters. This provides a certain measure of the stability of the simulations and also offers some insight into methods of improving the energy absorption behaviour of these composite tubes. Successful predictions contribute to the goal of safely and efficiently designing such parts for the ultimate incorporation into vehicles or other crashworthiness applications. 1.3 Outline of the Thesis This section presents a broad outline of the thesis. Chapter 2 presents a brief overview of the experimental and numerical studies that have been reported in the literature on the axial 3 Chapter 1 Introduction crushing o f composi te tubes. The ax ia l tubes serve to absorb the crash energy associated w i t h front-end impact situations i n vehicular structures. The exper imental studies indicate that these composi te structures offer a compet i t ive alternative to their meta l l i c counterparts. The numer ica l s imula t ion o f the progressive damage growth i n these structures is a diff icul t task because o f the compl ica ted nature o f mic romechan ica l damage i n composi te materials. There is a trade o f f between the accuracy p rov ided by mic romechan ica l models and the efficiency and robustness p rov ided by macromechanica l models . Chapter 3 concentrates o n the further development o f C O D A M , by expanding the appl icabi l i ty o f the consti tutive mode l to compress ively dominated failures i n composite materials. A b r i e f literature r ev iew o f the failure o f composi tes i n compress ion presents the general mic romechan ica l aspects o f this type o f failure. The consti tutive mode l is then refined based o n these observed failure mechanisms. Th i s chapter concludes by ver i fy ing the incorporat ion o f the m o d e l improvements v i a single element s imulat ions . Chapter 4 presents the numer ica l s imula t ion o f braided composi te panel failure due to tension and compress ion loads, as w e l l as load reversals. These s imulat ions serve to validate the mode l and calibrate the parameters for the ax ia l crushing o f tubes constructed o f identical material . Chapter 5 demonstrates C O D A M ' s abi l i ty to capture the progressive damage growth and the energy absorption behaviour i n the dynamic ax ia l crushing o f a variety o f braided composi te tubes. Success is measured i n terms o f overa l l predicted failure behaviour, force-displacement profi les and energy absorbed per unit mass o f crushed mater ia l . A parametric study indicates the important energy absorbing mechanisms and provides some insight into improv ing the specif ic energy absorbed by these structures. Chapter 6 presents the conclus ions drawn f rom the Current study and recommends future w o r k i n the area o f composi te damage mode l ing . 4 Chapter 2 Background and Literature Review C H A P T E R 2 B A C K G R O U N D AND L I T E R A T U R E REVIEW 2.1 Introduction The ability of a vehicle to absorb and dissipate impact energy during a crash is known as the crashworthiness of the structure. To ensure that the passenger compartment and thus the occupant(s) are not subjected to high decelerations, the impact energy should be absorbed in a controlled manner. Traditionally, the primary energy absorption devices for head on collisions have been side rail box sections typically made of steel or aluminium as shown in Figure 2.1. These sections absorb the crash energy by extensive plastic deformation leading to controlled buckling of the section. Figure 2.1 Principal scheme of the energy absorption devices in an automobile (Kerth and Maier, 1994). In the past two decades, a considerable amount of research has been done on using composite materials as the primary energy absorption system. Fibre reinforced composites offer a number of advantages over steel or aluminium. They generally have higher specific modulus and specific strength values, resulting in lighter structures. Additionally, they are highly tailorable to the designers needs. Composites, however, respond differently to loads than their metal counterparts. Most often, they respond in a brittle fashion, with resulting failure mechanisms that are quite different than those observed in metals. Composites fail through a 5 Chapter 2 Background and Literature Review sequence of fracture mechanisms involving fibre fracture, matrix cracking, fibre-matrix de-bonding, de-lamination and interply separation. The actual sequence of events and failure mechanisms are dependent on many properties including the geometry of the structure, lamina orientation, type of trigger, and crush speed. These properties can all by tailored to result in efficient energy absorbing systems. Extensive experimental and numerical research work has been done on the crushing behaviour and failure mechanism of axial loaded thin-walled polymer composite tubes, and many papers have been published on this topic. This chapter attempts to summarize the main features and results of these investigations. Experimental studies focused on the effects of varying design parameters on the specific energy absorption (energy absorbed per unit mass) and the failure modes. A limited number of numerical studies attempted to replicate the experimentally observed crushing characteristics, using various implicit and explicit finite element codes. 2.2 Composite Tube Crushing - Experimental Experimental studies were a necessary first step in the pursuit of understanding and ultimately predicting how composite tubes react to axial static and dynamic compressive loads. Due to the high tailorability of composites, a large compilation of papers that investigated the energy absorption trends of tubes with various geometry and dimensions, matrix and fibre types, trigger mechanisms, processing conditions, fibre volume fractions and testing speeds, were produced. Some of the more noteworthy conclusions are summarized below. 2 . 2 . 1 Test Methodologies Two main types of testing conditions, quasi-static and impact are typically used to quantify the energy absorption capacity of composite tubes. In quasi-static testing, the tube is progressively crushed at a constant low speed. This does not accurately represent the dynamic conditions that would exist if such structure were to be used in front end of a vehicle structure. Because many constituent materials are strain rate sensitive using quasi-static tests to draw conclusions on the effectiveness of tubes as energy 6 Chapter 2 Background and Literature Review absorbers can lead to false results. These types of tests, however, are commonly done because of their many advantages over impact tests. They are simple and easy to control, and the crushing process takes place over a longer period of time, making it easier to measure. , In impact testing, the speed of testing decreases from the initial impact speed to complete rest as the tube absorbs the energy. Most impact tests are done using drop towers. Arnaud and Hamelin (1998) used a different system called the block-bar, used to obtain reliable results for low to medium strain rates. Impact tests do represent the true simulation of a crash event, thus taking into consideration any strain-rate sensitivities of the system. However, because the test takes place over such a small time period, advanced equipment is necessary to obtain required results. 2.2.2 Crushing Modes and Mechanisms Tubes in axial compression fail either catastrophically or progressively. Catastrophic failures are characterized by an abrupt load increase to a peak force followed by a low post failure load. Catastrophic failures are common in long thin-walled tubes due to global Euler buckling or in brittle tubes due to unstable crack growth. Failure of this sort leads to low energy absorption. Progressive failures are characterized by an increase in load to a peak load followed by a slight decrease to a crushing load, which is sustained in a controlled and predictable manner throughout the remainder of the crushing process (Figure 2.2). The energy absorbed in progressive crushing is greater than that absorbed by tubes failing catastrophically. Incorporating a trigger mechanism, most often in the form of a chamfer at one end of the tube, is an effective way of initiating progressive failure. This trigger decreases the initial peak load, prolongs the triggering stage, reduces the probability of catastrophic failure and avoids unstable crushing. 7 Chapter 2 Background and Literature Review Load Triggering Progressive crushing O length Figure 2.2 Schematic showing the general load-displacement features characteristic of progressive axial tube crushing. Thornton (1990) investigated the energy absorption trends of pultruded glass fibre reinforced plastic tubes with either bevel or tulip triggers (Figure 2.3). The energy absorption was found to be consistently higher in both static and dynamic crushing tests in tubes treated with the tulip initiator. Sigalas et al. (1991) investigated the effect of various chamfer angles in glass cloth/epoxy tubes. The tubes had constant diameters and wall thicknesses with chamfer angles ranging from 10 to 90 degrees. The mechanisms of failure in the trigger were found to be dependent on the chamfer angle. Song et al. (2002) investigated the effect of using triggers at both ends of glass/epoxy circular tubes as shown in Figure 2.3. The initial slope and the peak load of a double-chamfered specimen were found to be less than that of either a single or non-chamfered tube. Additionally, the trigger length was longer. However, the overall energy absorbed was approximately the same in the single and double chamfered tubes. Recently, Mamalis et al. (2004) tested carbon fibre/epoxy square tubes that were not treated with any kind of trigger mechanism, and observed progressive crushing in 45% of the tubes. 8 Chapter 2 Background and Literature Review $55 , B E V E L . . T U L I P (a) (b) Figure 2.3 Two trigger mechanisms investigated by Thornton (1990) and double chamfered tubes studied by Song et al. (2002). Hull (1991) introduced the concept of two main modes of progressive failure, progressive folding and progressive crushing, as shown in Figure 2.4. The progressive crushing mode was further broken down to two extreme failure modes, splaying and fragmentation, as shown in Figure 2.5. A detailed study was done on the fundamental mechanisms of these two modes of progressive crushing. A unified approach to understanding progressive crushing was presented, which stated that the specific energy absorption was directly related to the behaviour and geometry of the crush zone. The crush zone was in turn dependent on the following variables: (1) Microfracture processes at the crush zone (2) Forces acting at the crush zone (3) Microstructural variables associated with the composite material (4) Shape and dimensions of the component (5) Testing variables 9 Chapter 2 Background and Literature Review Figure 2.4 Schematic showing the two main failure modes, (a) progressive folding and (b) progressive crushing (Hull, 1991). Load p laten (a) (b) Figure 2.5 Details of the two extreme modes ofprogressive crushing, (a) splaying and (b) fragmentation proposed by Hull (1991). Farley and Jones (1992a) later introduced and described three essential crushing modes exhibi ted by cont inuous fibre-reinforced composite materials. T h e y were transverse shearing, l amina bending and loca l buck l ing . A d d i t i o n a l l y they descr ibed a fourth mode cal led brittle fracturing, w h i c h was essentially a combina t ion o f transverse shearing and lamina bending modes. Far ley and Jones proposed that brittle fibre-reinforced composites cou ld fa i l b y any o f the 3 crushing modes, but most often exh ib i t bri t t le fracturing. D u c t i l e fibre-reinforced composi tes exhibi t loca l buck l ing crushing mode on ly . Tubes crushing i n the transverse shearing mode exhibi ted the highest c rushing eff ic iency. The three modes o f failure identif ied by H u l l are ident ical to the three modes o f fai lure descr ibed by Far ley and Jones. E a c h o f the modes can be summarized as fo l lows : 10 Chapter 2 Background and Literature Review Transverse Shearing or Fragmentation Mode • Brittle fibre reinforced tubes exhibit this crushing mode. • This mode is characterized by a wedge-shaped laminate cross-section with one or multiple short interlaminar and longitudinal cracks that form partial lamina bundles. The lengths of the cracks are less than the thickness of the lamina (Figure 2.6). • Mechanisms like the interlaminar crack growth and fracturing of lamina bundles control the crushing process. • The main energy absorption mechanism is fracturing of lamina bundles. /Shearing away of edges of lamina bundles Interlaminar cracks Scalloped end of crushed tube Longitudinal cracks Figure 2.6 Transverse shearing proposed by Farley and Jones (1992a). This failure mode is equivalent to Hull's fragmentation mode. Lamina Bending or Splaying Mode • Brittle fibre reinforced tubes exhibit this crushing mode. • This mode is characterized by very long interlaminar, intralaminar and parallel to fibre cracks but the lamina bundles do not fracture. The lengths of the cracks are greater than ten times the thickness of the lamina (Figure 2.7). • Mechanisms like the interlaminar, intralaminar and parallel to fibre crack growth control the crushing process. • The main energy absorption mechanism is matrix crack growth. Two secondary energy absorption mechanisms related to friction are also functioning. 11 Chapter 2 Background and Literature Review Extensive bending Figure 2.7 Lamina bending proposed by Farley and Jones (1992a). This failure mode is equivalent to Hull's splaying mode. Local Buckling or Progressive Folding • Both brittle and ductile fibre reinforced tubes exhibit this crushing mode. • This mode is characterized by the formation of local buckles (Figure 2.8). • Mechanisms like plastic yielding of the fibre and/or matrix control the crushing process. Formation of local buckle Interlaminar cracks at buckle Figure 2.8 Local buckling proposed by Farley and Jones (1992a). This failure mode is equivalent to Hull's progressive folding mode. Brittle Fracturing • Brittle fibre reinforced tubes exhibit this crushing mode. • This mode is a combination of transverse shearing and lamina bending (Figure 2.9). 12 Chapter 2 Background and Literature Review Figure 2.9 Fractured lamina bundles Interlaminar cracks Scalloped end of crushed tube Bending of lamina bundles Brittle fracturing proposed hy Farley and Jones (1992a). This failure mode is a mix of transverse shearing and lamina bending. Chiu et al. (1997) statically crushed 3D braided composite carbon/epoxy and Kevlar/epoxy square tubes and studied the crushing failure modes. The carbon/epoxy tubes displayed a splaying mode of failure, whereas the Kevlar/epoxy tubes exhibited a stable progressive folding failure. An inspection of the carbon/epoxy crush zone indicated that compressive shear failure, bending fracture, kink banding, shear and tensile fracture of the braiding bundles, matrix crush, interface debonding and crack propagation were the main failure processes (Figure 2.10). Chiu et al. concluded that in the 3D braided composite tubes, the axial yarns are responsible for absorbing most of the energy. 13 Chapter 2 Background and Literature Review Figure 2.10 Schematic offailure mechanisms observed in the crush zone of 3-D braided tubes failing in the splaying mode, showing (a) compressive fibre breakage, (b) kink banding, (c) bending failure, (d) shear slippage due to bending, (e) shear deformation and transverse cracking in braider yarns, (f) resin fracture and debonding (g) crack propagation (Chiu et al, 1997). 2.2.3 Energy Absorption Calculations The effectiveness of composite tubes as energy absorbers in crashworthy structures is most often indicated by the specific energy absorption (SEA). The SEA is defined as the energy absorbed per unit mass of material. In the case of tube crushing, the energy absorbed during crushing is generally taken as the area under the load-displacement curve from the beginning to the end of progressive crushing. This total energy absorption is then divided by the mass of the crushed portion of the tube. In some studies, the specific energy absorption was calculated as the mean crush load divided by the mass per unit length. Farley and Jones (Farley, 1983; Farley, 1986a; Farley, 1986b; Farley and Jones, 1992a; Farley and Jones, 1992b) typically used the specific sustained crushing stress, defined as the sustained crushing stress divided by the density to measure the crashworthiness of their test specimens. 2.2.4 The Effect of the Reinforcing Fibre Reinforcement material has a substantial effect on the crushing behaviour of composite tubes. Table 2.1 lists some of the fibre/matrix systems that have been investigated in studies 14 Chapter 2 Background and Literature Review on composite tube crushing. Thornton and Edwards (1982) and Farley (1983) observed that the type of reinforcement influenced the overall failure behaviour. Tubes reinforced with carbon or glass progressively crushed in either fragmentation or splaying modes. Kevlar reinforced tubes, however, failed by progressive folding. Farley believed the lower strain to failure of the glass and carbon fibres, as compared to Kevlar fibres, was the reason for the difference in crushing modes. The ability of composite systems to absorb energy was also observed to be dependent on the reinforcing fibre. A review conducted by Thornton et al. (1985) indicated that graphite/epoxy tubes, overall, had SEA values greater than that of Kevlar/epoxy and glass/epoxy tubes of similar configuration. This was credited to the lower density of carbon fibres. The same review indicated that the fibre modulus has only a small influence on the SEA. In a study conducted by Farley (1983), results indicated that tubes reinforced with longitudinally oriented graphite fibres had SEA values higher than tubes reinforced with longitudinally oriented glass or Kevlar fibres. Farley (1986a) also investigated the effect of maximum strain on the specific energy absorption. Tubes reinforced with higher strain-to-failure fibres had higher SEA values as long as the matrix had a greater strain-to-failure than the fibres. Karbhari et al. (1997) measured the energy absorbance of cylindrical composite tubes with triaxial [0/±45] braid architecture consisting of glass, Kevlar and graphite fibres either by themselves or as hybrids embedded in a vinyl ester resin. Tubes with glass fibre tows in the [±45] direction and carbon fibre tows in the axial direction had the highest values of both the peak and mean crush loads. They also displayed the best performance on the basis of both total energy absorbance, energy per unit deformation, and SEA. The failure mode of the glass-carbon triaxial hybrid was a mixture of splaying and extreme local fibre fracture. Chiu et al. (1999) also studied the influence of fibre type and hybrid structure of 2D triaxially braided composite tubes on the crush response, failure modes and specific energy absorption. Kevlar and carbon were used in the braided preform and set in an epoxy resin to create 6 different hybrid cylindrical structures. In each of the tubes, the braiding angle was 30°. All tubes failed by one of three different modes, splaying, folding or spiral-curling. The splaying mode was seen in tubes with biaxial carbon tows and either carbon or Kevlar axial tows. The folding mode was seen in tubes with biaxial Kevlar tows and either 15 Chapter 2 Background and Literature Review Kevlar or carbon axial tows. The spiral-curling mode of failure, which occurs when failure initially causes splaying but eventually, due to the fracture of some of the biaxial tows, the whole tube begins to rotate to form a spiral, was seen in tubes that had one biaxial carbon tow and one biaxial Kevlar tow and either a carbon or Kevlar axial tow. Based on these observations, Chiu et al. believed that the type of failure mode was determined by the braider yarns. The crush load-displacement curves of the tubes tested are shown in Figure 2.11. Tubes made entirely of carbon fibre tows exhibited the highest peak and mean crush load, and SEA. Additionally, it was observed that the braider yarns contributed significantly less to SEA than the axial tows. A carbon axial yarn supported higher crush loads and consumed higher crush energy than a Kevlar axial yarn. The authors concluded that tubes composed of one carbon biaxial yarn and one Kevlar biaxial yarn with a carbon axial yarn would fair the best in crashworthy structures because of its reasonably high value of SEA and its superior crush integrity over the pure carbon triaxial braided tubes. Figure 2.11 Typical crush load-displacement curves of braided composite tubes studied by Chiu et al. (1999). (K = Kevlar, C = carbon). 16 Chapter 2 Background and Literature Review Overall, two main conclusions can be made on the effect of the reinforcing fibres on the SEA capability of tubes. Firstly, graphite reinforced tubes generally lead to higher SEA values than glass or Kevlar reinforced tubes. This is attributed to the lower density of graphite fibres. Secondly, increasing the strain-to-failure of the fibre leads to a significant increase in SEA, whereas increasing the fibre stiffness will leads to only a slight increase in SEA. materials studied in composite tube crushing Table 2.1 List of various investigations. Glass-Epoxy Arnaud and Hamelin (1998) Fairfulland Hull (1987) Farley (1983) Gupta et al. (1997) Hamada a/. (1992a) Hul l (1988) Hull (1991) Schmueser and Wickliffe (1987) Sigalas et al. (1991) Song et al. (2002) Thornton (1979) Thornton and Edwards (1982) Thornton et al. (1985) Glass-Vinyl Ester Arnaud and Hamelin (1998) Karbhari et al. (1997) Karbhari and Haller (1998) Mamalis et al. (1997) Thornton (1990) Glass-Polyurethane Quek etal. (2001) Glass-Phenolic Thornton etal. (1985) Glass-Polyester Hull (1988) Hul l (1991) Mamalis etal. (1997) Thornton etal. (1985) Thornton (1990) Graphite-Epoxy Arnaud and Hamelin (1998) Chiu et al. (1997) Chiu et al. (1998) Chiu et al. (1999) Farley (1983) Farley (1986a) Farley (1986b) Farley (1991) Hamada e/a/. (1992b) Hull (1988) Hull (1991) Mamalis et al. (2004) Ramakrishna and Hull (1993) Schmueser and Wickliffe (1987) Schultze?a/. (2001) Thornton (1979) Thornton and Edwards (1982) Thornton et al. (1985) Graphite-Vinyl Ester Karbhari etal. (1997) Karbhari and Haller (1998) Graphite-PEEK Hamadaetal. (1992b) Hamada and Ramakrishna (1995) Hamada et al. (1995) Hamada etal. (1996) Kevlar-Epoxy Arnaud and Hamelin (1998) Chiu etal. (1997) Chiu etal. (1999) Farley (1983) Farley (1986b) Farley (1991) Schmueser and Wickliffe (1987) Thornton (1979) Thornton and Edwards (1982) Thornton etal. (1985) Kevlar-Vinyl Ester Karbhari etal. (1997) Karbhari and Haller (1998) Hybrid-Epoxy Arnaud and Hamelin (1998) Chiu etal. (1999) Farley (1983) Thornton and Edwards (1982) Hybrid-Vinyl Ester Castejon etal. (2001) Karbhari etal. (1997) Karbhari and Haller (1998) 17 Chapter 2 Background and Literature Review 2.2.5 The Effect of the Matrix Material For the most part investigators have looked at tubes constructed with a thermoset epoxy matrix material. Thornton et al. (1985) studied the effects of thermosetting resin and found that for glass reinforced tubes, SEA seemed to increase as you moved from phenolic to polyester to epoxy resin. It was proposed that this was due to the higher strength and stiffness of epoxy resin over the other two matrix materials. Based on his observations, Thornton proposed that as the stiffness of the thermosetting resin increases so will the specific energy absorption. Farley (1986a) studied the effects of the matrix strain-to-failure of thermosetting resin and found that increasing the strain-to-failure of the matrix to a value that is greater than the fibre strain-to-failure led to higher values of SEA. More recently, investigations have studied the advantages of using thermoplastic resins. Thornton (1990) investigated the crushing of tubes composed of glass fibre embedded in either a polyester (a thermoset) or a vinyl ester (a thermoplastic) resin. The tubes with vinyl ester resin had significantly higher values of SEA. Hamada et al. (1992b) also investigated the difference in energy absorption between tubes with a thermosetting epoxy resin and a thermoplastic PEEK resin. The superior energy absorption capability of carbon fibre/PEEK tubes (127 J/g for a ±30° architecture) was attributed to the high interlaminar fracture toughness of the thermoplastic matrix. The carbon/epoxy tubes, having much lower interlaminar fracture toughness values, absorbed less energy (53 J/g for a ±45° architecture). However, as these comparisons were made between two different lay-ups the results are not conclusive. Based on the observations, three general comments with respect to resin type are: 1) Increasing the matrix stiffness can increase the energy absorption capability of composite tubes reinforced with brittle fibres. 2) An increase in matrix strain-to-failure leads to an increase in SEA of brittle reinforced composite tubes. 3) Composite tubes with thermoplastic matrix materials seem to be better energy absorbers than thermosetting resin due to their higher fracture toughness, which provides resistance to crack growth. 18 Chapter 2 Background and Literature Review 2.2.6 The Effect of the Fibre Architecture Changes in fibre architecture, such as the fibre angle, the position of various plies within the wall thickness and the proportion of a given ply in the material significantly affect the crushing characteristics and the specific energy absorption capability of composite tubes. Additionally, the arrangement of the fibres within the laminate into forms such as woven, stitched or braided can be significant factors in the crushing response. Thornton and Edwards (1982) studied epoxy composite tubes with glass, graphite or Kevlar reinforcement in either a [0/90]n or a [45/-45]n lay-up. Changes in the lay-up that lead to higher values of the composite stiffness generally lead to higher energy absorption capabilities. Farley (1983) observed significant differences in the energy absorption trends of glass/epoxy, carbon/epoxy and Kevlar/epoxy composite tubes with fibre architectures of [0/±a] n , where the fibre angle varied from 15 to 90 degrees. The SEA of the glass and Kevlar reinforced tubes remained constant for fibre angles up to 45° and then increased afterwards, whereas SEA of the graphite/epoxy tubes decreased with increasing fibre angle up to 45° and then remained constant. These differences were attributed to the different crushing modes observed. Hull (1991) highlighted work that was done on the effect of fibre arrangement on progressive crushing. Research was completed on woven glass cloth/epoxy tubes and carbon fibre/epoxy unidirectional laminate tubes investigating the effects of the position, proportion and distribution of the hoop and axial layers had on the crushing modes and energy absorption of composite tubes. As the amount of hoop reinforcement in the woven glass cloth/epoxy tubes decreased, the crushing mode changed from a fragmentation mode to a splaying mode and, as a result, SEA decreased. The crushing modes and energy absorption in carbon fibre/epoxy unidirectional laminate tubes changed as the position, proportion and distribution of the hoop and axial layers changed. Unstable collapse occurred in tubes with all the hoop reinforcement on the inside, whereas stable collapse in the splaying mode occurred in tubes with the hoop reinforcement either all on the outside or split equally between the inside and outside. For tubes exhibiting fragmentation modes of failure, the SEA increased with a decrease in the ratio of hoop to axial reinforcement. Tubes with hoop reinforcement interleaved throughout the wall thickness crushed by fragmentation but absorbed significantly less energy than tubes with the same proportion of hoop reinforcement but 19 Chapter 2 Background and Literature Review placed equally on the inside and outside of the wall. Hull also summarized research that was done on filament-wound angle-ply glass fibre/polyester tubes, with winding angles varying from ±35 to 90 degrees. The specific energy absorption increased with increasing winding angle up to ±65 and then subsequently decreased. There was a systematic change in the load-displacement curves and crush zone morphology with the winding angle. Hamada et al. (1996) studied carbon fibre/PEEK composite tubes with fibre orientations of 0, ±5, ±10, ±15, ±20, ±25 and ±30 with respect to the axis of the tube. All of the tubes crushed in a progressive splaying manner except for the ±30° tubes, which failed catastrophically. The SEA and the crush mechanisms were found to be dependent on the fibre architecture. The SEA increased with increasing fibre orientations up to ±15° and thereafter decreased. This trend was believed to be due to the variation in the energy absorbed by fibre fracture, which increased initially with fibre orientation and then decreased. . Biaxial and triaxial braided composites of various fibre and matrix types have been investigated by Karbhari et al. (1997), Karbhari and Haller (1998), and Chiu et al. (1998, 1999). These types of designs proved to perform comparably as energy absorbers to tubes made of conventional unidirectional prepregs. Textile composite materials, achieved through weaving or braiding, offer an inexpensive and quick method of manufacturing composite tubes. As vehicle parts must be made quickly and cheaply in order to remain competitive, these low-cost materials are the most promising material for the automotive industry because of their low fabrication cost. Karbhari et al. investigated the effects of braiding parameters, such as braid yield, type of 2D braid pattern and number of layers used, on the crushing mechanisms and performance of cylindrical composite tubes. The crashworthiness of the tubes was measured using three characteristics, crush loads, energy absorption and damage and failure modes. Two types of braid architectures consisting of either glass, Kevlar and carbon fibres by themselves or as hybrids embedded in a vinyl ester resin were considered. The first, shown in Figure 2.12(a) was a biaxial braid with the braids oriented at [±45] and the second, shown in Figure 2.12(b) was a triaxial braid with biaxials at [±45] and a longitudinal [0] tow. A comparison of the specific energy absorption characteristics showed that overall, tubes with the triaxial braided architecture performed the best. As well, 20 Chapter 2 Background and Literature Review increasing the weight of the longitudinal tow and/or the number of layers led to increases in the specific energy absorption capability, but at the expense of higher weight. Figure 2.13 shows graphically the superior energy absorbance of the triaxial braided tubes over the biaxial braided tubes. Tubes failed in a splaying mode, an accordion type buckling mode or a progressive fragmentation into rings mode. The authors of this paper concluded that the effectiveness of braided composite tubes as an energy absorber is comparable to tubes made of conventional unidirectional prepregs. Figure 2.12 (a) Biaxial and (b) triaxial braid architectures (Karbhari et al., 1997). • •Jjg. >»' o> a> c LU U <D a. —• Biaxial (GB3) -----Triaxial.(GT3): 2 0 4 0 6 0 8 0 100 120 Deformation (mm) Figure 2.13 Comparison of specific energy absorbance traces for biaxial and triaxial braided circular tubes (Karbhari et al., 1997). Chiu et al. (1998) studied two different braiding parameters, the braiding angle and the axial yarn content, shown in Figure 2.14, to determine the effects these parameters have on the crush failure characteristics and the specific energy absorption capability. Two dimensional triaxially braided circular composite tubes composed of carbon biaxial yarns oriented at angles ranging from 20 to 45 degrees and carbon axial yarns with fibre content varying from 6k (6000 fibres/tow) to 18k embedded in epoxy resin were tested under uniaxial 21 Chapter 2 Background and Literature Review compression. The compression tests revealed that all of the tubes failed in a splaying mode. The width of the average frond increased with increasing braiding angle and decreasing axial yarn content. Tubes with a braider angle of 20°.showed the. narrowest fronds, and had the highest value of SEA, at 88:8 J/g. Figure 2.15 shows the changes observed in the crush stress and SEA when the two parameters, braiding angle and axial yarn content, were varied. A general trend of decreasing energy absorption with increasing braiding angle was not observed. This was thought to be due to the presence of additional energy consuming mechanisms in tubes with larger braiding angles. The value of SEA increased with increasing axial yarn content up to a limiting value of 6k. Thereafter the energy absorption remained relatively constant. The maximum energy absorption was seen in tubes with an axial content between 6 and 12k. Braking V«rn Figure 2.14 Sketch of a triaxially braided perform showing braiding angle and axial/braiding yarns (Chiu et al, 1998). 22 Chapter 2 Background and Literature Review 160 140 120 ==• 100 A 80 60 40 ^ 20 0 U fSE3 Mean crush stress : Specific energy,absorption'. ft • A20 A25 A30 A35 A40 Braiding angle type -(a) 100 80 60 40 20 A45 -A50 ••125-f r 75 50 H 25 S T ggg Mean crush stress ». Pss l^ Specific energy absortion 'If I m II 6603' - 6605 6509 - 6612 6615 . Axiahyarn content type -, (b) 80 •60-'. 40: •< : 20 Figure 2.15 Variation of crush stress and SEA of tubes with (a) varying braiding angle and (b) varying axial yarn content (Chiu et al, 1998). In general, an increase in SEA will result when the fibre architecture of a composite tube is modified such that either the number of fibre fractures, the material deformation, the axial stiffness of the composite material and/or the lateral support to the axial fibres increases. 2.2.7 The Effect of the Geometry of the Specimen Many different tube sections and geometry ratios have been studied to determine the effect such changes have on the crushing morphology and specific energy absorption. Table 2.2 lists the researchers that have studied circular and square/rectangular tube cross-sections. Thornton and Edwards (1982) statically crushed composite tubes having either circular, square or rectangular cross-sections with varying thickness to section dimension ratios (thickness/diameter (t/D) or thickness/side length (t/S)). There was a critical thickness to section dimension ratios, below which unstable collapse occurred. The cylindrical tubes demonstrated stable crush modes over the widest range of t/D ratios. Square and rectangular tubes required thicker wall sections to promote stable collapse and, in general, were less effective at absorbing energy. They concluded that for a given fibre lay-up and tube geometry circular cross-sections led to the most efficient crushing followed by square and then rectangular cross-sections. Thornton et al. (1985) made some additional observations 23 Chapter 2 Background and Literature Review on square, rectangular, and circular tubes. The SEA was found to be independent of the relative density, defined as the ratio of the volume of the tube to that of a solid of the same external dimensions. They also concluded that for composite structures that fail by fracture, the geometry dependence is negligible compared to that of structures that collapse by local buckling, have on the geometry. Farley (1986b) studied the effect of specimen geometry on the energy absorption capability in graphite and Kevlar reinforced epoxy tubes. Tubes with variable diameters and wall thickness to diameter ratios (t/D) were tested. The energy absorption was significantly affected by the t/D ratio and tended to increase nonlinearly as this ratio increased. Farley also determined that SEA remained relatively constant for Kevlar reinforced tubes with similar t/D ratios, but different diameters. This scalability was not seen in the graphite reinforced epoxy tubes. Fairfull and Hull (1987) also investigated glass cloth/epoxy composite tubes with various diameters and t/D ratios. The SEA increased as the diameter of the tubes decreased and increased with t/D ratios up to a value of about 0.2 and there after decreased. In this paper, Fairfull and Hull concluded that there could be "no universal relationship to predict the SEA of composite tubes". Hamada and Ramakrishna (1995) studied the effects of t/D ratios and absolute values of the wall thickness and diameter on the energy absorption behaviour of carbon fibre/PEEK composite tubes. Unstable collapse occurred in tubes with t/D ratios less than 0.015 whereas tubes with t/D ratios in the range of 0.015 to 0.25 crushed progressively. The specific energy absorption was dependent on the absolute value of t rather than the t/D ratio. It was concluded that in the case of thermoplastic composite tubes, the energy absorption characteristics are mainly influenced by the absolute value of t rather than the t/D ratio. Mamalis et al. (2004) studied the influence of the tube wall thickness on the peak load and SEA of tubes crushed with no triggering system. They noted that increasing the tube wall thickness generally led to an increase in the peak load, but did not produce any recognizable trend in the SEA. Based on the investigations of the effect of specimen geometry on SEA, the following conclusions can be made: 24 Chapter 2 Background and Literature Review 1) The crush zone fracture mechanisms, which determine SEA of a composite material, are greatly influenced by the tube dimensions and geometry. 2) For a given fibre lay-up and tube geometry, SEA generally increases as you move from rectangular to square to circular cross-section. 3) The specific energy absorption in thermosetting composites is a function of the t/D ratio, whereas in thermoplastic composites SEA is a function of the absolute value of the wall thickness. Table 2.2 List of researches that have studied circular, square and/or rectangular cross-sections. Circular Square and/or Rectangular Arnaud and Hamelin (1998) Chiu etal. (1997) Caliskan (2002) Kerth etal. (1996) Castejon et al. (2001) Mamalis etal. (1997) Chiu etal. (1998) Mamalis etal. (2004) Chiu etal. (1999) Schultz etal. (2001) Fairful land H u l l (1987) Thornton (1979) Farley (1983) Thornton and Edwards (1982) Farley (1986a) Thornton etal. (1985) Farley (1986b) Thornton (1990) Farley (1991) Gupta et al. (1997) Hamada etal. (1992a) Hamada <?f a/. (1992b) Hamada etal. (1995) Hamada and Ramakrishna (1995) Hamada etal. (1996) H u l l (1988) H u l l (1991) Karbhari etal. (1997) Karbhari and Haller(1998) Quek etal. (2001) Ramakrishna and Hu l l (1993) Schmueser and Wickliffe (1987) S c h u l t z ^ a / . (2001) Sigalas etal. (1991) Song et al. (2002) Thornton (1979) Thornton and Edwards (1982) Thornton etal. (1985) Thornton (1990) 25 Chapter 2 Background and Literature Review 2.2.8 The Effect of the Fibre Volume Fraction The amount of fibre within the matrix also affects the energy absorption characteristics of composite tubes. This aspect of composite tube design has not been studied as extensively as other areas and therefore results are not conclusive. Ramakrishna and Hull (1993) investigated the response of knitted-carbon-fibre-fabric/epoxy tubes under axial compressive loads with varying fibre content. The load displacement response, crush zone morphology and the SEA of tubes with fibre contents ranging from 5.25% to 31.5% volume were determined. Tubes with fibre content lower than 15% crushed in an irregular manner. The tubes with higher fibre content crushed progressively, in either the splaying or fragmentation mode, depending on the orientation of the inlay fibres. The value of SEA increased with increasing fibre content. The authors believed that increasing fibre contents to values as high as 40%o may further increase SEA. Other investigators have, however, reported a decrease in SEA in carbon fibre/epoxy tubes with an increase in fibre volume fraction. 2.2.9 The Effect of the Testing Speed A number of researchers have investigated the influence of strain-rate on the energy absorbing characteristics of composite tubes subjected to axial loading. Table 2.3 segregates the testing methodology used by various researchers into two categories, quasi-static testing and impact testing. Thornton (1979) studied the collapse of [0/90]n graphite/epoxy and glass/epoxy circular tubes at variable strain rates. The load-deformation curves became more jagged with an increase in crush rate, but there was no observable change in SEA. Farley (1983) tested [0/±a] n glass/epoxy, graphite/epoxy and Kevlar/epoxy tubes at both static and dynamic progressive crushing rates. The crushing load, failure modes, energy absorption mechanisms, and post crushing integrity observed in the dynamic tests were similar to results of the static tests. Schmueser and Wickliffe (1987) also studied graphite, glass or Kevlar/epoxy tubes, but contrary to Thornton and Farley's studies, a definite rate dependence was observed. They reported that although the failure modes were relatively consistent, there was up to 30% decrease in SEA in specimens tested dynamically. They concluded that these differences could be due to the sensitivity of interlaminar fracture toughness to loading 26 Chapter 2 Background and Literature Review rate. A study conducted by Schultz et al. (2001) on carbon fibre/epoxy tubes also showed lower SEA values in tubes tested dynamically compared to those statically tested. Farley (1991) reported opposite results in his study of graphite/epoxy and Kevlar/epoxy circular tubes with [0/±oc]2 and [±a]3 fibre architectures. Farley observed up to a 4 5 % increase in SEA, depending on whether the crushing mechanisms were functions of strain rate. Thornton (1990) tested pultruded glass fibre reinforced tubes in either a polyester or a vinyl ester resin at various crush rates. The glass/polyester resin SEA values generally increased (by up to 20%o) with increasing crush rate. Depending on the geometry, the glass/vinyl ester tubes showed either a positive or negative rate dependence, with up to a 10%> decrease in the SEA at higher testing rates. Thornton concluded that these differences were likely due to the different crush modes observed in specimens tested quasi-statically and dynamically. Karbhari and Haller (1998) reported increases in SEA as high as 2 8 % when braided tubes of carbon, glass and Kevlar fibres in a vinyl ester resin where tested dynamically, which contradicts Thornton's results. The damage mechanisms and failure modes were dependent on the crush rate. The peak and mean crush loads generally increased with the crushing rate with the largest recorded increase in the mean crush load. Triaxial braided architectures, in general had higher SEA values than the biaxially braided architectures. Karbhari and Haller concluded that the overall effects of changes in the crush rate depend on the strain rate sensitivity of the constituent fibre tows. Although the observed SEA trends with respect to testing rate are somewhat inconclusive, a general observation can be made. If the general failure characteristics, failure mechanisms and/or constituent material behaviour are strain rate sensitive, a change in the testing rate is likely to result in a change in the SEA value. 27 Chapter 2 Background and Literature Review Table 2.3 List of various testing rates investigated by researchers. Quasi-Static (rate < 1 mm/sec) Impact (rate > 1 mm/sec) Chiu etal. (1997) Arnaud and Hamelin (1998) Chiu etal. (1998) Caliskan (2002) Chiu etal. (1999) Farley (1983) Fairfulland Hul l (1987) Farley (1991) Farley (1983) Hull (1988) Farley (1986a) Hull (1991) Farley (1986b) Karbhari and Haller (1998) Farley (1991) Mamalis etal. (1997) Gupta etal. (1997) Schultz et al. (2001) Hamada etal. (1992a) Schmueser and Wickliffe (1987) Hamada et al. (1992b) Thornton (1979) Hamada et al. (1995) Thornton etal. (1985) Hamada and Ramakrishna (1995) Thornton (1990) Hamada et al. (1996) Hull (1988) Hull (1991) Karbhari etal. (1997) Karbhari and Haller(1998) Mamalis etal. (1997) Mamalis et al. (2004) Quek et al. (2001) Ramakrishna and Hull (1993) Schmueser and Wickliffe (1987) Schultz et al. (2001) Sigalase/a/. (1991) Song et al. (2002) Thornton (1979) Thornton and Edwards (1982) Thornton etal. (1985) Thornton (1990) 2.3 Composite Tube Crushing - Numerical Experimental studies confirmed that the specific energy absorption capacity of most composite structural parts under uniaxial crush load surpassed that of their metal or aluminium counterparts. The crushing process was dependent on many failure mechanisms, including fibre fracture, matrix cracking, fibre-matrix debonding, delamination, interply separation and friction. The cost and time necessary to develop and test composite crash elements or entire car front structures prompted the need to develop numerical methods to use as predictive tools in the design process. However, the complexity of the problem 28 Chapter 2 Background and Literature Review presents a considerable challenge to the goal of developing an accurate numerical model. The model must be able to handle both geometric and material nonlinearities, as well as large displacements, friction and contact between surfaces. Additionally, the material model must be able to represent the evolution of damage, failure initiation and post-failure behaviour of composite materials. Many researchers have evaluated the effectiveness and efficiency of various material models implemented into both implicit and explicit finite element programs, such as Engineering Analysis Language (EAL), PAM-CRASH, LS-DYNA, COMP-COLLAPSE, ABAQUS, Radioss, and MSC-Dytran at predicting the crushing behaviour of composite energy absorption devices. Some of the more relevant conclusions and results are discussed below. This area of study is relatively new and not very many numerical studies have been completed, as can be seen in the comparison of Table 2.4. 29 Chapter 2 Background and Literature Review Table 2.4 List of experimental and numerical studies on the axial crushing of composite tubes. Experimental Numerical Amaud and Hamelin (1998) Agaram et al. (1997) Caliskan (2002) Beard and Chang (2002) Castejon etal. (1998) Botkin et al. (1998) Chiu et al. (1997) Caliskan (2002) Chiu et al. (1998) Castejon et al. (1998) Chiu et al. (1999) Castejon etal. (2001) Fairfull and Hu l l (1987) Farley and Jones (1992b) Farley (1983) Haug et al. (1991) Farley (1986a) Kerth and Maier (1994) Farley (1986b) Kerth e ta l . (1996) Farley (1991) Mahmood and Wang (1997) Gupta et al. (1997) Morthorst and Horst (2004) Hamada et al. (1992a) Tay et al. (1998) Hamada et al. (1992b) Xiao et al (2003) Hamada and Ramakrishna (1995) Hamada et al. (1995) Hamada et al. (1996) H u l l (1988) H u l l (1991) Karbhari etal. (1997) Karbhari and Haller (1998) Kerth et al. (1996) Mamalis et al. (1997) Mamalis et al. (2004) Quek et al. (2001) Ramakrishna and H u l l (1993) Schultzetal. (2001) Schmueser and Wickliffe (1987) Sigalas etal. (1991) Song et al. (2002) Thornton (1979) Thornton and Edwards (1982) Thornton etal. (1985) Thornton (1990) 2.3.1 Finite Element Codes Numerical simulations of composite tube crushing have been done using both implicit and explicit finite element codes. An implicit code assembles large stiffness matrices at each time step and a considerable amount of memory is required, limiting the details or size of the problem. Implicit codes are useful when doing static or quasi-static analyses. 30 Chapter 2 Background and Literature Review Farley and Jones (1992b) used the implicit geometrically nonlinear finite element program, Engineering Analysis Language (EAL), to simulate a portion of the crushing process. The program was enhanced with the ability for material nonlinearities, crack initiation and propagation, and failure of elements. The following factors were considered by Farley and Jones when trying to successfully model the energy absorption capability of composite materials: (1) Crack initiation and growth (2) Material and geometrical nonlinearities (3) Fracture of individual lamina or laminae (4) Modeling around bifurcation loads (5) Element type and level of discretization Figure 2.16 illustrates the model used to predict the energy absorption capability of the tubes. The load was applied in increments until the maximum and minimum forces in one cycle of the progressive crushing phase were determined. The simulations were computationally more expensive than a corresponding linear-elastic model by approximately two orders of magnitude. The results from the finite element model were compared to test results obtained from an isotropic aluminium tube, and two composite tubes with varying lay-up of either Kevlar/epoxy or graphite/epoxy. Figure 2.17 shows that reasonable agreement was obtained between the analysis and the experiment, indicating that the important processes of tube crushing were included in the model. Improvements in the model were deemed necessary to improve the effectiveness and efficiency of the program. 31 Chapter 2 Background and Literature Review Zero-length springs for layer interconnection -v - \ ; - . Z e r o - l e n g t h springs have 3 axial and 3 rotational stiffness components Figure 2.16 Details of the finite model to predict the crushing response of circular tubes used by Farley and Jones (1992b). Specific sustained •: crushing-stress ' VN.rh/gi too •75' 50 •25 p Symbol size ofexperimental . data is large than data range O Analysis • Experiment O •'•"•' -" , ;"0'. ' " : n Q B. o. Q.: 15 30 45 60 ; 75 Ply angle, Q (a) o Specific;, sustained cushing •stress • . Nm/g 100 75 50 25 9( Symbol size; ofexperimental dataiis larger trian'data rarne 0 o Analysis , • Experiment • ,-, O f l . 15 : 30 • 45 S'J • Ply angle, (-) (b) •.;0-. 75. 90 Figure 2.17 Comparison of predicted and measured energy absorption for (a) Kevlar/epoxy and (b) graphite/epoxy Farley and Jones (1992b). Tay et al. (1998) used ABAQUS to perform an implicit analysis of the splaying mode of failure in a 20-ply carbon/PEEK composite tube. Due to the high computational energy required, the 3-dimensional model consisted of only one frond discretized into 4 layers of brick elements. Due to the lack of a post-failure constitutive relationship, all elements assumed to behave elastically. Delamination of the layers was simulated by the separation of node pairs once either a specified tensile or shearing break force was reached. Frictional effects and the relaxation of the hoop stresses as longitudinal cracks progressed were not included in this analysis. The initial peak force, the mean crush load and the overall shape of the force-displacement curve agreed well with the experimental results. However, in the simulation, the peak force was reached too early. A more computationally efficient axisymmetric model with 24 layers was also implicitly analyzed in ABAQUS. Both models 32 Chapter 2 Background and Literature Review captured the peak force and progressive crushing force reasonably well. The authors of this paper concluded that numerical modelling of the crushing behaviour of composite tubes is a challenging task because of the multiplicity of failure modes and their complex interactions. Additionally, some analyses are highly mesh sensitive and may provide unreliable results. Results obtained in this study showed fairly good agreement to the experimental force-displacement profile, but further investigation was deemed necessary. One major disadvantage of current finite element codes was the inability to model the post-failure constitutive relationships. tube wall plan view (a) (b) (c) Figure 2.18 (a) Typical failure behaviour in 20-ply carbon/PEEK composite tubes and response of (b) 3D model and (c) axisymmetric model (Fay et al. 1998). Beard and Chang (2002) successfully used ABAQUS implicit to simulate the complete crushing process of triaxially braided carbon fibre/epoxy-vinyl ester composite tubes using axisymmetrical and 3D elements. Experimental testing was completed on coupons to obtain information on the constitutive behaviour of the material necessary to develop an accurate material model. This material module calculates the overall effective material properties of a repeating unit cell (RUC) using the material properties of the individual constituents and the three-dimensional geometry of the composite. ABAQUS is used to calculate the stresses and strains, which are then passed to the module to predict the extent of damage and to update the material properties. Circular tubes with [0/±30]2, [0/±45]2, and [0/±60J2 braid architectures with plug type initiator were modelled using axisymmetric elements (Figure 2.19(a)). Figure 2.19(b) shows the comparison between the experimental and numerical tests for one of the circular tubes. Square tubes with a [0/±60]2 braid architecture were modelled using 3D elements (Figure 2.20(a)). Figure 2.20(b) shows the comparison of this simulation to 33 Chapter 2 Background and Literature Review experimental results. The numer ica l s imulations o f a l l tubes produced load-displacement profiles close to those measured experimentally. The mater ia l modu le implemented i n A B A Q U S was able to capture the load-displacement prof i le and predict the energy absorption properties fa i r ly accurately i n most tube s imulat ions . Rigid surface Axisymmetric < dements Coulomb friction (u.=0.3) between tube and plug 30 Load (kN)20 104 o R^ug ™ 7.94 mm [0ga/±60ijk]i (a) 10 20 Displacement (mm) (b) 30 Figure 2.19 Schematic showing (a) the finite element model and (b) comparison of predicted load-displacement profile to the experimental for circular tubes (Beard and Chang, 2002). 10 20 30 Displacement (ram) (b) Figure 2.20 Schematic showing (a) the finite element model and (b) comparison of predicted load-displacement profile to the experimental for square tubes (Beard and Chang, 2002). E x p l i c i t codes use the central difference method to integrate the equations o f mot ion for a system and do not fo rm or solve computat ional ly expensive stiffness matrices. A single iteration i n an exp l i c i t code is considerably less demanding numer i ca l ly than an i m p l i c i t iteration. E x p l i c i t codes are useful for dynamic problems, i n w h i c h the l oad is appl ied over a 34 Chapter 2 Background and Literature Review short time period, making them a suitable choice for the simulation of composite tube crushing. However, excessive element distortion and numerical instabilities in explicit analyses may cause results to be erroneous. Explicit codes used for the analysis of composite tube crushing include LS-DYNA (Kerth and Maier, 1994; Kerth et al, 1996; Agaram et al, 1997; Castejon et al, 1998; Botkin et al, 1998; and Xiao et al, 2003) ABAQUS Explicit (Castejon et al, 1998; Castejon et al, 2001; Tay et al, 1998), Radioss (Caliskan, 2002), and PAM-CRASH (Haug et al., 1991). Tay et al. (1998) explicitly analysed the same splaying mode of failure in a 20-ply carbon/PEEK composite tube in ABAQUS using both three-dimensional and axisymmetric elements. The 3D analysis significantly under-predicted the forces and had problems with severe element distortions. The axisymmetric analysis predicted a failure behaviour that was not observed experimentally and was also prone to severe element distortions. Tay et al. concluded that although the explicit analysis was computationally more efficient than the implicit analysis, the results may be unstable or unreliable unless element distortion is addressed. To summarize, explicit finite element codes seem to be the most suitable analysis method for the numerical simulation of progressive axial crushing of composite tubes because of their superior numerical efficiency. Excessive element distortion and numerical instabilities in explicit analyses can cause erroneous results and therefore must be understood and controlled. 2.3.2 Material Models The material models used to simulate damage progression in composite tube crushing events can be separated into two main categories: micro-mechanics models and phenomenological models. Micro-mechanics models consider the deformation and failure of the individual constituents in the composite. The stress-strain properties of the constituents and their interactions are modified as damage progresses and then used to determine the overall stress-strain properties of a repeating unit cell (RUC). The RUC properties are then used in the global model to 35 Chapter 2 Background and Literature Review define the behaviour of the elements that discretize the tube. Although these types of models are capable of capturing the damage properties of a composite quite accurately, they are quite complicated and lead to simulations with high computational energy requirements. Beard and Chang (2002) reported results of numerical braided tube simulations using a micro-mechanics based model. The model calculates the overall effective material properties of a repeating unit cell (RUC) using the material properties of the individual constituents and the three-dimensional geometry of the laminate. To simulate damage growth in the braided composite, the material properties of the individual constituents, namely the tows and the matrix, are degraded once stress levels reach a given maximum stress criterion. This model also captures braider tow scissoring and tow jamming. The material model was incorporated into the implicit finite element code ABAQUS, and was successfully used to predict the quasi-static crushing response of carbon fibre/epoxy-vinyl ester braided circular and square tubes. The model is limited by the fact that it is designed to simulate braided material of a particular architecture; each new architecture demands a new RUC. The elements used to discretize the tubes are equal in size to the RUC size. The authors did not comment on the applicability of the material model to elements of various sizes. Most of the current models available in commercial finite element programs are not micro-mechanical models. Rather, they are phenomenological models, which represent damage evolution in either the lamina or laminate by modifying the stress-strain relationships at the elemental level. The material properties are modified as damage increases in the element. The material behaviour is controlled by various parameters that are preferably experimentally determined, but most often are assumed and modified until the desired outcome is obtained. These assumed parameters limit the predictive capabilities of these types of models. Farley and Jones (1992b) produced one of the first phenomenological material models used to simulate composite tube crushing. The material model was used to describe the behaviour of shell elements that represented the layers in carbon reinforced and Kevlar reinforced circular tubes. The material was assumed to be nonlinear elastic with a piece-wise linear stress-strain relationship. Lamina failure, which was represented by a reduction in the elastic modulus, was determined by the maximum allowable strains in the principal material 36 Chapter 2 Background and Literature Review directions. The material model was used to calculate the maximum and minimum crushing forces in a single crushing cycle. Kerth and Maier (1994) presented the results of a numerical study on composite tube crushing using a composite damage material model implemented in the dynamic finite element code LS-DYNA. The material model was based on Classical Laminate Theory and continuum mechanics and models each lamina separately as a homogeneous and orthotropic medium. Four stress based failure modes, tensile fibre mode, compressive fibre mode, tensile matrix mode, and compressive matrix mode are evaluated to determine when damage initiates. Post-failure behaviour was characterized by the degradation of the material stiffness and strength properties. The model predicted the.failure behaviour, internal energy, and the force-displacement profile of the crushing of a filament wound graphite fibre reinforced tube quite accurately (Figure 2.21 and Figure 2.22). One of the shortcomings of the model was recognized as its inability to explicitly take into account delamination, resulting in slightly higher predicted peak forces and mean crush loads. Figure 2.21 Predicted failure behaviour offilament wound graphite fibre reinforced tube (Kerth and Maier, 1994). 37 Chapter 2 Background and Literature Review Figure 2.22 Comparison of numerical (a) force-displacement profile and (b) internal energy profile to experimentally measured trends (Kerth and Maier, 1994). Caliskan (2002) presented a simple model, in which the laminate was modelled by a single shell element. The thickness, orientation and position of plies defining the laminate were input into the model. Individual plies were deleted once an ultimate strength value, as determined by the Tsai-Wu failure criterion, was reached. This ply-by-ply fracturing leads to an overall laminate shell element response as shown in Figure 2.23 for a 7-ply laminate. Only the orthotropic material properties, tensile and compressive strengths, and the failure strains were required for this model, eliminating the need to define damage parameters. Although this results in a general-purpose model with easily determined parameters, the model does not accurately represent the crushing process because it does not accurately capture delamination or the post-failure response of individual plies. Additionally, the mesh sensitivity resulting from the strain softening elemental response is not accounted for. a Lamina 1 Lamina 2 Lamina 3 Lamina 4 Lamina 5 Lamina 6 Lamina 7 Figure 2.23 Representation of damage progression in a laminate shell element by successive failure of individual plies defining the laminate (Caliskan, 2002). 38 Chapter 2 Background and Literature Review Xiao et al. (2003) presented the results of a study done on braided composite tube crushing simulations using the phenomenological material model, MAT58, in LS-DYNA. The model represents the orthotropic material properties of composite laminae. The degradation of the elastic material moduli in each of the principal directions is dependent on the growth of damage parameters in these directions. The damage parameters are dependent on the elemental strain and defined for tensile and compressive failure in the longitudinal and transverse directions, as well as for shear. The elemental stress-strain profile is described using two independent parameters, £, and SLIMj, where i = IT, IC, 2T, 2C, S (TMension, C=compression, 5=shear). Et is defined as the strain at peak stress and SLIMj, is defined as the ratio of the minimum stress to the maximum stress. Figure 2.24 represents some stress-strain curves that result from different Et and SLIMj parameters. In the current study, the values of SLIMj were based on the plateau stress observed in many of the experimentally tested panel specimens, and the values of Et were chosen such that the peak stresses matched the experimentally observed values in the panel specimens. The material model captured the SEA and peak forces in braided tubes crushed with a plug initiator, but tended to under-predict the SEA of non-initiated tube crushes. The authors believed that the parameters could be experimentally determined, resulting in a physics based phenomenological model. One aspect of the material model that was not addressed was the issue of mesh sensitivity. 0 04 (1 06 Strain Figure 2.24 Examples of stress-strain profiles that result from various Et and SLIMj values (Xiao et al, 2003). 39 Chapter 2 Background and Literature Review Morthorst and Horst (2004) presented a material model incorporated into LS-DYNA to simulate failure of glass/epoxy conical shells under quasi-static crushing conditions. Their model focused on capturing the in-plane response of lamina plies as well as accurately modelling delamination between layers. A mixed model that uses shell elements to represent the plies and solid elements to represent the matrix represented the complete laminate as shown in Figure 2.25. Failure was either stress (ply failure) or strain (resin failure) based. Initiation of failure indicated the first stages of damage. As damage increased, as indicated by the degree that the failure criterion was exceeded, the elastic material modulus of the shells and solid elements were degraded. Although the model accurately predicted the crushing characteristics of glass/epoxy conical shells, the authors acknowledged that the model was fairly computational intensive and required a lot of parameters, most of which were determined using curve fitting methods. Additionally, although different mesh sizes were used in the simulations, no investigation into the effect of these changes was presented. Shell elements Solid elements (laminate plies) (resin pbase) Mixed model Figure 2.25 Mixed model used by Morthorst and Horst (2004) to capture the in-plane and delamination properties of a laminate. The majority of the phenomenological material models (Kerth and Maier, 1994; Xiao et al, 2003; Mothorst and Horst, 2004) use a strain-softening approach to represent the growth of damage in the element. A strain-softening response is characterized by the continual decrease of the elastic modulus with increasing strains once a certain damage initiation point has been reached (Figure 2.26). Strain-softening material models are susceptible to mesh 40 Chapter 2 Background and Literature Review sensitivity (Bazant and Planas, 1998). Mesh sensitivity is the tendency of the finite element analysis to produce significantly different results for different mesh densities, as shown in Figure 2.27. Mesh sensitivity occurs because the softening in the model is concentrated in one element. As the element size is reduced, the failure becomes localized in smaller volumes, which causes less energy to be dissipated by the softening. This can lead to instabilities or, at least, mesh-sensitive behaviour. The phenomenological material models reviewed above did not account for this mesh sensitivity incorporated in the model, limiting the robustness of the model and leading to debatable results. Figure 2.26 Schematic of general strain softening material model showing (A) damage initiation point and (B,C) unloading along a path defined by a reduced elastic modulus. (a) t f t t t t Zone of decreased A i strength A A A , (W : K mesh A (c) hi mesh B f t ' p. displacement Figure 2.27 Illustration of mesh sensitivity observed when a strain-softening material model is defined, showing (a) rectangular panel with a zone of decreased strength subjected to tensile load, localization offailure into a band of elements in (b) mesh A density and (c) mesh B density, and (d) dependency of the peak loads and post-peak response on the mesh density (Bazant and Planas, 1998). In order to accurately and efficiently model the crushing response of composite tubes, there is a definite need for a material model that is robust enough to represent the effect of damage growth in various types of composites. The model behaviour should be based on physically 41 Chapter 2 Background and Literature Review measurable model parameters and should be capable of simulating the response of systems discretized with various element sizes (i.e. mesh sensitivity should be limited). 2.3.3 Delamination Modeling Delamination growth is a fundamental energy absorbing mechanism in the crushing of composite tubes (Hull, 1991; Farley and Jones, 1992a). Therefore, it is necessary to account for delamination in order to accurately simulate tube crushing events. Most composite structures in the papers reviewed were modeled using shell elements, without the ability to model out-of-plane failures, such as delamination. Many researchers (Agaram et al. (1997), Kerth et al. (1996)) accepted this drawback and recognized it as a limitation to their model. Kerth and Maier (1994) accounted for delamination at the crash front by a global reduction of stiffness and strength values for the elements in the crash front. The authors concluded that the inaccurate predictions were due to this simplification. Authors that did confront the problem of delamination typically used a simple critical force spotweld technique, which ties nodes together using spring elements, rigid rods or other constraints on either side of an expected delamination front. The nodes are released, indicating delamination growth, once the forces exceed a defined criterion. Tay et al. (1998) used a spot weld technique to simulate the splaying mode of failure typical of carbon/PEEK tubes. Node pairs were defined at the interfaces and decoupled once a specified tensile or shearing break force was reached, resulting in delamination growth. Xiao et al. (2003) used a tiebreak contact definition between shell elements. The failure of the tied surfaces was based on interface strengths determined by a series of interlaminar fracture toughness tests and simulations. Although spotweld techniques are simple to incorporate, as only a force criterion need to be defined, they do not accurately capture the fracture phenomenon, which is dependent on more than just the forces transferred across the delamination front. Fleming (2001) compared the effectiveness of a simple critical force spotweld technique at capturing the growth of delamination in a double cantilever beam specimen to two other modeling techniques (cohesive fracture model and virtual crack closure technique). The cohesive fracture model uses spring element defined by a force-displacement response based on classical cohesive failure behaviour, as shown in Figure 2.28. This method assumes 42 Chapter 2 Background and Literature Review delamination growth is dependent on both the maximum force and the critical energy release rates, and allows energy to be absorbed via delamination, resulting in a more physical treatment of this failure mechanism than the spotweld technique. The virtual crack closure technique uses the forces and displacements near the delamination front to calculate the energy release rates for each mode of fracture (Mode I, II, and III), as shown in Equation 2.1 for Mode I, where Gj is the strain energy release rate, Fi is the force in the spring aligned with the mode I direction, AA is the interfacial area attributed to the spring, and u+ and u are the displacements at the nodes on either side of the delamination just ahead of the crack tip. If the combination of these energy release rates exceeds a mixed mode critical energy release rate the delamination grows. All three techniques were capable of capturing the growth of delamination in a double cantilever beam specimen. The parameters used in the spotweld technique were arbitrary, whereas the parameters in the other two methods were based on material properties. Fleming (2001) concluded that although the spotweld technique may not model delamination quite as accurately as the other methods (cohesive failure models and virtual crack closure techniques), at present it is the most efficient way of representing delamination growth. Figure 2.28 Typical stress-strain response defined for springs in cohesive failure delamination model (Fleming, 2001). 2.1 Other researchers used a more computationally intensive method of representing damage growth by defining solid delamination elements between plies. Botkin et al. (1998) simulated the crushing of a braided composite tube with and without the explicit use of solid delamination elements. The model with delamination element produced considerably better results. Morthorst and Horst (2004) defined solid elements between neighbouring shell 43 Chapter 2 Background and Literature Review elements (Figure 2.25) to represent delamination between plies. A strain based failure criterion and degradation model is used to adjust the properties of the solid element up to complete failure. Johnson et al. (2001) presented a method of modeling delamination between shell elements using a sliding interface. The failure of the interface, and therefore the growth of delamination, was governed by fracture mechanics criteria. 2.4 Summary Experimental research in the field of composite tube crushing has shown that the failure behaviour of composite tubes is affected by material parameters (fibre and matrix type, fibre architecture, and fibre volume fraction), specimen geometry, and external loading details. High SEA and stable progressive crushing can be obtained in brittle materials. The material is either sequentially damaged on a micro structural level (fragmentation mode) or deformed into a splaying type mode of failure. Tubes reinforced with ductile fibres, such as Kevlar, tend to progressively fold in a manner similar to that observed in metal tubes. At present, a universally accepted material model for predicting the crushing response of automotive parts has not been attained. Fleming (2000) stated, "It is a testament to the complexity of the crushing phenomenology that despite numerous efforts and years of study, no completely successful predictive model of composite crushing has been achieved". Composite crash modelling has proven to be a difficult area with many challenges. The complexity of the failure modes, the lack of reliable characterization of significant factors like friction and delamination, and numerical issues like mesh sensitivity all make the problem challenging. A successful material model should account for all of these issues, have parameters that are easily determined and physically based, and be general enough to apply to different dynamic crash problems. 44 Chapter 3 Numerical Implementation of a Compressive Damage Model C H A P T E R 3 N U M E R I C A L I M P L E M E N T A T I O N OF A COMPRESSIVE D A M A G E M O D E L 3.1 Introduction This chapter introduces, describes and verifies the material model, CODAM, used in the numerical simulations of the tube crushing experiments. The first section of the chapter presents a brief literature review of composite compressive failure in unidirectional, multidirectional and braided composites. Additionally, an overview of some of the existing analytical and numerical models illustrates some of the challenges faced in accurately modelling composites under severe compressive loads (such as in the tube crushing experiments). Next, the modifications incorporated into the existing material model, CODAM, are presented. Finally, CODAM is verified by comparing the numerical results of a single element under various loading and boundary conditions to analytical solutions. 3.2 Compressive Failure in Composite - Background Understanding how fibre reinforced composite materials react to severe compressive loads is of primary importance when attempting to design composite structures to absorb or dissipate energy through compressive damage propagation. Additionally, the ability to accurately simulate compressive failure in composites aids in using composites affordably and efficiently. Since the development of these materials, there has been a significant amount of research aimed at understanding how they respond to compressive loads. However, this area still remains to be one of the least understood in the field of composites today. This section provides a brief overview of some of the important literature, with emphasis on braided composites. 3.2.1 General Compressive Failure Mechanisms Compressive failure in laminated composites can occur at either the structural level (coupon geometry), macrostructural level (lamina level) or at the microstructural level (fibre-matrix 45 Chapter 3 Numerical Implementation of a Compressive Damage Model level). The compressive response is influenced by many factors, such as the specimen geometry, constituent properties, lamina (or fibre) orientation, fibre waviness, voids, and stress concentrations. The most notable failure mechanisms observed in composites are global Euler buckling, fibre or matrix compression failure, transverse tension failure, delamination, elastic microbuckling, and fibre kinking. Global Euler buckling is a structural level failure mode dependent on coupon geometry. Thin plates subjected to in-plane compressive loads are susceptible to buckling, as shown in Figure 3.1. Figure 3.1 Global Euler buckling of a thin plate under a remote applied stress. Fibre compression failure occurs when the fibres crush due to high axial strains. Typical crushing strains are 0.5% (for Kevlar and pitch-based carbon fibres) and 2.5% (for PAN-based carbon fibres) (Budiansky and Fleck, 1993). Laminates reinforced with fibres having a low crushing strength imbedded in a high yield strength epoxy may display this type of compressive failure. In most modern composites this type of failure is not significant due to the high crushing strength of fibres and the tendency to use tougher matrices (having lower yield strengths). Matrix compression failure is another failure mechanism in which the matrix microcracks or yields due to high compressive or shear stresses that build up in the matrix surrounding the fibres. In general, matrix damage occurs simultaneously with other fibre dominated failure mechanisms. Transverse tensile failure may precede other forms of failure in composites with low transverse tensile strength and high Poisson's ratio effects. Parallel splitting or 46 Chapter 3 Numerical Implementation of a Compressive Damage Model microcracking of the matrix due to the applied compressive load leads to ultimate tensile failure in the transverse direction, as shown in Figure 3.2. Figure 3.2 Transverse tensile failure of a unidirectional lamina. Delamination (interlaminar separation) occurs when the individual plies within a laminate partially or completely separate from each other (Figure 3.3). Delamination occurs due to high interlaminar (out-of-plane) stresses, which cause the bond between adjacent plies to fail. Initial defects due to either poor manufacturing or prior impact damage can significantly increase the risk of delamination. Splitting of the matrix Microcracking of the matrix a Figure 3.3 Delamination of outer plies in an idealized laminate comprised offour plies. Al Chapter 3 Numerical Implementation of a Compressive Damage Model Elastic Microbuckling is a failure mode that was first proposed by Rosen (1965) to describe the compressive failure of composites. In this failure mode, fibres buckle within the supporting flexible matrix in either a shear mode or in an extensional mode (Figure 3.4). For composites with fibre volume fraction's greater than 30% the shear mode dominates. SHEAR MODE EXTENSIONAL MODE Figure 3.4 Two modes of elastic microbuckling proposed by Rosen (1965). Fibre kinking is defined as the cooperative rotation of fibres within a narrow band, in conjunction with non-linear shear deformation of the supporting matrix. Figure 3.5 shows a typical kink band with characteristic fibre rotation angle, kink band width and kink band, rotation. Fibre kinking typically initiates in regions with local misalignments or discontinuities. Kinking can form either in-plane or out-of-plane (Figure 3.6), depending on the relative constraint in each direction. In typical laminates with thicknesses much less than the widths, out-of-plane kinking is the expected mode of failure unless some external boundaries prevent such deformation. There is a general consensus that kinking is the dominant failure mode in composites with fibres aligned or nearly aligned with the loading direction. Other failure mechanisms such as matrix failure, fibre failure and delamination tend to occur simultaneously with, or as a consequence of, kinking. 48 Chapter 3 Numerical Implementation of a Compressive Damage Model ( A <f> = fibre rotation 00 A V Fibres w = kink band width B = kink band inclination Thickness (out-of-plane) or width (in-plane) Figure 3.5 Details of kink banding with fibres aligned with remotely applied load. Thickness OUT-OF-PLANE KINKING IN-PLANE KINKING Figure 3.6 Graphical representation of differences between out-of-plane and in-plane kinking. 49 Chapter 3 Numerical Implementation of a Compressive Damage Model 3.2.2 Experimental Investigation of Failure Mechanisms 3.2.2.1 Unidirectional Composites Hahn and Williams (1986) tested unidirectional composites composed of graphite fibres in epoxy resin having various moduli in compression to study the observed failure modes. The experimental results indicated that when the fibres are strong enough to resist fibre compression failure and global Euler buckling is prevented, failure in unidirectional composites initiates as microbuckling. Microbuckling was observed in areas with low amounts of lateral support (near free boundaries, longitudinal splits, voids) or in areas with stress concentrations. In systems with soft resin, microbuckling proceeded until fibres ultimately failed in bending. In systems with an adequately stiff resin, microbuckling ultimately led to kinking as the final failure mode. Hahn and Williams concluded that microbuckling is the precursor of kinking in most composite systems. Increasing the resin modulus did not translate into a significant increase in the composite modulus, but did significantly increase the composite compressive strength. This strong influence on the strength was not observed in tension. A more recent study conducted by Sivashanker et al. (1996) focused on studying the initiation and progression of kink bands in unidirectional carbon fibre - epoxy matrix composites. The specimens tested were 145 x 50 x 3 mm (length, width, thickness) and had a 15 mm notch machined into one edge (Figure 3.7). To insure that damage would initiate at the notch, the root was dented before any compressive load was applied. Euler buckling was prevented by using anti-buckling guides in all the tests. In all of the specimens tested, out-of-plane kinking initiated at the notch tip and traversed the width of the specimen. The kink band was inclined at an angle of /? » 25° to the fibre direction. Examinations of the failed specimens indicated that the projected height of damage varied across the width of the specimen (from about 800 um at the notch tip to about 70 um at the final damage tip). A closer inspection revealed that the band of damage was composed of a series of individual kink bands, all with a relatively constant width of <y»100 um (20 fibre diameters). Significant fibre breakage was observed at the boundaries between successive kink bands. Strain gauges were placed above the predicted kink band path in elastic material that was not 50 Chapter 3 Numerical Implementation of a Compressive Damage Model damaged during the test (shown as location A and E). An estimate of the stresses bridging the damage zone at these locations was obtained by multiplying the measured gauge strains by the material elastic modulus. Although the strain gage readings may initially be affected by the large gradient in the applied stress, once the crack front passes the location, an estimate of the plateau stress can be obtained. Figure 3.8 illustrates how the stress (indicated by strain gauge A) increased to about 215 MPa just prior to the formation of the first kink band below the strain gauge. The stress then drops off to a relatively constant value of about 100 MPa. At this point the number of individual kinks bands (and thus the damage zone) continues to increase even though the stress at this location remains constant. This observation suggests that after the initial kink band has formed, less force (or energy) is required to propagate the damage zone lengthwise. The formation of additional kink bands, with fibre breakage at the boundaries (or in some cases an increasing kink band width without fibre breakage) has also been experimentally observed by others (Moran et al, 1995; Moran and Shih, 1998; Fleck, 1997; Vogler and Kyriakides, 1997; and Sivashanker, 1998) 51 Chapter 3 Numerical Implementation of a Compressive Damage Model ^ Applied displacement 50 mm •4 • 11 T t Figure 3.7 Specimen geometry, showing the location of the strain gauges and the distribution of damage along the width of the specimen after compressive loading (Sivashanker et al, 1996). 52 Chapter 3 Numerical Implementation of a Compressive Damage Model Growth of damage tip across specimen width (mm) Figure 3.8 Progression of stress (as estimated by measured product of gauge strain and elastic modulus) at location A as the damage tip approaches and passes. 3.2.2.2 Multidirectional Composites Sivashanker (2001) followed up his study on unidirectional laminates by a similar study on carbon fibre/epoxy multidirectional composites. In this study, three different lay-ups, [±45/04]2S, [±45/02]3s and [±45/0/90]3S, having specimen geometries similar to that shown in Figure 3.7, were loaded in compression. The progression of damage in all the laminates was studied using scanning electron microscope (SEM) and optical microscope. Strain gauges were again placed at 4 mm ahead of the notch tip in material that remained undamaged throughout the experiment to estimate the stress (stress = elastic modulus x strain gauge measurements) transferred across the damage zone as the damage tip progressed across the width of the specimen. All of the surface plies (45°) displayed similar damage patterns, characterized by regions of out-of-plane kink bands connected to one another by parallel diagonal splits between the fibres. The internal layers of all the laminates also displayed similar damage patterns characterized by fibre kinking in the 0° layers and off-axis damage in the form of matrix cracking and some kink banding oriented normal to the fibre direction. These damage mechanisms in each layer were connected to one another via inter-ply delaminations. Figure 3.9 shows the details of these failure mechanisms in the [±45/02]3S laminate. 53 Chapter 3 Numerical Implementation of a Compressive Damage Model Sivashanker observed that the overall damage height increased longitudinally in all specimens as the compressive load increased. This growth was represented by an increase in the width of kink bands, a growth of delaminations and an increase in the amount of off-axis damage. In all specimens tested, individual damage mechanisms were never seen in isolation, but rather damage seemed to initiate and grow in terms of a damage unit (comprised of kink bands, delaminations, and matrix damage). The measurements provided by the strain gauges were once again used to calculate an estimate of the stresses transferred across the damage zone as the damage initiates and progresses at location A. These calculated stresses again reached maximum values as the damage tip approached the location of the strain gauge and dropped off to a plateau level once the damage unit was well established and was in the process of growing longitudinally. The measured plateau stresses ranged from 120 to 200 MPa. The general trend of higher plateau stresses in the multidirectional laminates versus the unidirectional laminates was attributed to the progression of two mechanisms (kinking and delamination) versus one (kinking) in the unidirectional specimens. The energy (and thus the plateau stress) required to propagate the damage unit in the multidirectional case seems to be higher than in the unidirectional case. 54 Chapter 3 Numerical Implementation of a Compressive Damage Model Figure 3.9 Summary of the various damage mechanisms observed in multidirectional laminate (a) general through thickness out-of-plane failure mechanisms showing kinking, delaminations and off-axis matrix damage (b) region CI, showing out-of-plane kinking (c) close up of kink band (d) region C2, showing delamination and kink band formation in 45° layer (e) region C3, showing kinking in 0° layer extending into the 45° layer (Sivashanker, 2001). 55 Chapter 3 Numerical Implementation of a Compressive Damage Model 3.2.2.3 Braided Composites The failure mechanisms observed in braided composites under compression has also been extensively studied. The general trend of kink bands forming in the tows aligned (or nearly aligned) with the loading direction, leading to out-of-plane buckling of the tows combined with between-tow matrix damage and tow-matrix delamination (Falzon and Herszberg (1998), Quek et al. (2004a), Salvi et al. (2004)) is similar to the damage mechanisms observed in multidirectional laminates. One distinguishing feature between tape laminates (unidirectional or multidirectional) and braided composites is with respect to the waviness or initial misalignment of the fibres. The braided process inherently leads to high undulations or crimp in all of the tows. Due to the significant effect of initial fibre misalignment on the formation of kink bands, braided composites tend to have lower compressive strengths as compared to tape laminates. Quek et al. (2004a) studied the damage mechanisms in carbon 2D triaxial braided composites under both uniaxial and biaxial applied loads (Figure 3.10). The specimens tested in pure compression exhibited fibre kinking combined with matrix damage between the tows, which led to the ultimate out-of-plane buckling of the tows (Figure 3.11). These damage mechanisms localized into a band of about 8-12 mm oriented perpendicular to the direction of the applied load. In biaxially tested specimens, tow/matrix interfacial failure due to the high tensile failure precedes the kinking/buckling of the tows. 56 Chapter 3 Numerical Implementation of a Compressive Damage Model (a) 139.70 tarn 50,80 mm . H H 1 R = 12.70 ' 50.80 mm ^ Dimensions Applied Loads via Displacement Control I: * 4 • * • 4 er_ fb) - «• » 0- - • y »• Figure 3.10 Experimental study done by Quek et al. showing (a) details of the specimen dimension and (b) loading directions (Quek 2004a). Cross Section Front View Figure 3.11 Typical failure mechanism observed in braided composites due to compression showing global out-of-plane tow buckling and inter-tow matrix damage (Quek, 2004a). West and Adams (1999) studied the effect various degrees of crimp have on the compressive modulus, and strength of 2D triaxially braided composites under static compressive loading. To produce specimens with different degrees of axial crimp, panels were manufactured with tensioned and non-tensioned axial yarns (Figure 3.12). The stiffness of the panels was not affected by the degree of crimp, having an average value of 55 MPa. However, the axial strength of the crimped panels was 381 MPa, 30% lower than the average strength of 544 57 Chapter 3 Numerical Implementation of a Compressive Damage Model M P a for the low crimped panels. Minguet (1995) also observed up to a 40% reduction in the un-notched compressive strength of carbon/epoxy braided composites, as compared to tape laminates with similar elastic modulus and percentage of axial fibres. Figure 3.12 Degree of crimp observed in non-tensioned (a) and tensioned (b) braided panels (West and Adams, 1999). 3.2.3 Predicting the Compressive Response of Fibre-Reinforced Composites Researchers have proposed various theoretical and numerical models to predict the compressive response of composites. Rosen (1965) developed one of the first theoretical equations to estimate the compressive strength of unidirectional composites. His theory was based on the failure of initially straight fibres in a flexible matrix in either the shear or extensional mode (see Figure 3.4). Based on energy principles, Rosen developed the following equation for the ultimate compressive strength of unidirectional composites failing in the shear mode (limiting case for most composites), where <JC is the ultimate compressive strength, Gm is the shear modulus of the matrix material, and Vf is the volumetric fibre fraction. Rosen's formula typically over-predicted compressive strengths, prompting researchers to improve the microbuckling model by including such things as initial fibre waviness and misalignment (Davis, 1975; Hahn and (a) (b) 58 Chapter 3 Numerical Implementation of a Compressive Damage Model Williams, 1986), or material nonlinearity (Hahn and Williams, 1986; Schapery, 1995). Guynn et al., 1992 modeled the initial stages of kinking in unidirectional composites using a finite element analysis that incorporated non-linear material properties. Analysis was completed on initially straight and wavy fibres. These efforts resulted in predictions that were lower than Rosen's but still did not consistently predict the compressive strength of all composites. Other formulas were developed to predict the compressive strength of unidirectional composites assuming kinking as the dominant failure mechanism. Argon (1972) was perhaps the first to analyze the phenomenon of kinking. He suggested that a failure nucleus, formed in areas with local misalignment, led to the formation of a kink band at a stress less than that predicted by Rosen's elastic microbuckling formula. His analysis was based on the continuous shearing of the matrix between initially misaligned fibres, causing the fibres to further rotate. Argon developed the following formula for the ultimate compressive strength, where k is the interlaminar shear strength, and (/> is the initial misalignment of the fibres. Budiansky (1982) formulated an equation for the kinking stress, which accounted for an initial elastic portion in the composite shear stress-strain response. His formula, as shown below, reduces to Rosen's equation for an initial fibre misalignment of zero (if yy = r/G and G is assumed to equal the ratio of Gm over (l-Vj)). Additionally, if the yield strain is set to zero and the composite yield stress in shear (r ) is equal to the interlaminar shear strength (k), it reduces to Argon's formulation. 3.3 Although these analytical theories provide an estimate of the ultimate strength of a composite, they do not provide any insight into the response of composites under progressive deformation and failure. Additionally, they do not accurately represent the behaviour of 59 Chapter 3 Numerical Implementation of a Compressive Damage Model more complicated composites, such as multidirectional, woven or braided laminates. Therefore, they tend not to be very robust methods of analysis. Numerical models have the potential of predicting the pre- and post- failure response of a variety of composites under different loading situations up to complete failure. Most numerical models fall into one of two categories: micromechanical models and macromechanical models. Micromechanical models attempt to capture the overall stress-strain response of a repeating unit cell by modelling the failure mechanisms and material properties at the constituent level. This stress-strain response can then be used in a structural stress analysis code to describe the response of individual elements. Because micromechanical models incorporate the individual constituent properties and the interaction of these constituents, they are potentially capable of reproducing the different failure mechanisms observed in compression. They do, however, require a significant amount of computational energy and are fairly complicated, limiting their range of applicability. Fleck et al. (1996) discretized the fibres and the matrix in a finite element analysis to model kinking in notched unidirectional composites. The analysis assumed linear elastic material properties for the fibres and elastic-plastic material properties for the matrix. The initial growth of the damage was successfully captured, but damage growth and response past initiation was difficult to capture. Effendi et al. (1995) numerically analyzed compressive failure in unidirectional composites using a finite element model that included the effects of initial fibre waviness, matrix plasticity and large displacements. The model predicted that compressive failure shifted between fibre compressive failure (for low degrees of initial waviness) and fibre kinking (for higher values of initial misalignment). Although the model was successful at capturing the ultimate compressive load, the stress-strain response to failure was significantly stiffer and did not show the softening response observed experimentally. Micromechanical models have also been developed to simulate the compressive response of more complicated composites, such as woven and braided composites (Naik, 1995; Beard and Chang, 2002; Quek et al., 2004b). Naik developed a 3-dimensional model implemented in TEXCAD (Textile Composite Analysis for Design), capable of predicting the mechanical 60 Chapter 3 Numerical Implementation of a Compressive Damage Model properties, damage initiation/progression and strength of woven and braided composite materials. The model accurately represented the textile architecture, including the geometric nonlinearities associated with the straightening of the undulations in the braider tows. Additionally, the nonlinear shear response of both the impregnated yarns and the resin was included in the model. Constituent failure was based on the maximum stress and strain criteria. Naik believed the model could potentially improve the predictive capabilities of current structural analysis codes by incorporating it as a material model. As the model is fairly complicated, this would require significant computational effort. Quek et al. (2004b) developed an in-depth micromechanics based finite element model to simulate the compressive behaviour of braided composites, which accounted for the fibre/tow architecture, elastic properties of the braid, fibre volume fraction and non-linear stress-strain response of the matrix. A repeating unit cell (RUC) of the material was analyzed in ABAQUS to determine the compressive strength and the dependence of strength on tow misalignment. Figure 3.13 shows the results of the uniaxial compressive simulation for different imperfection magnitudes (representing the misalignment of the tows). The initial stiffness of the RUC agreed well with the experimental data. The peak load was dependent on the amount of imperfections and was generally over-predicted. As loading continued, the numerical response became progressively nonlinear, due to the geometric and matrix nonlinearities. A post-peak softening response that approached a plateau stress was also captured in the simulation. 61 Chapter 3 Numerical Implementation of a Compressive Damage Model — F E Mode! tm = 0.05 —<•— FE Model im = 0.10 —.- FE Model im - 0.20 Experimental 4-0.00 0.01 0.02 0 03 0,04 0.05 0 06 0.07 0,08 0,09 strain, Figure 3.13 Uniaxial numerical compressive stress-strain response compared to the experimental data for varying degrees of imperfection. Initial modulus is well predicted, but the peak stress is over-predicted (Quek et al. 2004b). Macromechanical models, on the other hand, attempt to capture the global system response by defining an elemental stress-strain response, which theoretically embodies the overall effect of all the local failure mechanisms. As these types of models are very computationally efficient, they tend to be the most commonly used methods to analyze large-scale problems (e.g. material models 22, 54-55, 58, 59, 161 in the commercial Finite Element code, LS-DYNA). Most of these models incorporate orthotropic material properties and define element failure based on some maximum stress or strain criterion. The problem with many of these models is that they tend not to be physically based. Instead, curve fitting or tweaking methods are often used to obtain acceptable results for a given simulation. This creates non-robust models, which cannot be used with confidence for predictive purposes. In order for these types of material models to become trusted as efficient and effective modelling methods, they must be physically based. A physically based macromechanical material model must be based on physically measurable model parameters. To the best of the author's knowledge, a purely predictive physically based macro-mechanical model capable of simulating composites under a variety of loading conditions does not exist. 62 Chapter 3 Numerical Implementation of a Compressive Damage Model 3.3 Compressive damage modelling using CODAM 3.3.1 Overview of the CODAM material model The CODAM (COmposite DAMage) model is a physically-based macromechanical material model that attempts to represent the constitutive behaviour of polymer composite laminates through both the elastic (pre-failure) and post-initial failure regimes. CODAM was originally developed at the University of British Columbia by Williams (Williams, 1998; Williams et al., 2003). It was further extended by Floyd (Floyd, 2004) to include three dimensional brick elements and a scaling law based on energy principles that ultimately results in objective predictions that exhibit mesh-independence. CODAM has been incorporated into the finite element code, LS-DYNA, and has accurately modelled tensile failure in composite systems under many loading conditions (Floyd, 2004). The ability of CODAM to capture the response of laminated composites under compressive loads has not yet been thoroughly studied. The most important aspects of the CODAM material model are summarized below. For a complete description of the model formulation see Williams et al. (2003) and Floyd (2004). • The development of the model is based on a Continuum Damage Mechanics formulation. • The model represents damage growth in a Representative Volume Element (RVE). The dimensions of the RVE have significant meaning. The height and width are equal to the characteristic height of damage in each of these directions. The characteristic height of a quasi-brittle material is defined as the actual height the damage would grow to if experimentally tested under conditions that lead to stable damage growth. The thickness corresponds to the thickness of the sub-laminate or laminate. Defining the RVE thickness in such a way allows layer interactions to be captured as well as making the model more computationally efficient. • Damage growth in the RVE is represented by an overall decrease in the normalized secant modulus (E - E/En, where E0 is the initial undamaged modulus). The modulus reduces with increasing amounts of damage, as indicated by damage 63 Chapter 3 Numerical Implementation of a Compressive Damage Model variables (co) as shown generally in Figure 3.14. A normalized modulus reduction curve is defined for each principal local direction (x, y, and z) for both tension and compression. co Figure 3.14 Generalized normalized modulus curve. • The damage variables vary with the RVE strain state, as defined by an effective strain function (F). The effective strain is a function of the various strain components as shown below, which allows strain interactions to be incorporated into the model. A damage growth curve (co vs. F), as shown in Figure 3.15, is defined for each principal direction (x, y, and z) for both tension and compression. F=. f _ \ K, 2 /„ \ 2 M, i = x, y, and z M, 2 f \ l 3.4 64 Chapter 3 Numerical Implementation of a Compressive Damage Model • The overall one-dimensional RVE stress-strain curve (assuming F is equal to the strain in that direction) that results from the damage growth and normalized modulus curves is shown in Figure 3.16. 1.2 0.05 0.06 Average strain Figure 3.16 General one-dimensional stress-strain response, assuming previous damage growth and normalized modulus curves and that F equals the strain in the direction of interest. • Two scaling laws, crack band method (previously referred to as equivalent work method by Floyd, 2004) and the modified crack band method (previously referred to as crack band method by Floyd, 2004) are incorporated into the material model to 65 Chapter 3 Numerical Implementation of a Compressive Damage Model account for element sizes that are different from the RVE size. The scaling laws are designed to overcome mesh sensitivity problems, often encountered when trying to model strain-softening materials (Floyd, 2004). These methods are based on the idea that the fracture energy, Gf, of a given element with height, he, should be equivalent to the fracture energy of the characteristic RVE element with defined characteristic height, hc. To maintain the fracture energy in an element that is not equal to the characteristic size, the stress-strain curve is typically modified by a scaling factor. In general, the requirement for constant fracture energy can be written in equation form as shown below. (Gf)c=K(rf)c=K(rf)e 3.5 where Gf is the fracture energy, yf is the fracture energy density, hc is the characteristic height, he is the element height and subscripts c and e refer to characteristic and element, respectively. • Recent investigations by McClennan (2004) indicate that in addition to the element height, the element width is an important factor in strain-softening material models when modeling crack growth near a notch tip. Elements must be small enough to capture any notch tip stress concentrations. In cases when element width is larger than this (for improved efficiency), the average stress of the element (underintegrated shell elements) is much lower than the actual stress at the notch tip. McClennan observed that scaling the stress of elements with the width can help to capture the proper growth of damage from a notch tip. 3.3.2 Overview of the Analog Model Preliminary tests indicated that the previous version of CODAM, as described above, did not adequately capture the complete compressive response of laminated composites. This is reasonable, as the model was originally developed to study tensile dominated failures, which are characterized by failure mechanisms that are different than those observed in compression. In tension, failure is generally characterized by matrix cracking and fibre 66 Chapter 3 Numerical Implementation of a Compressive Damage Model breakage, in combination with delaminations (Floyd, 2004). In compression, experimental evidence suggests that the main failure mechanisms in a laminate consist of fibre kink banding in laminae with fibres parallel to the loading direction. Within these kink bands, matrix is cracking (or plastically deforming) and fibres may or may not be breaking (Sivashanker, 2001; Moran and Shih, 1998). In addition to the matrix and fibre damage in the kink band, there are delaminations between dissimilar laminae and matrix cracking in off-axis laminae (reinforced with fibres perpendicular (or near perpendicular) to the loading direction). As the current CODAM does allow for individual representation of damage growth in each of the constituents, capturing the response due to the degradation of the constituents is possible with the current CODAM model. The problem arises when one attempts to capture the constant plateau (or bridging stress) observed by Sivashanker et al., (1996), Sivashanker, (2001), Moran and Shih, (1998), and Vogler and Kyriakides, (1997). This plateau stress arises from the continued ability of the bulk damaged material to support compressive load even though the material is completed damaged, in comparison to the lack of tensile load carrying capacity of similarly damaged material. This occurs because although the matrix and fibres may be significantly damaged in the original kink band, as long as the bulk material is still present and the material above and below the kink band have not completely separated and slid past each other, they can still carry a significant amount of load. An analogy can be made between this scenario and a bucket of sand. As long as the bucket (or confinement) is provided, the sand can support high compressive loads. If the bucket were removed, the sand would easily collapse under the compressive load. This is completely opposite to the tensile response, where once the fibre or matrix becomes damaged, all load carrying capacity is lost. Another important difference between tensile and compressive failures is that the damage height in tension remains relatively constant throughout the loading history, whereas in compression, the damage height tends to increase as the damage zone propagates longitudinally (into previously undamaged areas) under increasing compressive deformation (Sivashanker et al; 1996; Moran and Shih, 1998; and Vogler and Kyriakides, 1997). The fact that the damage height is growing is important when considering the energy absorbed and the scaling laws to be incorporated into the constituent model. 67 Chapter 3 Numerical Implementation of a Compressive Damage Model To aid in the incorporation of compressive damage modelling capabilities into CODAM, an analog model was developed jointly with Navid Zobiery (2004) to describe the response of fibre reinforced polymer composites at the RVE level in the force-displacement domain. A brief outline presenting the main aspects of the model follows. For a more detailed explanation of the model refer to Zobiery (2004). The model was created to enhance our understanding of the mechanisms of failure in composites under tension, compression and load reversals. The analog model consists of two main boxes, the rubble box and the laminate box, as shown in Figure 3.17. In addition to these two boxes, there is a slider element and a lock-spring element. The laminate box represents the one-dimensional softening response of the laminate due to increasing damage growth in both the fibre and matrix, and actually equates to the constitutive model (in the force-displacement domain) currently implemented in CODAM. The rubble box represents the ability of the bulk damaged material to support compressive loads, even when severely damaged. The rubble box is therefore only active under compressive displacements, whereas the laminate box is active regardless of the loading direction. The slider element represents the broadening of the damage unit longitudinally into undamaged material under a constant load, Fy. The lock-spring element leads to a permanent compressive deformation upon unloading. Ultimate failure of the analog model is represented by complete failure of the laminate box. 68 Chapter 3 Numerical Implementation of a Compressive Damage Model Rubble Box Represents the response of the bulk material once completely damaged Slider-Represents the response of the laminate during band broadening / / / / / / / / / / Laminate Box Represents the response of the remaining damaged laminate Lock-spring element Represents the permanent compressive deformation upon unloading Figure 3.17 Schematic of overall analog model showing the two main boxes, the slider and the lock-spring element. The laminate box is composed of a series of infinite internal repeating units (RU) in parallel, as shown in Figure 3.18 (a). The RU's are in turn composed of two types of elements, T/C spring and T/C fuse. A T/C spring is simply a spring with different tensile and compressive stiffness values. A T/C fuse similarly has different tensile and compressive strengths. The fuse strengths incrementally increase with each successive RU. As displacement is applied in either tension or compression, fuses begin to fail (representing the increasing damage in either the matrix or fibre). Failure of all of the fuses represents complete laminate failure. Figure 3.18 (c) shows the overall response of the laminate box to either tensile or compressive displacements. Note that this response is completely captured (in the one dimensional stress-strain domain) by the current CODAM material model. The laminate box responds to load reversals by unloading or reloading along linear paths with a stiffness value defined by the remaining intact springs. If tensile displacement were applied (sufficient enough to cause damage), upon removal of the displacement, unloading would follow a linear path. This linear path would have a stiffness value less than the original "undamaged" stiffness of the laminate box. If the analog model were then subjected to compressive displacements, the force would increase linearly (defined by a reduced stiffness) until it reaches a displacement value that causes the next fuse to fail (i.e. intersects the original undamaged force-displacement profile). 69 Chapter 3 Numerical Implementation of a Compressive Damage Model Repeating unit-f fails + 1 Km Kf \FT\]FC / • • / r i / v t m LJ m J LJ ]> T/C Spring T/C Fuse s ; ; v ; ; ; ; ; / ,y , ; 1 2 n (a) A Single repeating unit Forcef rf m fails f fails (b) Infinite number of repeating units Force Softening due to fibre damage Softening due to matrix damage Total response Figure 3.18 (a) Schematic of Laminate Box, showing the arrangement of the T/C springs and fuses in a parallel infinite series. Repeating units (RU) have the same stiffness values, but each successive fuse has a higher failure load (b) Force-displacement response of a single RU. (c) Force-displacement response of laminate box consisting of an infinite number of RU. 70 Chapter 3 Numerical Implementation of a Compressive Damage Model The rubble box is composed of gap and spring elements, which are arranged in an infinite parallel series (Figure 3.19). The gaps do not transfer any tensile forces and only transfer compressive forces once the gap has closed. Each successive RU has an incrementally larger defined gap length. All of the springs have equivalent stiffness values. Basically, as compressive displacements increase the gaps begin to close sequentially, thus increasing the overall stiffness (and thus the load carrying capacity) of the rubble box. The closure of the first gap occurs concurrently with the failure of the first compressive matrix fuse in the laminate box. At some point, all of the gaps have closed (which occurs concurrently with the failure of all compressive matrix fuses in the laminate box) giving the rubble box a stiffness value equal to n • ke. A release of the compressive displacement results in the successive unloading of the springs and opening of gaps. Unloading in the rubble box occurs along a path identical to but offset from the loading path. 71 Chapter 3 Numerical Implementation of a Compressive Damage Model 1 2 n (a) A Single repeating unit Force Infinite number of repeating units Force Figure 3.19 (a) Schematic of Rubble Box, showing the arrangement of gap and spring elements (parallel infinite series) and the slider element below. Repeating units (RU) have the same stiffness values, but each successive RU has a higher gap size (b) Force-displacement response of a single repeating unit showing inactivation before gap closure, (c) Force-displacement response of rubble box. 72 Chapter 3 Numerical Implementation of a Compressive Damage Model The overall force-displacement response predicted by the analog model is shown in Figure 3.20. Additionally, dashed lines represent the unloading (or reloading with established damage) loading paths. Figure 3.21 and Figure 3.22 demonstrate the failure mechanisms that the analog model represents in both tension and compression for an idealized [90/0/90 [. RVE. Figure 3.20 Response of individual boxes and analog model in force-displacement domain showing elastic and damaging regions. Possible unloading paths (or reloading paths with prior damage) are also shown as dashed lines. 73 Chapter 3 Numerical Implementation of a Compressive Damage Model R V E 1 I it?*" | 5 2 Force Elastic Matrix cracks, delamination Complete and some fibre cracks failure (a) (b) Figure 3.21 Idealized tensile damage growth in a representative volume of a [90/O/90] laminate and force-displacement response predicted by the analog model. I Some kink bands, Damage in all plies Damage propagating Elastic delamination and off-axis (kink bands, lengthwise. (Kink band matrix damage delaminations, off-axis width increasing, matrix damage) delaminations growing). Force Figure 3.22 Idealized compressive damage growth in a [90/0/90js RVE and resulting force-displacement response predicted by the analog model. 7 4 Chapter 3 Numerical Implementation of a Compressive Damage Model 3.3.3 Incorporation of the Analog Model into CODAM To incorporate the features of the analog model into the stress-strain constitutive model, CODAM, requires a few minor changes, which are summarized below. As mentioned previously, CODAM is already capable of capturing the complete response of the laminate box in the stress-strain domain, so no changes were made regarding this aspect. • The modified Bazant scaling laws were further expanded to enable different scaling factors in compression. Previously, one scaling factor (k) was defined for each principal direction (x, y, and z) and applied to both tension and compression. As k is only a function of n (ratio of the characteristic height and the element height) and the element height is fixed for any given finite element model, this in effect forces the characteristic damage height to be identical in tension and compression, which may not be the case. The scaling factors were increased to 6 to allow different scaling factors to be defined for tension and compression. See Appendix A for more details. See Floyd (2004) and Zobiery (2004) for further details on how to calculate the scaling factors. Currently a trial and error method (incorporated into an Excel spreadsheet) easily leads to an accurate elemental scaling factor by ensuring that the characteristic fracture energy equals the elemental fracture energy. • The aspects of the rubble box were incorporated into CODAM. The response of the rubble box is closely associated with the parameters describing the laminate box. The rubble box begins to carry load at the point of matrix damage initiation, indicating the failure of the first matrix compression fuse (location A in Figure 3.20). The load plateaus to a constant force at the saturation of matrix damage, indicating the failure of all of the matrix compression fuses (location B in Figure 3.20). Finally, the plateau load drops to zero at complete matrix and fibre damage saturation, indicating the failure of all of the matrix and fibre compression fuses (location C in Figure 3.20). The rubble box in effect becomes completely dependent on the growth of the parameters already in CODAM, with exception of the magnitude of the actual plateau stress. Therefore, the incorporation of the plateau stress into CODAM becomes quite 75 Chapter 3 Numerical Implementation of a Compressive Damage Model simple, only requiring 3 additional parameters, namely the plateau stress in each principal direction (x, y, and z). See Appendix A for more details. • Two minor corrections/modifications were made to the existing implementation of CODAM. When scaling was active, the elemental unloading or reloading stress-strain profile was bi-linear, with an apparent change in the modulus as the strains crossed the strain-at-peak-stress. This unloading or reloading path should theoretically be linear, as the amount of damage (and thus the modulus) is not changing. Further investigation into this matter indicated that the definition of stress needed to be redefined for cases with changing strains but a constant damage state. Figure 3.23 illustrates the modifications incorporated into CODAM for a one-dimensional case with modified Bazant scaling, where ee is the elemental strain, sc is the corresponding characteristic strain at the same stress level, k is the scaling factor, e k is the strain at peak stress, E0 is the initial undamaged modulus, E is the instantaneous normalized damaged modulus, and £"e(max) is the maximum elemental strain. Figure 3.23 Schematic showing the erroneous bi-linear unloading path and corrected . linear unloading path for a one-dimensional case assuming modified Bazant scaling. 76 Chapter 3 Numerical Implementation of a Compressive Damage Model A minor error in the Bazant scaling sub-routine was also corrected. In this subroutine a table is created relating the post-peak characteristic strains to the post-peak elemental strains. Previously, the coded formula incorrectly calculated this relationship. This error resulted in predicted elemental fracture energies that were not equal to the characteristic fracture energy. A simple correction to the calculation of the post-peak elemental strains resulted in accurate fracture energy predictions. Figure 3.24 compares the calculation of the post-peak elemental strains for the case of ri>\ before and after the correction, where he is the element height, hc is the characteristic height, n is the ratio of hjhe, cr'is a given stress value, s'e is the elemental strain at a', s'c is the characteristic strain at a', s'L is the loading strain at a', ande'y is the unloading strain at cr'. 3.3.3.1 Physically Based CODAM parameters Assuming plane stress conditions, CODAM involves several parameters, which can be organized into 5 main categories (see below). Floyd (2004) presented some general guidelines for approximating some of these parameters in tension. This section briefly 77 Chapter 3 Numerical Implementation of a Compressive Damage Model reviews some of these assumptions and discusses their applicability in compressive damage modelling. . Elastic Constants ((E0 )x , (E0 )y, vxy, (G0 )xy ) Given that the fibre and matrix properties are known, a simple rule-of-mixtures theory, followed by laminate plate theory (LPT) can provide a good approximation for.these parameters. • Effective Strain Interaction Constants (K'f, L], Mj, Sj, T,c, Uf,mdKf, i f , Mf , Sf, 7 f , Uf in Equation 3.4 / = x and y) These parameters remain to be supplied by the user. A simple approach is to choose the parameters such that the effective strain is equal to the strain in that direction. For yr y, y* y* y, example in the principal x-direction (tension), setting Lx , Mx , Sx , Tx , and Ux to very large values and K7X equal to one results in Fx = £[ . . Effective Strain Parameters (F T m i,Fl, F j ' , FTfi,,F^, Fj, F%) These parameters define the effective strain values that cause damage to initiate and saturate in the sublaminate constituents in both tension and compression. These values remain as the most subjective parameters in the model and are highly dependent on the lay-up and constituent properties. A thought experiment was presented by Floyd (2004) to estimate these parameters in tension and some of the general guidelines are presented below. 1. Matrix cracking (in off-axis layers) is generally assumed to be the first manifestation of damage in a sub-laminate under tension. The failure strain of the matrix material is generally used as the matrix damage initiation strain. 2. Fibre damage in tension typically occurs shortly after matrix damage initiation. The failure strain of a unidirectional lamina of the same material is generally used as the fibre damage initiation strain (given that fibres are oriented parallel to the loading direction). 78 Chapter 3 Numerical Implementation of a Compressive Damage Model 3. The matrix and fibre saturation strains in tension are much more difficult to predict. They may saturate together or independently of each other. For simplicity, they are typically assumed to saturate at the same strain state. Saturation strains in the range of 2 to 4% in tension are suggested (Floyd, 2004). In compression, the mechanisms of failure and behaviour are different than in tension. A similar thought process can be used to estimate the initiation and saturation effective strain parameters. 1. Matrix damage is again likely to be the first manifestation of damage. Matrix damage will likely form in 0° (or near 0 degree) plies as cracking or plastic deformation of the matrix between the fibres due to high shear stresses that develop as the fibres continue to bend under the applied compressive load (Figure 3.25). The initiation of matrix damage is highly dependent on the misalignment of the fibres, with earlier initiation in laminates with higher degrees of misalignment. Assuming a value of strain equal to or slightly below (depending on the degree of misalignment) the failure strain of the matrix material. At the same time, matrix damage is likely to form in the off-axis plies as matrix cracking and/or delamination. 2. Fibre damage in compression is typically assumed to initiate shortly after matrix damage initiation. This fibre damage occurs at the boundaries of the kink band due to significant bending. As the initiation of damage in the matrix typically means that a kink band has started to form, assuming fibre damage initiation values of strain slightly larger than or equal to the initiation strain of the matrix seems reasonable. In laminates without significant kinking (i.e. with little or no fibres aligned with the loading direction), fibre damage may initiate at a later strain value. 3. The fibre and matrix saturation strains are once again rather difficult to estimate. Matrix damage is typically assumed to saturate once the initial damage unit has formed across the entire RVE and before longitudinal 79 Chapter 3 Numerical Implementation of a Compressive Damage Model broadening of the kink band (stage 3 in Figure 3.22). The matrix damage saturation strain depends on the ease with which the damage unit (composed of kinking, matrix cracking, and delaminations) progresses through the RVE. Matrix saturation strains in the range of 1 to 3% seem to be reasonable. As the compressive strains continue to increase, the damage unit spreads longitudinally into previously undamaged material. With this progression, the two sections of the RVE begin to slide past one another (stage 4 in Figure 3.22). Fibre damage saturation in compression corresponds to the point at which all the fibres have failed and the two sections of the RVE slide past one another. This fibre failure is of the form of breaking due to high shear or bending strains at the boundaries of the kink bands. As fibres can bend significantly before breaking, the fibre saturation strain can be quite large (DeTeresa et al, 2001). Assuming values at least twice the matrix saturation strain in compression seems to be reasonable. Figure 3.25 Schematic showing the progression of failure mechanisms in compression. Matrix damage is assumed to initiate first, followed by fibre damage in the kink band. • Damage Parameters ( co'm, co7fs, cocmx, coLfli) These parameters are simply based on the area projections of each constituent in the failure plane. For example if 45% of the failure plane area consisted of fibres, co'is would 80 Chapter 3 Numerical Implementation of a Compressive Damage Model equal 0.45. A simple approximation of these parameters can be made by considering a sub-laminate cross-section (taken perpendicular to the loading direction). Approximating the area of the cross-section represented by each of the constituents easily leads to an approximation of the damage parameters. Additionally, assuming the compressive parameters are equivalent to the tensile parameters is reasonable, as the proportion of fibres and matrix within the failure plane is not likely to vary significantly. • Modulus Reduction Parameters (ETX, E%s) These parameters represent the remaining modulus (in terms of a normalized secant modulus) in the laminate assuming complete matrix damage. To estimate normalized modulus at matrix damage saturation in tension, Floyd (2004) suggested using the ply discount method and laminate plate theory (LPT). At matrix damage saturation, 90° plies were assumed to retain 0% of their original lamina stiffness, 0° plies were assumed to retain 100% of their original lamina stiffness, and angle plies retained some percentage (between 0% and 100%) of their original lamina stiffness. In addition to these guidelines, the shear stiffness of all the plies is typically reduced to zero. He suggested two ways of representing matrix damage in the laminate. The first method, removes the damaged 90° and angle plies (i.e. 90° or 45°) completely, effectively reducing the lamina stiffness to zero in all directions. In the second method the stiffness in the matrix direction of all plies is reduced. In this case, a 90° lamina would have no stiffness in the matrix direction, but full stiffness in the fibre direction. Floyd presented an example of the calculations for a [0/±45/90/±45/0] carbon/epoxy laminate with an initial elastic modulus of 50 MPa. In the first method, at matrix saturation the 90° and ±45° layers are completely removed. Only the 0° layers remain, and the reduced laminate modulus is 33 MPa. For this case, the normalized sub-laminate modulus at matrix saturation (ETS) would be 33/50=0.66. In the second method, only the matrix direction stiffness values of the plies are set to zero. The laminate modulus reduces to 42 MPa, resulting in a normalized sub-laminate modulus at matrix saturation ( E^ ) equal to 0.84. The same ply-discount method used by Floyd can be used to predict the compressive modulus reduction parameters. However, there is one key difference. The assumption 81 Chapter 3 Numerical Implementation of a Compressive Damage Model that 0° plies maintain 100% of their lamina stiffness at matrix saturation is no longer valid. In tension, this was assumed because the tensile lamina stiffness of a 0° ply with no matrix is approximately the same as that of a 0° ply with matrix, as a majority of the stiffness comes from the fibres. In compression, this is not the case. The lamina stiffness of a 0° ply without the matrix is significantly reduced because the damaged matrix no longer provides confinement to the embedded fibres leading to instability and buckling under compressive load. This buckling of the fibres leads to an overall reduction in the lamina stiffness. In conclusion, the lamina stiffness (at matrix damage saturation) is generally assumed to be 0 to 10% of the original in the fibre direction and zero in the direction perpendicular to the fibre direction. Using the same [0/±45/90/±45/0] laminate with all lamina stiffness values reduced to 10% in the fibre direction and 0% in the direction perpendicular to the fibre direction results in normalized sub-laminate at matrix damage (E„s) equal to 0.066. Table 3.1 summarizes the guidelines for the reduction of lamina stiffness values in tension and compression for use in LPT. Table 3.1 Guidelines for remaining lamina stiffness values (at matrix damage saturation) for use in LPT Lamina Direction TENSION COMPRESSION (% of original) (% of original) Parallel to fibre direction 100 0-10 Perpendicular to fibre direction 0 0 Shear stiffness 0 0 3.4 Verification of CODAM To verify that the compressive model has accurately been incorporated and coded into the existing CODAM material model, a series of single element simulations were performed on both a shell and a brick element. Strains were introduced into the elements by applying displacements to the nodes. The resulting stress-strain profile was compared to the analytical solution. Due to the importance of energy absorption in tube crushing simulations, a thorough investigation of the various LS-DYNA outputs for energy was also completed. The 82 Chapter 3 Numerical Implementation of a Compressive Damage Model material properties and CODAM specific parameters used in all of the single element simulations are shown in Table 3.2. Table 3.2 CODAM parameters for the single element simulations General Laminate Properties Ex 50GPa Gyy 9 GPa Ey 50 GPa 0.0 tlaminate 2.40 mm CODAM Specific Parameters (x = y) Tension Compression K 0.008 FTfl 0.008 Fc mi 0.008 Fl 0.010 K 0.025 F'l 0.025 Fc 0.018 Fl 0.080 coTms 0.55 0)Tfi 0.45 ©I ms 0.55 a>% 0.45 E~l 0.87 Ec 0.01 hi' 8 mm hcc 7 mm GTf 54 kJ/m2 Gcf 30 kJ/m2 kT (he= 5mm) 2.04 &c(he=5mm) 1.71 100 MPa A few observations can be made with respect to the CODAM parameters. The parameters are assumed to be equal in both local x and -^directions. As the parameters are based on a fictitious material and are meant for verification purposes only, details relating to how the values were chosen are not discussed. Rather, some of the general similarities and differences between the tensile and compressive parameters are presented. In both tension and compression, the damage due to matrix saturation (col,col) w a s assumed to be 55%. This is equivalent to saying that if all of the matrix in a given cross-section was completely damaged and there was no fibre damage, 55% of that cross-sectional area would be damaged. It is reasonable to assume that the damage parameters are independent of loading direction because regardless of how the damage is introduced into the laminate (tension or compression), the same proportion of the cross-section will be damaged at matrix damage saturation. Unlike the damage parameters, the modulus reduction parameters in tension are considerably different from those in compression. In tension, complete matrix damage 83 Chapter 3 Numerical Implementation of a Compressive Damage Model reduces the laminate modulus to 87% of its original value, whereas in compression, complete matrix damage reduces the laminate modulus to 1%> of its original value. This illustrates the high dependence of the compressive behaviour on matrix damage growth. The characteristic height of damage is assumed to be equal to 8 mm in tension and 7 mm in compression. 2 2 These parametric assumptions result in estimated fracture energies of 54 kJ/m and 30 kJ/m , in tension and compression, respectively. Scaling factors greater than unity account for the fact that the element size of 5 mm used in all of the simulations is smaller than the characteristic height of damage (he < hc). Using element sizes that are smaller than the characteristic height (effectively going into the damage zone) is numerically admissible. The overall response of an element with characteristic height is equivalent to the response of n arbitrary elements in sequence (where n = characteristic height/element height). In each case, the scaling factor modifies the post-peak relationship between the effective strain and secant modulus such that the overall fracture energy is constant. 3.4.1 Analytical Solution 3.4.1.1 Tension The analytical tensile normalized modulus-effective strain and stress-strain responses are shown in Figure 3.26. The normalized modulus curve (for simplicity, only shown for the characteristic element height) graphically shows how the secant modulus reduces as the effective strain increases. In all of the verification runs, parameters were chosen such that the effective strain fF. ) equals the strain in the direction of interest (i.e. Fx - sTx). This figure additionally illustrates the significant dependence of the normalized laminate modulus on the fibre behaviour. The stress-strain profile is shown for both the characteristic height (un-scaled stress-strain profile) and the element height of 5 mm (scaled stress-strain profile). In tension the stress reaches a peak of about 450 MPa at a strain of about 0.012. 84 Chapter 3 Numerical Implementation of a Compressive Damage Model 1.05 Effective Strain, tension (mm/mm) Tensile Strain (mm/mm) Figure 3.26 Analytical Tensile Response (x=y): (a) Normalized modulus curve and (b) Stress-strain curve for the single element simulations. Analytically, energy is either stored, to be recovered upon unloading, or dissipated as strains increase. Figure 3.27 illustrates how the total energy density of a given element subjected to a particular level of tensile strain (sTx) can be separated into three main types, pre-peak energy dissipation, post-peak energy dissipation and recoverable elastic strain energy. These energy densities are represented as areas under the stress-strain profile (having units of kJ/m3). Pre-peak energy dissipation, y'y, represents the energy dissipated by the element before the peak stress is reached. It corresponds to the energy dissipated by the bulk material before damage localizes into a band. Post-peak energy dissipation, {yTf) , represents the energy dissipated by the element after the peak stress has been surpassed. It corresponds to the energy dissipated by the material within the damage zone after damage localization into a band. Recoverable elastic strain energy, (y'e \ , represents the energy that is released as the element unloads from a given strain (less than the ultimate failure strain) to a strain and stress free state. 85 Chapter 3 Numerical Implementation of a Compressive Damage Model 0-012 T 0.038 Tensile Strain (mm/mm) Figure 3.27 Graphical representation of various energy densities for a given element subjected to a given pre-failure tensile strain (x=y). Figure 3.28 illustrates how the values of the recoverable and dissipative energy densities vary as the strain in the element increases from zero to its ultimate value of 0.038. The recoverable elastic strain energy increases from zero (at zero strain) to a maximum of 4500 kJ/m3 and then back to zero at the ultimate strain of 0.038. The pre-peak energy reaches and maintains its maximum of approximately 900 kJ/m at a strain equal to 0.012 (corresponding to the strain at peak stress). The post-peak energy increases from zero (at a strain value of 0.012) to its maximum value (equal to 10900 kJ/m3) at final failure. Figure 3.28 also shows the total energy profile of the element. 86 Chapter 3 Numerical Implementation of a Compressive Damage Model 12 0.00 0.01 0.02 0.03 0.04 Tensile Strain (mm/mm) Figure 3.28 Breakdown of various energies (from stress-strain profile) for an element, showing three components of the total energy density. 3.4.1.2 Compression Figure 3.29 shows the analytical compressive response, illustrating the variation of both the normalized modulus and the stress with increasing values of compressive strain. As noted previously, in all of the simulations, the effective strain equals the primary strain. The normalized modulus-effective strain profile demonstrates how, .unlike in tension, the secant modulus of the system is almost independent of the fibre damage. The stress-strain profile is shown for both the characteristic height and the element height. The peak stress of about 400 MPa occurs at a strain of about 0.009. 87 Chapter 3 Numerical Implementation of a Compressive Damage Model c Effective Strain, compression (mm/mm) Compressive Strain (mm/mm) Figure 3.29 Analytical Compressive Response (x=y): (a) Normalized modulus curve and (b) Stress-strain curve for the single element simulations. As was done in tension, the analytical stress-strain profile for the element can be used to illustrate the various forms of stored and dissipated energy. For simplicity, two loading cases with different maximum compressive strains will be considered. The first loading case considers an element with a maximum compressive strain (£%) that places the elemental stress in the pre-plateau stress range. Figure 3.30 shows three energy density components, pre-peak energy dissipation, post-peak initial fracture energy dissipation and recoverable elastic energy. The pre-peak energy dissipation,^, once again represents the energy dissipated by the bulk material in the element before the peak stress (and damage localization) is reached. The post-peak initial fracture energy dissipation, {ycf ) , represents the energy dissipated by the element as the first kink band (or more generally, the damage zone) forms. And finally, the recoverable elastic strain energy, ( y f )e, represents the energy that would be released if the element was unloaded to a strain and stress free state. 88 Chapter 3 Numerical Implementation of a Compressive Damage Model 0.009 £x 0.132 Compressive Strain (mm/mm) Figure 3.30 Graphical representation of energy development in compression (x—y). The second loading situation considers the case where the maximum strain (scx ) places the elemental stress within the plateau range as shown in Figure 3.31. The relative stored and dissipated energies are shown on this figure as well. In addition to the three previously mentioned energy types (pre-peak energy dissipation, post-peak initial fracture energy dissipation and recoverable elastic energy), there is a fourth type, band broadening energy dissipation, {y^b) . The post-peak band broadening fracture energy arises from the propagation of the initial fracture zone lengthwise (i.e. along the direction of loading), therefore causing an overall growth in the damage band height. This can be visualized as the successive addition of new bands of damage onto the original one. 89 Chapter 3 Numerical Implementation of a Compressive Damage Model Post-peak energy dissipation (energy dissipation due to initial fracture) 0.009 sc 0.132 X Compressive Strain (mm/mm) Figure 3.31 Graphical representation of energy development in compression (x=y). As was done in tension, the various forms of dissipative and recoverable energies can be plotted as the compressive strain increases from zero to the ultimate value of about 0.132 (Figure 3.32). The recoverable strain energy again begins to increase upon loading and reaches a maximum value of 2600 kJ/m3 at a strain of approximately 0.014. It then proceeds to drop off to a relatively constant value that is maintained until it drops to zero at final failure. The pre-peak energy dissipation reaches its maximum at 0.009 strain (corresponding to the strain at peak stress) and maintains this value to failure. The initial fracture energy is zero up to the strain at peak stress. It then grows up to a maximum value of 5700 kJ/m3, which is maintained up to final failure. The band broadening fracture energy dissipation begins to grow at a strain of 0.018 and continues to increase linearly up to a maximum value of 11700 kJ/m . Figure 3.32 also shows how the various energies sum up to the total energy. 90 Chapter 3 Numerical Implementation of a Compressive Damage Model 0.00 0.02 0.04 0.06 0.08 0.10 0.12 0.14 Compressive Strain (mm/mm) Figure 3.32 Breakdown of various energies (from stress-strain profile). 3.4.2 Numerical Solution The details of the shell and brick elements used in the verification runs including the local and global co-ordinate systems are shown in Figure 3.33. Both elements were 5x5 mm and had a defined thickness of 2.40 mm. The elements were oriented such that the assumed laminate x-direction was aligned with the global Z-direction. Displacements were applied to the nodes in the global Z-direction (local x-direction) and the development of stresses and energy in the element was studied. 91 Chapter 3 Numerical Implementation of a Compressive Damage Model 3isplacement applied to nodes 3,4,7 and i in global Z direction (local x-direction) Displacement applied to nodes 3 and 4 in global Z direction (local x-direction) \ r I i^ . i l oiwirjiniitc i\stem Global co-ordinate system z X 5 x 5 mm SHELL ELEMENT (a) 5 x 5 x 2.4 mm BRICK ELEMENT (b) Figure 3.33 Loading applied to a single (a) shell element and (b) brick element in the local and global co-ordinate systems. 3.4.2.1 Tensile Response Simulation Figure 3.34 illustrates the tensile stress-strain response of both the shell and solid element under a ramped tensile displacement. The stresses reach a peak stress of approximately 450 MPa at a strain of about 0.012, which agree with the analytical solution. Additionally, the stresses drop to zero, indicating complete failure in this direction at the analytical failure strain (0.038). 92 Chapter 3 Numerical Implementation of a Compressive Damage Model 500 0 0.01 0.02 0.03 0.04 Tensile Strain (mm/mm) Figure 3.34 Numerical stress-strain profde for 5 mm shell and solid element (Tension). Results agree with the analytical solution. LS-DYNA provides several means of measuring the energy in a particular simulation. The global data (GLSTAT) can be used to output the kinetic, internal and total energies of the entire model. The material energies (MATSUM) can also be used to provide the same energy outputs on a material basis. Finally, the force-displacement profile can be manually integrated to represent the energy in the system. All of these methods provide the energy in units of mJ. Dividing them by the volume of the element, a comparison can be made between the LS-DYNA energy outputs and the analytical solution in Figure 3.28. Figure 3.35 shows that all three numerical outputs result in profiles for the single shell and solid element simulations that are identical to the analytical total energy profile. Therefore, LS-DYNA outputs provide the total energy in the system (including any elastic recoverable energy). 93 Chapter 3 Numerical Implementation of a Compressive Damage Model 0.00 0.01 0.02 0.03 0.04 Tensile Strain (mm/mm) Figure 3.35 Numerical output of total energy density for shell and solid element (Tension). Agrees with the analytical solution. 3.4.2.2 Compressive Simulations The numerical stress-strain profile (Figure 3.36) also matches well with the analytical solution, reaching a peak stress of 400 MPa at a strain of about 0.009. Additionally, the plateau stress of 100 MPa is observed and maintained up to failure at a strain of 0.132. CS cn t in <u > in co cu s-I O U 450 400 350 300 250 0.02 0.04 0.06 0.08 0.1 Compressive Strain (mm/mm) 0.12 0.14 Figure 3.36 Numerical stress-strain profile for 5 mm shell and solid element (Compression). Agrees with the analytical solution. 94 Chapter 3 Numerical Implementation of a Compressive Damage Model LS-DYNA again provides numerical energy density outputs for the total energy of the system that correspond to the analytically calculated profile as shown in Figure 3.37. 0.00 0.02 0.04 0.06 0.08 0.10 0.12 0.14 Compressive Strain (mm/mm) Figure 3.37 Numerical output of total energy density (Compression). Agrees with the analytical solution. 3.4.2.3 Strain Reversal Simulations In addition to the simple tension and compression simulations, the shell and solid elements were also subjected to a strain reversal simulation. Ramped displacements were again applied to the nodes in a manner resulting in initial compressive strains followed by unloading and then tensile strains to failure. Figure 3.38 shows the stress-strain profiles of both the shell and solid element. The profiles follow the compressive stress-strain profile in Figure 3.36 up to a compressive strain of about 0.063. At this point, the strain was reversed causing the stress to follow the unloading path back to a strain and stress-free state. As the strains moved into the tensile domain, the stresses increased linearly with a slope defined by the reduced secant modulus. Eventually, tensile strains reached a value that causes damage to grow (and thus the secant modulus to further decrease) resulting in a softening curve that followed the stress-strain profile in the pure tensile simulations (Figure 3.34) to ultimate failure. 95 Chapter 3 Numerical Implementation of a Compressive Damage Model 300 Loading - no damage growth Shell Solid 200H iooH 0.01 0.02/0.03 0.0^ Loading - / damage growing -500 Strain (mm/mm) Figure 3.38 Stress-strain profde for strain reversal simulation. Loaded in compression, then unloaded and loaded to failure in tension. The numerical output of the total energy as strains progress from compression to tension is shown in Figure 3.39. The initial growth in the total energy up to compressive strains of 0.063 follow the total energy profile in Figure 3.37 reaching a maximum of 10800 kJ/m . Upon unloading the stored elastic strain energy is released, resulting in a drop in the total energy to a value of about 10000 kJ/m . As the strains increase into the tensile domain, the total energy increases. Initially, this increase is attributable only to the growth in the recoverable strain energy (because damage is not growing and unloading would occur along the same loading path). Once damage begins to grow, the total energy is comprised of a combination of the recoverable strain energy and the post-peak fracture energy. The total energy in the strain reversal simulation reaches 15000 kJ/m , a value larger than was absorbed in the pure tensile simulations (11800 kJ/m ) but smaller than the value in the pure compression simulations (17500 kJ/m ). 96 Chapter 3 Numerical Implementation of a Compressive Damage Model -0.07 -0.06 -0.05 -0.04 -0.03 -0.02 -0.01 0.00 0.01 0.02 0.03 0.04 Strain (mm/mm) Figure 3.39 Numerical output of total energy density (Compression - Tension). 3.5 Summary In tension, the main failure mechanisms in fibre reinforced polymer composites are fibre and matrix cracking, and delamination. In compression, kinking is the dominant damage mechanism observed in composites with fibres aligned with the loading direction. In multidirectional or braided composites, this kinking is joined by off-axis matrix damage and delamination. Current predictive methods available to capture the compressive response of composites are limited. Analytical methods are generally not very robust and do not provide any information regarding the post-peak behaviour of the composite. Micro-mechanical models, although capable of accurately modelling the response of composites under various loading conditions, are very complex and thus become very inefficient when modelling structure of considerable size. Macro-mechanical models offer a sufficiently accurate method of modelling the response of composites that is considerably more efficient and robust than either of the two previous methods. CODAM is a physically based macro-mechanical model that captures the laminate or sub-laminate constitutive behaviour at the level of a Representative Volume Element (RVE). 97 Chapter 3 Numerical Implementation of a Compressive Damage Model The dimensions of the RVE have significant meaning. The height and width are equal to the characteristic height of damage growth in each perspective direction. The thickness of the RVE is equal to the laminate or sub-laminate thickness. CODAM predicts the overall stress-strain response of the RVE by defining relationships between the effective strain and the modulus in each local direction. In tension, the growth of damage in the RVE results in an overall softening stress-strain. In compression, the overall response of the RVE is characterized by initial strain softening followed by a relatively constant plateau stress. The material model has been incorporated into the commercial finite element code, LS-DYNA, to allow efficient modelling of composite structures under a variety of loading conditions. 98 Chapter 4 Calibration and Validation of Compressive Model C H A P T E R 4 CALIBRATION AND VALIDATION OF COMPRESSIVE M O D E L 4.1 Introduction This chapter presents the results of tensile and compressive simulations on braided panel specimens. These simulations serve two purposes, to validate CODAM, and to determine the parameters to be used in the tube crushing simulations. The first section reviews the experimental details of the braided panel specimens tested by Lawrence Livermore National Laboratory (LLNL) (DeTeresa et al., 2001). The next section presents and describes the CODAM specific parameters used in panel simulations. The last three sections present the simulation results and compare the predicted peak stresses to experimentally observed data. 4.2 Panel Material Properties and Experimental Details The panels were made of Hetron 922 resin, (an epoxy vinyl ester) and two layers of a [0/±30] braid of Fortafil #556 80k carbon fibre tows, with an overall fibre volume fraction of 45% and a thickness of 3 mm. All panels were of very poor quality, with warped surfaces, resin-rich regions, poor tow distribution/alignment, and dry and poorly bonded tows. The average mechanical properties of the resin and fibre are listed in Table 4.1 below. Table 4.1 Average mechanical properties of the Fortafil #566/Hetron 922 composite constituents (DeTeresa etal, 2001) Constituent Modulus (GPa) Shear Modulus (GPa) Tensile Strength (MPa) Compressive Strength (MPa) Shear Strength (MPa) Tensile Elongation at Break (%) Resin 3.2 1.38 83 114 36 3.8 Fibre 231 - 3800 - - 1.6 Experimental tension and compression tests were performed along the local x and y-directions, as well as at 30° and 60° (from the local x-direction), as shown in Figure 4.1. All specimens were un-tabbed and loads were introduced via hydraulic grips and surfalloy 99 Chapter 4 Calibration and Validation of Compressive Model wedges. The tension specimens were 163 mm x 38 mm, as shown in Figure 4.2, and were loaded in an Instron test machine at a rate of lxlO"4 s"1. The compression specimens were 22mm x 38mm (Figure 4.2) and were loaded in a servohydraulic test machine at a rate of lxlO"3 s"1. Note that in the off-axis (i.e. 30° and 60°) tests, the specimen is subjected to shear strains in addition to the applied axial strains because the restrictions of the fixed supports prevent shear coupling. If supported on rollers, the anisotropy of the composite would cause the panel to shear in-plane as the displacements are applied. Local x-direction 30 °from Braider Figure 4.1 Schematic of basic 2-dimensional braid architecture, showing the four loading directions in the panel tests. 100 Chapter 4 Calibration and Validation of Compressive Model TENSION COMPRESSION Figure 4.2 Schematic showing the testing configuration for tension and compression. In tension, the general g loba l failure mechanisms o f the panel specimens were a combinat ion o f failure between the tows, fibre fracture, and de laminat ion as s h o w n i n F igure 4.3. Damage typ ica l ly ini t ia ted at the edges o f the specimens as interfacial fai lure between the tows not a l igned w i t h the load ing direct ion (off-axis tows) and the mat r ix . F igure 4.4 shows the orientation o f the a l igned and off-axis tows i n each o f the four l o a d i n g direct ions. A l l o f the off-axis tows are oriented such that they intersect at least one edge. It was at these intersection points that c rack ing was observed dur ing the in i t i a l stages o f loading. A s loading continued, a larger crack generally grew f rom one (or two) o f these in i t ia l edge cracks across the spec imen along a path paral lel to the off-axis tows. F i n a l fracture was characterized by fibre fai lure i n tows al igned wi th the load ing d i rec t ion and spli t t ing (wi th very litt le or no fibre fracture) a long off-axis tows. The tensile strength o f matr ix dominated directions ( local y and 6 0 ° direction) was consistently less (by up to a factor o f 3) than the tensile strength o f the mat r ix alone, supporting the s ignif icant effect the poor tow/matr ix bond has o n the tensile behaviour o f these specimens. In compression, fai lure was characterized by b u c k l i n g o f the braider and ax i a l tows and fibre failure (Figure 4.3). D a m a g e was often observed over the entire gauge length. In a l l cases, 101 Chapter 4 Calibration and Validation of Compressive Model the peak stress measured in compression was 18 to 52% less than the compressive strength of the resin alone, indicative of the poor quality of the panels. A post-peak plateau stress, equal to approximately half the peak stress value was sustained by all specimens for a considerable amount of time before final failure (Figure 4.5). Table 4.2 summarizes the initial modulus and the peak stresses measured for all the specimens in both tension and compression. Braider tows remained relatively intact, with little or no fibre fracture f surface « and axial tows Figure 4.3 Typical failure characteristics in tension (a) and compression (b). Local x-direction (DeTeresa etal, 2001). Figure 4.4 Schematic showing the orientation of the tows in the four loading directions (shown for one axial-braider-braider tow group). The off-axis tows (tows not aligned with the loading direction) intersect at least one edge. 102 Chapter 4 Calibration and Validation of Compressive Model 0 0 . 0 5 0 .1 0 . 1 5 Apparent Compression Strain (in./in.) Figure 4.5 Example of the plateau stress observed in compression (local x-direction). Apparent strain calculated using machine displacement (DeTeresa et al, 2001). Table 4.2 Experimentally measured initial modulus and peak stresses (DeTeresa et al., 2001) INITIAL MODULUS PEAK STRESS PEAK STRESS COMPRESSION (MPa) DIRECTION (GPa) TENSION (MPa) LB* UB* Average LB* UB* Average LB* UB* Average Local x 33.1 50.2 41.4 204 363 300 71 110 94 30° from x 35.2 44.0 41.7 184 234 214 83 110 100 60° from x 13.1 14.5 13.8 28 35 32 41 64 54 Local y 5.7 8.0 7.0 18 30 24 44 81 57 * LB = lower bound, UB = upper bound 4.3 CODAM Parameters The CODAM specific parameters reflect the complete orthotropic response of a representative volume element (RVE) of the braided material in tension and compression. Both a lower and upper bound material model was defined to capture the large variation in the experimental results. Parameters were only defined for the local x and _y-directions, and the response in the 30° and 60° directions were predicted by rotating the material axes appropriately. The parameters were based on experimental data, when available. In other 103 Chapter 4 Calibration and Validation of Compressive Model cases, supplementary experiments on braided composites or previously successful simulation parameters were used to approximate the values. Table 4.3 presents the upper and lower bound CODAM parameters. Following the table is a discussion on each of the parameters. Table 4.3 CODAM parameters in both the local x and ^ -directions showing upper and lower bound assumptions CODAM PARAMETERS Local jc-direction Local j-direction Tens ion Compres s ion Tens ion Compress ion Initial elastic modulus Ea (GPa) 3 3 - 5 0 5 . 7 - 8 . 0 Matrix damage initiation F 7,, F^ 0.006 - 0.008 0.002 - 0.003 0.006-0.008 0.006 - 0.008 Matrix damage saturation F T F c m.v' ms 0.021 -0.025 0.010-0.011 0.021 -0.025 0.021 -0.025 Fibre damage initiation F Jfi , Ff: 0.006 - 0.008 0.002 - 0.003 0.011 -0.012 0.007 - 0.009 Fibre damage saturation FJS, F^s 0.021 - 0.025 0.120-0.150 0.050-0.060 0.120-0.150 Damage due to matrix saturation 0.55 0.55 Damage due to fibre saturation 0.45 0.45 Normalized modulus at matrix damage saturation ETms, E%s 0.87 0.01 0.48 0.02 Plateau stress <JPLATEAU (MPa) - 5 5 - 6 0 - 4 5 - 5 0 Damage height hTc ,h^ (mm) 12 - 15 11 - 14 1 0 - 1 3 9 - 12 Fracture energy G Tf , GJ (kJ/m2) 36 -101 7 . 3 - 2 0 1 6 - 4 2 5 - 1 3 Modified Bazant scaling factor for element size o f 2 mm kT, kc 9.75-12.25 6 .58-8 .25 6.3 -8 .2 5 .9 -7 .96 INITIAL ELASTIC MODULUS The elastic moduli were based on the upper and lower bounds measured in the experimental panel tests using 1-inch (25.4 mm) foil strain gauges. Additionally, strain gauges placed perpendicular to one another provided an estimate of the major and minor Poisson ratios. The major Poisson's ratio ranged from 1.06 - 1.25, with an average of 1.15, and the minor Poisson's ratio ranged from 0.198 - 0.215, with an average of 0.206. 104 Chapter 4 Calibration and Validation of Compressive Model MATRIX DAMAGE INITIATION Damage in composites typically initiates in the matrix due to high stress levels that develop near discontinuities. As the quality of the panels was extremely poor with significant discontinuities, damage is expected to initiate sooner than in panels without this degree of variability. Matrix damage initiation values were assumed to be equal in both directions. In tension, matrix damage initiation values were assumed to be equal in both local directions. Floyd (2004) suggested that this parameter could be related to the elongation at break of the resin, which was reported to be 3.8%, but varied from 1.1% to 6.25% for the vinyl ester epoxy. The damage initiation strains in tension were assumed to be much lower, ranging from 0.6% to 0.8%. This is reasonable as the poor tow/matrix bonds will likely lead to interfacial failure (along the off-axis tows) well before the total potential of the resin is reached. Nonlinearities in the experimental tensile stress-strain profiles of the panels also indicate an early onset of damage. Studies conducted by Burr (1996) on similar braided composite panels noted that at strain values between 0.3 and 0.4%>, audible popping sounds were heard in most tensile tests, corresponding to the initiation of cracks oriented parallel to the off-axis tows. In compression, matrix damage initiation values are not assumed to be equal for both local directions. In the local x-direction, the initiation of damage in the matrix is related to the initiation of kinking in the axial tows. As a kink band forms, the shear strains between the successive fibres in the kink band have increased substantially. These high shear strains cause matrix damage, either in the form of micro-cracking or interfacial bond failure. Due to the significant waviness of the axial tows (caused by the braided architecture) and the poor interfacial bond, kinking (and thus matrix damage initiation) is assumed to initiate fairly early. Initiation values between 0.2 and 0.3%> are thus assumed. In the local ^-direction, kinking is unlikely to form, due to the lack of fibres aligned with the loading direction. Alternately, matrix damage will likely be caused by a simple increase in the compressive strains, rather than the high shear strains characteristic in kink bands. It is assumed that matrix damage will initiate later than in the local x-direction, at values of strain between 0.6 to 0.8%. 105 Chapter 4 Calibration and Validation of Compressive Model FIBRE DAMAGE INITIATION In tension, the orientation of fibres significantly affects the damage initiation strains. In the local x-direction, the axial fibres are aligned with the loading direction. As the matrix begins to micro-crack, the fibres adjacent to these cracks experience an almost instantaneous increase in their local tensile strains. These strains are assumed to initiate some form of fibre damage, and thus fibre damage is assumed to initiate at the same strain as matrix damage (0.6 to 0.8%). In the local ^ -direction, there are no fibres directly aligned with the loading direction (i.e. axial tows are now perpendicular and the braider tows are off by 60°). A micro-crack in the matrix will not initially cause an increase in the tensile fibre strains. Alternately, this micro-crack (which will likely form parallel to and possibly within the off-axis tows) will lead to the rotation of any fibres bridging the crack as they slowly align themselves with the loading direction (Figure 4.6). The tensile and bending strains in these fibres will continue to increase as the crack grows until the crack has opened enough to cause fibre failure. The exact initiation of damage to the fibres is hard to predict without detailed experimental testing, and thus arbitrary values ranging from 1.1 to 1.2% strain are assumed. These values, which are approximately twice the matrix damage initiation strain, attempt to represent the delayed onset of fibre damage. Crack growing within tow Original fibre direction Figure 4.6 Fibre bridging observed when a crack propagates parallel to a fibre tow, forming a large damage (or process) zone (Salvi et al., 2003). 106 Chapter 4 Calibration and Validation of Compressive Model In compression, the initiation values are again assumed to be different for each of the loading directions. In the local x-direction, kink band formation is assumed to be the first indication of damage. It is reasonable to assume that some fibre damage will occur simultaneously with kink band formation (Sivashanker et al, 1996, Sutcliffe and Fleck, 1994). This fibre damage is likely in the form of cracking at the boundaries of the kink band due to significant rotation of the fibres. Fibre initiation values in this direction are thus equivalent to the matrix damage initiation values (0.2 to 0.3%). Fibre damage in the local _y-direction is not likely to follow the same trend, as significant kinking is unlikely to form because of the lack of aligned fibres. Rather matrix damage will eventually lead to either transverse crushing of the fibres or failure due to significant out-of-plane buckling of the surface braider tows. It is assumed that fibre damage initiates slightly after matrix damage, and thus initiation values between 0.7 and 0.9% strain are assumed. MATRIX DAMAGE SATURATION The saturation strains are the hardest parameters to predict. Floyd (2004) used a value of 3% in tension for a quasi-isotropic carbon/epoxy [±45/0/90] s lay-up. In tension, values between 2.1 and 2.5% are assumed for matrix damage saturation for both the local x and -^directions. These values, which are smaller than assumed by Floyd, represent the poor quality and expected premature failure of the panels. In compression, matrix damage saturation values were not equal in both local directions. In the local x-direction, matrix damage saturation occurs when the first kink band (or damage unit) completely traverses the representative volume element (RVE). This is assumed to occur shortly after matrix damage initiation, as the kink bands easily traverse across the axial tows and individual kink bands connect together via matrix cracks and delaminations between the tows. Values between 1.0 and 1.1% are thus assumed. In the local y-direction, matrix damage is more distributed, as it is not localized into kink bands. Thus, saturation of this matrix damage across the entire projected cross-section of the RVE is likely to occur at a higher strain value. Values between 2.1 and 2.5% are assumed. 107 Chapter 4 Calibration and Validation of Compressive Model FIBRE DAMAGE SATURATION In all cases, saturation of damage in the fibres indicates final failure. In tension, fibre damage is assumed to saturate at the same strains as matrix damage saturation (2.1 to 2.5%) in the local x-direction. As most of the fibres are aligned or near aligned with the applied load, once matrix damage completely saturates the cross-section of the RVE, the strains in the fibres bridging these cracks will quickly cause failure. In the local ^ -direction, the fibres bridging the matrix cracks are not aligned with the loading direction, and thus can rotate significantly before strains reach failure levels. Higher fibre saturation values (between 5.0 and 6.0%) account for this. Experimentally, a significant post-peak stress observed in this direction was attributed to the load borne by the braider tows that bridged the fracture zone. In compression, fibre damage is assumed to saturate well after saturation of the damage in the matrix. This accounts for the ability of fibres to bend significantly (either in the kink band or in the buckled tows) before the strains reach a level that causes fibre failure. In all of the panel tests, ultimate failure occurred at displacements well after the peak stress was reached. Values of 12 to 15% strain are thus assumed in the local x and .y-directions. DAMAGE DUE TO MATRIX/FIBRE SATURATION These parameters represent the proportion of each of the constituents within the projected area of the fracture zone. An estimate of these parameters for a braided composite is based on the volume fractions. Given that the fibre volume fraction is 45%, it is reasonable to assume that roughly 45% of any fracture zone will be fibres, leaving 55% to the matrix. NORMALIZED MODULUS AT MATRIX DAMAGE SATURATION As matrix and fibre damage progresses in the RVE, the normalized modulus appropriately decreases. These parameters, which describe how the modulus is affected by damage growth in each of the constituents, are defined using the LPT and the ply discount method outlined in Section 3.3.3.1. For the LPT analysis, the braided composite is idealized as a [0/±30J2 tape laminate. This results in laminate moduli estimations (58.8 and 8 GPa in the local x and y directions) that are close to the measured values. The ply discount method can then be used, assuming complete matrix damage and no fibre damage. The individual plies stiffness values are reduced according to the suggestions outlined in Table 3.1. 108 Chapter 4 Calibration and Validation of Compressive Model In tension, the "fibre direction" ply stiffness of the individual plies remains at its initial value and the "matrix direction" ply stiffness reduces to zero. Additionally, the shear stiffness of the individual plies is assumed to reduce to zero. These modifications lead to laminate moduli estimations of 51.05 and 3.81 GPa. The normalized moduli at matrix damage saturation are thus equal to 0.87 (51.05/58.8) for the local x-direction, and 0.48 (3.81/8.0) for the local ^ /-direction. In compression, the "fibre direction" ply stiffness of the individual plies reduces to 1% of its original value and the "matrix direction" ply stiffness reduces to zero. The ply shear stiffness is again assumed to reduce to zero. These modifications lead to laminate moduli estimations of 0.84 and 0.142 GPa. The normalized moduli at matrix damage saturation are thus equal to 0.01 (0.84/58.8) for the local x-direction, and 0.02 (0.142/8.0) for the locals-direction. PLATEAU STRESS The plateau stress was estimated from the compressive force-displacement profiles of the panels. All of the tests indicated that the panel specimens could support a load equal to approximately half of the peak load. Based on this experimental data, values between 55 to 60 MPa and 45 to 50 MPa for the local x and y-directions, respectively, were assumed. CHARACTERISTIC DAMAGE HEIGHT This parameter describes the actual height the damage (fracture process zone) would grow to in a specimen under the applied tensile or compressive loads. Floyd (2004) assumed characteristic heights in tension between 5 and 10 mm for the quasi-isotropic laminate mentioned previously. For the case of the braided specimens, the tensile damage height was assumed to vary between 10 and 15 mm. This larger damage height is attributed to the coarse architecture (i.e. large 80k tows) of the braided panels. The post-failure zone of damage in all the specimens was observed to be fairly large. In compression, the characteristic damage height is the height to which the damage grows before it begins to propagate under the constant plateau stress. In unidirectional composites, Fleck et al. (1995) observed kink band heights between 10 and 15 fibre diameters (between 70 and 105 pm for typical carbon fibres). Sivashanker et al. (1996) also observed kink band heights in UD composites of the same magnitude. In multidirectional laminates, damage 109 Chapter 4 Calibration and Validation of Compressive Model progresses in terms of a damage unit, which includes matrix damage and delamination. This damage unit can extend over heights much larger than a single kink band, as previously illustrated in Figure 3.9. The height of damage before the self-similar progression of the damage unit lengthwise in the braided panel specimens was not thoroughly discussed, and further investigation into this matter is necessary to pinpoint the exact value. Figure 4.7 shows some of the damage heights observed in compression for various braided panels at different levels of strain. Damage heights are assumed to be slightly smaller than the corresponding tensile values. Values between 9 and 14 mm are assumed. Figure 4.7 Examples of compressive damage heights observed in (a) the local y-direction of [0/±30J2 carbon/vinyl ester panels(DeTeresa et al, 2001) and (b) the local x-direction, (c) the local y-direction and (d) 45° from the local x-direction of [0/±60]2 carbon/vinyl ester panels(Xiao, 2003b). FRACTURE ENERGY The fracture energy per unit area is the energy absorbed by the material during fracture in either tension or compression. With respect to CODAM, its value is a direct consequence of the parametric values assumed thus far. As the fracture energy for the braided panels was not experimentally measured, engineering judgement is used to decide whether the resulting values are admissible. In tension, the fracture energy in the local x-direction varies between 36 and 101 kJ/m2. In the local v-direction this reduces to between 16 and 42 kJ/m2. The higher relative fracture energy in the x-direction is a direct consequence of the larger percentage of fibre fracture in 110 Chapter 4 Calibration and Validation of Compressive Model this direction. These tensile fracture energies are in the same range as the values assumed by Floyd (2004) for a quasi-isotropic [±45/0/90] s laminate (50 and 100 kJ/m2). In compression, the fracture energy is the energy associated with the formation and complete growth of the first kink band (or damage unit) across the RVE (Zobiery, 2004). It does not include the additional energy absorbed by the RVE during band broadening. The fracture energy in the local x-direction varies between 7.3 and 20 kJ/m2. In the local v-direction this reduces to between 5 and 13 kJ/m2. The consistently lower fracture energy in compression is a consequence of the different failure mechanisms and the lack of complete fibre fracture at boundaries of the kink bands. The complete energy, including the energy required to broaden the fracture zone lengthwise into previously undamaged material (band broadening) may very well surpass the fracture energies reported in tension. Refer to Figure 3.31 and Figure 3.32 for a comparison between the fracture and band broadening energy density for a theoretical element. MODIFIED BAZANT SCALING FACTOR Scaling factors are defined to appropriately scale the post-peak stress-strain response for an element size of 2 mm. These scaling factors ensure that the fracture energy in the element is equivalent to the fracture energy of the characteristic element size (height hrc , hcc ). EFFECTIVE STRAIN FUNCTION The effective strain function, F„ (Equation 4.1) determines the extent that individual strain components affect the overall behaviour of the laminate in the local x and v-directions. An individual strain function is defined for tensile and compressive loading in each direction. In pure axial tension and compression tests, it is reasonable to equate the effective strain function, to the strain in the corresponding direction (i.e. Fx = ex and Fy = Sy). However, as was stated earlier, due to the fixed constraints, both the 30° and 60° tests were under high shear strains in addition.to the axial strains. These shear strains understandably lead to damage growth in both the local directions, which should be captured in the numerical model. This can be accomplished by defining the shear constants, Sx and Sy. Preliminary investigations indicated that the effective strain functions shown in Equation 4.2 accurately capture the effect both axial and shear strain have on the behaviour in the local x and y-111 Chapter 4 Calibration and Validation of Compressive Model directions. These strain functions were assumed to be equivalent in tension and compression. v 2 K, v 2 L, F, = Jl — I + i = x, y, and z M. 2 \ 2 4.1 F, =. V 1 J I 3 J 2 (£ \ y + I1 J { 3 / 4.2 4.4 Tensile Simulation Results The experimental peak stress values are compared to simulation predictions using CODAM. The simulation details are shown in Figure 4.8. By defining different rotation angles in the material description, the various loading directions (x, y, 30° from x, and 60° from x) are simulated. Loading was applied via ramped displacements applied to the top nodes. To capture the edge effects that are suspected of causing premature damage initiation at the edges, some weak elements with identical damage parameter but a 50% reduction in the elastic moduli have been defined along both edges. A typical failure pattern and resulting panel stress-strain profile is shown in Figure 4.9, for the upper bound local x simulation. Damage initiated in one of the weak elements and proceeded until all the elements in that row failed. The predicted initial modulus of 51'GPa and peak strength of 386 MPa are similar to experimentally observed values of 50 GPa and 363 MPa. 112 Chapter 4 Calibration and Validation of Compressive Model Global Coordinate System ' ' Figure 4.8 Tensile simulation details. Overall dimensions are 38 x 162 mm. Ramped displacement is applied to top nodes, while bottom nodes are fixed. 113 Chapter 4 Calibration and Validation of Compressive Model 450 7 G l o b a l C o o r d i n a t e S y s t e m (a) Figure 4.9 Example of the failure behaviour (a) and nominal stress-average strain response (b) of the panel in tensile simulations, shown here for the upper bound local x-direction. Predicted initial modulus and peak stress are similar to upper bound experimentally measured values. Figure 4.10 shows the comparison of the experimental and numerically predicted peak stress values. Peak stresses reported by Beard (2001) on identical braided material are also included in the figure. The figure demonstrates that peak stresses are accurately predicted using CODAM as a material model. The slight over-prediction of the peak stress values in the matrix dominated directions (local y and 60°) is attributed to the extreme edge effects in these directions, which may not be completely captured by the numerical mesh. The edge effects are significant in these directions because of the shallow orientation of the tows, as shown in Figure 4.4, which cause significant fracture initiation points distributed along the edges. In many cases, these off-axis tows pull out instead of fracture, resulting in a measured strength that is less than the true strength of the braided specimen. Defining the specimen with elements designed to capture the true response of the material is thus not accurate. 114 Chapter 4 Calibration and Validation of Compressive Model Degrading the elastic properties of a four elements along each side may not be the ideal approach to take to capture this effect. 400 o r x -•—Experimental (UB and LB) -±— Numerical (UB and LB) 0 Beard (2001) 30° from x 60° from x y Loading Direction Figure 4.10 Comparison of upper and lower bound predictions of tensile peak stress with experimentally measured values. 4.5 Compressive Simulation Results The experimental compressive peak stress values are compared to the numerical predictions. Figure 4.11 shows the details of the numerical model. As the edge effects are not considered to be significant in compression (mainly due to the small gauge length of 22 mm), no weak elements were incorporated into the mesh. Figure 4.12 shows the failure behaviour and resulting panel stress-strain profile, including a post-peak plateau stress, predicted for the local x-direction (upper bound). The plateau stress was maintained until an average strain of about 6.5%. The predicted initial modulus of 8 GPa and peak stress value of 75 MPa is close to experimentally measured values of 8 GPa and 81 MPa in this direction. 115 Chapter 4 Calibration and Validation of Compressive Model Nodes are fixed in X and Y z x Global Coordinate System Time Bottom nodes are fixed in X, Y, and Z 2 mm x 2 mm elements Ramped Displacement Figure 4.11 Compressive simulation details. Overall dimensions are 38 x 22 mm. Ramped displacement is applied to top nodes, while bottom nodes are fixed. lit Ill /, k •'>/ 11 Wt. \ %. M. ft. •A. Experimental U B_strength _ _ Numerical UB strength Numerical UB initial modulus Experimental LB strength Experimental LB and UB initial modulus Failure of one complete row of elements (a) 0.01 0.02 0.03 0.04 0.05 Average Strain (mm/mm) (b) 0.06 0.07 Figure 4.12 Example of the failure behaviour (a) and nominal stress-average strain response (b) of the panel in compressive simulations, shown here for the upper bound local x-direction. Note that compressive stresses and strains are shown as positive values. Predicted values are similar to the experimental values. Figure 4.13 compares the predicted peak stress values to the values measured experimental ly. Peak stresses reported by B e a r d (2001) on ident ical braided mater ia l are also inc luded i n the figure. The C O D A M material mode l predicts the peak stresses quite w e l l over the four different load ing direct ions. The slight over-predict ion o f the peak stress values i n the loca l x-di rect ion c o u l d be due to deficiencies o f the effective strain funct ion presented i n Equa t ion 4.2 at capturing the complete response. F o r example , due to the h i g h major Po isson ' s ratio o f this braided mater ial (reported to be near 1.1), h i g h transverse tensile strains cou ld develop and contribute to the softening o f the ax ia l response. T h i s was not captured by the effective 116 Chapter 4 Calibration and Validation of Compressive Model strain function. Additionally, inadequate modelling of through-thickness delamination could be the reason for this over-prediction. 200 • Experimental (UB and LB) -^ —Numerical (UBandLB) ® Beard (2001) 1 1 I x 30° from x 60° from x Loading Direction y Figure 4.13 Comparison of upper and lower bounded predictions of compressive peak stress with experimentally measured values. 4.6 Strain Reversal Simulation Results In addition to the multiple tensile and compressive tests, tension-after-compression experimental tests were carried out in the local x-direction using the compressive test geometry. These tests quantified the amount of damage induced by large compressive strains. Uniaxial compression was first applied to a maximum machine displacement of 1.27 mm, and then the load was reversed and the specimen was pulled in tension to failure. Figure 4.14 shows the typical stress-displacement response and failure characteristics of the specimen. Note that the displacements reported were those of the machine, which inevitably contain errors and typically over-estimate the actual specimen displacement. The panels were able to support quite a substantial tensile load upon the reversal of strain. These results indicate that the damage induced by the high compression strains did not lead to complete fracture of the fibres. 117 Chapter 4 Calibration and Validation of Compressive Model -0.06 -0.04 -0X12 0 0.02 0.W 0.06 Displacement (In.) Figure 4.14 Stress-displacement profile of specimen first compressed before being pulled in tension to failure. Panel still has significant strength in tension. Note that the reported displacement is measured from the machine. (DeTeresa et al. 2001) To measure the effectiveness of the CODAM material model at capturing these strain reversal characteristics, panel simulations (with geometry shown in Figure 4.11) were subjected to a ramped compressive displacement and then unloaded and subsequently loaded in tension to failure. A maximum compressive displacement of 1.27 mm was applied in the simulation (thus ignoring the discrepancy between the machine and actual displacement). Figure 4.15 displays the predicted nominal stress-average strain profiles. Failure occurred along a row of elements near the middle of the specimen. The simulation accurately captured the peak compressive stress, as well as the unloading behaviour and the subsequent tensile strength. Additionally, the overall energy absorbed by the load cycle (as indicated by the area within the stress-strain profile) is accurately captured. The simulation does seem to predict a more rapid drop to the plateau stress and perhaps over-estimates the remaining tensile strength slightly, but these details are not as important as capturing the overall general behaviour and energy dissipation characteristics. 118 Chapter 4 Calibration and Validation of Compressive Model 100 -200 A : , , : ! : , 1 -0.06 -0.04 -0.02 0.00 0.02 0.04 Average Strain (mm/mm) Figure 4.15 Numerical prediction of the nominal stress-average strain response during a strain reversal (compression-unload-tension). The peak stresses, unloading behaviour, and energy dissipation are reasonably predicted. 4.7 Summary The CODAM parameters required to define the complete orthotropic response of the braided material were estimated from the constituent properties, observed failure behaviour and braid architecture. The large variation in the quality and measured material properties of the panel specimens warranted the use of upper and lower bound parameters. Tensile, compressive and strain reversal simulations using CODAM were compared to the experimental results. The predicted peak stresses in the two local directions agreed well with the experimentally measured values. Additionally, the material model successfully captured 30° and 60° off-axis strength. These simulations validated the compressive model incorporated into CODAM as well as determined the upper and lower bound parameters to be used in the simulations of the axial crushing of the braided tubes, as presented in the next chapter. 119 Chapter 5 Tube Crushing Simulations C H A P T E R 5 T U B E CRUSHING SIMULATIONS 5.1 Introduction This chapter presents the results of the numerical study of the axial crushing of braided composite tubes. The first section reviews the overall test geometry and loading conditions of the experimental tube crushing program carried out by the Automotive Composites Consortium (ACC) on the braided composite tubes of interest in this study. The energy balance of the tubes is also reviewed in this section. The subsequent sections deal with the results of tube simulations with and without plug initiators. Numerical and experimental results are compared in terms of the force-displacement profile, peak force, general failure behaviour, and Specific Energy Absorption (SEA). Parametric studies help to understand the energy absorbing mechanisms and provide some insight into designing such structures to efficiently absorb the energy in front end low-velocity vehicle impacts. 5.2 Experimental Details and Energy Balance of Tubes The braided composite tubes tested in the experimental study were approximately 360 mm in length and had nominal outside dimensions of 55 x 55 mm. The wall thickness varied in all the tubes, but nominal thickness values were reported to be 2.3, 4.9 and 8.1 mm for the single ply, double ply, and four ply tubes, respectively. The tubes had a braided [0/±30] lay-up consisting of carbon Fortafil fibre tows in a vinyl ester matrix. This is the same material used in the panel specimens studied in Chapter 4. The density of the material was calculated from the reported specimen weight and geometry, and was found to be roughly equal to 1.23, 1.25, and 1.36 g/cm , for the single ply, double ply, and four ply tubes respectively. This agrees well with the predicted density of approximately 1.4 g/cm3, as calculated from the constituent densities and volume fractions. The leading edges of all the tubes were chamfered at a 45° angle. In tubes tested with a plug initiator, a metal plug was inserted into the bottom of the tube prior to the weight drop. The specimens were attached to a weight of 140 kg and dropped from a height of 2.54 m. This drop height and weight resulted in an 120 Chapter 5 Tube Crushing Simulations initial velocity just prior to impact of about 7.1 m/s. The force-displacement profile was recorded for each of the tests. Figure 5.1 and Figure 5.2 show the general failure characteristics of tubes crushed without and with a plug initiator by ACC. The tubes split into very indistinct fronds with tearing concentrated at the corners. Transverse tensile fracture was observed between and within the tows. Away from the corners, the main failure mechanism was observed to be delamination (in multiply lay-ups) followed by bending of the separated material either inwards or outwards. (a) (b) (c) Figure 5.1 Typical failure characteristics of braided tubes without a plug initiator from the front (a and b) and from the bottom (c). A significant amount of loose material gathers within the tube (DeTeresa et al, 2001). 121 Chapter 5 Tube Crushing Simulations Figure 5.2 Typical failure characteristics of braided tubes crushed with a plug initiator. Tube tears at the corner and the four fronds bend away from the plug (Xiao, 2003c). Energy absorption is an important aspect of all the tube crushing tests. The law of conservation of energy states that the total energy within a system must be conserved; energy cannot be created or destroyed. In the case of tube crushing, the total energy of the system at the initiation of crushing is in the form of kinetic and potential. As crushing proceeds, this energy is absorbed or dissipated through various processes, such as material yielding, fracture, and friction, either internally within the damaged material or externally between the tube and the plug. Additionally, heat may be generated and lost from the system due to friction. This energy loss is assumed to be negligible. Upon termination, the total energy of the system is in the form of kinetic energy (final velocity at termination), absorbed energy and energy lost to heat generation. This energy balance can be represented as shown in the equation below. 2 EP„< = mglsh Eahs = energy absorbed via material straining or fracture, and energy dissipated via internal or external friction Eheat = energy lost due to heat generation 5.1 where: 122 Chapter 5 Tube Crushing Simulations Assuming that the mass (w=140 kg) and initial velocity (v = 7.1 m/s) are the same for all tubes and using the measured crush distance (Ah), the total energy input into the system is calculated and summarized in Table 5.1. The energy absorbed by the system is also presented in the table, as estimated by integrating the measured force-displacement profile. Any remaining energy in the system is attributed to a final velocity or to energy lost via heat generation. The energy absorbed by the system can be normalized with respect to the mass of the material crushed to obtain the specific energy absorbed (SEA). This provides a useful measure of the energy absorption efficiency of the system. Table 5.1 shows that the SEA increases as the number of plies increases. Additionally, in the single ply tubes, the SEA is about 50% higher in tubes crushed without the initiator plug. The table also summarizes the experimentally measured peak forces and final displacements. Table 5.1 Summary of experimentally measured tube crushing properties Tube Total initial Peak force Crush distance (mm) Energy absorbed (J) Residual SEA configuration energy (J) (kN) energy (J) (J/g) lply F[0/±30] with plug 3858 16.7 240 2153 1705 15 lply F[0/±30] 3880 50.7 256 3390 490 22 without plug 3882 38.7 258 3720 162 24 2ply F[0/±30] 3632 74.1 75 2617 1015 28 without plug 3641 75.7 82 2808 833 28 4ply F[0/±30] without plug 3588 107.8 43 2724 864 30 5.3 General Simulation Details The simulations consisted of 3 main components: a drop weight, a tube and either a plug initiator or a rigid plate as shown in Figure 5.3. In all of the numerical simulations, the drop weight was modelled at a rigid body with a mass of 140 kg. The tube was modelled with the CODAM material model and was rigidly attached to the bottom of the drop weight. The plug initiator was modelled as a rigid body and was fixed in space. An initial velocity of 7.1 123 Chapter 5 Tube Crushing Simulations m/s was assigned to the drop weight causing the tube to impact the plug (or in case of simulation without the plug initiator, the rigid plate). All of the simulations used the system of consistent units shown in Table 5.2. Drop weight modeled as a rigid body with mass = 140 kg Local x-direction of braided material (aligned with the longitudinal axis of the tube) Composite tube modeled with either 5 mm x 5mm or 2.5 mm x 2.5 mm shell elements (5 mm mesh discretization shown) Global coordinate system Local y-direction of braided material (aligned at right angles to the longitudinal axis of the tube) Plug initiator modeled as rigid body Figure 5.3 Rigid plate used in simulations without a plug initiator Schematic of the basic components in the tube simulation. The tube contacts either the plug or the rigid plate. Table 5.2 Consistent System of Units Mass Length Time Force Stress Energy g mm ms N (kg-m/s2) MPa mJ 124 Chapter 5 Tube Crushing Simulations To ensure the mesh refinement of the tube was adequate enough to capture the response, separate elastic simulations on tubes discretized with either 5 mm x 5 mm elements or 2.5 mm x 2.5 mm elements were compared. In both cases, the global system elastic response and the local elemental response were identical, indicating that both of these refinements are adequate. Elements larger than 5 mm were deemed too large to accurately represent the curved corners of the tube. As the tubes were made of the same material as the panel specimens discussed in Chapter 4, the CODAM material parameters outlined in Table 4.3 for the local x and v-directions were used in all the tube crushing simulations. The local x-direction is aligned with the longitudinal axis of the tube and the local v-direction is aligned at right angles to this axis as shown in Figure 5.3. Both lower and upper bound parameters were used in the tube simulations. The rigid components (drop weight and plug initiator) were modeled with solid elements and the tubes were discretized with shell elements. Comparisons were done using under-integrated (Ul) and fully-integrated (FI) shell elements. FI cases represent a benchmark to which Ul cases are measured against. Ul elements dramatically decrease the computational burden, often decreasing the run time by up to 70%. This time saving characteristic is the main reason for the popularity of Ul elements in explicit finite element codes. However, because these elements have only one integration point, hourglassing into non-physical, zero-energy shapes can occur. To prevent these deformation modes, hourglassing (HG) control is applied. Two main types of HG control are available in LS-DYNA: stiffness based HG control and viscosity based HG control. Stiffness based HG control is generally more effective for structural parts and is thus used in all of the tube crushing simulations. The energy required to prevent HG, as a rule of thumb, should be less than 10% of the internal energy of the part. This requirement is monitored in all simulations. The option to erode elements, thus removing them completely from the analysis, is included in CODAM. Erosion is helpful when a strain softening material model is defined because degraded elements can distort easily, leading to high hourglass energy requirements and visually disturbing results. An erosion criterion can remove these low-strength elements before they become a burden on the simulation. In all the simulations, erosion was based on 125 Chapter 5 Tube Crushing Simulations a total strain criterion, which erodes the element once the sum of the strain component reaches or surpasses a defined limit (Erode), as shown in Equation 5.2. In the tube simulations, Erode was set to a high value of 300%, which ensures that the element is completely damaged upon removal. Erode = sx+£y+£z+£xy+£yz+£xz 5.2 Modelling contact between the tube and the plug (or rigid wall), as well as self-contact of the tube, is a crucial aspect in accurately capturing the response of the tube crushing simulations. LS-DYNA offers a variety of contact types. In general, the penetration of a "slave" node or surface into a "master" surface is checked at each time step. If a penetration is found, a force is applied to eliminate the penetration. In contact situations involving rigid bodies it is important that the mesh refinement of the rigid body be as fine, or finer, than the deformable body coming into contact with it. The plug initiator is thus defined with fairly small elements, ranging from 1 to 2 mm. In crash analyses, it is important to use contact definitions capable of handling unpredictable and complicated contact interfaces. An Automatic Single Surface contact type (type 13) is used in all of the simulations. This contact definition is the most popular in LS-DYNA for crash analyses due to its reliability and accuracy (Bala, 2001a and Bala, 2001b). Contact between all of the parts defined in the slave set, including self-contact is monitored. Additionally, contact on either side of shell elements is considered. A few general comments follow on some of the contact parameters used in the simulations. • Static and dynamic coefficients were set equal to a value ranging from 0.22 to 0.35. These values were based on experimental frictional tests done on panels made of the same material (DeTeresa et al., 2001). Setting the static and dynamic coefficients equal to each other can help to decrease the numerical noise. • Viscous damping in the contact definition is set to values between 20 and 40 (corresponding to 20-40%) of critical damping) to damp out oscillation normal to the contact surfaces typical in crash or impact simulations. • Contact forces in LS-DYNA are determined using linear contact springs between the slave nodes and nearest master segments. Both penalty-based and soft-constraint-126 Chapter 5 Tube Crushing Simulations based contact approaches are available to define the stiffness of these springs. The penalty-based approach uses the size of the contact segment and its material properties to determine the contact forces, whereas the soft-constraint-based approach uses the nodal masses that come into contact and the global time step size to determine the contact forces. The second approach is more effective for contact between materials having dissimilar material properties. As the composite tube has a much smaller through-thickness modulus (~8 GPa) than the metallic plug (£'«200 GPa), this soft-constraint-based approach was used. In addition, a segment-based contact algorithm, which treats contact between two segments (instead of the typical nodes-to-segments), was used to produce a more realistic distribution of the contact forces. • All other parameters in the contact definition were left at the default values. The success of the numerical model at capturing the observed tube crushing response was based on the comparisons of the force-displacement profile, peak force, overall failure behaviour and final predicted Specific Energy Absorption (SEA). Two methods were used to obtain the forces that develop longitudinally as the crushing progressed in the simulation. The first tracks the acceleration of the top mass and uses Newton's Second Law of Motion (F = ma) to calculate the longitudinal force. The second method uses the longitudinal forces that develop at the contact interface. In all the simulations, these two methods resulted in the same force-displacement profile. The energy development was numerically monitored in LS-DYNA using many global and material based energy outputs. Time histories of various forms of energy (kinetic energy, strain energy, total energy) were monitored as a consistency check on energy calculations. 127 Chapter 5 Tube Crushing Simulations 5.4 Simulations with Plug Initiator 5.4.1 Single Ply 5.4.1.1 Elastic Orthotropic Analysis An initial elastic simulation with orthotropic material properties illustrates how strains develop and eventually cause damage in the tube. Figure 5.4 shows the distribution of the longitudinal and circumferential strains along the length of the tube under an applied vertical displacement of 7 mm. The longitudinal strains are relatively constant along the length with a slight increase at the bottom. The circumferential strains are also relatively constant along the length until approximately 20 mm from the bottom, where the strains dramatically increase. The results of this analysis indicate that damage will likely initiate at the bottom of the tube due to the high strains at this location. Applied displacement = 7 mm (-0.02 mm/mm) 0.049 (a) (b) (c) Figure 5.4 (a) Displacement applied to the orthotropic tube, (b) Distribution of longitudinal and (c) circumferential strains along the length of the tube due to the applied displacement. Figure 5.5 illustrates the relative variation of the longitudinal, circumferential and shear strains predicted by the elastic analysis at the bottom of the tube. The elements in the tube 128 Chapter 5 Tube Crushing Simulations are under a biaxial state of stress, characterized by compressive longitudinal strains and tensile circumferential strains. The shear strains are negligible compared to either the longitudinal or circumferential strains. The higher longitudinal and circumferential strains at the corners of the tube suggest that damage is likely to initiate at the corners. Figure 5.5 Variation of (a) longitudinal and (b) circumferential strains at the bottom of the tube. Terms in brackets refer to strains normalized with respect to the applied global strain of-0.02. 5.4.1.2 Predicted Failure Behaviour and Energy Absorption This section presents the results of the numerical crushing simulation on the single ply braided tube crushed with a plug initiator. The tube was modeled with 5 mm x 5 mm shell elements and the material was defined by the upper and lower bound parameters discussed in Chapter 4. The resulting force-displacement, in general reached a peak force followed by an oscillations about a force slightly smaller than the peak, as shown in Figure 5.6 for the upper and lower bound parameters. The force-displacement profile from the experimental tube crush is also shown on this figure. Compared to the experimental results, the numerical force-displacement matches well, but tends to show a higher frequency of oscillation, which is likely caused by the progressive failure of the elements that discretize the tube. The numerical simulation predicts damage localization in the corners of the tube, ultimately leading to progressive tearing at the four corners. Eventually, the tube separates into four distinct fronds as shown in Figure 5.7. The numerically predicted bounds on the peak force 129 Chapter 5 Tube Crushing Simulations (11.8 to 17.9 kN) straddle the experimentally measured value of 16.7 kN. Estimates of the total energy absorbed, obtained by integrating the force-displacement profile, vary between 1430 and 2850 J, comparing well to the experimental result of 2153 J. Dividing this total energy by the mass of the crushed tube (calculated as the product of the mass/unit length of tube and the crush distance) leads to estimations of SEA between 10 and 20 J/g. Figure 5.8 compares these predicted peak force and SEA values to experimentally measured values. Displacement (mm) Figure 5.6 Comparison ofpredicted upper and lower bound force-displacement profdes to the experimentally measured trace for single ply tube with plug. Experimental results provided by (Xiao, 2003a). 130 Chapter 5 Tube Crushing Simulations Figure 5.7 General failure behaviour predicted for single ply braided tubes crushed with a plug initiator. Figure 5.8 Comparison ofpredicted (a) peak forces and (b) Specific Energy Absorption (SEA) to experimental values. Error bars refer to upper and lower bound values. The var ia t ion o f the longi tud ina l , circumferential and shear strains around the perimeter o f the tube at a arbitrary crush distance o f 165 m m is shown i n F igure 5.9, F igure 5.10, and Figure 5.11. D u e to the self -s imilar progression o f c rushing a long the length o f the tube, these strain distr ibutions are representative o f any cross-sect ion as it approaches the crush 131 Chapter 5 Tube Crushing Simulations zone. Figure 5.12 represents the relative strain distribution of this cross-section. From this figure, it is immediately evident that in-plane shear strains are insignificant. The longitudinal strains are slightly concentrated at the sides and the circumferential strains are severely localized to the corners. Fringe Levels -4.470e-08 Figure 5.9 Fringe plot of longitudinal strains in single ply tube crushed with plug. Absolute minimum = -0.001, absolute maximum = -0.70. Figure 5.10 Fringe plot of circumferential strains in single ply tube crushed with plug. Absolute minimum = 0.0, absolute maximum = 0.74. 132 Chapter 5 Tube Crushing Simulations Figure 5.11 Fringe plot of in-plane shear strains in single ply tube crushed with plug. Absolute minimum = -0.15, absolute maximum = 0.026. Max • -0.70 Max = 0.74 Min = 0.0 -0.64 Max = 0.0: (b) (c) Figure 5.12 Variation of (a) longitudinal strains and (b) circumferential strain at a cross section 165 mm from the bottom of the tube. The profiles of various forms of energy provided from LS-DYNA are shown for the lower bound simulation in Figure 5.13. The total energy in the system remains relatively constant, as the kinetic energy decreases and the internal and factional energy increases. Figure 5.13(b) shows various forms of energy absorption in the system. The sliding energy is the energy absorbed at the contact interface via friction. The material internal energy is the energy absorbed by the composite tube, but subtracts the energy absorbed by elements once they erode. The global internal energy represents the energy absorbed by the composite tube, and does not subtract the energy absorbed by elements once they fail. The energy indicated 133 Chapter 5 Tube Crushing Simulations by the integral of the force-displacement profile, which represents the total energy absorbed by the system is also shown in Figure 5.13(b). Note that the sum of the global internal energy and the sliding energy approximately equals the energy indicated by the force-displacement profile. For this lower bound case, approximately 27% of the total energy was absorbed via friction and 73%) was absorbed via material straining and damage growth. The energy required to control hourglassing was negligible (typically less than 1 % of the internal energy). 4 0 0 0 3 5 0 0 3 0 0 0 2 5 0 0 2 0 0 0 1 5 0 0 1 0 0 0 5 0 0 0 Total energys Global kinetic energy' Global internal energy ^ Global sliding (frictional) energy.. D i s p l a c e m e n t ( m m ) (a) 1 6 0 0 1 4 0 0 -1 2 0 0 -_ 1 0 0 0 | 2 8 0 0 " 6 0 0 4 0 0 2 0 0 0 Integral of force-displacement „ Global internal energy . Material internal energy . Global sliding (frictional) energy 8 0 1 2 0 1 6 0 D i s p l a c e m e n t ( m m ) (b) Figure 5.13 Growth and decay of various (a) global and (b) dissipative energy measures predicted by LS-DYNA as the tube crushes with lower bound parameters. One method of measuring the energy absorbing efficiency of the tube is to evaluate the amount of damage induced during crushing. This can be accomplished via the elemental damage parameters, cox and coy, which quantify the amount of damage in the elements in each of the local directions. Damage parameters range from zero (indicating no damage) to one (indicating complete damage). The most efficient tube crush would fully damage all of the elements in both local directions (i.e. a>x = co =Y). Figure 5.14 shows the fringe plots of the damage parameters at the initial stages of the simulation in the local x-direction and the local v-direction. All of the elements are equally damaged longitudinally, but local y-direction damage is localized to the corner elements. 134 Chapter 5 Tube Crushing Simulations (a) fr inge l e v e l s 1.0000*00 9.0006 in BJMOa-01 7.000B-01 t, nolle id 5 MOe-01 •1.0006 01 3.000e 01 Z.OOOe-Ol 1.0006-01 0.<HH3e*0B 1 (b) Figure 5.14 Fringe plot of the (a) local x-direction and (b) local y-direction damage parameters during initial stages of tube crushing with the plug initiator. Monitoring these damage parameters in the elements defining the tube can provide some insight into the mechanisms responsible for dissipating the crash energy. The values of these damage parameters vary depending on the location within the cross-section (i.e. local y-direction damage is more severe at the corners than in the sides). However, the damage distribution will not vary significantly along the length of the tube, because crushing proceeds in a self-similar manner. Thus, given the distribution of the damage parameters within a cross-section, the overall proportion of energy absorbed via damage growth in the two local directions can be approximated. The damage parameters in eighteen elements (ten elements from the side, eight elements from the corner) located about 75 mm from the bottom of the tube were monitored to estimate the amount of damage typical of a given cross-section. The average value of the local x-direction and ^-direction damage parameters of the elements are presented in Table 5.3. The corner elements are, on average, damaged significantly in both local directions, whereas the side elements are generally highly damaged in the local x-direction but only moderately damaged in the local _v-direction. These results indicate that longitudinal damage growth in the sides of the tube and local y-direction damage growth in the corners of the tube are important energy dissipation mechanisms. 135 Chapter 5 Tube Crushing Simulations Table 5.3 Average damage parameters of side and corner elements. Location of Elements (longitudinal) coy (circumferential) Side 0.89 0.22 Corner 0.50 0.78 5.4.1.3 Parametric Study Several factors affect the overall behaviour of the tube crushing experiments. Among the most important are the material properties and the structural geometric details. This section presents the results of a parametric study that determines how some of these properties affect the overall behaviour of the tube crush. The sensitivity of the tube crush to changes in the material properties was evaluated by changing the fracture energy associated with failure in each of the local directions. The effect of changing the magnitude of the plateau stress was also investigated. Additionally, the effect of changing the frictional coefficient between the tube and the plug was studied. The sensitivity of the simulations to these changes was evaluated in terms of the percent change in the predicted Specific Energy Absorption. As longitudinal compression (local x-direction) and local ^ -direction tension are considered to be the dominant failure mechanisms, only the fracture energies in these local directions were manipulated. The fracture energy was increased or decreased by changing the fibre and matrix saturation points (Fms and Ffs). Figure 5.15 demonstrates how the SEA did not vary significantly with changes in the local v-direction tensile fracture energy. Doubling of the fracture energy resulted in only a 6% increase in the predicted SEA. 136 Chapter 5 Tube Crushing Simulations The SEA was similarly only loosely dependent on the compressive fracture energy in the local x-direction as shown in Figure 5.16. A 50% increase in the compressive fracture energy resulted in only a 6% increase in the SEA. In both cases, increasing the fracture energies did not result in a significant change in the failure behaviour of the tubes. Failure was dominated by tensile failure of the corner elements and compressive failure of the side elements. The damage induced into elements that did not fail remained relatively constant, as indicated by average damage variables that remained close to those reported in Table 5.3. Thus, increasing the fracture energy associated with local v-direction fracture did not promote more elements to damage, but rather increased the fracture energy associated with the failure of the corner elements. 137 Chapter 5 Tube Crushing Simulations • G\ =20 kJ/m2 20 i HI GCF decreased 25% g QC increased 25% ^ QC increased 50% 15 < UJ 0 3 5 o J. i i i 11111 m 11111 ii — i w w w v i , Figure 5.16 Variation of the SEA with compressive fracture energies (local x-direction) for single ply with plug tube simulations. At first, this rather insensitivity of the SEA to increases in the fracture energy seems counterintuitive. However, a simple analysis that follows indicates that these predictions are reasonable. Given that the total energy absorbed by the tube, as predicted by upper bound parameters, is approximately 2750 J , a few approximations can be made to divide this energy into components. These components are the energy absorbed by the complete failure of the elements, energy absorbed by the partial failure of the elements and energy absorbed via friction. By estimating the number of elements that completely fail and knowing the fracture energy attributed to this failure, one can approximate the energy absorbed via complete failure of the elements. The following assumptions are made with respect to the number of elements that completely fail in either tension or compression, noting that there are a total of 12 corner elements and 24 side elements per row and that these elements can absorb energy by partial failure in other directions in addition to the energy absorbed by complete failure in one direction: • Corner elements (per row): 8 fail in tension, 6 fail in compression. • Side elements (per row): 6 fail in tension, 18 fail in compression. The total crush distance of approximately 240 mm means crushing involves 48 rows of elements. The energy attributed to complete tensile and compressive failure in both the 138 Chapter 5 Tube Crushing Simulations corner and side elements can now be calculated as the product of the number of elements that completely fail, the elemental cross-sectional area, and the fracture energy as shown in Table 5.4. Table 5 . 4 Predicted variation of energy absorbed by the complete failure of corner and side elements for different fracture energies. Descr ip t ion Local y G", (kJ/m 2) Local x G/ C (kJ/m 2) Energy a )sorbed via complete fracture Corner -tension (J) Corner -compression (J) S i d e -tension (J) S ide -compression (J) Base l ine 42 20 190 410 140 1240 GTf increased 5 0 % 63 20 280 410 210 1240 Gcf increased 5 0 % 42 30 190 450 140 1340 The frictional energy reported by LS-DYNA was consistently equal to approximately 600 J. Knowing the total energy (2750 J), energy absorbed by friction (600 J), and energy absorbed by complete failure (1980 J), the energy absorbed by the partial failure of elements can be calculated (as the difference between the total energy and energy absorbed via either friction or failure). The baseline prediction of the energy absorbed by the partial failure of elements is equal to approximately 200 J. This equates to an average of 0.12 J of energy absorbed per element or 10 kJ/m , which seems reasonable as it is less than either the tensile or compressive fracture energies. Based on the fact that the average damage state of the elements remained approximately constant between runs with different fracture energies, the energy absorbed by this partial failure is assumed to also remain constant at 200 J. Figure 5.17 shows the relative distribution of the various components of energy for the baseline (no modifications to the fracture energies) tube simulation with upper bound parameters. The majority of the energy is absorbed by the complete failure of elements. The predicted increases in the energy absorbed by increasing the fracture energy, as shown in Table 5.4 can now be used to approximate how these modifications affect the SEA. Increasing the local y-direction tensile fracture energy by 50% causes an overall increase in the energy absorbed by about 160 J, corresponding to a 6%> increase in the SEA. Increasing the local x-direction compressive fracture energy by 50% causes an overall increase in the energy absorbed by 139 Chapter 5 Tube Crushing Simulations about 140 J (5% increase in the SEA). These percent increases in the SEA, as estimated using a simple analysis, agree well with the simulation results of a 6% increase in the SEA for both cases. 3000 Total energy absorbed Energy absorbed by partial failure of elements Energy absorbed by friction Energy absorbed by complete — < C failure of elements Energy absorbed by complete failure 1240 J4Q_ Corner compression Corner tension Side compression Side tension Figure 5.17 Assumed distribution of the energy absorbed in the tube simulation with upper bound parameters, showing the components of energy absorbed by complete failure of the elements. In addition to the investigation on the effect of the fracture energy, the sensitivity of the SEA to changes in the plateau stress, translating to an increase in the band broadening energy, was numerically investigated. The plateau stress was increased from 60 MPa to either 75 MPa (25% increase) or 90 MPa (50% increase), and the resulting SEA predicted by the simulation was compared to the baseline prediction. Once again, these modifications resulted in insignificant changes in the predicted SEA values as shown in Figure 5.18, with only a 15% increase for a 50%) increase in the plateau stress. The average damage parameters in the side and corner elements remained close to those reported in Table 5.3, indicating that changes in the plateau stress did not translate into changes in the global response of the tube. Rather, they simply increased the amount of energy absorbed locally by failing elements. The same simple energy analysis presented previously can be used to understand the insignificant increase in the SEA. Increasing the plateau stress increases the band broadening energy absorbed by the compressive failure of the corner elements by 115 J and side elements by 140 Chapter 5 Tube Crushing Simulations 350 J. This additional energy absorption translates into an overall increase in the SEA of 17%, which is very close to the 15%> increase predicted by the numerical simulation. = 60 MPa increased 25% increased 50% Figure 5.18 Variation of SEA with changes in the longitudinal plateau stress for single ply with plug tube simulations. The final parameter to be investigated with respect to the sensitivity of the SEA is the coefficient of friction between the plug initiator and the tube. Understandably, increasing the coefficient will lead to a direct increase in the energy absorbed via frictional dissipation. However, increasing the resistance to sliding also increased the longitudinal compressive damage induced in the corner elements, as indicated by a slight increase in the average local x-direction damage parameter from cox = 0.50 to about cox = 0.95 for a 50%> increase in the coefficient of friction. All other average damage parameters were relatively unaffected. Figure 5.19 shows the relative effect that increasing the coefficient of friction has on the SEA. Numerical simulations predicted a 15%> increase in the SEA for a 50%> increase in the coefficient of friction. The increase in the frictional dissipation was directly measured from the numerical simulation to be approximately 125 J and 300 J for increases in the coefficient of friction of 25% and 50%, respectively. Increasing the number of corner elements that fail in compression from 6 to 10 per row accounts for the reported increase in the average damage variable. This leads to an estimated increase of the energy absorbed by the corner elements from 410 J to 680 J. The predicted increase in the SEA using the simple analysis is 141 Chapter 5 Tube Crushing Simulations 20% (for a 50%) increase in the coefficient of friction), which is close to the 15%) increase predicted by the numerical simulation. 25 20 H 2 5 , wa 10 H • n = 0.22 5 H increased 25% EI ti increased 50% Figure 5.19 Variation of SEA with changes in the coefficient offriction for single ply with plug tube simulations. This parametric study indicates that that energy absorbed by the crushing of the tube is relatively insensitive to variations in individual parameters, such as the fracture energy, band broadening energy, or the coefficient of friction. The SEA was slightly more sensitive to changes in the plateau stress (band broadening energy) or coefficient of friction. In all cases, the global failure of the tube was independent of the parametric modifications studied above; main failure characteristics were local v-direction tensile failure at the corners combined with local x-direction compression failure in the elements around the perimeter of the tube. Additionally, except for the coefficient of friction, the parametric modifications did not significantly affect the average amount of damage induced into either the side or corner elements. Analytically dividing the total absorbed energy into components and predicting how each component was affected by the parametric changes helped to describe and understand the numerical predictions. 5.4.1.4 Mesh Insensitivity To illustrate that the model predictions are independent of the mesh size used to discretize the tube, predictions using two element sizes were compared. The upper bound prediction 142 Chapter 5 Tube Crushing Simulations from the tube discretized with 5 mm x 5 mm elements is compared to the response predicted using elements half this size. The only material parameters that differ in these two cases are the scaling factors, k, used to maintain the fracture energy associated with the characteristic height of damage. Figure 5.20 demonstrates that the overall predicted failure behaviour was similar with each of the mesh densities. Additionally, Figure 5.21 compares the individual force-displacement profiles, which aside from the difference in the frequency of the oscillations, fluctuate between similar minimum and maximum values with comparable average force levels. The higher number of oscillation in the simulation modeled with 2.5 mm elements is a direct result of the smaller element size. The drops coincide with element erosion. This discrepancy is not considered to be as important as capturing the overall response in terms of failure behaviour, peak force and SEA. The comparison between the peak force and SEA predictions, as shown in Figure 5.22, demonstrate that there is good agreement between the two simulations. Figure 5.20 Comparison of the overall failure behaviour predicted by different mesh densities for single ply with plug tube simulations. 143 Chapter 5 Tube Crushing Simulations 20 i 0 -P 1 1 — — i — — i — i — i 0 40 80 120 160 200 240 Displacement (mm) Figure 5.21 Comparison of the force-displacement profdes predicted for different mesh densities. (a) (b) Figure 5.22 Comparison of (a) peak forces and (b) SEA values predicted for different mesh densities for single ply with plug tube simulations. 5.5 Simulations without Plug Initiator 5.5.1 Single Ply 5.5.1.1 Predicted Failure Behaviour and Energy Absorption The next challenge of this study was to accurately model the response of the braided tubes when crushed without a plug initiator. Experimentally, the crushing of single ply braided tubes without a plug initiator was about 50% more efficient at absorbing energy, and thus offers an obvious advantage over plug-initiated tube with respect to absorbing automotive crash energy. Initial numerical simulations using 5 mm x 5 mm shell elements to discretize 144 Chapter 5 Tube Crushing Simulations the tube exhibited severe hourglassing, leading to unstable and inaccurate results. For this reason, all tubes were modeled with 2.5 mm x 2.5 mm shell elements, which did not display any instability under reasonable hourglassing control. Two experimental force-displacement profiles are shown in Figure 5.23. The predicted upper and lower bound force-displacement profiles are compared to Experimental A results in Figure 5.24 and Figure 5.25. Figure 5.26 illustrates the predicted failure morphology at two different crushing distances. Failure was characterized by the compressive failure of the elements surrounding the perimeter. This failure initiated at the bottom and progressed lengthwise up the tube as the displacement increased. Unlike the experimental tube crushes and the LB numerical simulation, which proceeded to maximum displacements of 254 mm, the UB numerical simulation became unstable at about 150 mm. SEA calculations for the UB simulation will thus be calculated using a maximum crush distance of 150 mm. The predicted force-displacement profiles match the experimental profile quite well, with reasonable peak force and post-peak predictions. Integrating the numerical and experimental force-displacement profiles leads to similar estimates of total energy absorption (2300 J numerically (LB) and 3390 J to 3720 experimentally). Energy data files from LS-DYNA indicate that approximately 15% of this total energy was absorbed via friction and 85%> absorbed via damage accumulation in the material. The numerically predicted SEA and peak force match the experimental measured values fairly well as shown in Figure 5.27. The simulations without a plug initiator display some significant differences to the results from the plug-initiated simulations. Firstly, the overall failure behaviour is different. In the plug-initiated tubes, failure was mostly confined to the four corners of the tube, whereas failure was evenly distributed around the perimeter in the non-initiated simulations. Secondly, the predicted force-displacement profile has some distinct differences. The predicted peak force was up to 3 times higher and the post-peak force oscillated about a slightly higher average force level in the non-initiated case. These two aspects result in a SEA prediction that is up to 50%> higher for the non-initiated case. 145 Chapter 5 Tube Crushing Simulations 0 -f 1 1 1 1 1 1— 0 40 80 120 160 200 240 Displacement (mm) Figure 5.23 Comparison of two experimental force-displacement profiles for single ply tube without plug. Experimental results provided by (Xiao, 2003a). Displacement (mm) Figure 5.24 Comparison ofpredicted upper bound force-displacement profiles to the experimentally measured trace for single ply tube without plug. Experimental results provided by (Xiao, 2003a). 146 Chapter 5 Tube Crushing Simulations 40 80 120 160 Displacement (mm) Experimental A Numerical L B 200 240 Figure 5.25 Comparison ofpredicted lower bound force-displacement profiles to the experimentally measured trace for single ply tube without plug. Experimental results provided by (Xiao, 2003a). Elements erode along the perimeter due to high strains Global coordinate system (a) (b) (c) (d) Figure 5.26 General failure behaviour as predicted for single ply braided tubes crushed without a plug initiation showing morphology from the (a, c) side and (b, d) bottom for 15 mm and 130 mm crush distances. 147 Chapter 5 Tube Crushing Simulations 60 -, Figure 5.27 Comparison ofpredicted (a) peak forces and (b) Specific Energy Absorption (SEA) to experimental values. Error bars refer to upper and lower bound values. The predicted damage parameters in each of the local directions, as determined by the average of about 24 side and corner elements, are presented in Table 5.5. The damage parameters indicate the dominance of longitudinal damage growth over damage growth in the local ^ -direction, with almost 100% longitudinal damage growth in all the elements. In comparison to the plug-initiated case, the predicted longitudinal damage parameters in the corner elements and the local y-direction damage parameters in the side elements have increased, from 0.50 to 0.96 and 0.22 to 0.55, respectively. In effect, more energy is being absorbed via the growth of damage in the material. The higher SEA prediction in the non-initiated case is a direct consequence of the increase in these damage parameters. Table 5.5 Average damage parameters of side and corner elements. Location of Elements (longitudinal) (circumferential) Side 0.90 0.55 Corner 0.96 0.80 5.5.1.2 Parametric Study In the plug-initiated simulations, the SEA was observed to be rather insensitive to changes in the material parameters. This section examines the sensitivity of the non-initiated single ply simulation results to similar modifications. As indicated by the high longitudinal damage 148 Chapter 5 Tube Crushing Simulations parameters, the compressive material properties are the most significant factor in the overall SEA response. For this reason, only the sensitivity to the longitudinal compressive fracture energy and the plateau stress (or band broadening energy) is investigated. The UB simulation was used as a baseline for comparison and sensitivity was based on changes in the SEA values. As Figure 5.28 indicates, increasing the compressive fracture energy by 50% (from 20 kJ/m to 30 kJ/m2) results in an 18% increase in the predicted SEA. Increasing the plateau stress by 50%> (from 60 MPa to 90 MPa) resulting in a more dramatic 28% increase in the SEA prediction (Figure 5.29). In both cases, the modifications did not dramatically alter the overall failure behaviour or the average damage parameters presented in Table 5.5. These predicted increases in the SEA for the non-initiated single ply tubes is much more dramatic than the corresponding increase predicted in the initiated tube simulations (6% SEA increase for a 50%> increase in the compressive fracture energy and 15%> SEA increase for a 50%) increase in the plateau stress). This demonstrates the importance of the compressive properties in the non-initiated tube simulations. Figure 5.28 Variation of SEA with changes in the longitudinal plateau stress for single ply without plug tube simulations. 149 Chapter 5 Tube Crushing Simulations ESj ^plateau = 60 MPa E2 ® plateau increased 25% ID °'plateau increased 50% Figure 5.29 Variation of SEA with changes in the longitudinal plateau stress for single ply without plug tube simulations. 5.5.2 Double Ply 5.5.2.1 Predicted Failure Behaviour and Energy Absorption Double ply braided tubes were also crushed without a plug initiator. Experimental force-displacement profiles for two tubes are shown in Figure 5.30. Forces increased to maximum values of approximately 75 kN before dropping off to levels between 25 and 45 kN. Crushing proceeded to maximum displacements of approximately 80 mm. To model the tube in LS-DYNA, shell elements of 2.5 mm x 2.5 mm were again used due to their resistance to hourglassing. The numerically predicted force-displacement profiles for the upper and lower bound parameters match the experimental results well, as shown in Figure 5.31 and Figure 5.32. The predicted failure morphology was characterized by compressive crushing of the elements surrounding the perimeter, as shown in Figure 5.33. Integrating the numerical and experimental force-displacement profiles leads to similar estimates of total energy absorption (2395 to 3160 J numerically and 2617 to 2808 J experimentally). Approximately 13% of this total energy was absorbed via friction and 87%> absorbed via damage accumulation in the material. The numerically predicted SEA and peak forces match the experimental measured values well as shown in Figure 5.34. 150 Chapter 5 Tube Crushing Simulations Experimental A Experimental B 30 40 50 60 Displacement (nun) 90 Figure 5.30 Comparison of two experimental force-displacement profdes for double ply tube without plug. Experimental results provided by (Xiao, 2003a). Numerical Upper Bound Experimental A Numerical UB 10 20 30 40 50 60 Displacement (mm) 70 80 Figure 5.31 Comparison ofpredicted upper bound force-displacement profdes to the experimentally measured trace for double ply tube without plug. Experimental results provided by (Xiao, 2003a). 151 Chapter 5 Tube Crushing Simulations 10 Experimental 20 30 40 50 60 Displacement (mm) Experimental A Numerical L B Figure 5.32 Comparison ofpredicted lower bound force-displacement profdes to the experimentally measured trace for double ply tube without plug. Experimental results provided by (Xiao, 2003a). Elements erode along the perimeter due to high strains Global coordinate system (a) (b) Figure 5.33 General failure behaviour as predicted for double ply braided tubes crushed without a plug initiator showing (a) side and (b) bottom morphology. 152 Chapter 5 Tube Crushing Simulations Experimental Numerical Figure 5.34 Comparison ofpredicted (a) peak forces and (b) Specific Energy Absorption (SEA) to experimental values. Error bars refer to upper and lower bound values. Table 5.6 presents the average local damage parameters predicted by the numerical simulation. The damage parameters are very similar to those reported for the single ply simulation crushed without a plug initiator, indicating similar damage patterns. The numerical simulations indicate major damage growth longitudinally, caused by the high compressive strains in this direction. The local ^ -direction damage is somewhat localized to the corner elements, but does not play as significant of a role in absorbing the crash energy as the longitudinal damage growth. Table 5.6 Average damage parameters of side and corner elements. Location of Elements (longitudinal) (circumferential) Side 0.98 0.54 Corner 1.00 0.90 5.5.3 Four Ply 5.5.3.1 Predicted Failure Behaviour and Energy Absorption The last axial tube non-initiated crushes to be simulated using the CODAM material model are tubes constructed of four plies. As shown in Figure 5.35, the longitudinal force 153 Chapter 5 Tube Crushing Simulations experimentally progressed to a maximum value of 108 kN before dropping off to levels around 70 kN up to a final displacement of about 43 mm. The numerical model was constructed using 2.5 mm x 2.5 mm shell elements to discretize the tube. The numerically predicted force-displacement profiles for the upper and lower bound parameters match the experimental results well, as shown in Figure 5.35. The total energy absorbed is obtained by integrating the numerical and experimental force-displacement profiles. Numerically predicted results of 1993 to 2840 J agree well with the experimental result of 2724 J. As indicated by the energy output files in LS-DYNA, approximately 15% of this total energy was absorbed via friction and 85%> absorbed via damage accumulation in the material. Figure 5.36 displays the failure morphology predicted by the simulation, which was characterized by compressive crushing of the perimeter elements. The numerically predicted SEA and peak force match the experimental measured values nicely as shown in Figure 5.37. Displacement (mm) Figure 5.35 Comparison ofpredicted upper and lower bound force-displacement profiles to the experimentally measured trace for four ply tube without plug. Experimental results provided by (Xiao, 2003a). 154 Chapter 5 Tube Crushing Simulations Elements erode along the perimeter due to high strains Global coordinate system Figure 5.36 General failure behaviour as predicted for four ply braided tubes crushed without a plug initiation showing (a) side and (b) bottom morphology. o u u o U-« OH 140 -| 35 -j 120 - i 30 -100 -DJJ 25 -80 - 20 -60 - SEA 15 -40 - 10 20 - 5 0 •rl 0 -Experimenta 1 Numerical (a) Figure 5.37 Comparison ofpredicted (a) peak forces and (b) Specific Energy Absorption (SEA) to experimental values. Error bars refer to upper and lower bound values. The average local damage parameters predicted by the numerical simulation are shown in Table 5.7. The damage parameters are very similar to those reported for both the single ply and double ply simulation, indicating similar damage patterns. These damage parameters indicate that the majority of the energy is absorbed via longitudinal damage growth caused 155 Chapter 5 Tube Crushing Simulations by the applied compressive strains in this direction. The most effective method of increasing the SEA of these tube crushes would be to increase the energy absorbed via the axial crushing of the material in the local x-direction, as long as these material alterations do not influence the relative values of the damage parameters in each of the local directions. Table 5.7 Average damage parameters of side and corner elements. Location of Elements ®* (longitudinal) <°y (circumferential) Side 0.98 0.52 Corner 0.98 0.80 5.6 Summary This chapter presented the results of the numerical study on the axial crushing of braided composite tubes. The aim was to accurately predict the response of initiated and non-initiated tubes. Upper and lower bound CODAM material parameters were used to obtain a window of prediction due to the high variability of the material. The effectiveness of the simulations at capturing the overall response of the tube crushing was based on the predicted force-displacement profile, general failure behaviour and the Specific Energy Absorption. As shown in Figure 5.38 and Figure 5.39, the numerical upper and lower bound predictions captured the general trend in the peak force and SEA values accurately. 156 Chapter 5 Tube Crushing Simulations 140 120 3 100 CU o 80 u o ta 60 2 40 PM H Experimental 0 Numerical l p l y with plug l p l y without plug 2 ply without plug 4 ply without plug Figure 5.38 Comparison of measured and predicted peak force values of single ply, double ply and four ply tubes. Error bars refer to upper and lower bound values. 35 - i a 30 -o • ^ Q . 25 -O fAb 20 -li u 15 -a> e td 10 -W '3 5 -a 0 -I Experimental 0 Numerical l p l y with plug l p l y without plug 2 ply without plug 4 ply without plug Figure 5.39 Comparison of measured and predicted SEA values of single ply, double ply and four ply tubes. Error bars refer to upper and lower bound values. Failure was characterized by tearing at the four corners in conjunction with moderate compressive longitudinal damage in the elements around the perimeter in the plug-initiated case. In the non-initiated simulations, failure was characterized by severe compressive failure of the perimeter elements in conjunction with slight localization of local v-direction tensile strains at the corners. This difference in the overall failure behaviour between initiated and non-initiated cases was also depicted by the difference in the average damage parameters in the side and corner elements. 157 Chapter 5 Tube Crushing Simulations A parametric study on the initiated single ply tubes indicates that the response, as evaluated by the SEA, is not sensitive to changes in the material fracture energies. The response is only slightly dependent on the band broadening energy and the friction coefficient between the plug and the tube. The response of the non-initiated tubes was more sensitive to changes in the longitudinal compressive fracture energy and the band broadening energy. In summary, the CODAM material model was capable of capturing the general features and trends of the tube crushing experiments on braided material. 158 Chapter 6 Conclusions and Future Work C H A P T E R 6 CONCLUSIONS AND F U T U R E W O R K 6.1 Summary and Conclusions The main objective of this thesis was to enhance the capabilities of CODAM by incorporating a compressive damage model based on experimentally observed compressive failure mechanisms in composites. The main conclusions drawn from this investigation are summarized. A thorough review of the literature on composite tube crushing indicated that these structures could be designed to efficiently and safely absorb crash energy in a progressive manner. The design flexibility, environmental durability, and superior specific energy absorption (SEA) as compared to traditional materials are the most attractive features of composites for use in the automotive industry. The main disadvantages holding back the widespread use of composite materials in the automotive industry include uncertainties about performance, high material cost, and apprehension to design with composites. The manufacturing of low cost composite materials, such as woven or braided materials, is helping to resolve the issue of high cost. Continued experimental research on composite materials is helping to alleviate the uncertainties and apprehension to design with these materials. Compressive damage growth is an important mechanism responsible for absorbing the front-end crash energy. Research over the past forty years in the area of compressive composite failure has led to the general consensus that cooperative fibre kinking is the dominant failure mode in unidirectional composites (Sivashanker et al., 1996; Moran et al., 1995). In multidirectional composites kinking in aligned (or nearly aligned) plies is combined with matrix cracking (or yielding) in the off-axis plies and delamination between dissimilar plies to form a damage unit (Sivashanker, 2001). Additionally, research indicated that after initial formation, the kink band (or damage unit) broadened longitudinally with further applied displacement. Fibre waviness and matrix shear stiffness have been shown to greatly influence the initiation of kinking and the resulting compression strength. Due to the significant fibre waviness introduced during the braiding process, braided composites are 159 Chapter 6 Conclusions and Future Work extremely susceptible to premature kinking ultimately reducing the compressive strength (West and Adams, 1999). A previously developed plane-stress continuum damage mechanics based model (CODAM), which was suitable for simulating composite tensile damage growth and failure, has been extended to include compressive predictive capabilities. The improvements were based on the observed compressive failure mechanisms and behaviour. The model parameters, which include the elastic material properties, damage initiation and saturation in each of the constituents and the mechanical effect of damage growth are related to the observed composite behaviour and have physical meaning. The incorporated modifications were verified via single element simulations. The improved material model successfully predicted the tensile and compressive strengths, as well as the load reversal response of flat carbon/vinyl ester braided panels. The model was further evaluated based on its predictions of the damage propagation and energy absorption in the axial crushing of braided composite tubes. Tubes with and without plug initiator were numerically simulated and, in all cases, the damage propagation, force-displacement profile, and specific energy absorption trends agreed well with the experiments. In plug-initiated tubes the dominant global failure mechanisms were tearing at the corners due to high circumferential strains, combined with moderate longitudinal compressive damage along the entire perimeter. A parametric study indicated that the SEA of plug-initiated tubes was relatively insensitive to changes in the tensile and compressive fracture energies, band broadening energy, or the coefficients of friction between the plug and initiator. In tubes simulated without plug, the dominant global failure mechanisms were mild tearing at the corners, combined with significant longitudinal compressive damage along the entire perimeter. The SEA of the un-initiated tubes was approximately twice as sensitive as the plug-initiated cases to changes in the longitudinal compressive fracture energy and the band broadening energy. The improvements incorporated into CODAM result in a material model capable of representing the energy absorbed via significant damage growth in the axial crushing of the braided composite tubes studied in this thesis. This success suggests that CODAM could offer an attractive design aid for future incorporation of lightweight composite energy 160 Chapter 6 Conclusions and Future Work absorbers into crashworthy structures. Additionally, the success of this study improves the confidence of using CODAM as a predictive tool to simulate the effect of damage in composites subjected to either tensile or compressive loads. 6.2 Future Work Although the model has been shown to be successful at simulating the crushing response of braided composite tubes, simulation of other compressive dominated failures is needed to further validate and assess the robustness of the material model. With respect to the tube crushing, a thorough experimental investigation on the damaging behaviour of the braided composite material would help to improve the confidence in the model parameters and provide some insight into improving the energy absorbing capacity of these types of materials. A broader numerical sensitivity analysis, which considers other aspects such as geometry and trigger characteristics, is another area that needs addressing. Plug-initiated circular tubes should theoretically absorb more energy because circumferential damage is not localized to the corners. Additionally, the sensitivity of circular tubes to changes in the circumferential model parameters may be more significant. Issues regarding CODAM that remain to be addressed include an investigation into the three-dimensional interactions and effects these interactions have on damage growth. For example, transverse tensile strains could potentially influence the longitudinal compressive behaviour because of the reduction in fibre support leading to premature kinking as cracks grow at the fibre matrix interface due to the applied tensile strains. This interaction is captured in CODAM via the effective strain function. An investigation into the appropriate effective strain function constants is required to accurately define the damage surface. Additionally, the treatment of shear damage in CODAM needs to be addressed. The current linear degradation of the shear modulus as defined by the growth in the local damage variables may not be an accurate representation of the shear behaviour. Another deficiency in the model is the lack of permanent strain capabilities, which may be important during load reversals. Finally, a further investigation into physically determining the model parameters, 161 Chapter 6 Conclusions and Future Work especially the saturation points, which remain the hardest parameters to estimate, is necessary to maintain the physical basis and predictive capabilities of the model. 162 References R E F E R E N C E S Agaram, V., Bilkhu, S., Fidan, S., Botkin, M., and Johnson, N., (1997). "Simulation of Controlled Failure of Automotive Composite Structures with LS-DYNA3D", Proceedings of the 12th ESD/S AE Advanced Composites Conference, The United New Generation Vehicle Conference and Exposition, Detroit, MI, USA, 7-10 April 1997, pp. 181-190. Argon, A., (1972). "Fracture of Composites", Treatise on Materials Science and Technology, Vol. l,pp. 79-114. Arnaud, P. and Hamelin, P., (1998). "Dynamic Characterization of Structures: A Study of Energy Absorption in Composite Tubes", Composite Science and Technology, Vol. 58, No. 5, pp. 709-715. Bala, S., (2001a). "Contact Modelling in LS-DYNA - Some Recommendations - Part 1", FEA Information International News, August 2001, pp. 8-16. Bala, S., (2001b). "Contact Modelling in LS-DYNA - Part 2 - Contact Parameters - Default and Recommended Values", FEA Information International News, September 2001, pp. 6-16. Bazant, Z.P. and Planas, J., (1998)."Fracture and Size Effect". CRC Press LLC. Beard, S.J., (2001). "Energy Absorption of Braided Composite Tubes", Ph.D. Thesis, Department of Aeronautics and Astronautics, Stanford University. Beard, S.J., and Chang, F., (2002). "Energy Absorption of Braided Composite Tubes", International Journal of Crashworthiness, Vol. 7, No. 2, pp. 191-206. Botkin, M., Johnson, N., Zywicz, E., and Simunovic, S., (1998). "Crashworthiness Simulation of Composite Automotive Structures", Proceedings of the 13th Annual Engineering Society of Detroit Advanced Composites Technology Conference and Exposition, Detroit, MI, USA, 28-29 September 1998. Budiansky, B., (1982). "Micromechanics", Computers and Structures, Vol. 16, No. 1-4, pp. 3-12. Budiansky, B., and Fleck, N.A., (1993). "Compressive Failure of Fibre Composites", Journal of the Mechanics and Physics of Solids, Vol. 41, No. 1, pp. 183-211. Burr, S.T., (1996). "Characterization of Damage Mechanics and Behavior in Two-Dimensionally Braided Composites Subjected to Static and Fatigue Loading", Ph.D. Thesis, Department of Engineering Mechanics, Virginia Polytechnic and State University. 163 References Caliskan, A.G., (2002). "Design & Analysis of Composite Impact Structures for Formula One Using Explicit FEA Techniques", Proceedings of the 2002 SAE Motorsports Engineering Conference and Exhibition, Indianapolis, Indiana, 2-5 December, 2002. Castejon, L., Cuartero, J., Clemente, R., and Larrode, E., (1998). "Energy Absorption Capability of Composite Materials Applied to Automotive Crash Absorbers Design", from Polymer Composites and Polymeric Materials, Proceedings of the 1998 SAE International Congress & Exposition, Detroit, MI, USA, 23-26 February 1998, pp. 37-46. Castejon, L., Miravete, A., and Cuartero, J., (2001). "Analytical Formulation of (0°, ±ot°) Braided Composites and its Application in Crashworthiness Simulations", Mechanics of Composite Materials and Structures, Vol.8, No. 3, pp. 219-229. Chiu, C.H., Lu, C.K., and Wu, C M . , (1997). "Crushing Characteristics of 3-D Braided Composite Square Tubes", Journal of Composite Materials, Vol. 31, No. 22, 1997, pp. 2309-2327. Chiu, C.H., Tsai, K.-H., and Huang, W.J., (1998). "Effects of Braiding Parameters on Energy Absorption Capability of Triaxially Braided Composite Tubes", Journal of Composite Materials, Vol. 32, No. 21, pp. 1964-1983. Chiu, C.H., Tsai, K.-H., and Huang, W.J., (1999). "Crush-Failure Modes of 2D Triaxially Braided Hybrid Composite Tubes", Composites Science and Technology, Vol. 59, No. 11, pp. 1713-1723. Davis, J.R., (1975). "Compressive Strength of Fibre Reinforced Composites", in Composite Reliability ASTM STP 580(ed.), ASTM, Philadelphia, pp. 364-377. DeTeresa, S.J., Allison, L.M., Cunningham, B.J., Freeman, D.C., Saculla, M.D., Sanchez, R.J., and Winchester, S.W., (2001). "Experimental Results in Support of Simulating Progressive Crush in Carbon-Fiber Textile Composites", Lawrence Livermore National Laboratory, UCRL-ID-143287, March 12, 2001. Effendi, R.R., Barrau, J.-J., and Guedra-Degeorges, D., (1995). "Failure Mechanism Analysis Under Compression Loading of Unidirectional Carbon/Epoxy Composites Using Micromechanical Modelling", Composite Structures, Vol. 31, No. 2, pp. 97-98. Fairfull, A.H. and Hull, D., (1987). "Effects of Specimen Dimensions on the Specific Energy Absorption of Fibre Composite Tubes", Proceedings of Sixth International Conference on Composite Materials (ICCM VI) and Second European Conference on Composite Materials (ECCM 2), Vol. 3, 20-24 July 1987, pp. 3.36-3.45. Falzon, P.J. and Herszberg, I., (1998). "Mechanical Performance of 2-D Braided Carbon/Epoxy Composites", Composites Science and Technology, Vol. 58, No. 2, pp. 253-265. 164 References Farley, G.L., (1983). "Energy Absorption of Composite Materials", Journal of Composite Materials, Vol. 17, No. 3, pp.267-279. Farley, G.L., (1986a). "Effect of Fiber and Matrix Maximum Strain on the Energy Absorption of Composite Materials", Journal of Composite Materials, Vol. 20, No. 4, pp. 322-334. Farley, G.L., (1986b). "Effect of Specimen Geometry on the Energy Absorption Capabilities of Composite Materials", Journal of Composite Materials, Vol. 20, No. 4, pp. 390-400. Farley, G.L., (1991). "The Effects of Crushing Speed on the Energy-absorption Capability of Composite Tubes", Journal of Composite Materials, Vol. 25, No. 10, pp. 1314-1329. Farley, G.L. and Jones, R.M., (1992a). "Crushing Characteristics of Continuous Fiber-reinforced Composite Tubes", Journal of Composite Materials, Vol. 26, No. 1, pp. 37-50. Farley, G.L. and Jones, R.M., (1992b). "Prediction of the Energy-Absorption Capability of Composite Tubes", Journal of Composite Materials, Vol. 26, No. 3, pp. 388-404. Fleck, N.A., (1997). "Compressive Failure of Fibre Composites", Advances in Applied Mechanics, Vol. 33, pp. 43-113. Fleck, N.A., Deng, L., Budiansky, B., (1995). "Prediction of Kink Width in Compressed Fibre Composites", Journal of Applied Mechanics, Transactions ASME, Vol. 62, No. 2, pp. 329-337. Fleck, N.A., Sutcliffe, M.P.F., Sivashanker, S., and Xin, X.J., (1996). "Compressive R-curve of a Carbon Fibre-Epoxy Matrix Composite", Composites Part B: Engineering, Vol. 27, No. 6, pp. 531-541. Fleming, D.C, (2001). "Delamination Modeling of Composites for Improved Crash Analysis", Journal of Composite Materials, Vol. 35, No. 19, pp. 1777-1792. Fleming, D.C, (2000). "Modeling Delamination Growth in Composites using MSC.Dytran", Proceedings of the 2nd Worldwide Automotive Conference, MSC Software Corporation, Dearborn, MI, USA, 9-11 October 2000, pp. 1-15. Floyd, A.M., (2004). "An Engineering Approach to the Simulation of Gross Damage Development in Composite Laminates", Ph.D. Thesis, Department of Civil Engineering, The University of British Columbia. Gupta, N.K., Velmurugan, R., and Gupta, S.K., (1997). "Analysis of Axial Crushing of Composite Tubes", Journal of Composite Materials, Vol. 31, No. 13, pp. 1262-1286. 165 References Guynn, E.G., Ochoa, O.O., and Bradley, W.L., (1992). "Analysis of Fiber Microbuckling in Thermoplastic Composites", International Journal of Non-Linear Mechanics, Vol. 27, No. 6, pp. 1039-1047. Hahn, H.T., Williams, J.G., (1986). "Compression Failure Mechanisms in Unidirectional Composites", ASTM Special Technical Publication, pp. 115-139. Hamada, H., Coppola, J.C., and Hull, D., (1992a). "Effect of Surface Treatment on Crushing Behaviour of Glass Cloth/Epoxy Composite Tubes", Composites, Vol. 23, No. 2, pp. 93-99. Hamada H., Coppola J. C , Hull D., Maekawa Z., and Satoh H., (1992b). "Comparison of Energy Absorption of Carbon/Epoxy and Carbon/PEEK Composite Tubes", Composites, Vol. 23, No. 4, pp. 245-252. Hamada, H. and Ramakrishna, S., (1995). "Scaling Effects in the Energy Absorption of Carbon-fiber/PEEK Composite Tubes", Composites Science and Technology, Vol. 55, No. 3, pp. 211-221. Hamada, H., Ramakrishna, S., and Satoh, H., (1995). "Crushing Mechanism of Carbon Fibre/PEEK Composite Tubes", Composites, Vol. 26, No. 11, pp. 749-755. Hamada, H., Ramakrishna, S., and Satoh, H., (1996). "Effect of Fiber Orientation on the Energy Absorption Capability of Carbon Fiber/PEEK Composite Tubes", Journal of Composite Materials, Vol. 30, No. 8, pp. 947-963. Haug, E., Fort, O., and Tramecon, A., (1991). "Numerical Crashworthiness Simulation of Automotive Structures and Components Made of Continuous Fiber Reinforced Composite and Sandwich Assemblies", SAE Technical Paper Series 910152. Hull, D., (1988). "Energy Absorbing Composite Structures", Science and Technology Review, Vol. 3, pp. 22-30. Hull, D., (1991). "A Unified Approach to Progressive Crushing of Fibre-reinforced Composite Tubes", Composites Science and Technology, Vol. 40, No. 4, pp. 377-421. Hyer, M.W., (1998). "Stress Analysis of Fiber-Reinforced Composite Materials", WCB/McGraw-Hill. Jacob, G.C., Fellers, J.F., Simunovic, S., and Starbuck, J.M., (2002). "Energy Absorption in Polymer Composites for Automotive Crashworthiness", Journal of Composite Materials, Vol. 36, No. 7, pp. 813-850. Johnson, A., Pickett, A., and Rozycki, P., (2001). "Computational Methods for Predicting Impact Damage in Composite Structures", Composites Science and Technology, Vol. 61, No. 15, pp. 2183-2192. 166 References Karbhari, V.M. and Haller, J.E., (1998). "Rate and Architecture Effects on Progressive Crush-of Braided Tubes", Composite Structures, Vol. 43, No. 2, pp. 93-108. Karbhari, V.M., Falzon, P.J., and Herzberg, I., (1997). "Energy Absorption Characteristics of Hybrid Braided Composite Tubes", Journal of Composite Materials, Vol. 31, No. 12, pp. 1164-1186. Kerth, S., Dehn, A., Ostgathe, M., and Maier, M., (1996). "Experimental Investigation and Numerical Simulation of the Crush Behavior of Composite Structural Parts", Proceedings of the 1996 41st International SAMPE Symposium and Exhibition, Vol. 41, No. 2, pp.1397-1408. Kerth, S. and Maier, M., (1994). "Numerical Simulation of the Crushing of Composite Tubes", Proceedings of the 4th International Conference on Computer Aided Design in Composite Material Technology, Southampton, UK, June 1994, pp. 141-148. Mahmood, H. and Wang, R., (1997). "Finite Element Analysis of Automotive Fiber Reinforced Composite Structures Subjected to Crush Load", Proceeding of the 12th ESD Advanced Composites Conference and Exhibition, Detroit, MI, 7-10 April 1997, pp. 167-180. Mamalis, A.G., Manolakos, D.E., Demosthenous, G.A., and Ioannidis, M.B., (1997). "The Static and Dynamic Axial Crumbling of Thin-walled Fibreglass Composite Square Tubes", Composites Part B: Engineering, Vol. 28, No. 4, pp. 439-451. Mamalis, A.G., Manolakos, D.E., Ioannidis, M.B., Papapostolou, D.P., (2004). "Crashworthy Characteristics of Axially Statically Compressed Thin-walled Square CFRP Composite Tubes: Experimental", Composite Structures, Vol. 63, No. 3-4, pp. 347-360. McClennan, S.A., (2004). "Crack growth and damage modeling of fibre reinforced polymer composites", MASc. Thesis, Department of Materials Engineering, The University of British Columbia. Minguet, P.J.A., (1995). "Comparison of Graphite/Epoxy Tape Laminates and 2-D Braided Composites Mechanical Properties", Proceedings of the 36th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference, 10-13 April 1995, New Orleans, LA, USA, pp. 17-26. Moran, P.M., Liu, X.H., Shih, C.F., (1995). "Kink Band Formation and Band Broadening in Fibre Composites Under Compressive Loading", Acta Metallurgica et Materialia, Vol. 43, No. 8, pp. 2943-2958. . Moran, P.M., and Shih, C.F., (1998). "Kink Band Propagation and Broadening in Ductile Matrix Fibre Composites: Experiments and Analysis", International Journal of Solids and Structures, Vol. 35, No. 15, pp. 1709-1722. 167 References Morthorst, M. and Horst, P., (2004). "Failure Model for Composite Materials Under Quasi-static Crushing Conditions", Journal of Strain Analysis for Engineering Design, Vol. 39, No. 5, pp. 411-421. Naik, R.A., (1995). "Failure Analysis of Woven and Braided Fabric Reinforced Composites", Journal of Composite Materials, Vol. 29, No. 17, pp. 2334-2363. Quek, S.C., Waas, A., Hoffman, J., and Agaram, V., (2001). "The Crushing Response of Braided and CSM Glass Reinforced Composite Tubes", Composite Structures, Vol. 52, No. 1, pp. 103-112. Quek, S.C., Waas, A., Shahwan, K.W., Agaram, V., (2004a). "Compressive Response and Failure of Braided Textile Composites: Part 1 - Experiments", International Journal of Non-Linear Mechanics, Vol. 39, No. 4, pp. 635-648. Quek, S.C., Waas, A., Shahwan, K.W., Agaram, V., (2004b). "Compressive Response and Failure of Braided Textile Composites: Part 2 - Computations", International Journal of Non-Linear Mechanics, Vol. 39, No. 4, pp. 649-663. Ramakrishna, S. and Hull, D., (1993). "Energy Absorption Capability of Epoxy Composite Tubes with Knitted Carbon Fibre Fabric Reinforcement", Composites Science and Technology, Vol. 49, No. 4, pp. 349-356. Rosen, V.W., (1965). "Mechanics of Composite Strengthening", Fibre Composite Materials. American Society of Metals, Metals Park, Ohio, pp. 37-75. Salvi, A.G., Chung, J., Waas, A.M., and Caliskan, A., (2003). "Strain-Rate Effects on Unidirectional Carbon-Fiber Composites", AIAA Journal, Vol. 41, No. 10, pp. 2020-2028. Salvi, A.G., Waas, A.M., and Caliskan, A., (2004). "Specimen Size Effects in the Off-Axis Compression Test of Unidirectional Carbon Fiber Tow Composites", Composites Science and Technology, Vol. 64, No. 1, pp. 83-97. Schapery, R.A., (1995). "Prediction of Composite Strength and Kink Bands in Composites Using a Work Potential", International Journal of Solids and Structures, Vol. 32, No. 6-7, pp. 739-765. Schmueser, D.W. and Wickliffe, L.E., (1987). "Impact Energy Absorption of Continuous Fiber Composite Tubes", Journal of Engineering Materials and Technology, Transactions of the ASME, Vol. 109, No. 1, pp. 72-77. Schultz, M.R., Hyer, M.W., and Fuchs, H.P., (2001). "Static and Dynamic Energy-absorption Capacity of Graphite-Epoxy Tubular Specimens", Mechanics of Composite Material and Structures, Vol.8, No. 3, pp. 231-247. 168 References Sigalas, I., Kumosa, M., and Hull, D., (1991). "Trigger Mechanisms in Energy-absorbing Glass Cloth/Epoxy Tubes", Composites Science and Technology, Vol. 40, No. 3, pp^ 265-287. Sivashanker, S., (1998). "Damage Growth in Carbon Fibre-PEEK Unidirectional Composites Under Compression", Materials Science & Engineering A: Structural Materials: Properties, Microstructure and Processing, Vol. A249, No. 1-2, pp. 259-276. Sivashanker, S., (2001). "Damage Propagation in Multidirectional Composites Subjected to Compressive Loading", Metallurgical and Materials Transactions A: Physical Metallurgy and Materials Science Physical Metallurgy and Materials Science, Vol. 32, No. l,pp. 171-182. Sivashanker, S., Fleck, N.A., Sutcliffe, M.P.F., (1996). "Microbuckle Propagation in a Unidirectional Carbon Fibre-Epoxy Matrix Composite", Acta Materialia, Vol. 44, No. 7, pp. 2581-2590. Song, H.W., Du, X., and Zhao, G., (2002). "Energy Absorption Behavior of Double-chamfer Triggered Glass/Epoxy Circular Tubes", Journal of Composite Materials, Vol. 36, No. 18, pp. 2183-2198. Sutcliffe, M.P.F., Fleck, N.A., (1994). "Microbuckle Propagation in Carbon Fibre-Epoxy Composites", Acta Metallurgica et Materialia, Vol. 42, No. 7, pp. 2219-2231. Tay, T., Lee, K.H., Ramakrishna, S., and Shen, F., (1998). "Modelling the Crushing Behaviour of Composite Tubes", Key Engineering Material, Vol. 141-143, pp. 777-790. Thornton, P.H., (1979). "Energy Absorption in Composite Structures", Journal of Composite Materials, Vol. 13, pp. 247-262. Thornton, P.H., (1990). "The Crush Behavior of Pultruded Tubes at High Strain Rates", Journal of Composite Materials, Vol. 24, pp. 594-615. Thornton, P.H. and Edwards, P.J., (1982). "Energy Absorption in Composite Tubes", Journal of Composite Materials, Vol. 16, No. 6, pp. 521-545. Thornton, P.H., Hardwood, J.J., and Beardmore, P., (1985). "Fiber-reinforced Plastic Composites for Energy Absorption Purposes", Composites Science and Technology, Vol. 24, No. 4, pp. 275-298. Vogler, T.J., and Kyriakides, S., (1997). "Initiation and Axial Propagation of Kink Bands in Fibre Composites", Acta Materialia, Vol. 45, No. 6, pp. 2443-2454. West, A.C. and Adams, D.O., (1999). "Axial Yarn Crimping Effects in Braided Composite Materials", Journal of Composite Materials, Vol. 33, No. 5, pp. 402-419. 169 References Williams, K.V., (1998). "A Physically-Based Continuum Damage Mechanics Model for Numerical Prediction of Damage Growth", Ph.D. Thesis, Department of Materials Engineering, The University of British Columbia. Williams, K.V., Vaziri, R., Poursartip, A., (2003). "A Physically-Based Continuum Damage Mechanics Model for Thin Laminated Composite Structures", International Journal of Solids and Structures, Vol. 40, pp. 2267-2300. Xiao, X., Johnson, N., and Botkin, M., (2003). "Challenges in Composite Tube Crash Simulations", Proceedings of the American Society for Composites 18th Technical Conference, ASC 154, Gainesville, FL, USA, 19-22 October 2003. Xiao, X. (xinran.xiao@gm.com), (2003a). Personal correspondence. Experimental axial tube crushing force-displacement profiles, 21 July, 2003 and 1 April, 2004. Xiao, X., (2003b). Personal correspondence. File name: 2003-07 compressive_test results conducted at GM.pdf, 7 July, 2003. Xiao, X., (2003c). "Axial Crash Simulations of Braided Carbon Composite Tubes using LS-DYNA", Presentation, University of British Columbia, 3 April, 2003. Zobiery, N., (2004). "Progressive Damage Modeling of Composite Materials under Compressive Loads", M.A.Sc. Thesis, Department of Civil Engineering, The University of British Columbia. 170 Appendix A: Implementation of Compressive Damage Model APPENDIX A I M P L E M E N T A T I O N OF COMPRESSIVE D A M A G E M O D E L This appendix presents the details of the two main modifications made to CODAM in order to incorporate compressive damage modeling capabilities. These modifications were incorporated into a new version of CODAM (version 8.1). The first deals with increasing the number of scaling factors (k) from 3 to 6. This allows scaling factors in tension and compression to be independently defined in the three principal directions. The following steps were taken to accomplish this: • The scaling factors were moved from the material (cm) vector (defined in *MAT_USER_ DEFLNEDMATERIALJVIODEL) to the end of the ssparam.dat-auxiliary file because of lack of space in the cm vector. The following three lines were added to the end of the ssparam. dat file: c Modified Bazant s c a l i n g factors . c Kscalext Kscalexc Kscaleyt Kscaleyc Kscalezt Kscalezc 2.0 1.2 1.0 1.7 1.0 1.0 A flag (cm(2S) = BSFLAG) was defined in the material model input deck to define whether modified Bazant scaling is active. If this parameter is set to 1.0, scaling is active, otherwise it is not. The k-scaling section in the user material model code was modified to include the calculation of strain at maximum stress in compression (xcmaxeps, ycmaxepx, and zcmaxeps). This parameter determines when and how the post-peak response is modified by the scaling factors in compression. The second modification involves incorporating the plateau stress into CODAM. The following steps were taken to accomplish this: • The plateau stress in each of the three principal directions is dependent on the compressive matrix damage initiation, matrix damage saturation, and fibre damage saturation. The plateau stress (platx, platy, platz), compressive matrix damage. 171 Appendix A: Implementation of Compressive Damage Model initiation effective strains (Fxci, Fyci, Fzci), and compressive matrix damage saturation effective strains (Fxcs, Fycs, Fzcs) must be defined for each principal direction. This is done in the s spar am. dat file as shown below. c Plat e a u s t r e s s c p l a t x F x c i Fxcs -50 0 . 005 0 . 010 c p l a t y F y c i Fzcs -30 0 . 005 0 . 012 c p l a t z F y c i Fzcs 0 . 0 The plateau stress is calculated in the model when the strains are such that the element is being sufficiently loaded in compression (damage is induced in the matrix) using the following equations. If the strains are such that the element is completely damaged (fibre damage saturation), the plateau stress is set to zero. Unloading and reloading follows a path identical, but offset, to the initial loading path: [° plateau \ = PlatX " (wXtmpf where: wxtmp - min 1.0, Fx - Fxci Fxcs - Fxci Fx = instantaneous effective strain function A.l A.2 Fibre damage saturation Figure A. 1 Variation of plateau stress with effective strain (local x-direction). • At each time step, the plateau stress is added on the level of stress already calculated by the material model. sigx(i) = c \ \ * s x + c \ 2 * e y + c n * s z + ( a p l a l e a u )x sigy(i) = c\2*sx+c22*ey + c23 * sz + (aplaleau \ A.3 sigz(i) = c\3*£x + c23 *sy+c33 * sz + (crplateau \ 172 Appendix B: LPT and Ply Discount Calculations for Braided Material APPENDIX B L P T AND PLY DISCOUNT C A L C U L A T I O N S FOR BRAIDED M A T E R I A L This appendix presents the method used to estimate the modulus reduction due to the presence of damage in the braided composite laminate. To apply this method, the [0/±30]2 braid architecture is first approximated as a symmetric laminate made of 6 individual plies in a [±30/0] s lay-up. The individual ply properties are determined from the experimentally reported volume fractions and constituent properties (DeTeresa et al., 2001). Laminate Plate Theory (LPT) (see Hyer, 1998 for details) is then applied to estimate the undamaged elastic properties of the laminate. The next step is to estimate how the elastic laminate properties are reduced by damage growth. To do this, the ply discount method in conjunction with LPT is used. In summary, the approach is used to estimate the reduction in the elastic laminate properties at 100% matrix damage. Complete matrix damage is represented in the individual plies by a complete loss in the 2-direction (matrix direction) and shear stiffness. The 1-direction stiffness (fibre direction) is assumed to be unaffected by the matrix damage. The appropriate reduction in. the ply stiffness values, as outlined in Table 3.1, is then incorporated into a new LPT analysis to attain an approximation of the damaged laminate elastic properties. 1. DETERMINE VOLUME FRACTIONS For this analysis, the volume fractions of the three individual laminae are assumed to be equal to the volume fractions of the overall braided lamina reported by DeTeresa etal., 2001. Overall braided lamina properties = individual lamina properties Vm = 0.513 Vf = 0.450 Vv = 0.037 2. DETERMINE INDIVIDUAL PLY PROPERTIES Given the following constituent properties: 173 Appendix B: LPT and Ply Discount Calculations for Braided Material Determined experimentally (DeTeresa et al., 2001): Ef = 2 3 1 , 0 0 0 M P a Em = 3,250 M P a (average o f compressive and tensile modu l i ) Gm = 1,380 M P a Derived: v m = 0 . 1 8 (from Gm = Em/2{l + vm)) Assumed (from typical carbon fibre properties (Hyer, 1998)): Gf 8,960 M P a vf = 0.20 The i n d i v i d u a l p l y elastic properties are determined us ing a rule-of-mixtures approach (Hyer , 1998): £, = VfEf + VmEm = (0 .45)231000+ (0.514)3250 = 105,620 MPa E. + (0.514) > E7 =6,246 MPa 231000 3250 ' 2 ^f Vn =Vfvf+Vmvm = (0 .45)0 .20+ (0 .514)0 .18= 0.18 MPa — = Vf — + Vm — = ( 0 . 4 5 ) - ^ — + ( 0 . 5 1 4 ) — > G 1 2 = 2,366 MPa G 1 2 / Gf Gm 8960 1380 E2 6246 2 1 £ , 1 2 105620 (0.18) =0.0106 MPa 3. DETERMINE THE INDIVIDUAL PLY STIFFNESS MATRIX [Q] 'Qn Qn 0 " '*1 Qn Q22 0 £2 0 0 & 6 _ 7l2. ( B . l ) where: a,- E' l - v 1 2 v 2 1 0 , - r ^ ~ l - v 1 2 v 2 1 Q 2 2 = T J ^ ~ l - v , 2 v 2 , 066 = G \ 2 The [Q] matr ix o f the i nd iv idua l unidirect ional pl ies is equal to: (B.2) (B.3) (B.4) (B.5) 174 Appendix B: LPT and Ply Discount Calculations for Braided Material \Q} = 105822 1126 0 1126 6258 0 0 0 2366 MPa 4 TRANSFORM THE PLY STIFFNESS MATRIX FROM THE LAMINA TO THE LAMINATE CO-ORDINATE SYSTEM [Q ] (1,2,3 TO x,y,z) (B.6) "fin fi.2 fi16" '<*> Qn fi22 fi26 • X xy fil6 fi26 fi66 7 *y, where: m Qn =Qum4 +2(QU +2Q66)n2m2+Q22n\ Qn = (fin + fi22 - 4fi 6 6 )n2m2 + Qu («4 + m4) fi,6 =(fin 'Qn-IQaW +(Qn ~ f i 2 2 + 2 f i 6 6 ) « Q22 =Quni+2(Qu +2Q66)n2m2 +Q22m" Qie =(Qu~Qn- 2fi 6 6 Vm + (Qn - Q22 + 2Q66 )nm3 Q66 = (fin + fi22 - IQn ~ 2fi 6 6 )n2m2+ Q66(n4 + mA) m = cos(d) n = sin(0) (B.7) The [Q ] matrix of the individual unidirectional plies is equal to: 0=[Q} = 105822 1126 0 1126 6258 0 0 0 2366 MPa -30 621113 19944 32421 1126 12331 10692 32421 10692 21184 MPa 1-30 621113 19944 -32421 1126 12331 -10692 -32421 -10692 21184 MPa 175 Appendix B: LPT and Ply Discount Calculations for Braided Material 5. DETERMINE THE AVERAGE LAMINATE STIFFNESS MATRIX [ Q ] For symmetric laminates: Q n (B.8) The average stiffness matrix for the braided laminate is equal to: Q 76683 13672 0 13672 10307 0 0 0 14912 MPa 6. DETERMINE THE UNDAMAGED LAMINATE ELASTIC CONSTANTS QnQn-Q 12 Ex= _ Q: p _E&n Qn Gxy = 066 22 xy Q: 22 (B.9) The undamaged laminate elastic constants for the braided laminate are: (Ex\=58.5 GPa (Ey)o=1.9 GPa fej0=14.9 GPa 7. DETERMINE THE DAMAGED LAMINATE ELASTIC CONSTANTS The damaged laminate elastic constants at matrix damage saturation are determined by reducing the appropriate ply stiffness values. Following the guideline presented in Table 3.1, the individual ply properties are reduced appropriately and then steps 3 through 6 are followed to estimate the damaged laminate elastic constants. A t matrix damage saturation in tension, the 2-direction and shear stiffness values are reduced to 176 Appendix B: LPT and Ply Discount Calculations for Braided Material zero. In compression, in addition to these reductions, the 1-direction stiffness is decreased to 1% of its undamaged value. The resulting damaged laminate properties at matrix damage saturation in tension and compression are as follows: The damaged laminate elastic constants for the braided laminate in tension are: Gxy =13.2 GPa vxy = 3.00 GPa The damaged laminate elastic constants for the braided laminate in compression are: Ex=0.59GPa -> ( l £ ) x = ° - % 5 = 0.01 Ey = 0-16 GPa -> [Ecm\ = 0 . 1 ^ g = 0.02 Gxy =0 .13 GPa v =3 .00 GPa 50.9 GPa £ =3 .8 GPa xy 111
The Open Collections website will be unavailable July 27 from 2100-2200 PST ahead of planned usability and performance enhancements on July 28. More information here.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Simulation of progressive damage development in braided...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Simulation of progressive damage development in braided composite tubes under axial compression McGregor, Carla Jane 2005
pdf
Notice for Google Chrome users:
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
Page Metadata
Item Metadata
Title | Simulation of progressive damage development in braided composite tubes under axial compression |
Creator |
McGregor, Carla Jane |
Publisher | University of British Columbia |
Date Issued | 2005 |
Description | Composite materials have recently been investigated as an attractive alternative to traditional materials used as energy absorbers in the automotive industry. One of the factors inhibiting the implementation of these lightweight energy absorbers is the lack of accurate, robust, and physically meaningful material models capable of simulating damage growth in composite materials under a variety of applied loads. A plane-stress continuum damage mechanics based model for composite materials, CODAM (COmposite DAMage), recently developed at the University of British Columbia (Williams, 1998; Williams et al., 2003; Floyd, 2004) has proved to be capable of predicting the complete tensile behaviour of composites from initiation of damage to complete failure. The current study focuses on extending the model to capture the complete compressive response of composite materials as well as their behaviour under load reversals. Refinements to the model are based on the experimentally observed compressive failure mechanisms presented in the literature. In particular, the mechanical consequences of kinking and kink band broadening, in conjunction with matrix cracking and delamination, are represented in the model. The model parameters defining the compressive response are related to experimentally observed behaviour, thus maintaining the physical basis of CODAM. These improvements contribute to the development of CODAM as a predictive tool in the analysis of gross damage growth in composite materials. To demonstrate the capability of the model, the predicted response of braided composite panels under tensile loads, compressive loads and load reversals are compared to experimentally observed behaviours. The model is also evaluated on its ability to predict the damage propagation and energy absorption in axial crushing of braided composite tubes. Model predictions for both panel and tube specimens correlate well with the experimental results. The successful tube crushing predictions suggests that CODAM could offer an attractive design aid for the future incorporation of lightweight composite energy absorbers into crashworthiness structures. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2009-12-10 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0063295 |
URI | http://hdl.handle.net/2429/16352 |
Degree |
Master of Applied Science - MASc |
Program |
Civil Engineering |
Affiliation |
Applied Science, Faculty of Civil Engineering, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 2005-05 |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-ubc_2005-0268.pdf [ 23.56MB ]
- Metadata
- JSON: 831-1.0063295.json
- JSON-LD: 831-1.0063295-ld.json
- RDF/XML (Pretty): 831-1.0063295-rdf.xml
- RDF/JSON: 831-1.0063295-rdf.json
- Turtle: 831-1.0063295-turtle.txt
- N-Triples: 831-1.0063295-rdf-ntriples.txt
- Original Record: 831-1.0063295-source.json
- Full Text
- 831-1.0063295-fulltext.txt
- Citation
- 831-1.0063295.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
data-media="{[{embed.selectedMedia}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0063295/manifest