Thesis

LOW-ENERGY IMPACT BEHAVIOUR OF LAMINATED COMPOSITE Thesis submitted for the degree of bachelor of aeronautical engineering by Student Nguyen Duc Anh

Advisor: Dr.Nguyen The Hoang & Eng.Nguyen Song Thanh Thao

Ho Chi Minh, july 2010 i

Abstract Impact behaviour of fibre-reinforced composite materials brings together for the first time the most comprehensive and up-to-date work undertaken in leading research centres worldwide. Impact response, damage tolerance and failure are extensively investigated from a number of viewpoints. The low-energy impact on laminated composite is the research topic in this thesis. The main objective is to simulate the phenomenon of low-energy impact on laminated composite. An integrated modelling for carbon/epoxy composite plates has been developed to evaluate the impact damage initiation and propagation. Three failure modes: resin crazing, delamination and fibre failure have been investigated and implemented for the damage simulation. Resin crazing and fibre failure are studied based on Hashin criterion. And delamination is studied based on cohesive zone model. Three-dimensional dynamic finite element method was used to simulate the real impact events. The simulation is based on Abaqus CAE package. This thesis is divided into six chapters and appendix Chapter 1: Introduction; Chapter 2: An overview of the impact behaviour of fibre-reinforced composites; Chapter 3: Contact mechanics; Chapter 4: Modelling with Abaqus/CAE 6.9; Chapter 5: Results; Chapter 6: Conclusion and recommendation; Appendix

Keywords Composites; Impact damage; Failure modes; Delamination; Cohesive zone model (CZM); Finite element analysis

ii

Acknowledgements I would like to express my gratitude to my advisor, Dr. Hoang, for his guidance. Since I have been studing Aviation Engineering, I have become interested in structural numerical computing. It was him who gave me the chance to approach it. I also want to thank engineer Thao who has given me a great deal of advice and help at every stage of my research. I may not finish this thesis without their kindness. I wish to thank all teachers who was giving me instructions during my study at university, especially the teachers in Aeronautical Department. Finally, special thanks to my family and friends who have supported and encouraged me throughout this work. I love all.

iii

List of Figures Figure 2.1 Composite material ............................................................................................ 5 Figure 2.2 Fracture mechanisms observed in laminates ..................................................... 9 Figure 2.3 Fracture propagation in the case of poor fibre-matrix bonding ......................... 9 Figure 2.4 Fracture surface associated with poor fibre/matrix bonding in the case of a carbon fibre composite (ONERA document) ................................................................... 11 Figure 2.5 Fracture surface associated with high fibre/matrix bonding in the case of a carbon fibre composite (ONERA document) ................................................................... 12 Figure 2.6 X-ray observation of the fracture state of carbon fibre composites after fatigue (105 cycles) in the case of impacted specimens (ONERA document) .............................. 12 Figure 2.7 The acoustic emission process......................................................................... 13 Figure 2.8 Acoustic emission signals recorded during bending tests on unidirectional carbon fibre-epoxide composites ...................................................................................... 13 Figure 2.9 Delamination ................................................................................................... 17 Figure 2.10 Ultrasonic scanning images of the dlemination in cross-ply laminates after impact: (a) after 1J impact; (b)after 2J impact; (c) after 3J impact (ELSEVIER document) ........................................................................................................................................... 17

Figure 3.1 impact mechanism ........................................................................................... 18 Figure 3.2 contact mechanical models .............................................................................. 20 Figure 3.3 definition of parameter R................................................................................. 20 iv

Figure 3.4 indentation of glass-epoxy laminate from Yang and Sum (1982). (a) Loading; (b) Unloading; (c) Reloading ............................................................................................ 22 Figure 3.5 (a) two-degree-of-freedom model. (b) single-degree-of-freedom model. ....... 25 Figure 3.6 Linear two-degree-of-freedom spring-mass impact model ............................. 26 Figure 3.7 Free vibrations of two-degree-of-freedom sustem: (a) first non-dimensional frequency; (b) second non-dimensional frequency (solid line: exact; dashed line; approximate); (c) Curve veering phenomenon for k2/k1=10 (solid line: exact; dashed line; approximate). .................................................................................................................... 29 Figure 3.8 Apprpximate mode shapes of linear two-degree-of-freedom system.............. 29

Figure 4.1 The contronlled drop-weight rig ...................................................................... 34 Figure 4.2 modelling of composite plate and impactor .................................................... 35 Figure 4.3 Material orientations of the composite plate ................................................... 35 Figure 4.4 Boundary conditions ........................................................................................ 37 Figure 4.5 The mesh ......................................................................................................... 38 Figure 4.6 pressure-overclosure relationship of hard contact ........................................... 39 Figure 4.7 cohesive layer mesh ......................................................................................... 41 Figure 4.8 model with cohesive elements for delaminations ............................................ 41

v

Figure 5.1 Impact process at t=0.00ms; 0.08ms; 0.16ms; 0.24ms; 0.32ms; 0.40ms; 0.48ms; .............................................................................................................................. 45 Figure 5.2 Velocity of impactor during and after impact process .................................... 45 Figure 5.3 Impact force at right impact point during and after impact process. ............... 46 Figure 5.4 Contact stress at right impact point during and after impact process .............. 46 Figure 5.5 Relation between contact force and indentation during impact process including loading and unloading, t=0 to t=0.28ms. .......................................................... 47 Figure 5.6 Delamination process at t=0.28ms; 0.52ms;.................................................... 48 Figure 5.7 Stress field in the surface of impact at t=0.36ms............................................. 48 Figure 5.8 Mises stress distributes along thickness at t=0.08ms, distance from centre of impact is normalized by percent ....................................................................................... 49 Figure 5.9 Shear stress in-plane 1-0-2 at first interface, t=0.04ms; distance from centre of impact to distance of 10mm is normalized by percentage. ............................................... 49 Figure 5.10 Shear stress in-plane 1-0-2 at second interface, t=0.04ms; distance from centre of impact to distance of 10mm is normalized by percentage. ................................ 50 Figure 5.11 Status of elements after impact, t=0.8ms....................................................... 51 Figure 5.12 Tensile fibre damage initiation after impact, t=0.8ms. .................................. 52 Figure 5.13 Compressive fibre damage initiation after impact, t=0.8ms.......................... 52 Figure 5.14 Compressive fibre damage evolution after impact, t=0.8ms. ........................ 53 Figure 5.15 Tensile matrix damage initiation after impact, t=0.8ms. ............................... 53 Figure 5.16 Tensile matrix damage evolution after impact, t=0.8ms. .............................. 54

vi

Figure 5.17 Compressive matrix damage initiation after impact, t=0.8ms. ...................... 54 Figure 5.18 compressive matrix damage evolution after impact, t=0.8ms. ...................... 55 Figure 5.19 Delamination ................................................................................................. 55 Figure 5.20 Shear damage evolution after impact, t=0.8ms. ............................................ 56 Figure 5.21 Delamination at both two interlaminars at t=0.52ms .................................... 57 Figure 5.22 State of first interlaminar (left) and second interlaminar (right) after impact, t=0.8ms ............................................................................................................................. 58 Figure 5.23 Traction damage initiation at first interlaminar (left) and second interlaminar (right) after impact, t=0.8ms ............................................................................................. 58 Figure 5.24 Interlaminar delamination shaped a peanut at second interlaminar, t=0.04ms ........................................................................................................................................... 59 Figure 5.25 Interlaminar delamination at second interlaminar after impact, t=0.8ms ...... 60

vii

Contents Abstract

ii

Acknowledgements

iii

List of Figures

iv

Chapter 1: Introduction

1

1.1

Background ...................................................................................................................... 1

1.2

Objective and scope of this thesis .................................................................................... 2

1.3

Thesis layout .................................................................................................................... 2

Chapter 2 An overview of the impact behaviour of fibre-reinforced composites

4

2.1

Introduction ...................................................................................................................... 4

2.2

Composite materials ......................................................................................................... 4

2.2.1

Properties of fibres and resins ................................................................................... 5

2.2.2

Composite properties ................................................................................................ 6

2.3

Failure processes .............................................................................................................. 8

2.3.1

Why are composites prone to impact damage?....................................................... 10

2.3.2

Estimation of the stresses acting in an impact ........................................................ 10

2.4

Observation of fracture mechanisms .............................................................................. 11

2.4.1

Observation by Microscopy .................................................................................... 11

2.4.2

Radiography analysis .............................................................................................. 12

2.4.3

Acoustic emission analysis ..................................................................................... 13

2.5

Failure criteria ................................................................................................................ 14

Chapter 3 Contact mechanics

18

Abstract ..................................................................................................................................... 18 3.1

Introduction .................................................................................................................... 18

3.2

Contact mechanics.......................................................................................................... 19

viii

3.3

Energy-balance models .................................................................................................. 22

3.4

Impact on infinite composite plates ............................................................................... 23

3.5

Spring-mass models ....................................................................................................... 24

3.5.1

General case ............................................................................................................ 24

3.5.2

Linear two-degree-of-freedom model ..................................................................... 25

3.5.3

Linear single-degree-of-freedom (SDOF) spring-mass model ............................... 30

3.5.4

Nonlinear SDOF spring-mass models .................................................................... 30

3.6

Complete models ............................................................................................................ 30

3.7

Selection of an impact dynamics model ......................................................................... 31

3.8

Conclusions .................................................................................................................... 32

Chapter 4 Modelling with Abaqus/CAE 6.9

33

4.1

Introduction .................................................................................................................... 33

4.2

Simulation procedure of the impact model .................................................................... 34

4.2.1

Model ...................................................................................................................... 34

4.2.2

Material properties of composite ............................................................................ 35

4.2.3

Boundary conditions ............................................................................................... 37

4.2.4

Mesh........................................................................................................................ 37

4.2.5

Technique for simulation the impact ...................................................................... 38

4.2.6

Jobs of simulation ................................................................................................... 39

4.3

Simulating the delamination using cohesive zone model .............................................. 39

4.3.1

Mesh the surrounding bulk material ....................................................................... 40

4.3.2

The cohesive layer properties ................................................................................. 40

Chapter 5 Results

43

5.1

Overview ........................................................................................................................ 43

5.2

Observation of impact process ....................................................................................... 43

5.2.1. Contact properties ...................................................................................................... 46 5.2.2. Delamination process ................................................................................................ 47

ix

5.3

Investigation of stress field ............................................................................................ 48

5.4

Resin crazing, fibre failure and prediction of delamination ........................................... 50

5.4.1 Fibre failure ................................................................................................................ 51 5.4.2 Matrix failure .............................................................................................................. 53 5.4.3 5.5

Prediction of interlaminar delamination ................................................................. 55

Delamination .................................................................................................................. 57

Chapter 6 Conclusion and recommendation

61

6.1

Conclusion...................................................................................................................... 61

6.2

Recommendation ............................................................................................................ 62

References

63

Appendix

65

x

Chapter 1: Introduction 1.1

Background The last four decades have witnessed a steady increase in the use of light-weight polymetric composites for structural applications. In addition to their excellent quasistatic mechanical properties such as high specific stiffness and strength, it has become essential for these composite structures to perform well under various types of impact loading. For example, fibre-reinforced composite materials are used in aircraft, modern vehicles and light-weight structures. In the aerospace industry, the residual compressive strength of an impact-damaged composite struture has become the design-limiting factor. Major aspects of current damage-tolerant design philosophies in other industrial sectors are similar to that of the aerospace industry. There has been a significant increase in the application of composites in the aerospace industry. Composite are being selected for their high stiffness, high strength to weight ratio and energy absorbing properties. Composite offers significant weight savings compared to metal. Reduced weight results in opportunities to improve fuel economy, which subsequently helps to reduce emissions. However, their application remains limited because of not only high material and production cost but also long manufacturing cycle times. Nowadays, finite element simulations are use extensively in airplane design. These tools have allowed engineers to reduce design cycles and avoid the high costs associated with experimental testing. Simulations of metallic structures are now well established; however, advanced finite element methods for fibre-reinforced composites are just emerging. In particular, the study of impact damage progression and failure in composite structures has become one of the most researched areas over the last two decades. However, in contrast to monolithic composite structures, composite sandwich structures have received limited attention. The level of structural and material complexity of the composite sandwich is greater than that of metallic and monolithic composite strutures due to the inhomogeneous and anisotropic nature of the sandwich construction; in addition to the more complex failure mechanisms of the composite skins and core material. Advanced finite element composite material models can simulate damage progression and fracture are emerging; however the development of a predictive modelling methodology for composite sandwich structure is still in its infancy.

1

1.2

Objective and scope of this thesis The overall objective of this thesis is to develop a predictive modelling of impact damage in fibre-reinforced composite plates. Specifically, the focus of this study will be on lowenergy impact on sandwich structure composite plate. Composite structures have a high strength-to-weight and stiffness-to-weight ratio; however, because they are laminar systems with weak interfaces, they are very susceptible to impact loading. This may cause invisible cracks and delaminatios to occur in the material, which are often difficult to predict and detect. The damage of composite structures to impact events is one of the most important aspects to behaviour that inhibits more widespread application of composite materials. Understanding the deformation and damage mechanisms involved in the impact of composite targets is important in the effective design of a composite structure. A prerequisite for increasing the use of composite materials is the need to predict damage from low-energy impact. Under lowenergy impact loading where the elastic energy absorbing capability of the structure is important. Experiments show that low-energy impact on sandwich plates introduces significant reductions in tensile, comprehensive, shear and bending strengths. They aslo show that damage due to low-energy impact are mainly matrix cracks and delaminations, not fibre failures. The main aims are to: 



1.3

Develop an appropriate finite element modelling for the dynamic impact loading of composite sandwich structures which predicts the damage and deformation behaviour of these structures. Determine the failure modes of composite sandwich structures under low-energy dynamic impact.

Thesis layout This thesis is divided into six chapter. Following this introduction (chapter 1), the thesis is composed as follows: Chapter 2 presents a comprehensive literature review of the composite material and its impact behaviour. Chapter 3 presents an overview of contact mechanics for understanding the principle of the impact event. Chapter 4 presents modelling using Abaqus for simulation the dynamic impact. This modelling uses the Hashin criterion for prediction of damages; and uses the cohesive zone model for simulation of interlaminar delamination. 2

Chapter 5 presents the results of simulation with discussion. Chapter 6 presents the conclusion of this thesis. Future research themes are also outlined.

3

Chapter 2 An overview of the impact behaviour of fibrereinforced composites 2.1

Introduction Impact may be defined as the relatively sudden application of an impulsive force, to a limited volume of material or part of a structure. The proviso is that relatively and limited are capable of an extraordinarily wide range of interpretations. The effects of impact are widely known and yet analysing the phenomenon and relating effects to the forces acting and the materials’ properties, in order to predict the outcome of a particular event , can be very difficult. The results of an impact can be largely elastic, with some energy dissipated as heat, sound, internally in the material, ect. Alternatively there may be deformation, permanent damage, complete penetration of the body struck or fragmentation of the impacting or impacted body, or both. For fibre composite materials it is permanent damage, possibly subsurface and barely visible, penetration and fragmentation, that are of interest. There are various ways of analysing the impact process; in terms of the energy deposited and gross damage produced, micro energy dissipation or by considering the stresses acting on flaws in the material and the effects that are generated. The latter method, which is known as fracture mechanics.

2.2

Composite materials What is composite material? As the term indicates, composite material reveals a material that is different from common heterogeneous materials. Currently composite material refers to materials having strong fibers—continuous or noncontinuous—surrounded by a weaker matrix material. The matrix serves to distribute the fibers and also to transmit the load to the fibers.

4

Figure 2.1 Composite material To discuss the impact behaviour of composites it is helpful initially to consider the nature of the materials. Modern polymer composites based on glass, carbon, aramid, ceramic or polymer fibres in a polymer matrix are heterogeneous and anisotropic materials. They have a low density, high strength and stiffness and hence excellent specific properties (the ratio of the property to the density or specific gravity) in the fibre direction. They can be relatively easily formed into components with complex shapes by vacuum pressure compaction or pressing at temperatures of no more than ~200oC for 1h or more, though post curing may be needed. Because of the need to have long fibres for optimum properties and, in many cases, easier and effective processing, composite structures are often shells. They have excellent properties within the plane of the shell, in the fibre direction in tension, and to a lesser extent, in compression. The performance tranversely, through the thickness and in shear and impact, is poorer. Generally, properties which are matrix dependent are much lower than those which are governed by the fibre. Nevertheless the role of the matrix is vital in composite behavior; protecting the fibre, transferring stresses and, in some cases, alleviating brittle failure by providing alternative paths for crack growth. 2.2.1 Properties of fibres and resins [1] Reinforcing fibres usually have a diameter of 5- 10 m , through it may be up to ~100

m for certain types ( e,g boron), and a modulus and strength of the order of 70-800GPa and 1000-7000MPa, respectively. The failure strain is ~0.27-5.0%. Carbon and aramid fibres may exhibit an increasing modulus with increasing strain due to changes in their internal structure during stressing. The resins, whether thermosets or thermoplastics, have a modulus and strength of the order of 2-5GPa and 50-100 MPa, respectively, and a strain to failure of ~1% upwards. It should be noted that resin failure strain is usually well in excess of the minimum figure. Properties of resins are usually strain rate dependent 5

and markedly influenced by temperature. This occurs because the glass transition temperature, Tg, which is indicative of the upper working temperature, is ~80-200oC for most resins. This is much closer to room temperature than is the case for the melting or degradation point of most structural materials. Hence creep can be a problem with polymers at relatively low temperatures. Some fibre and resin properties are summarised in table 2.1. It should be noted that these are indicative of the property range for the materials. Individual grades/types of material, the way the property is measured, etc., might cause variations outside the limits given. Density 3 ( Mgm ) 2.49-2.55 1.7-2.0 1.39-1.45 2.0-3.97 1.0-1.35 1.1-1.23 1.2-1.2

Tensile modulus Tensile strength Material (Gpa) (Mpa) Glass fibre 73-86 3400-4500 Carbon fibre 160-827 1400-7070 Aramid fibre 73-160 2400-3400 Inorganic fibre 152-462 1720-3900 Phenolic resin 3.0-4.0 60-80 Polyester resin 3.1-4.6 50-75 Epoxy resin 2.6-3.8 60-85 Bismaleimide resin 1.2-1.32 3.2-5.0 48-110 Table 2.1 Some fibre and resin properties [1]

Strain to failure (%) 3.5-5.4 0.27-1.9 1.4-4.6 ~1.8 1.0-6.5 1.5-8.0 1.5-3.3

2.2.2 Composite properties [1] Table 2.2 lists some composite properties. These are again indicative of the performance of 60v/o (fibre volume fraction), unidirectional fibre composites. Values will vary with the type of fibre (e.g. grade of carbon fibre), method of test, matrix and way the specimens were fabricated. The important points to appreciate are very significant differences between longitudinal and transverse properties, the lesser difference between longitudinal tensile and compressive strengths and the disparity between the interlaminar chear strength (ILSS) and tensile/compressive strengths. Shear modulus is in the range 2.5-5.0GPa for all of the materials. Longitudinal tensile and compressive strain to failure is essentially that of the fibre. Transversely the figure may be lower, ~0.3-0.5%, and is governed by the properties of the interface and resin matrix. Laminates, containing layers of fibre with different orientations within a plane, may have higher impact, shear and transverse properties, but lower longitudial performance. Charpy impact energy for unnotched materials, stressed in the longitudinal of fibre direction, varies widely depending on the interface, matrix and experimental details. There may also be errors due to the scaling up of results for different cross section specimens. Charpy impact strength in the transverse direction (i.e. at right angles to the fibre) can be very low.

6

Table 2.2 Some properties of unidirectional composite materials [1] Material Glass Carbon Aramid Material Glass Carbon Aramid

ρ (Mg m-3) 1.8-2.0 1.54-1.66 1.36-1.4

Elt (GPa) 37-53 125-330 66-107

G1c (J m-2) 125-3300 60-3300 1100

Ett (GPa) 9.0-13.6 5.9-10 4.1 G2c (J m-2) 600-1200 154-2033

σlt (GPa) 1.2-1.9 1.76-2.9 1.29-1.5

σtt (GPa)

σlc (GPa)

30-80 27

0.78-1.6 0.19-0.28

Imp st. l, (kJm-2) Up to 609 8-150 20-40

τ (MPa) 53-83 70-127 38-69

Imp st. t, (kJm-2) Up to 28 0.2-4.0

Where ρ is the density, Elt and Ett are the longitudinal tensile and transverse tensile moduli, respectively, σlt, σtt, σlc, are the longitudinal tensile, transverse tensile and longitudinal compressive strengths, respectively, τ is the shear strength (also known as the interlaminar shear strength or ILSS), G1c, G2c are the critical works of fracture for the 1 and 2 deformation modes and Imp st. l, and Imp st, t, are the longitudinal and transverse Charpy impact strengths. It is common to use Gc, the critical strain energy release rate, as a measure of toughness, though some workers distinguish between initiation and propagation energies for this quantity. Often Gc is equated to 2γ, where γ is the work of fracture or the energy required to produce a new surface in a material. It is important to remember that both γ and Gc refer to a specific model of deformation (a tensile type opening for G1c) and that the direction of the fibres in relation to the sense of deformation must be clearly defined. In table 2.2, G1c refers to the extension of a crack parallel to the fibre direction and is thus analogous to the transverse Charpy impact energy. G1c defined in this way is also known as the mode 1 interlaminar shear work of fracture. More information on the fibre, resin and composite properties is available, for instance, as is a simple description of the various deformation modes associated with G1c, G2c and G3c. A general account of composites fracture mechanics and toughening mechanisms is given in Hancox and Mayer. More details and descriptions of the experimental difficulties of determining the critical strain energy release rates and the influence of materials’ and testing parameters on the results have been recorded.

7

2.3

Failure processes New situations and materials tend to be judged initially by comparison with existing ones. Since the principle class of construction materials used in high performance applications is metals, composites are usually compared with them. Generally, metals tend to be dense, stiff, strong, good conductors of heat and electricity, isotropic and inexpensive. In particular many show plastic behaviour before failure and it may be possible to remove the effects of deformation/damage by annealing and/or reworking the material. Unfortunately, at least on a significant scale, polymer matrix fibre composites do not exhibit this type of behaviour-elastic deformation is followed by irreversuble damage and failure. Another important influence on composite performance, especially in failure, is the heterogeneous nature of the material on a micro/meso scale. Composites are generally strong, and have a reasonable impact resistance, if the applied stress is in the fibre direction. In other directions they tend to be weak and to have a low impact resistance. Small secondary stresses in the transverse direction or unexpected stresses due to an impact, in a weak direction, can easily cause damage. Bacause of the ease with which localised permanent damage can be produced and the low values of γ or G1c for failure parallel to the fibres, or between plies in a laminate, G1c is sometimes regarded as critical and the focus of most attention in impact work. Fracture mechanics methods of determining G1c are favoured over traditional measurements of impact strength, from which γ and arguably G1c can be determined, because of their greater accuracy and virtual elimination of extraneous energy absorption. However, the assumptions implicit in fracture energy determinations may not be true for fibre composites, where crack propagation may be dominated by crack blunting, fibre pull-out and fibre/matrix delamination. It has been suggested that in such circumstances conventional fracture toughness. Impact tests can be more useful, lathough no so readily analysed as fracture toughness measurements, and can be conducted at strain rates more representative of the deformation rates seen in service. There are five basic machanical failure modes that can occur in a composite after initial elastic deformation. These are:     

Fibre failure, fracture, and, for aramids. Defibrillation. Resin crazing, microcracking and gross fracture. Debonding between the fibr and matrix. Dlamination of adjacent plies in a laminate. Fibre pull out from the matrix and stress relaxation.

In practice, debonding and delamination are not always distinguished.

8

Figure 2.2 Fracture mechanisms observed in laminates [3]

Figure 2.3 Fracture propagation in the case of poor fibre-matrix bonding [3]

9

2.3.1 Why are composites prone to impact damage? [1] The lack of plastic deformation in composites means that once a certain stress level is exceeded permanent damage, resulting in local or structural weakening, occurs. Unlike a metal, which may undergo plastic deformation by can retain its integrity (e.g. water tightness), composites stressed above a certain level, though possibly retaining some structural properties, are permanently damaged. A blow with an energy of ~1J or less at ~2m/s can cause irreversible damage in a realistic composite plate. To summarise, the reasons for a low impact strength are:   

Low transverse and interlaminar shear strength. Laminar construction, which is required if the reinforcing fibres are to be used efficiently and anisotropy reduced. No plastic deformation.

2.3.2 Estimation of the stresses acting in an impact [1] It is informative to try to estimate the stress generated by an impact. The situation is very complex and our simple approach merely indicates the order of magnitude of the stress. Consider a steel ball of radius ~25mm and mass 0.5kg, with an energy of 1J impacting a surface at a velocity, V, of 2m/s. The compression stress generated, σc, is, given by

 c  V ( Ect  )0.5 where Ect is the transverse compressive modulus and ρ the density. Typical values of these quantities are 10GPa and 2x103kg/m3, respectively, giving σc=9MPa. In pratice, compressive loading could cause localised bending and transverse failure in suitably oriented plies and generate interlaminar shear stresses. The compressive stress will continue to act until the impacting body rebounds. The stress wave will propagate into the impacted material, with a decreasing amplitude, and be reflected at the back surface as a tension wave. It is possible to apply the theory of Hertzian contact stresses to the estimation of the impact stress. This indicates that if the diameter of the contact area is 2a the maximum shear stress appears below the point of contact at a depth of ~a/2. The maximum contact stress, q0, for indentation on a plate surface is

q0  0.7aE / R This assumes a common Poisson’s ratio of 0.3 and that the controling modulus is that of the composite through its thickness ~10GPa. If a=0.25mm and R, the radius of the impactor is 25mm as before, the maximum compressive stress is found to be ~70MPa and the maximum shear stress, ~0.31q0, is 22MPa. Both approaches neglet, among other 10

points, the differences in the properties of the impacting and impacted bodies, geometrical effects and the heterogeneous, anisotropic nature of the composite. In addition, whether a static thoery of contact stresses can be applied to an impact situation is open to debate. In view of all these factors the results are approximate but they do indicate the potential for damage.

2.4

Observation of fracture mechanisms The observation of the fracture mechanisms in laminates can be carried out by different techniques, We give hereafter some basic elements on these techniques. 2.4.1 Observation by Microscopy [3] Optical observation with a microscope is a very simple technique to carry out for the continuous observation of fracture mechanisms during tests. However, this technique is restricted to local observation and the depth field is limited. Scanning electronic mocroscopy increases the depth field, allowing high magnifications to be obtained. Figures and show the microsgraphs obtained in the case of transverse fracture of composites with poor fibre-matrix bonding (Figure 2.4) and high bonding (Figure 2.5).

Figure 2.4 Fracture surface associated with poor fibre/matrix bonding in the case of a carbon fibre composite (ONERA document) [3]

11

Figure 2.5 Fracture surface associated with high fibre/matrix bonding in the case of a carbon fibre composite (ONERA document) [3]

2.4.2 Radiography analysis [3] The technique of analysis by X-radiography consists in impregnating the test specimens by means of opacifying agent (as zinc iodide) and then taking an X-radiograph of the test specimens. Radiography gives a two-dimensional image of the fracture stage (Figure 2.6). It is, however, easy to localize the damage in the body of laminates when one knows the orientations of the layers. Raidography allows a very fine observation of the cracks, and of the cracks transverse to the thickness of the laminates in particular. It should be noted that it is necessary to demount the test specimen for each radiography, and then to remount it in the testing machine in order to carry on the test. This makes the tests considerably time consuming.

Figure 2.6 X-ray observation of the fracture state of carbon fibre composites after fatigue (105 cycles) in the case of impacted specimens (ONERA document) [3] It is also possible to observe the fracture state of test specimens by radiography with a medical scanner. The analysis of the density variations allows us to obtain information in three dimensions. 12

2.4.3 Acoustic emission analysis [3] The preceding techniques permit observations at different times. They are also time consuming to carry on because of the mounting and demounting of the test specimens, necessary for the observations on the fracture state. In contrast, acoustic emission is a physical process which allows us to access, in real time, information about the fracture mechanisms at they happen. When a fracture mechanism is induced inside a material, it creates a local discontinuity, called an event, generates a strain wave which propagates through the material. At the surface of material, an adapted transducer converts the wave received into an electric signal (the acoustic emission signal) which is next amplified, then analysed. The transducers are piezoelectric transducers, developed specifically for acoustic emission so that they have a high sensitivity. The frequency domain studied generally extends from 50kHz to 1MHz

Figure 2.7 The acoustic emission process [3]

Figure 2.8 Acoustic emission signals recorded during bending tests on unidirectional carbon fibre-epoxide composites [3]

13

The processes of recording and analysis of the acoustic emission signals are greatly improved by the numerical equipment which is available today for the engineer.

2.5

Failure criteria The objective of the failure criteria is to allow the designer to have an evaluation of the mechanical strength of laminates. Quite generally, the mechanical resistance of a material is associated to an irreversible degradation. In fact, the definition of failure may change from one application to another. In case of composite materials, the end of the elastic domian is generally associated with the development of microcracking: matrix microcracking, fibre-matrix debonding, etc,. In the initial stage of fracture process, the initiated microcracks do not propagate, and their development changes the stiffness of the material very gradually. Failure criteria have been established in the case of a single layer of a laminate and may be classified as:   

The criterion of maximum stresses, The criterion of maximum strains, The interactive criteria, usually called as energy criteria.

In some applications, the maximum stress and strain criteria do not allow us to describe the experimental results observed. Moreover, these criteria do not consider interactions between the modes of fracture: longitudinal, transverse and shear fracture. So, the different fracture mechanisms are assumed to occur independently. Thus, the interactive criteria have been investigated, extending to the case of orthotropicmaterials the Von Mises’ criterion introduced for the isotropic materials. Von Mises’ criterion is related to the deformation energy stored per unit volume of the strained material. These interactive criteria are    

Hill’s criterion, The Tsai-Hill criterion, Hoffman’s criterion, Tsai-Wu criterion.

In this paper, we introduces the Hashin criterion which predicts anisotropic damage in elastic-brittle materials. It is primarily intended for use with fiber-reinforced composite materials and takes into account four different failure modes: fiber tension, fiber compression, matrix tension, and matrix compression. According to the Hashin criterion [14]: The initiation criteria have the following general forms: Fiber tension ( ˆ11  0 ): 14

 ˆ   ˆ  Fft   11T     12L  X  S  2

2

Fiber compression ( ˆ11  0 ):

 ˆ  F   11C  X 

2

c f

Matrix tension ( ˆ 22  0 ):

 ˆ   ˆ12  Fmt   22  L  T  Y  S  2

2

Matrix compression ( ˆ 22  0 ): 2 2 2 C  ˆ  ˆ 22   Y   ˆ12  F   T    T   1 22   L C  2S   2S  S   Y c m

In the above equations: XT denotes the longitudinal tensile strength (fiber tensile strength); XC denotes the longitudinal compressive strength (fiber compressive strength); YT denotes the transverse tensile strength (matrix tensile strength); YC denotes the transverse compressive strength (matrix compressive strength); SL denotes the longitudinal shear strength (longitudinal shear strength); ST denotes the transverse shear strength (transverse shear strength); α is a coefficient that determines the contribution of the shear stress to the fiber tensile initiation criterion. The initiation criteria presented above can be specialized to obtain the model proposed in Hashin and Rotem (1973) by setting   0.0 and ST  Y C / 2 or the model proposed in Hashin (1980) by setting   1.0 .

ˆ11 , ˆ 22 , ˆ12 are components of the effective stress tensor, ˆ , that is used to evaluate the initiation criteria and which is computed from: ˆ  M  , where  is the nominal stress and M is the damage operator:

 1 1  d f   M  0    0 

0 1 1  dm 0

 0    0   1   1  ds 

15

d f , d m , and d s are internal (damage) variables that characterize fiber, matrix, and shear t

c

damage, which are derived from damage variables d f , d f

t

c

, d m , and d m ,

corresponding to the four modes previously discussed, as follows:

d tf if ˆ11  0  df   c  d f if ˆ11  0  d mt if ˆ 22  0  dm   c  d m if ˆ 22  0  d s  1  (1  d tf )(1  d cf )(1  d mt )(1  d mc ) Prior to any damage initiation and evolution the damage operator, M, is equal to the identity matrix, so ˆ   . Once damage initiation and evolution has occurred for at least one mode, the damage operator becomes significant in the criteria for damage initiation of other modes. The effective stress, ˆ , is intended to represent the stress acting over the damaged area that effectively resists the internal forces. An output variable is associated with each initiation criterion (fiber tension, fiber compression, matrix tension, matrix compression) to indicate whether the criterion has been met. A value of 1.0 or higher indicates that the initiation criterion has been met. In addition, we introduce an interlaminar delamination criterion to predict whether the delamination occurred in the impact event. It is described below [6]: Interlaminar delamination appears when the following criterion is reached: 2

2

2

  33    13    23           1 ( 33  0)   33           This failure criterion includes all stress components ( 33 , 13 and  23 ) on a potential fracture plane (1-0-2). Each component produces its own damage fraction in this failure mode. Therefore there is an interaction between the tensile stress  33 normal to a layer plane and the in_plane shear stresses  13 and  23 . The tensile stress  33 causes two adjacent layers to move realtively along the thickness direction. The in-plane shear stress

13 and  23 cause the two adjacent layers to move relatively (parallel to the layer plane). This interaction can lead to interlaminar delamination. This above equation describes an ellipsoidal failure envelope (give a value of the right-hand side equal to unity). Any point inside the envelope shows no failure in the material. The composite structure does not

16

support any load locally along these three directions when the failure occurs. Hence the local stiffness matrix is modified to simulate the situation.

Figure 2.9 Delamination [7]

Figure 2.10 Ultrasonic scanning images of the dlemination in cross-ply laminates after impact: (a) after 1J impact; (b)after 2J impact; (c) after 3J impact (ELSEVIER document) [8]

17

Chapter 3 Contact mechanics [5] Abstract Impacts of foreign objects on composite structures can create internal damage that reduces the strength of the structure significantly. The study of such impacts requires understanding the dynamics of the event, predicting the extent of the induced damage, and estimating the residual properties of the structure. The impact event involves the motion of the target, the motion of the projectile, and the local indentation in the contact zone. A large number of parameters affect the impact dynamics, many types of responses can be obtained, and many models have been proposed in the literature. These models can be classified into three categories: (1) energy-balance models that assume a quasistatic behavior of the structure; (2) spring-mass models that account for the dynamics of the structure in a simplified manner; (3) complete models in which the dynamic behavior of the structure is fully modeled. Simple models can bring insight into the problem and be efficient but have limited applicability. Complex models may have wider applicability but require significantly higher modeling and computational effort. There is a need for a general understanding of the impact dynamics and for a method for developing efficient and accurate models.

Figure 3.1 impact mechanism Keywords: Impact dynamics; Models; Laminated composite materials.

18

3.1

Introduction A first step towards understanding the effect of impacts is to develop a model for predicting the contact force history and the overall response of the structure. It involves modeling the motion of the projectile, the dynamics of the structure, and the local indentation of the structure by the projectile. Experimental modal analysis showed that low velocity impact damage has only minor effects on the dynamic properties of laminated plates. Small shifts in the natural frequencies of higher-order bending modes are observed which confirms that damage needs not modeled in the impact dynamic analysis. Some impacts produce deformations in a small zone surrounding the point of impact while others involve deformations of the entire structure. In some cases, a major portion of the impact energy is transferred to the plate and in other cases most the energy is restituted to the projectile. For some problems, the indentation absorbs a significant portion of the impact energy so that it must be modeled adequately in the analysis. In other cases, the effect of indentation are negligible. Sorting out these different types of behavior is necessary for the interpretation of experimental results and for the selection of an appropriate mathematical model. The objectives of this chapter are to study the various models available for analyzing the impact dynamics and to present an approach for selecting an appropriate model for each particular case. The many models used to study the impact dynamics are classified here according to how the structure is modeled: spring-mass models, energy balance models, complete models, and a model for impact on infinite plates.

3.2

Contact mechanics Local elastic deformation properties have been considered as early as 1880 with the Hertzian Theory of Elastic Deformation. This theory relates the circular contact area of a sphere with a plane (or more general between two spheres) to the elastic deformation properties of the materials. In the theory any surface interactions such as near contact Van der Waals interactions, or contact Adhesive interactions are neglected. An improvement over the Hertzian theory was provided by Johnson et al. (around 1970) with the JKR (Johnson, Kendall, Roberts) Theory. In the JKR-Theory the contact is considered to be adhesive. Hence the theory correlates the contact area to the elastic material properties plus the interfacial interaction strength. Due to the adhesive contact, contacts can be formed during the unloading cycle also in the negative loading (pulling) regime. Such as the Hertzian theory, the JKR solution is also restricted to elastic spheresphere contacts. A more involved theory (the DMT theory) also considers Van der Waals interactions outside the elastic contact regime, which give rise to an additional load. The theory simplifies to Bradley's Van der Waals model if the two surfaces are separated and 19

significantly appart. In Bradley's model any elastic material deformations due to the effect of attractive interaction forces are neglected. Bradley's non-contact model and the JKR contact model are very special limits explained by the Tabor coefficient.

Figure 3.2 contact mechanical models [4]    

Hertz: fully elastic model, JKR: fully elastic model considering adhesionin the contact zone, Bradley: purely Van der Waals model with rigid spheres, DMT: fully elastic, adhesive and Van der Waals model.

Local deformation in the contact zone are not modeled with beam, plate or shell theories since those theories usually assume that the structure is inextensible in the transverse direction. However, in many cases, local indentation has a significant effect on the contactforce history and mush be accounte for in the analysis. The contact phenomenon is recognized as being rate independent for most laminated formost laminated composite materials and statically determined contact laws are used by most investigators. During the loading phase of the impact, the contact force P is related to the indentation α by

P  k 3/2 .

Figure 3.3 definition of parameter R (1)

The contact stiffness is given by k

4 1/2 ER , 3

(2)

20

where the parameters R and E are defined as

1 1 1 1 1  12 1  22   and  ,  R R1 R2 E E1 E2

(3)

where R1 and R2 are the radii of curvature of the two bodies. The Young’s moduli and Poisson’s ratios of the two bodies are E1, ν1 and E2, ν2, respectively. Subscripts 1 denotes properties of the indentor, while subscript 2 identifies properties of the target. Eq. (1) is usually referred to as the Hertzian law of contact. Permanent indentation occur even at relatively low loading levels, and the unloading phase of the process is significantly different from the loading phase. During unloading, the contact law is

P  Pm    0  /  m   0  , 2.5

(4)

where Pm is the maximum force reached before unloading, αm the maximum indentation, and α0 is the permanent indentation. α0 is zero then the maximum indentation remains below a critical value αcr. When αm>αcr,  0   m 1  ( cr /  m )2/5  .

(5)

Typically , during impact the contact force increases to maximum value and then decreases back to zero. In some cases, multiple impacts and reloading occurs. During subsequent reloading, the reloading curve is distinct from the unloading curve but always returns to the point where unloading began. The unloading curve is modelled by

P  Pm    0  /  m   0 

3/2

(6)

21

Figure 3.4 indentation of glass-epoxy laminate from Yang and Sum (1982). (a) Loading; (b) Unloading; (c) Reloading In some cases, only a small fraction of the impact energy is used in the local indentation process and therefore it it not necessary to distinguish between loading ang unloading branches, and Eq.(1) is used throughout the indentation process. Sometimes the indentation preocess need not be modeled at all.

3.3

Energy-balance models One approach for analyzing the impact dynamics is to consider the balance of energy in the system. The initial kinetic energy of the projectile is used to deform the structure during impact. Assuming that the structure bahaves quasi-statically, when the structure reaches its maximum deflection, the velocity of the projectile becomes zero and all the initial kinetic energy has been used to deform the structure. Therefore, the energy balance equation can be written as 1 MV 2  Eb  Es  Em  Ec 2

(7)

where the subscrifts b, s, m refer to the bending, shear, and membrane components of the overall structural deformation and Ec is the energy stored in the contact region during indentation.

22

When the overall deflections of the structure are neglibible compared to the local indentation, the problem is reduced to that of an impact on a half-space and the maximum contact force, and the contact duration are given by 5 P  4

3/5

1/5

 M 3V 6 k 2 

(8a)

1/5

M2  Tc  3.2145  2  Vk 

(8b)

These simple expressions show the effect of the projectile mass and velocity and the contact stiffness on the contact force.

3.4

Impact on infinite composite plates Here to present an approximate solution for impacts on infinite plates. With an Hertzian contact law (Eq.(1)), non-dimensional indentation and time variables can be defined as



 TV

, t 

t , T

(9)

where T   M / (kV 1/2 ) 

2/5

.

(10)

The non-dimensional indentation is governed by the single, nonlinear. Ordinary differential equation

d 2 3 d    1/2   3/2  0 2 dt 2 dt

(11)

which depends on a single non-dimensional parameter

  kc2/5V 1/5 M 3/5 / 8 mD*  

(12)



Called the inelasticity parameter. The equivalent bending rigidity in Eq. (12) is defined as A 1 ( D11 D22 )1/2 2 A=(D12  2 D66 ) / ( D11D22 )1/2 D* 

(13)

Eq. (11), with the initial conditions

 (0)  0,

d (0) 1 dt

(14)

23

must be solved numerically. The contact force is given by 1/5

P   k 2 M 3V 6   3/2

(15)

After contact ceases, the deflection at the point of impact remains constant while the deformation propagates outward. During an impact on an infinite plate, the wavefront is nearly elliptical. The distance a between the wavefront and the impact point in the xdirection, can be estimated using 1/4

D  a  2   11   m 

 2( A  1)

1/8

t

(16)

In the y-direction, the wavefront is approximately located at a distance

b  a / ( D11 / D22 ) 1/4

(17)

Eqs. (16) and (17) provide an estimate of the size of the deformed zone during impact.

3.5

Spring-mass models 3.5.1 General case Spring-mass models are simple and provide accurate solutions for some types of impacts often encountered during tests on small size specimens. The most complete model consists of one spring representing the linear stiffness of the structure (Kbs), another spring Km for the nonlinear membrane stiffness, a mass M2 representing the effective mass of the structure, the nonlinear contact stiffness, and M1 the mass of the projectile. If the effect of shear deformation is negligible, the spring constant Kbs is replaced by Kb which account for bending deformations only. From the free body diagrams of the two masses M1 and M2, the equations of motion of the system can be written as ..

M 1 x1  P  0

(18a)

..

M 2 x2  Kbs x2  K m x23  P  0

(18b)

where P is the contact force which is a highly nonlinear function of the indentation x1 – x2.

24

Figure 3.5 (a) two-degree-of-freedom model. (b) single-degree-of-freedom model. .

.

The initial conditions x1  V , x 2  0, x1 (0)  x2 (0)  0 (19) can be studied numerically. The stiffness used in the spring-mass models can be determined from formulas available in many handbooks or, numerically, using the finite element method for example [Shivakumar K N. Prediction of impact force and duration due to low-velocity impact on circular composite laminates.]. For a completely clamped, isotropic, circular plate, the bending and membrane stiffnesses are given by 4 Eh3 (353  191 ) Eh Kb  and K m  2 2 3(1  )a 648(1  )a 2

(20&21)

where E is the elastic modulus, m the Poisson's ratio, h the thickness, and a is the radius of the plate. The effective mass of the plate is taken as one-fourth of the total mass of the plate [Olsson R. Impact response of orthotropic composite plates predicted from a oneparameter differential equation]. 3.5.2 Linear two-degree-of-freedom model In many cases, the transverse deflections are small and membrane-stiffening effects are negligible and the stiffness of the structure can be represented by the linearspring k2. In order to understand the dynamics of the impact, we also assume that the local indentation can be represented by a linear spring k1 (see Figure 3.6). The motion of the linear twodegree-of-freedom system is governed by

25

..

m1 x1  k1 ( x1  x2 )  0

(22)

..

m2 x 2  k2 x2  k1 ( x2  x1 )  0

Figure 3.6 Linear two-degree-of-freedom spring-mass impact model with the initial conditions .

.

x1 (0)  x2 (0)  0, x1 (0)  V , x2 (0)  0

(23)

Introducing the non-dimensional variables k1 and yi  xi m1

 t

k1 / m1 V

(24)

Eqs. (22) can be written as y1''  y1  y2  0 y''2 

(25)

m1 m  k  y1 + 1 1  2  y 2  0 m2 m2  k1 

And the initial conditions become

y1 (0)  y2 (0)  0, y1' (0)  1, y'2 (0)  0

(26)

The behavior of the system depends on the two non-densional parameters m1/m2 and k2/k1. The non-dimensional natural frequencies of the system are obtained by solving the biquadratic equation

 m  k  m k  4   2 1  1 1  2    1 2  0  m2  k1   m2 k1

(27)

26

The natural frequencies can be approximated as indicated in Table in which we distinguish four special cases. Figure 3.7(a) shows that as k2/k1 becomes large, the first natural frequency tends to the limiting value of one. The convergence rate is strongly affected by the massratio when m1/m2<<1. When k2/k1 is larger than one, the second nondimensional frequency increases with k2 / k1 and the validity of the approximations given in Table is shown in Figure 3.7(b). Figure 3.7(c) shows a ―curve veering‖ phenomenon for the two natural frequencies as m1/m2 increases for a fixed k2/k1 ratio. Table 3.1 approximate formulas for natural frequencies of two-degree-of-freedom system

For case I ((m1/m2)<<1, (k2/k1)<<1), the natural frequencies are approximated by k k 1  2 , 2  1 (28) m2 m1 and the free vibration modes are as shown in Figure 3.8(a). With mode I, the two masses move together but m1 is negligible compared to m2. For mode II, the heavy mass remains stationary and oscillates as a single-degree-of-freedomsystem with stiffness k1. For case II ((m1/m2)>>1, (k2/k1)<<1), the mass of the target is small compared to the mass of the impactor and the stiffness of the target is small compared to the contact stiffness and the natural frequencies are approximated by k k 1  2 , 2  1 (29) m1 m2 In mode I, the indentationis negligible and the projectile and the target form a SDOF system (Figure 3.8(b)). In model II, the projectile remains stationary, the stiffness af the target is begligible, and the target oscillates as a SDOF system. For case III ((m1/m2)<<1, (k2/k1)>>1), the stiffness of the target is much larger than the contact stiffness, the mass of the projectile is small compared to that of the target. Then k2 k1 1  , 2  (30) m2 m1 27

The mode shapes are shown in Figure 3.8(c). For case IV ((m1/m2)>>1, (k2/k1)>>1), the stiffness of the target is much larger than the contact stiffness and

1 

k1 k2 , 2  m1 m2

(31)

28

Figure 3.7 Free vibrations of two-degree-of-freedom sustem: (a) first nondimensional frequency; (b) second non-dimensional frequency (solid line: exact; dashed line; approximate); (c) Curve veering phenomenon for k 2/k1=10 (solid line: exact; dashed line; approximate). [5]

Figure 3.8 Apprpximate mode shapes of linear two-degree-of-freedom system. [5] For mode I,the target remains stationary and the projectile oscillates (Figure 3.8(d)). For mode II, the motion of the plate is unaffected by the presence of the projectile due to the weak coupling provided by the contact stiffness.

29

3.5.3

Linear single-degree-of-freedom (SDOF) spring-mass model

A significant simplification occurs when membrane-stiffening effects are negligible and the indentation is small compared with the overall deformation of the structure. In that case, a linear SDOF system (Fig. 3.5(b)) can be used and the contact force is then given by P  V K bs M sin





(32) Eq. (32) predicts that the maximum contact force is directly proportional to the initial velocity of the projectile. Similarly, the maximum contact force is proportional to the square root of the kinetic energy. The contact force increases with the square root of the stiffness of the structure and the square root of the mass of the impactor. A stiffer structure will cause a harder impact and, for the same initial velocity, a larger mass will have a larger kinetic energy which will also increase the contact force. K bs M t

3.5.4 Nonlinear SDOF spring-mass models There are two situations for which a nonlinear SDOF can provide accurate predictions of the contact force history. In the first instance, the overall deflection of the structure is negligible compared to the local indentation. In that case, the spring in Fig. 3.5(b) represents the contact stiffness, and the equation of motion is ..

(33) M 1 x kx13 / 2  0 The second situation for which nonlinear a SDOF model can yield accurate predictions of the contact force history is when the local indentation is negligible but the deflections of the structure becomes large and membrane stiffening is significant. The equation of motion of the single-degree-of-freedom model is ..

M x kb x  km x3  0

(34) A numerical solution of these nonlinear equations of motion (Eqs. (33) and (34)) yields the dynamic response and contact force histories.

3.6

Complete models With a complete model, the dynamic behavior of the structure is described accurately. This means that the appropriate structural theory is used. For example, in many cases the classical plate theory can be used but, in some cases, transverse shear deformations become significant and higher-order theories must be used. If the initial velocity of the projectile is sufficiently large, damage can be introduced plate motion is established. In that case, a three-dimensional analysis is required. Once an appropriate theory is selected, all the vibration modes participating in the response have to be predicted accurately and must be retained in the model. For a simply supported plate, for example, an analytical solution can be found for the natural frequencies and mode shapes. The transient response is then expressed in terms of these mode shapes and all participating modes can be included. For other geometries or boundary conditions, variational or finite element models must be used. With such approximate methods, a sufficiently large number of degrees of freedom must be selected so that the participating modes are predicted accurately. 30

If N equations are needed to describe the motion of the structure and one equation for the projectile, the N + 1 differential equations can be written in matrix form as

M  X   K X   F  ..

(35)

 

and integrated using Newmark's step-by-step time integration method. The contact force is unknown and is a nonlinear function of the indentation. Therefore, the force vector in Eq. (35) is assumed to be known at the end to the nth time step but its value at the end of step n + 1 is unknown. In order to determine the displacements at the end of step n . 1, we start assuming that Pn1  Pn and solve equation (35) for a first estimate of X n1 . A new estimate of Pn1 can be calculated from these displacements and a new iteration can be performed. After several iterations, the solution converges so that both the equations of motion and the contact law are satisfied and the process is repeated for the following time steps.

3.7

Selection of an impact dynamics model The selection of an appropriate impact dynamics model starts by neglecting the deflections of the plate and using the energy-balance model (Eqs. 8) to obtain a first estimate of the maximum contact force and of the contact duration. The infinite plate model is used to determine the effect of plate deflections and the type of response of the target. The inelasticity parameter  determines how much of the impact energy is absorbed by the deformation of the plate. When  << 1, the deformation of the plate is negligible and the energy balance approach is satisfactory. Otherwise, the dynamics of the plate play an important role. If the deformation front has not reached the boundary of the plate, the infiniteplate model is a very efficient and accurate way of analyzing the impact. When the deformation front propagates many times more than the distance from the impact point to the boundary, the finite size of the target must be accounted for. In the examples discussed here as well as in many examples in the literature, the contact deformation is modeled using Hertz's law for both loading and unloading. However, contact laws have been investigated at length and some authors account for different behavior in the unloading and reloading phases (Eqs. 4, 6). When  is small, when the contact force reaches its maximum, most of the energy is used to indent the target so that energy losses due to indentation can be significant. Therefore, when  is large, the indentation uses only a small fraction of the impact energy so that it is not necessary to model the contact behavior carefully.

31

3.8

Conclusions This chapter presents an overview of mathematical models used for the analysis of the dynamics of impacts between a foreign object and a composite structure. Currently available models are classified into four categories: Spring-mass models, energy-balance models, complete models, and an impact on infinite plate model. Simple models are easy to use and efficient but have limitations due to the simplifying assumptions on which they are based. Selecting an appropriate model requires on understanding of the effects the many factors affecting the impact dynamics. In this chapter, a procedure is presented in order to determine the type of impact to be expected and to select an appropriate model. The process starts by assuming that the plate is infinite. In this case, a simple approximate solution requiring to solve a single nonlinear ordinary differential equation is available. This model predicts the contact force history and the overall deformation of the plate. Then we can determine whether the deformation reached the plate boundaries during the duration of the impact. If bending waves travel from the impact point to the edge of the plate and back many times during the predicted contact duration, we have a boundary controlled impact and a spring mass model or an energy balance approach might be adequate because in that case, the plate behaves in a quasi-static manner. If the deformation never reaches the edges of the plate, we have a wave-controlled impact and the approximate solution provides very good results. For intermediate cases, the infinite plate model might be adequate initially but reflected waves will affect the contact force history. Then, a complete model taking into account the full dynamic behavior of the plate and the boundary conditions will be necessary.

32

Chapter 4 Modelling with Abaqus/CAE 6.9 4.1

Introduction The numerical simulation of impact damage to the composite plate was performed using the finite element package Abaqus/CAE 6.9. There are two models in this simulation. They are named as impact model and czm. The impact model is used to investigate the stress field during and after impact process for prediction damage. The czm model is used to simulate interlaminar delamination due to matrix cracking, which is the most common damage of the laminated composite. The czm model is more exact because it is implemented cohesive zone model to take into account the interaction between failure modes. Both of them simulate the low-energy impact on laminated composite. Velocity of impactor is 6.26m/s with mass of the impactor is 0.1kg, so the energy impact is 1.96J. Experimental set-up For the composite tests, 16-ply symmetric crossply composite plates (04/904)S, were used. The laminated plates used were 2mm thick with circular plate specimens having an effective diameter of 100mm could be cut from each fabricated plate. A controlled dropweight rig shown in Fig. 4.1 was employed to perform the impact test. A steel ring was bolted on the platform of the test rig. The composite specimen was clamped using eight equally spaced bolts through an upper ring into the fixed lower steel ring in order to simulate fixed boundary conditions. The maximum height of the rig was 3.75 m which can produce approximately 8.57 m/s initial velocity on impacting the target. Dropping heigh (m) Impact velocity (m/s) Impact energy (J) 1.00 4.43 0.98 1.25 4.95 1.26 1.50 5.42 1.47 1.75 5.86 1.72 2.00 6.26 1.96 2.25 6.64 2.20 2.50 7.00 2.45 2.75 7.34 2.69 3.00 7.67 2.94 3.25 7.98 3.18 3.50 8.28 3.43 3.75 8.57 3.67 Table 4.1 Impact energy of projectile 0.1kg of mass

33

Figure 4.1 The contronlled drop-weight rig [6]

4.2

Simulation procedure of the impact model 4.2.1 Model For this simulation, 16-ply symmetric crossply composite plate (04/904)s was used. The laminated plates were 2mm thick with circular shape having an diameter of 100mm. Because of the symmetry of the circular plates, it was necessary to model only a quarter of the plates. In addition, the appropriate symmetrical boundary conditions is needed to apply. The impactor shaped ball having an diameter of 10mm was modelled as a rigid body. Its potential contact area was taken as a ―Master surface‖, and the potential contact area of the composite plate was taken as a ―Slave surface‖. Under these definitions, the individual nodes in the composite plate were constrained not to penetrate the surface of the impactor whereas the nodes on the surface of the impactor could penetrate segments of the surface of the composite plate.

34

Figure 4.2 modelling of composite plate and impactor

Figure 4.3 Material orientations of the composite plate The figure shows that the fibre direction is the 1-axis and the in-plane direction that is perpendicular the fibre direction is the 2-axis. 4.2.2 Material properties of composite The plate was made of composite CFRP of density 1600kg/m3 which includes 16 plies [04/904]s, each ply has a thickness of 0.125mm.

35

Mechanical properties of laminates E11 E22 E33 G12 G13 G23 ν12 ν13 ν23

145Gpa 9.2Gpa Assuming = E22 4.6GPa 5.2Gpa 3.0Gpa 0.3 Assuming = ν12 Assuming = ν13

To model the initiation and the propagation of the damage in the composite plate, we use the Hashin criterion by identifying the strength of composite and the damage variables in type fracture energy. Strength of composite

Longitudinal tensile strength Longitudinal compressive strength Transverse tensile strength Transverse compressive strength Longitudinal shear strength Transverse shear strength

3600MPa 1500MPa 72MPa 102MPa 41MPa 41MPa

Damage evolution Longitudinal tensile fracture energy Longitudinal compressive fracture energy Transverse tensile fracture energy Transverse compressive fracture energy

12.5N/mm 12.5N/mm 2.5N/mm 2.5N/mm

For analysis with Hashin criterion, we must assign shell section for the part composite plate. Especially, to investige the shear stress field, Abaqus need to define thickness integration points. As a result, in this analysis we set 5 integration points along the thickness of shell element.

36

4.2.3 Boundary conditions Since the model in this simulation is only a quarter of the whole model, we need to apply appropriate boundary conditions to achieve the reasonable response as the whole model. We shall apply fixed condition at curved boundary of composite plate, that contraints six DOFs of this curved boundary not to permit all the displacements there. Besides, we apply four symmetrical boundary conditions to the impactor and the composite plate. We also assign a condition to contraint vertical motion of the impactor.

Figure 4.4 Boundary conditions 4.2.4 Mesh For the impactor, it is modeled as rigid body, so we mesh it with discrete rigid element. For this mesh, it need not to be fine compared with the composite plate. We discriminate the impactor to 180 triangle elements. With the composite plate, whose mesh is important because it affect to the consistency and the convergence of model. If the mesh is not good enough, it may even cause the divergence of the result as some elements are distorted during calculated process. As a result, we need using an appropriate mesh. In this analysis, we mesh the composite plate regularly with total 6720 elements, including 6400 8-node quadrilateral continuum shell elements and 320 6-node triangular in-plane continuum shell elements. We must notice the mesh stack orientation in according to Abaqus’s recommendation so that it complies with the thickness-direction of the composite plate.

37

Figure 4.5 The mesh The mesh of composite plate includes total 6720 elements which generated by 7514 nodes. These elements are 6400 SC8R and 320 SC6R. The mesh of impactor includes 180 elements which generated by 92 nodes. These elements are R3D3. 4.2.5 Technique for simulation the impact This model is used to simulate the impact between the impactor and the composite plate according to the drop-weight experiment. To simplify, we do not simulate the drop process of the impactor. Instead of this, we assign the impact velocity and set position exactly at the impact point to the impactor. The impact velocity is 6.26m/s in this model. To model the impact between the impactor and the composite plate, we need to define the interaction properties between them. Here is an interaction with behavior of hard contact in normal direction and damping coefficient of 0.0152. The damping coefficient of the composite plate is needed to determine the system damping matrix in the impact simulation. In experiment, this damping coefficient can be measured using a free vibration decay method. And hard contact implies that:  



the surfaces transmit no contact pressure unless the nodes of the slave surface contact the master surface; penetration is not allowed at each constraint location (depending on the constraint enforcement method used, this condition will either be strictly satisfied or approximated); there is no limit to the magnitude of contact pressure that can be transmitted when the surfaces are in contact.

The hard contact is the most common contact pressure-overclosure relationship which is shown in below figure. When surfaces are in contact, any contact pressures can be 38

transmitted between them. The surfaces separate if the contact pressure reduces to zero. Separated surfaces come into contact when the clearance between them reduces to zero.

Figure 4.6 pressure-overclosure relationship of hard contact [14] 4.2.6 Jobs of simulation We create a job for this simulation. In this model, the time of impact is about 0.6ms, so we set the time for simulation as 1ms. This time includes both impact process and afterimpact process.

4.3

Simulating the delamination using cohesive zone model We can create a model using cohesive elements to model the following cases:  

Adhesive joints – two components connected by a glue which is like material having a finite thiskness. Fracture at bonded interfaces – crack propagation in glue material that is very thin and for all practical purposes may be considered to be of zero thickness.

The two main approaches that we can use to include cohesive elements in our model are:  

Embedding one or more layers of cohesive elements in the mesh of an existing model; or Creating the analysis model using the geometry and mesh tools.

We can model the connection at the interface between the cohesive layer and the surrounding bulk material by sharing nodes or by defining a tie contraint. The tiecontraint approach allows you to model the cohesive layer using a finer discretization than that of the bulk material and may be more desirable in certain modeling situations. So in this analysis, we use creation of a model with cohesive elements using geometry and mesh tools, and define tie contraints between the cohesive layer and the surrounding bulk material. 39

4.3.1 Mesh the surrounding bulk material We can use any of the mesh tools to mesh the surrounding bulk material. In this model, we create a regular mesh. Especially, since the stress gradient at the vicinity of impact point is high, so we choose the fine mesh at this region. At the far region from the impact point, we choose the coarse mesh. 4.3.2 The cohesive layer properties In this model, since the prediction interlaminar delaminations from above model that is quite reasonable with the results of experiments, we create the cohesive layers with the length of 10mm. Therefore, we determine a fine mesh to these cohesive layers. We can just observe the delamination when damages occur in the cohesive elements, at that time these damaged cohesive elements will be removed. So the mesh should be as smooth as possible to observe. The cohesive layers are zero-thickness ones, so we select traction-separation for defining the constitutive behavior of the cohesive layer. With this behavior, Abaqus offers an elasticity definition that can be written directly in terms of the nominal tractions and the nominal strains. Both uncoupled and coupled behaviors are supported. For uncoupled behavior, each traction component depends only on its conjugate nominal strain, while for coupled behavior the response is more general. In the local element directions the stress-strain relations for uncoupled behavior are as follows:

tn   K nn    ts    t   t 

K ss

  n      s  Ktt   t 

The quantities tn, ts, and tt respectively represent the nominal tractions in the normal and the two local shear directions; while the quantities εn, εs, and εt represent the corresponding nominal strains. For coupled traction separation, behavior the stress-strain relations are as follows:

tn   K nn    ts    K ns t   K  t   nt

K ns K ss K st

K nt   n    K st   s  Ktt   t 

The stability criterion for uncoupled behavior requires that Knn>0, Kss>0, and Ktt>0. For coupled behavior, the stability criterion requires that:

40

K nn  0, K ss  0, K tt  0; K ns  K nn K ss ; K st  K ss K tt ; K nt  K nn K tt ;  K nn det  K ns  K nt

K ns K ss K st

K nt  K st   0 K tt 

In this model, we use uncoupled behavior

Figure 4.7 cohesive layer mesh

Figure 4.8 model with cohesive elements for delaminations 41

The mesh of each cohesive layer includes 420 elements wich generated by 884 nodes. These elements are 400 COH3D8 and 20 COH3D6. To obtain a successful FEM simulation using CZM, two conditions must be met: (a) The cohesive contribution to the global compliance before crack propagation should be small enough to avoid the induction of a fictitious compliance to the model, and (b) the element size must be less than the cohesive zone length. In this simulation, the parameters of cohesive properties are set as follows (as per experience in the simulation process): Elastic 106N/mm3 106N/mm3 106N/mm3

E/Knn E/Kss E/Ktt Quads damage Nominal stress normal-only mode Nominal stress first direction Nominal stress second direction

61MPa 68MPa 68MPa

Damage evolution Normal mode fracture energy Shear mode fracture energy first direction Shear mode fracture energy second direction

42

0.02J/mm 0.2J/mm 0.2J/mm

Chapter 5 Results 5.1

Overview Composites are generally strong, and have a reasonable impact resistance if the applied stress is in the fibre direction. In other directions they tend to be weak and to have a low impact resistance. Small secondary stresses in the transverse direction or unexpected stresses due to an impact, in a weak direction, can easily cause damage. There are five basic mechanical failure modes that can occur in a composite after initial elastic deformation. They are:     

Fibre failure, fracture. Resin crazing, microcracking and gross fracture. Debonding between the fibre and matrix. Delamination of adjacent plies in a laminate. Fibre pull out from the matrix and stress relaxation.

In this analysis, low-velocity impact, the damages are mainly resin crazing and interlaminar delamination. So this analysis focuses on investigating the stress field and uses the Hashin criterion for prediction the resin crazing that may cause delamination damage.

5.2

Observation of impact process Time of impact is about 0.28ms. After 0.28ms, the impactor separates from the composite plate. The observation shows that the impact process is separated into two periods, loading and unloading. In the loading, the impactor transfers energy to the composite plate by contact, and silmutaneously also reduces the velocity to zero. The unloading is right after the loading. During this period, the impactor receives energy from the composite plate, increases the velocity and reflects backward.

43

44

Figure 5.1 Impact process at t=0.00ms; 0.08ms; 0.16ms; 0.24ms; 0.32ms; 0.40ms; 0.48ms;

Figure 5.2 Velocity of impactor during and after impact process Process from beginning to approximated 0.18ms is loading process; Impactor reduces velocity to zero. Then after process is unloading process to 0.28ms when impact process is over.

45

5.2.1. Contact properties

Figure 5.3 Impact force at right impact point during and after impact process.

Figure 5.4 Contact stress at right impact point during and after impact process

46

Relation between contact force and indentation is showed in Fig. 5.5. The maximum indentation is about 0.55mm. The maximum contact force is 64.81N reached in loading process. This relation one difference with contact mechanics theory which reveals that contact force reaches the maximum value when maximum indentation occurs.

Figure 5.5 Relation between contact force and indentation during impact process including loading and unloading, t=0 to t=0.28ms.

5.2.2. Delamination process The impactor was beginning to separate from the composite plate at t=0.28ms, three layers of the plate started to separate. Though delamination has been occurred before.

47

Figure 5.6 Delamination process at t=0.28ms; 0.52ms;

5.3

Investigation of stress field The result shows that stress concentrates on the centre of impact and the vicinity. There is also a region at the bind inflicting the stress, the other regions infict the much smaller stress. Since the anisotropic property of composite, the stress at the layers that have the different fibre directions distributes at the different regions.

Figure 5.7 Stress field in the surface of impact at t=0.36ms

48

Figure 5.8 Mises stress distributes along thickness at t=0.08ms, distance from centre of impact is normalized by percent The Mises stress reached the maximum value at t=0.08ms, right after beginning of impact event. The maximum value occurs at a point which is below impact point, not right impact point. It is about 2.826GPa. We investigates shear stress in-plane 1-0-2 at interfaces. The region investigated is circular shape having 10mm of distance from centre of impact.

Figure 5.9 Shear stress in-plane 1-0-2 at first interface, t=0.04ms; distance from centre of impact to distance of 10mm is normalized by percentage.

49

Figure 5.10 Shear stress in-plane 1-0-2 at second interface, t=0.04ms; distance from centre of impact to distance of 10mm is normalized by percentage.

5.4

Resin crazing, fibre failure and prediction of delamination This modeling uses the Hashin criteria to predict damage for the composite plate. The parameters of the damage stabilization is set as follows: Viscosity Viscosity coefficient coefficient in the in the longitudinal longitudinal compressive tensile direction direction 10-3 10-3 Table 5.1 Viscosity cofficients

Viscosity coefficient in the transverse tensile direction 10-3

Viscosity coefficient in the transverse compressive direction 10-3

The results of the model using Hashin criteria shows that status of all the elements of the composite plate are active. In Abaqus, the status of an element is 1.0 if the element is active, 0.0 if the element is not, the value of this variable is set to 0.0 only if damage has occurred in all the damage modes.

50

Figure 5.11 Status of elements after impact, t=0.8ms. The experiments with the low-velocity impact shows that the damage has been identified as interlaminar delamination coupled with matrix failure. There is no fibre failure in 0.59J impact [6]. In this analysis, the energy is 1.96J impact, this paper investigates whether the fibre failure has occurred. For the fibre and matrix failure, the results of the Hashin criteria demonstrated that 5.4.1 Fibre failure The lower composite plate inflicts tension but there is no tensile fiber fracture initiation. The variable HSNFTCRT is of maximum of around 0.717 that is lower than 1.0 of limited value. In contrast, the upper composite plate inflicts so much compression caused by the impact. There is small region at the impact point which has begun to initiate fibre failure, the variable HSNFCCRT has reached 1 of limited value. After the initiation, the value 0.225 of DAMAGEFC demonstrates that the fibre failure has grown significantly. In conclusion, with 1.96J impact, the fiber occurs damage as compressive failure in the 1mm radius of circular region right at impact point.

51

Figure 5.12 Tensile fibre damage initiation after impact, t=0.8ms.

Figure 5.13 Compressive fibre damage initiation after impact, t=0.8ms.

52

Figure 5.14 Compressive fibre damage evolution after impact, t=0.8ms. 5.4.2 Matrix failure There are so much matrix damages occurred. State of the matrix damage is significant. The tensile matrix damage has begun at lowermost surface and propagated along thickness direction. While the compressive matrix damage has begun right at the impact surface and propagated along thickness direction. Region of matrix damage is observed as about 5.5mm from the centre of impact.

Figure 5.15 Tensile matrix damage initiation after impact, t=0.8ms. 53

Figure 5.16 Tensile matrix damage evolution after impact, t=0.8ms.

Figure 5.17 Compressive matrix damage initiation after impact, t=0.8ms.

54

Figure 5.18 compressive matrix damage evolution after impact, t=0.8ms. 5.4.3

Prediction of interlaminar delamination

Figure 5.19 Delamination For the composite material, damages are usually resin crazings and shear damages. The shear damages include interlaminar delamination, debonding between fibre and matrix. The fibre damage always occurs last. In particular, since the interlaminar fracture toughness of the laminate is inferior to its in-plane one, interlaminar delamination between the plies with different fiber orientations is easily induced by transverse loading 55

to the laminate. In this analysis, with 1.96J impact, there is small region that damaged of fibre failure. Therefore, the interlaminar delamination is sure to have been occurred. For this demonstration, this analysis investigates shear damage and the shear stress field in the interlaminars.

Figure 5.20 Shear damage evolution after impact, t=0.8ms. There is no normal stress σ33 in the interlaminars, so the stress components causing interlaminar delaminations are shear stress in the interlaminars. This analysis uses following interlaminar delamination criterion [6]: Interlaminar delamination appears when the following criterion is reached: 2

2

2

  33    13    23           1 ( 33  0)          33       This failure criterion includes all stress components ( 33 , 13 and  23 ) on a potential fracture plane (1-0-2). Each component produces its own damage fraction in this failure mode. Therefore there is an interaction between the tensile stress  33 , which is normal to a layer plane, and the in-plane shear stresses  13 and  23 . The tensile stress  33 causes two adjacent layers to move relatively along the thickness direction. The in-plane shear stress 13 and  23 cause the two adjacent layers to move relatively (parallel to the layer plane). This interaction can lead to interlaminar delamination. This above equation describes an ellipsoidal failure envelope (give a value of the right-hand side equal to unity). Any points inside the envelope show no failure in the material. In the first interface, the shear stress field shows that influence of the shear stress τ13 is larger than shear stress τ23. However, in the second interface, it is reverse. The graphics in 56

Fig. 5.5 and Fig. 5.6 show that state of the shear stress field in the vicinity region of the centre of impact overcomes the limited state of delamination, with 41MPa of interlaminar shear strength (the same value with the transverse shear strength). In conclusion, the damage has been identified as interlaminar delamination coupled with matrix failure from the simulation. The damages have occurred dramatically. It also has been revealed from the simulation that there is very small fiber failure at the impact point.

5.5

Delamination According to the above conclusion, the interlaminar delaminations have occurred. So here we use cohesive zone model to simulate these delaminations. Here is results of the cohesive zone model. There are delaminations at both interlaminars. State of the cohesive layer after impact shows shape of the interlaminar delamination. With the fine mesh of the cohesive layer, we can see this shape and evaluate the length of delamination quite accurately.

Figure 5.21 Delamination at both two interlaminars at t=0.52ms

57

Figure 5.22 State of first interlaminar (left) and second interlaminar (right) after impact, t=0.8ms

Figure 5.23 Traction damage initiation at first interlaminar (left) and second interlaminar (right) after impact, t=0.8ms

58

Figure 5.24 Interlaminar delamination shaped a peanut at second interlaminar, t=0.04ms

59

Figure 5.25 Interlaminar delamination at second interlaminar after impact, t=0.8ms

The results show that shape of delamination at the second interlaminar is a peanut. This shape of peanut is orientated along fibre-direction of the lower material layer, that is validated by the above results of simulation. Besides, shape of delamination at the first interlaminar is quite circular. This shape may be caused by the mesh, that is not fine enough to demonstrate the peanut shape of the interlaminar delamination. In addition, the delamination at the second interlaminar is stronger than the first one. The length of delamination at the second interlaminar is about 9.34x2mm, and at the first one is about 6.40x2mm. These length-delaminations are accepted with the result of experiments which have done by the others.

60

Chapter 6 Conclusion and recommendation 6.1

Conclusion The main objective of this thesis is to simulate the phenomenon of low-energy impact on laminated composite. An integrated modelling for carbon/epoxy composite plates has been developed to evaluate the impact damage initiation and propagation. Three failure modes: resin crazing, delamination and fibre failure have been investigated and implemented for the damage simulation. Resin crazing and fibre failure are studied based on Hashin criterion. And delamination is studied based on cohesive zone model. Threedimensional dynamic finite element method was used to simulate the real impact events. The simulation is based on Abaqus CAE package. The objectives obtained of this thesis can be summarised as follows:     

Understanding of impact behaviour of fibre-reinforced composites; Overview of contact mechanics; Modelling with Abaqus/CAE package; Introduction Hashin criteria for prediction matrix and fibre damage; Introduction cohesive zone model for simulation impact behaviour with interlaminar delamination;

The results obtained from simulation is not validated by experiment. However, it is compared with results of other similar researches. The comparison shows that the results of this thesis is satisfied about physical nature of the impact phenomenon. For modelling, works focus on simulation low-energy impact behaviour which includes failure modes. Using Hashin criteria for prediction of matrix and fibre failure of composite materials, we recommend that the mesh and element type need to be chosen appropriately. In addition, for cohesive zone model to simulate delamination, we note that the cohesive contribution to the global compliance before crack propagation should be small enough to avoid the introduction of a fictitious compliance to the model. We recognize that main problem of the simulation with cohesive zone model is to set the behaviour of cohesive zone.

61

6.2

Recommendation The present research highlights a number of issues which require further investigation. The following recommendations are therefore made: 





There is a need to establish drop-weight experiment for characterising laminated composite under low-energy impact loading. This would allow for better comparison between results obtained from simulation and experimental results. For simulation of impact problem, the modelling of low-energy impact in this research is completed. Its results can be used for prediction of damages and residual after impact. Therefore, next development of this research is to establish a model for high-energy impact. Since limit of time, the research do not study effect of mesh and increment of time of calculation. These may not be so much serious in the low-energy impact simulation since the gradient of the variable fields are not extremely. However, for high-energy impact simulation, these troubles are important. Therefore we need to care all these troubles in high-energy impact simulation. So it will be a general and complex problem which takes a lot of time of calculation and needs support of many processors for parallel calculation.

62

References

1.

S R Reid and G Zhou. Impact behaviour of fibre-reinforced composite materials and structure. 2000. p. 1-29.

2.

Nguyen Tuan Anh. Graduate Thesis Dap ung cua tam composite hang khong duoi tac dong va cham nang luong thap. 2007.

3.

Jean-Marie Berthelot. Mechanical Behaviour of Composite Materials and Structures. 1999. Chapter 12.

4.

J.A.Greenwood. Adhesion of Elastic Spheres. 1997. 1277-1297.

5.

Serge Abrate. Modeling of impacts on composite structures. Composite structures 51 (2001) 129-138.

6.

R.K.Lou, E.R.Green, C.J.Morrison. Impact damage analysis of composite plates. International journal of impact engineering 22 (1999) 435-447.

7.

C.Soutis & P.T.Curtis. Prediction of the post-impact compressive strength of CFRP laminated composites. Composites Science and Technology 56 (1996) 677684.

8.

Masaaki Nishikawa, Tomonaga Okabe, Nobuo Takeda. Numerical simulation of interlaminar damage propagation in CFRP cross-ply laminates under transverse loading. International Journal of Solids and structures 44 (2007) 3101-3113.

9.

M.F.S.F.de Moura, J.P.M. Goncalves. Modelling the interaction between matrix cracking and delamination in carbon-epoxy laminates under low velocity impact. Composites Science and Technology 64 (2004) 1021-1027.

10.

J.P.Hou, N.Petrinic, C.Ruiz. A delamination criterion for laminated composites under low-velocity impact. Composites Science and Technology 61 (2001) 20692074.

11.

Rikard Borg, Larsgunnar Nilsson, Kjell Simonsson. Modeling of delamination using a discretized cohesive zone and damage formulation. Composites Science and Technology 62 (2002) 1299-1314.

63

12.

A.Turon, C.G.Davila, P.P.Camanho, J.Costa. An engineering solution for mesh size effects in the simulation of delamination using cohesive zone models. Engineering Fracture Mechanics 74 (2007) 1665-1682.

13.

C.F.Li, N.Hu, Y.J.Yin, H.Sekine, H. Fukunage. Low-velocity impact-induced damage of continuous fiber-reinforced composite laminates. Part I. An FEM numerical model. Composite: Part A 33 (2002) 1055-1062.

14.

Abaqus documents: Abaqus/CAE User’s Manual; Abaqus Analysis User’s Manual; Abaqus Example Problems Manual.

15.

Eric Qiuli Sun .Shear Locking and Hourglassing in MSC Nastran, Abaqus, and Ansys.

64

Appendix A. Flowchart of FEM program [13]

65

B. Shear locking and Hourglassing [15] Abstract A solid beam and a composite beam were used to compare how MSC Nastran, Abaqus, and Ansys handled the numerical difficulties of shear locking and hourglassing. Their tip displacements and first modes were computed, normalized. It was found that fully integrated first order solid elements in these three finite element codes exhibited similar shear locking. It is thus recommended that one should avoid using this type of element in bending applications and model analysis. In FEA, shear locking and hourglassing are two major numerical problems because they may cause spurious solutions in certain situations. It is, therefore, necessary to have a solid understanding of how major commercial FEA codes handle these numerical difficulties before performing any serious engineering analysis. Shear locking Fully integrated linear brick elements of Nastran, Abaqus, and Ansys are overly stiff in bending applications and modal analysis. This numerical problem is called shear locking. Fully integrated first order elements, such as solid elements, Tomoshenko beam elements, Mindlin plate elements, may suffer from the locking (Prathap, 2005). FEA codes could therefore give false results when this type of element is used. In an ideal situation, a block of material under a pure bending moment experiences a curved shape change (Fig. 1). Suppose there are straight dotted lines on the surface of the block. Under the bending moment, horizontal dotted lines and edges bend to curves while vertical dotted lines and edges remain straight. The angle A remains at 90 degrees after bending, just as predicted by classical beam theory (Gere and Timoshenko, 1997).

Figure 1. Shape change of the Material Block under the Moment in the Ideal Situation To correctly model the ideal shape change, an element should have the ability to assume the curved shape. The edges of the fully integrated first order element are, however, not 66

able to bend to curves. The linear element will develop a shape show in Fig. 2 under a pure bending moment. The top surface experiences tensile stress, and the lower surface experiences compressive stress. All dotted lines remain straight. But the angle A can no longer stay at 90 degrees.

Figure 2. Shape change of the Fully Integrated First Order Element under the Moment

To cause the angle A to change under the pure moment, an incorrect artificial shear stress has been introduced. This also means that the strain energy of the element is generating shear deformation instead of bending deformation. The overall effect is that the linear fully integrated element becomes locked or overly stiff under the bending moment. Wrong displacements, false stresses, and spurious natural frequencies may be reported bacause of the locking. The fully integrated second order element behaves differently since its edges are able to bend to curves. Under a bending moment, the shape change of the element will correctly assume that of the material block (Fig. 3). The angle A continues to remain at 90 degrees after the bending. No artificial shear stress is introduced and the element can correctly simulate the behaviours of the material block. There is no shear locking associated with this type of element,

Figure 3. Shape change of the Fully Integrated Second Order Element under the Moment

67

Hourglassing To address the shear locking and to increase computational efficiency, a reduced integration scheme is proposed and widely implement in FEA codes. For example, for the reduced-integration first order 8-node brick element, a single integration point scheme is used while its fully integrated version employs eight integration points. For the reducedintegration, second-order 20-node brick element, an 8-integration-point scheme is used while its fully integrated version in Abaqus needs 27 integration points. In addition, the reduced integration element is tolerant of shape distortions, which is significantly beneficial in finite element modeling. However, nothing is perfect. The reduced integration first order element suffers from its own numerical difficulty called hourglassing since it tend to be excessively flexible. The hourglassing has to be properly controlled. If not, the results from this type of element are often not usable. Figure 4 demonstrates the deformation of such an element under a bending moment. To visualize the deformation, notice that vertical and horizontal dotted lines and the angle A remain unchanged. This means that normal stresses and shear stresses are zero at the integration point and that there is no strain anergy generated by the deformation. This zero-energy mode is a nonphysical response, which may propagate when a coarse mesh is used. The propagation of such a mode may therefore produce meaningless results. The results often indicate that the structure is excessively flexible. In order to make the reduced integration elements useful, the FEA codes provide default hourglassing control internally. The user may be able to adjust control parameters.

Figure 4. Shape change og the Reduced Integration Element under the Moment For the second order solid element with reduced integration, it may also suffer from hourglassing when only one layer of elements is used. But this rarely causes numerical problems because it virtually vanishes with two layers of elements. No special technique is needed to control it.

68