Failure Modeling of Curved Composite Beams

(1)

DEGREE PROJECT IN MECHANICAL ENGINEERING, SECOND CYCLE, 30 CREDITS

STOCKHOLM, SWEDEN 2018

Failure Modeling of

Curved Composite Beams

Numerical Modeling of Failure Onset and Propagation in L-Profile Beams

SUHAS GURURAJ SHETTY

KTH ROYAL INSTITUTE OF TECHNOLOGY SCHOOL OF ENGINEERING SCIENCES

(2)

i

This thesis work was carried out as a part of degree project of the Master’s programme in Aerospace Engineering – Lightweight Structures at KTH Royal Institute of Technology.

The author takes full responsibility for the Master Thesis work presented in this report.

GKN Supervisors - Olofsson Niklas and Tsampas Spyros KTH Supervisor/Examiner - Stefan Hallström

(3)

ii

Abstract

The high strength/stiffness-to-weight ratio that composite materials exhibit has led to the utilization of composites as alternative to traditional materials in weight-critical applications.

However, the highly anisotropic nature of composites renders the strength prediction under complex loading challenging. To efficiently predict the failure of composite structures especially in cases where out-of-plane stresses are dominant, the modeling of damage onset and propagation plays an essential role in accurate strength predictions.

Firstly, in this Thesis work the analysis of a composite L-profile, which is loaded such that significant out-of-plane stresses are generated in the curved region, is conducted. However, the inherent heterogeneity at the micro/meso scale is not modeled for the stress analysis.

Secondly, in this project the target was to accurately predict the initiation of failure at the ply level, modal based Puck’s matrix failure criteria have been implemented to the failure analysis. Maximum stress failure criteria were however retained to check the possible fiber- based failure which is not directly captured with in Puck’s failure criterion.

Thirdly, Cohesive Zone Material Model has also been employed to model the growth of interlaminar damage (delamination). The delamination study is based on the Inter Fibre Fracture crack initiation and doesn’t include other causes like edge effects, voids, manufacturing defects etc.

Finally, the attempt to validate the analysis results with the available test results was made.

Further development of the existing model and several tests are required to be carried out for material characterization and complete validation of the developed damage model for composite structure.

(4)

iii

Acknowledgments

This master Thesis work was written as a part of the master’s degree programme in Aerospace Engineering under Lightweight Structures track at Kungliga Tekniska Högskolan. The work has been carried out at GKN Aerospace Engine Systems in Trollhättan, Sweden during the spring term of 2018. Firstly, I would like to thank my supervisors Olofsson Niklas, Tsampas Spyros and Stefan Hallström for their time and continuous support during the project. I would substantially like to thank them for their valuable feedbacks which motivated me to understand further down the line in the world of composites. Special gratitude to all the employees at the R&T- department for their warm welcome and nice work environment who kept me propelled throughout the project especially due to their helpful nature and patience with my questions.

I would like to thank Erik Marklund, Thomas Bru and Fredrik Ahlqvist at Swerea SICOMP who guided me with the initial input for the project, despite their busy schedule. I would also like to thank ANSYS and Hyperworks technical support team for their valuable response during this Thesis work. This Thesis work would have been impossible without the inputs from the Journals, Text Books and other technical papers published in the related topics and hence I would also like to thank all the authors for publishing their articles which indirectly helped me to carry out my Master Thesis work. I would like to dedicate this master Thesis work to all teachers and professors who had upskilled me till date. Last but not the least I would like to thank my family, friends and relatives for continuously reinforcing me during the Thesis work.

- Suhas Gururaj Shetty Trollhättan, July 2018

(5)

iv

List of Symbols

€ : Euros (Currency)

ρ : Density

E : Elastic Modulus

𝜎̂ : Tensile Strength

b : Width of Specimen

t : Thickness of Specimen

δ : Applied Displacement

σ n(θ) : Normal Stress Component on Action Plane Ͳ nt(θ) : Component of shear stress on Action Plane Ͳ n1(θ) : Other component of shear stress on Action Plane

θ : Action plane search angle

θfp : Fracture angle

fE(θ) : Local Stress Exposure

σmax : Maximum Normal Traction

δnc : Normal displacement jump at the completion of debonding

Tmax : Maximum Tangential Traction

δtc : Tangential displacement jump at the completion of debonding

(6)

v

List of Abbreviations

UD : Unidirectional

NCF : Non-Crimp Fabrics

RTM : Resin Transfer Moulding

LR : Long Range

GE : General Electric Company

OMC : Organic-Matrix Composites

MMC : Metal-Matrix Composites

CMC : Ceramic-Matrix Composites

PMC : Polymer-Matrix Composites

WWFE : Worldwide Failure Exercise

FE : Finite Element

VCCT : Virtual Crack Closure Technology

CTE : Crack Tip Element Method

CZM : Cohesive Zone Modeling

M1 : Module 1

M2 : Module 2

M3 : Module 3

CTE : Coefficient of Thermal Expansion

T : Tension

C : Compression

LC : Loading Case

IFF : Inter Fiber Fracture

(7)

vi

List of Figures

Figure 1: Locations with complex geometry [6]. ... 2

Figure 2: Typical T-profile cross Section [7]. ... 3

Figure 3: GEnx-1B Engine [4]. ... 4

Figure 4 : Gantt Chart of Thesis work. ... 6

Figure 5: Generic workflow of the thesis work and modeling process (marked in the blue frame). ... 7

Figure 6: Example of a composite material in nature [9]. ... 8

Figure 7: Classification based on reinforcement [10]. ... 9

Figure 8: Chipped stone tool in paleolithic age [9]. ... 10

Figure 9: Three modes of Fracture [9]. ... 11

Figure 10: Illustration showing different damage in composites [13]. ... 12

Figure 11: Example of L-shaped laminate showing section loads and delamination [16] ... 13

Figure 12: Three major analysis modules. ... 16

Figure 13: GKN Workflow. ... 17

Figure 14: Relevant scale for composite analysis [25]... 18

Figure 15: L-Shaped specimen used in this thesis work [7]. ... 19

Figure 16: Layered Solid Model. ... 20

Figure 17: Mesh quality check. ... 21

Figure 18: General Hooke’s Law [25]. ... 21

Figure 19: Laminate and transformation of co-ordinate system [26]. ... 22

Figure 20: Material Properties at the ply level. ... 22

Figure 21: Global co-ordinate system and Local co-ordinate system [26]. ... 23

Figure 22: Tool setup used for testing the specimen in Tension and Compression [7]. ... 24

Figure 23: Out-of-plane stress components. ... 25

Figure 24: Out-of-plane stress components at δ = 0.6 mm. ... 26

Figure 25: Out-of-plane stress components at δ = 0.7 mm. ... 26

Figure 26: Stresses in UD-Laminae [27]. ... 28

Figure 27: Inter Fiber Fracture [28]. ... 29

Figure 28: Stresses on the action plane [28]. ... 29

Figure 29: Searching for fracture plane angle [28]. ... 30

Figure 30: Initiation of failure for LC Tension and LC Compression. ... 31

Figure 31: Schematic showing crack bridging tractions in cohesive zone [29]. ... 32

(8)

vii

Figure 32: Various typical cohesive traction-displacement curves: (a) Triangular, (b)

exponential, (c) trapezoidal, (d) perfectly plastic and (e) linear/Polynomial [29]. ... 33

Figure 33: Schematic diagrams of different formulations of interface elements [29]. ... 33

Figure 34: Mode 1 and mode 2 dominated Bilinear CZM law [24]. ... 34

Figure 35: Typical CZM mesh introduced into the existing model. ... 35

Figure 36: Example illustrating bilinear behavior check. ... 36

Figure 37: At δ = 0.9 mm. ... 37

Figure 38:At δ = 1.45 mm. ... 37

Figure 39: At δ = 2 mm. ... 37

Figure 40: At δ = 0.9 mm. ... 38

Figure 41: At δ = 1.45 mm. ... 38

Figure 42: At δ = 0.9 mm. ... 38

Figure 43: At δ = 1.45 mm. ... 39

Figure 44: Experimental test rig [7]. ... 40

Figure 45: Failure strains in Tensile tests. ... 41

Figure 46: Failure strains in Compressive tests. ... 41

Figure 47: Failure displacement in tensile tests ... 42

Figure 48: Failure displacement in compressive tests. ... 42

(9)

viii

List of Tables

Table 1: Lifetime data for LR A330 [5]. ... 3

Table 2: Estimation for composite fan case [5]. ... 4

Table 3: Mechanical properties of some UD fiber composites [8]. ... 5

Table 4: Pre-study on modeling approaches in composite structures. ... 15

Table 5: Dimensions of the four test specimens... 19

Table 6: Defined Boundary Conditions [7]. ... 24

Table 7: Results based on max stress failure criteria, using small deformation analysis... 27

Table 8: Results based on max stress failure criteria, using large deformation analysis. ... 27

Table 9: IFF failure criteria for LC Tension and LC Compression. ... 30

Table 10: Mixed mode bilinear material model input [24]. ... 35

Table 11: Max stress failure criteria for all load steps ... 51

(10)

ix

1 Introduction

Failure of fiber-reinforced composite laminates is a complicated process and because of that, failure prediction models become overly simplified since accounting for all the different physical phenomena occurring in the failure process would be very complex [1]. In this Thesis, composite failure modeling has been studied with the major focus in out-of-plane direction. It is vital to know how the structure responses to the applied load, before it is further developed and implemented in the design and manufacturing. Interpreting the failure onset and propagation in composite structures is very important. Hence, it is important to understand and appreciate the difficulties prior so that the limitations of these can be understood properly [1].

Non-planar complex composite structures are prone to out-of-plane loading which might lead to catastrophic failure due to for example delamination must be avoided when designing such lightweight structures. Accordingly, it is necessary to precisely predict the initiation by using suitable failure criteria. Finally, the model was developed to study how delamination grows by introducing cohesive elements. Hence, the major objective of Thesis work is to develop a suitable composite failure modeling approach.

1.1 Background and Origin of Thesis

GKN Aerospace Sweden has for several years been exploring the potential of using composite materials in aero-engines to reduce weight without compromising the overall performance.

For complex thermo-mechanical applications like aero engines, composites with high temperature capability can play an important role in enabling the use of composites in such applications. To optimize such composites for aero-engine applications and predict their structural performance, numerical simulations are performed. In this project, failure modeling has been carried out on curved L-profile geometries with the aim to validate the analysis results with the test results from a previous development program. However, to accurately predict with accounting for all the different phenomena further work is needed in all the major modules of the developed composite failure modeling.

1.2 GKN Aerospace Sweden Presentation

GKN Aerospace is the aerospace operation of GKN plc, serving a global customer base. With sales of £3.6 billion in 2017, the business is focused around 3 major product areas – aerostructures, engine products and transparencies and several specialist products like electro- thermal ice protection, fuel and floating systems and bullet resistant glasses. The business has

(14)

2

significant participation on most major civil and military programs. GKN Aerospace is a major supplier of integrated airframe and aero-engine composite structures. GKN also offers one of the most comprehensive capabilities in high performance metallics processing. GKN is also the world leading supplier of cockpit transparencies and passenger cabin windows [2].

1.3 Aims and Objectives of Thesis Work

In the General Electric GEnx engine, fan blades and the fan case are made of fibre-reinforced composites that enable a significant reduction in the overall weight [3] [4]. Such lightweight fan blades and fan case can effectively lead to reductions in the operational cost by a notable amount [5]. Similarly, the replacement of frame structures in the compressor module by composites as shown in the Figure 1 can further reduce the weight, and thus cost and CO2

emissions. The complex geometries of the modules result in new challenges where the effects due to out-of-plane loading cannot be neglected. L and T-joints are typical examples of such complexity (Figure 2). To predict the behavior of composites in such complex geometries due to out-of-plane loading there is a need for accurate composite failure modeling.

Figure 1: Locations with complex geometry [6].

The aims of this Thesis are briefly listed below:

• Numerical modeling of composite structure, capable of handling the stress analysis of the composite laminates at the ply level

• To predict onset of inter-laminar damage by using suitable failure criteria

(15)

3

• To study how delamination grows and develops by introducing cohesive elements

• To validate the analysis with the help of available test results

Figure 2: Typical T-profile cross Section [7].

This Thesis however is limited to prediction of the damage of simple L-profile composite geometries instead of T-profiles (to reduce the computational time due to geometrical complexity) and to develop the suitable modeling techniques. The proposed modeling approach is also capable of handling multidirectional laminates provided the inputs are carefully defined. However, in this report the major focus is on failure modeling for Unidirectional (UD) non-crimp fabrics (NCF), Resin Transfer Moulding (RTM6) composite material system [7]. The model is in-between the meso and macro scale, thus in this study the effect of glass yarns and the fiber waviness is not included.

1.4 Why do we need Composites?

Despite several complexity involved in handling composites over other conventional material system it is important to understand why there is a need for alternate material system. To understand this further let us consider lifetime data for LR A330 as given by SAS airlines, tabulated in Table 1 [5].

Table 1: Lifetime data for LR A330 [5].

Total fuel consumption One billion liters of jet fuel per engine Average fuel prize € 0,40/L

(16)

4

The presented data is statistical and thus have variation based on the flight cycle which is dependent on the airliner and the fuel prize which is also varying. However, from the above data, 1 kg of structural weight reduction should save around 4000 – 5000 L of fuel which could be further estimated to € 1,000 – 2,000 direct cost saving per kg of structural weight [5].

GEnx engines have replaced fan blade with composite which weighs 10-15% less than a hollow-core titanium blade [3]. A typical GEnx engine is presented in Figure 3 for reference.

Since, the blades are lighter, they also explored the possibility of building a lighter fan case.

The resulting fan case which is made of composites is in total 158 kg less weight per engine [3]. So, from the statistics mentioned above one could quickly estimate lifetime fuel and cost saving as tabulated in Table 2:

Table 2: Estimation for composite fan case [5].

Fuel saved 700,000 – 800,000 L per engine Cost saved € 200,000 – 320,000 per engine

Firstly, from the above study one could infer that the run time cost plays a vital role when the cost modeling is done for structures involving composites. Also, from the above statistics it is

Figure 3: GEnx-1B Engine [4].

also clear (and a well-known fact) that composites could be one of the possible solutions to reduce operating costs for an airliner.

Secondly, composites can be tailored in such a way that better mechanical properties or any other desired properties can be achieved, which is further discussed in Section 2.1. At the material selection stage, an engineer working to reduce weight must investigate specific

(17)

5

stiffness and strength properties. Typical mechanical properties of some UD fiber composites with other conventional materials is given in Table 3 [8]:

Table 3: Mechanical properties of some UD fiber composites [8].

Material ρ (kg/m³) E (GPa)¹ 𝝈̂ (MPa)²

Mild Steel 7800 206 250-500

Stainless Steel 7900 196 200

Aluminum alloy 2024 2700 73 300

Titanium alloy 4500 108 980

UD³ - Carbon/Epoxy 1600 180/10 1500/40

UD³ – Glass/Epoxy 1800 39/8 1060/30

UD³ – Kevlar/Epoxy 1300 76/6 1400/12

Note: 1 - Elastic modulus in/perpendicular to fiber direction, 2 – Tensile strength in/perpendicular to fiber direction and 3 – Prepregs with high fiber volume fraction

From the table one should note that it is possible to make a material system with higher specific mechanical properties. Thus, with the major goal of reducing weight in aero-engine one must overcome the barriers created due to non-planar complex composite structures as mentioned in Section 1.3. Hence there is a need for accurate numerical modeling of composite structures which will be discussed in the upcoming chapters.

1.5 Flow Chart for the Thesis Work

Composite Failure Modeling is a vast topic by itself. Hence, there was a need for the plan before going ahead with the Thesis work. A detailed work plan of the thesis work is illustrated in Figure 4 and overall workflow of the thesis work is also shown in Figure 5.

(18)

6

Figure 4 : Gantt Chart of Thesis work.

(19)

7

Figure 5: Generic workflow of the thesis work and modeling process (marked in the blue frame).

(20)

8

2 Literature Review

The purpose of this chapter is to summarize the literature relevant to this work investigated in a survey that was carried out in an early stage of this thesis work. One must note that the sections included in this chapter are brief and the main aim is to revisit some basic concepts.

2.1 Brief Introduction to Composite Materials

Composite material is defined as ‘‘A macroscopic combination of two or more distinct materials into one with the intent of suppressing undesirable properties of the constituent materials in favor of desirable properties’’. In a material science perspective, a composite material is thus composed of several different distinct materials [1].The application of composite materials is not a unique invention by mankind. Other species excel in making composites, and they have benefited from it for millions of years. For example, the nest of the Chinese bird shown in Figure 6 uses a similar concept as advanced carbon fiber reinforced composites from a mechanical point of view. Clay basically plays the role of a matrix that holds intact the reinforcements and protect them from being affected by the environment [9].

Figure 6: Example of a composite material in nature [9].

Composites are used not only for their structural properties, but also for their electrical, thermal, tribological and environmental properties. They are usually optimized to achieve a balance of properties based on applications and are commonly classified in two major levels.

The first level of classification is usually made with respect to the matrix constituent. This level of classification includes: organic-matrix composites (OMC), metal-matrix composites (MMC) and ceramic-matrix-composites (CMC). OMCs are further classified to two classes of composites: polymer-matrix composites (PMC) and carbon matrix composites. The second

(21)

9

level of classification refers to the reinforcement form which includes: particulate reinforcements, whisker reinforcements, continuous fiber composites and woven composites as depicted in Figure 7 [10].

Figure 7: Classification based on reinforcement [10].

The final category of fiber architecture is formed either by weaving, braiding, knitting the fiber bundles or also known as ‘tows’ to create interlocking fibers that often have orientations slightly or fully in an orientation orthogonal to the primary structural plane. This approach is taken for variety of reasons, including the ability to have structural, thermal, electrical properties etc. in the out-of-plane direction [10].

Thus, it becomes evident that there are infinite ways by which the composite part could be built. This is one of the major tasks of any composite engineer to optimize the design with the right selection of constituent materials.

2.2 Application of Composite Materials

The aim of lightweight construction is to preserve or even expand a product’s functionality while reducing the overall weight of the product. Some of the existing approaches for reducing mass include the use of less dense materials such as metal foams and honeycombs, composite materials etc. or decrease the material volume by reducing wall thickness in structural components. The main reasons for the application of lightweight composites are weight savings and possible cost savings. If there are significant weight reductions with

(22)

10

improved performance, it will also mean that there is less fuel consumption and CO2

emissions. In addition, there are several other advantages like noise and vibration reduction, impact resistance and energy absorption capability. Composites can also be tailored to meet specific design requirements in ways that are not possible for most conventional materials.

This could be done by correctly choosing the constituent materials and the orientation of the reinforcement fibers. This is of primary importance for performance optimization and hence lightweight construction could play a vital role in such applications [11].

2.3 Brief Introduction to Fracture Mechanics

Human ancestors made use of fracture phenomena more than 2 million years ago. Brittle solids, such as flint stones, usually crack in terms of cleavage when they are tapped, and sharp edges are formed on the stones as shown in Figure 8 which can be used as tools for cutting food or hunting [9].

Fracture Mechanics is the study of mechanical failures which can be of many different kinds.

A failure is a sudden loss of functionality of a mechanical component or structure by exhaustion of its load bearing capacity. The failure mechanism is thus the mode by which this occurs like for instance buckling, fracture etc. This is to be distinguished from the concept of damage mechanism. Damage is under most circumstances a non-favourable change of the material properties and the material behavior which in general develops over some time span or possibly with increasing loading. Failure is often preceded and promoted by damage formation and sometimes the difference between the two concepts is difficult to perceive [12].

Figure 8: Chipped stone tool in paleolithic age [9].

(23)

11

The goal of fracture mechanics is to enable predictions of initiation and propagation of growth of existing or postulated cracks of given configurations in structures of arbitrary shape. In general, for three-dimensional elastic crack problems, three stress intensity factors are enough to fully characterize the state at a point along the crack front. To obtain a visual impression of the three different modes idealized illustration are provided in Figure 9 [12].

Figure 9: Three modes of Fracture [9].

2.4 Failure Mechanisms in Composite Structures

Composite failure is the result of competition between different damage mechanisms. Failure mechanisms or modes of failure of a laminated composite can generally be divided into three types:

1. Translaminar: Through the thickness in which fibers have been broken (Fiber- dominated failure)

2. Intralaminar: Through the thickness in which only matrix, or fibre/matrix interface have failed (matrix and fiber/matrix interface-dominated failure)

3. Interlaminar: In the laminate plane, in which the layers (or plies) have been separated (matrix and fiber/matrix interface-dominated failure). This damage of laminate is also known as delamination [13]

(24)

12

Illustrations of these 3 failure modes are given in Figure 10.

Figure 10: Illustration showing different damage in composites [13].

2.5 Failure Criteria in Composite Structures

In 1991, an ‘Expert Meeting’ was held at St Albans (UK) on the subject of ‘Failure of Polymeric Composites and Structures’. Two key findings emerged during the meeting are as follows [14]:

1. There is lack of faith in the developed failure criteria

2. There is no Universal definition of what constitutes ‘failure’ of composite

This meeting later led to the Worldwide failure exercise also known as WWFE [14].

Composite plies and laminates have directionally-dependent strength and they exhibit several distinct failure modes. These anisotropic materials display more complex interaction of multiaxial stresses and strains, making the development of reliable failure theories much more difficult [15]. Failure criteria in a broader sense could be classified into two main categories as listed below:

1. Interactive Failure Criteria 2. Modal Failure Criteria

Depending on the load case and stress state at laminae, the failure predicted by these criteria can be fiber-dominated, matrix-dominated and fiber/matrix interface-dominated. In this Thesis work, Max stress and Pucks matrix failure criteria has been implemented to the failure analysis.

(25)

13

2.6 Brief introduction to Composite Curved Beams

As mentioned earlier in Section 1.3, Composite laminates in a wide variety of shapes start to replace metallic counter-parts and L-shaped geometry is frequently encountered composite curved beam. Interlaminar normal stresses are induced at the interfaces between the plies in addition to the well-known interlaminar shear stresses due to the geometry. Mixed-mode delamination failure occurs in the curved region of the L-shaped composite laminates under sectional forces as illustrated in Figure 11 [16].

Figure 11: Example of L-shaped laminate showing section loads and delamination [16]

The origins of such failures can often be associated simply with transverse strength limitations. When the loading or environmental condition is such that these interlaminar stresses are tensile, failure may occur at load levels much less than predictions based on the in-plane strength properties would indicate [17]. The analytical methods to determine load bearing capabilities of curved beams are under constant development. For example; an analytical technique suitable for calculating the stresses, strains and maximum load for a symmetric UD layup is developed by Fu-Kuo Chang and George S. Springer in 1985 [18].

Several Numerical methods are under constant development for the stress evaluation of curved structures. For example; stress variation along the curved beam’s width in bending using a 3D finite element analysis was investigated. It was observed that the assumption of plane strain for the analysis model resulted in a close solution to the 3D analysis in case of large specimen’s width, while significant errors were obtained with the assumption of plane stress [19].

(26)

14

Delamination research has mostly dealt with the initiation and growth of delamination.

Initiation can be predicted using stress-based criteria with some characteristic lengths.

Methods using fracture mechanics were developed for simulating delamination growth successfully. Another approach for the numerical simulation of delamination is the cohesive zone method, in which the framework of damage mechanics and softening is employed.

Delamination is interpreted as the creation of a cohesive damage zone in front of the delamination front, separating the adjacent plies. This method can handle both delamination onset and growth [20]. To predict accurately the damage of curved composite structures damage modeling is thus important.

2.7 Modeling Approaches in Composite Structures

The use of classical (continuum) methods of stress analysis has been developed over many decades to give techniques that can be applied satisfactorily to a vast range of situations.

Classical methods are however very limited to simple geometries and ‘real structural features’

for example the details of attachment of a stringer to a skin panel, cannot be analyzed. In such cases one must resort to Finite Element (FE) methods. FE analysis is merely an alternative approach to solving the governing equations of a structural problem [21]. In this Thesis, a pre- study was carried out to decide on a suitable numerical modeling for composite L-profile specimen and is briefly presented in Table 4.

(27)

15

Table 4: Pre-study on modeling approaches in composite structures.

Type – Stress Analysis

Application Advantages Disadvantages

Homogenized solid model

Useful when the geometry is complex and large

Easy to model and not heavy modeling files

Only smeared stresses on laminate level

Shell model Works for any layup orientation

Suitable for in-plane composite analysis

Extrapolated results for out of plane stresses

Layered solid model Works for any layup orientation

3D stress and strains are obtained at the laminae level

It might be heavier model when

implemented on the large geometry Type – Failure and

Growth Analysis

Application Advantages Disadvantages

Virtual crack closure technology (VCCT)

Works for any layup orientation

Based only on critical fracture toughness [9]

Crack tip needs to be defined

Crack tip element method (CTE)

Based only on critical fracture toughness [9]

Crack tip needs to be defined

Cohesive Zone Modeling (CZM)

Crack tip is not needed

Extensive input and mesh size influence

(28)

16

3 Methodology

The composite failure modeling carried out in this thesis can be divided into three major modules as shown in Figure 12 (for detailed workflow of the thesis work refer to Section 1.4).

Module 1 [M1] is the first part of composite numerical modeling where the stress analysis is carried out to obtain 3D stressing in the ply. The stresses obtained are later input to Module 2 [M2] to perform the required composite failure analysis. Once the initiation of damage is identified then progression of failure can be simplified in Module 3 [M3] to simulate the damage growth. It is important to note that all the three major modules are highly inter- dependent. Therefore, the accuracy of the analysis carried out in M2 and M3 is very much dependent on the precision achieved in M1 and M2 respectively.

Figure 12: Three major analysis modules.

3.1 FEM Software used

The software used during this thesis as per the GKN specific workflow is given in Figure 13 for the better understanding of subsequent work carried out in the thesis work.

3.1.1 Hyperworks

• HyperMesh: In this thesis HyperMesh was used as the pre-processor tool, due to its ability to quickly generate quality meshes. The advanced geometry and meshing capabilities provide an environment for rapid model generation [22].

Stress Analysis [M1]

Failure Analysis

[M2]

Progression Analysis

[M3]

(29)

17

Figure 13: GKN Workflow.

• HyperView: HyperView was used as the post processing tool to visualize results interactively [22].

3.1.2 MATLAB

MATLAB is a tool which combines a desktop environment tuned for iterative analysis and design processes with a programming language that expresses matrix and array mathematics directly [23]. Puck’s failure criteria are implemented in this tool for swift computation.

3.1.3 ANSYS

ANSYS structural analysis software is capable to solve complex structural engineering problems. In this thesis work, ANSYS was used as the solver in M1 and M3 modules. CZM pre and post processing were also performed using ANSYS [24].

3.2 Material and Component Testing at Swerea SICOMP

Test results used for comparison were obtained from the tests which were performed at Swerea SICOMP in Piteå using an Instron 8800 testing machine with a 100 kN load cell (calibration date 2012-11-26). The fiber preform was a Sigmatex UD weave consisting of alternating E-glass yarns (13.5 g/m²) and 12K Toho Tenax HTS40 F13 carbon fiber roving (242 g/m²) in the warp direction and a combi-yarn (areal weight 8.5 g/m²) in the weft direction. The resin material used is the monocomponent epoxy system RTM6 supplied by Hexcel corporation [7].

Hyperworks and ANSYS

[M1]

Matlab [M2]

ANSYS and Hyperworks

[M3]

(30)

18

4 Numerical Modeling of Composite Structures

As stated in Chapter 3, this Chapter briefly introduces to all the three modules developed in this thesis. Numerical Modeling of Composites can be broadly classified as follows;

• Structural analysis of the behavior of a fully consolidated composite structure

• Process modeling; the analysis of the manufacturing and forming of composite materials and parts [25]

The thesis focuses on the structural analysis of the composite structure and process modeling is not considered.

4.1 Module 1 – Stress Analysis

Stress Analysis of composites can be categorized as follows;

• Macromechanical approach: This approach involves constructing models strictly at the global scale. This approach is straightforward in the linear elastic regime

• Micromechanical approach: In the nonlinear regime and when trying to predict damage and failure, the macromechanical approach becomes problematic. This approach explicitly considers the constituent materials and how they are arranged [25]

Illustration of the relevant levels of scale for composite analysis is given in Figure 14. Note that the stress analysis carried out in this work is slightly inclined towards the mesoscale.

With the current model it is possible to obtain the 3D stress field at the ply level.

Figure 14: Relevant scale for composite analysis [25].

(31)

19

4.1.1 L-profile Geometry Modeling

As mentioned earlier in Section 1.3, the geometry under consideration is L- shaped. Also, one should note that L/T – profiles are the typical cross-section of the complex geometry as stated in Section 1.3. The 2D geometry was created in Hypermesh as per the cross-section dimensions given in Figure 15. Variation in the dimensions of the four test specimens manufactured as per this cross section is also listed in Table 5 for reference.

Figure 15: L-Shaped specimen used in this thesis work [7].

Table 5: Dimensions of the four test specimens.

Width of specimen (b in mm) Thickness of specimen (t in mm)

20.20 3.28

20.39 3.43

20.27 3.30

15.68 3.43

(32)

20

4.1.2 Pre-Processing

In this step, the composite model is setup for stress analysis. Accuracy of stress output is controlled at this stage for any typical composite structural analysis.

4.1.2.1 Meshing

Due to simple geometric shape, the results are not sensitive to meshing. Hence, a quadrilateral mapped plane strain element [22] was created initially by meshing the geometry mentioned in the Section 4.1.1. This mesh was later extruded in the z-direction as per the width of the geometry which is illustrated in Figure 16. There was no significant variation in stress results based on the element type chosen. Finally, the stress analysis was carried out using Solid 185 elements. The element selection should be made in such a way that there are no compatibility issues when the interface or contact elements is implemented to the existing FE model.

SOLSH190 was also tested for the model which could be alternative to the element type chosen [24].

Figure 16: Layered Solid Model.

As a general routine in the process of meshing, the quality of the mesh was also checked as shown in Figure 17 [22] [24]. It was observed that the minimum Jacobian was above 0.98.

Mesh was checked for default hypermesh parameters. The mesh was further refined based on the quality check. The quality of the mesh could affect the stress and strain outputs which are the major inputs for failure analysis. However, this check will play a vital role for complex shaped real structures.

(33)

21

Figure 17: Mesh quality check.

4.1.2.2 Material and Properties

Hooke’s Law for a fully anisotropic material (such as fiber-reinforced composites) is shown in Figure 18. In its most general form, it has 21 independent elastic constants [1] and 6 independent coefficients of coefficient of thermal expansion (CTE) [25].

Figure 18: General Hooke’s Law [25].

For computational purposes, general anisotropic laminae can be simplified using orthotropic laminae or a transversely isotropic material. As mentioned earlier in Section 3.2 the material system under consideration is a weave-based lamina. For computational simplicity UD laminae can be characterized as transversely isotropic materials. But weave-based laminae are better treated as orthotropic layered (materials with three orthogonal planes of symmetry for its material properties) [1] and the model is thus developed to run with orthotropic material input and is characterized by nine independent elastic engineering constants. As of now thermal aspects are not included in the analysis [25]. With the stable set of inputs (thermal coefficient) it is possible to include them to the existing model. For the typical laminate, all the plies may not have the same ply orientation as shown in Figure 19 [26]. In such cases, the user must be careful while defining the material properties by controlling the necessary transformation of the properties to the laminate axes or controlling the element orientation.

(34)

22

Figure 19: Laminate and transformation of co-ordinate system [26].

The model developed in this thesis is capable of handling multidirectional laminates by defining the transformed material properties for each ply. A typical illustration after defining material properties for composite numerical modeling is shown in Figure 20.

Figure 20: Material Properties at the ply level.

However, in this project the material system under consideration is UD as mentioned in Section 3.2. Thus, laminate and laminae axes are in the same direction and hence the laminae axes coincide with the laminate axes. But as mentioned earlier if the model has to be updated with additional layers, change of orientation of the laminae or in general change the stacking sequence, then the illustration in Figure 21 is helpful to define the orthotropic material properties. The tiny blue arrows shown in the image is the out-of-plane direction for the composite.

(35)

23

Figure 21: Global co-ordinate system and Local co-ordinate system [26].

4.1.2.3 Loads and Boundary Conditions

The interaction of tool on the specimen is defined with the help of boundary conditions. An attempt has been made to reach closer to reality with the help of suitable boundary conditions as listed in Table 6. Also, as shown in Figure 22, the tool setup is such that the specimen is ideally restricted to move in the z-direction. However, during testing it was observed that the displacement in the z-direction is miniscule for the 4 specimens. The load is controlled by prescribing the displacement in tension (T) or compression (C) respectively as illustrated in Figure 22 and table 6. It is important to note that tooling is not modeled in the current model.

However, for more accurate stress analysis modeling including the tool is needed.

(36)

24

Figure 22: Tool setup used for testing the specimen in Tension and Compression [7].

Table 6: Defined Boundary Conditions [7].

Global directions Horizontal (H) part of L-profile Vertical (V) part of L-profile

X free to roll gripped

Y gripped +/- prescribed δ

Z gripped gripped

4.1.3 Processing

Once the composite FE model is setup in the pre-processing stage, it needs to be processed with help of a suitable solver. In this project ANSYS was used as the solver as per the GKN workflow discussed in Chapter 3.

4.1.3.1 ANSYS Solver

ANSYS allows the user to control the analysis settings as per the structure under consideration. For simple linear static analysis, it is not needed to change any of the solver settings [24]. But if one needs to consider the incremental load step of 0.12 mm/min which was applied during the testing [7] then it might be needed to consider large deflections.

Activating large deflection will consider stiffness changes resulting from changes in element

T C

V

H

(37)

25

shape and orientation. While small deflection and small strain analysis assume that displacements are small and the resulting stiffness changes are insignificant. However, a rule of thumb is that one should use large deflection if the transverse displacements in a slender structure are more than 10% of the thickness. Due to this uncertainty, analysis was carried out for both and the significant variation in stress and strain outputs was observed. However, the computation time is increased when the large displacement setting is switched on [24].

4.1.4 Post-Processing

Post-processing of the stress analysis was done with the help of HyperView. It is obvious to obtain the stresses or strains to use the suitable failure criteria in Module 2. Post processing in composite analysis is rather substantial when one considers all the three modules. For better clarity and to explain the methodology mentioned in the previous Chapter, the results for only one case is considered (when the specimen is tested in tension). The results for the compression case are presented in Appendix 4.

4.1.4.1 3D Stress Field

The stress output from linear small displacement stress analysis is illustrated in Figure 23.

Stress (XX) represents the out-of-plane normal stresses and Stress (XY) represents out-of- plane shear stresses for the curved part of the L-profile when the prescribed displacement is1 mm. The other stress components are also given in Appendix 1 for reference.

Figure 23: Out-of-plane stress components.

(38)

26

The stress analysis was carried out with incremental load steps. The most critical sub-step based on Max stress failure criteria are illustrated in Figure 24 and Figure 25. Some additional stress output closer to the critical load for both the load cases is given in Appendix 3 and Appendix 4 for reference.

• Sub step 5: δ = 0.6 mm

Figure 24: Out-of-plane stress components at δ = 0.6 mm.

• Sub step 6: δ = 0.7 mm

Figure 25: Out-of-plane stress components at δ = 0.7 mm.

(39)

27

4.2 Module 2 – Failure Analysis

The application of realistic failure criteria for the UD composite layers in laminate design is a precondition for the successful use of laminated components in lightweight structures. Also, use of laminated composites as a primary/secondary structural component implies an accurate assessment of damage initiation, damage mechanisms, failure and post failure behavior. An attempt has been made to understand and implement suitable composite failure analysis to the model.

4.2.1 Max Stress Failure Criteria

Based on the stress output obtained in Module 1, Max stress criteria was initially used to predict the first ply failure due to out-of-plane stresses. Based on linear stress analysis the failure δ for different Loading Cases (LC) are listed in Table 7. All δ values are given in mm and all stresses and failure strengths are given in MPa

Table 7: Results based on max stress failure criteria, using small deformation analysis.

Loading case δ^c [mm]

δ^t [mm]

δ^s [mm]

Critical Ply number

LC Tension 3.69 0.36 0.35 8 or 9

LC Compression 2.97 0.44 0.35 9

Out of plane strengths 218 26.3 65

When the large displacement setting is switched on, for the incremental load step the failure δ comes out as in Table 8. Detailed failure analysis using Max Stress criteria for every sub-step is given in Appendix 2 for reference.

Table 8: Results based on max stress failure criteria, using large deformation analysis.

Loading case Failure δ [mm] Critical Ply number

LC Tension 0.6 – 0.7 8

LC Compression 0.3 – 0.4 9

As mentioned in the above table it is possible to predict the initiation of failure (the critical initiation points are marked with the red). However, Max Stress and Strain failure criteria

(40)

28

have several limitations and can be used only for the quick estimate of initiation of failure.

Thus, there is need for better failure criteria which account for the mutual interaction of stresses and for the damage mechanism involved.

4.2.2 Puck’s Failure Criteria

Puck’s failure criteria are modal based matrix failure criteria which have been implemented in the failure analysis. Puck’s failure criteria constitute one of the competing theories involved in WWFE-II. These criteria are based on 3D stress formulations and have already proven their capability under two-dimensional stresses in the WWFE-I [27].

4.2.2.1 Introduction to Puck’s action plane fracture criteria

As mentioned in Section 2.4, delamination is defined as the separation of layers from each other. This separation is caused by tensile stresses acting in the thickness direction and/or shear stresses acting in planes which are parallel to the layer interfaces. For better understanding, nine stresses which could possible lead to failure of the UD-laminae is also shown in Figure 26. Interlaminar stresses exist close to geometric discontinuities such as free edges and can be caused both by mechanical and hygrothermal loading. Even more important for the development of delamination zones are stress concentrations at inner defects such as tips of Inter Fiber Fracture (IFF) cracks. A figure illustrating a typical IFF crack tip is also shown in Figure 27. Higher local stresses occur and cause local delamination at each IFF- crack tip. Also, intensive testing has been carried out by Puck and his colleagues and their experimental and theoretical investigations even suggest that no delamination can occur in the absence of impact if no IFF-cracks have been developed in the laminate [28].

Figure 26: Stresses in UD-Laminae [27].

(41)

29

Figure 27: Inter Fiber Fracture [28].

For IFF analysis it is reasonable to use an adapted coordinate system. This is done because IFF can take place on an inclined fracture plane. It is clear from the UD laminae stress state that shear stresses τ12 andτ13 never lead to inter fiber failure. The interaction of all the stresses is accounted for in the Puck’s failure criteria. Puck simplifies the UD laminae stress state on the plane inclined by the angle θ as shown in the Figure 28. In the resolved plane only one normal stress (σ n(θ)) and two shear stresses (Ͳ nt(θ) and Ͳ n1(θ)) are acting. These three stresses potentially provoke IFF on their common action plane inclined by the angle θ.

Figure 28: Stresses on the action plane [28].

If fracture occurs on a plane inclined by a certain angle θ this plane is called fracture plane and the corresponding angle θ is called fracture angle (θfp). A search scheme for the fracture plane is illustrated in Figure 29 [27].

(42)

30

Figure 29: Searching for fracture plane angle [28].

4.2.2.2 Algorithm for Puck’s failure criteria

The failure analysis used in this Thesis is based on universal 3D-formulation of the action plane related to Puck’s IFF-criteria. To determine the stresses at fracture, it is necessary first to determine the fracture plane angle (θfp). This can be obtained by carrying out a numerical search of fE(θ). The fracture plane is characterized as the action plane with the maximum local stress exposure (fE(θ)) [27] [28].

4.2.2.3 Results from IFF failure analysis

The important results after the numerical search is listed in Table 9 and for better visualization the results are also illustrated in Figure 30. For further detailed outputs after failure analysis refer to Appendix 5.

Table 9: IFF failure criteria for LC Tension and LC Compression.

Loading case Critical Failure element ID θfp fE(θfp) Critical Ply number

LC Tension 103988 -83 3.9614 9

LC Compression 90362 85 2.6375 8

(43)

31

Figure 30: Initiation of failure for LC Tension and LC Compression.

From the IFF failure analysis the critical failure location could be obtained without any difficulties which are observed in max stress criteria. From the above results for LC Tension the initiation of failure happens at ply 9 interface while for the compression the initiation of failure should happen at ply 8 interface. Prediction of IFF crack initiation gives enough information needed to carry out the progression analysis in the laminated structure.

4.3 Module 3 – Progression Analysis

The prediction of progression of failure in composite materials is of great importance in the design of composite structures. However, the prediction of growth of failure in composite materials is not a trivial matter, because modes and mechanisms of failure are complex and varied, occurring at multiple-length scales and often interacting with each other to lead to global failure. As mentioned in the earlier chapters and sections, delamination is widely acknowledged as one of the most important failure modes and involves both opening and sliding modes. It can occur at relatively low load levels compared to ultimate failure by fiber fracture but still with significant consequences for the structural load-bearing capability. It is also caused by high interlaminar stress levels, which lead to through-thickness debonding of the individual plies. Within this area there are now several methods for predictions, as mentioned in Section 2.6 [29].

(44)

32

4.3.1 Introduction to delamination modeling with cohesive interface elements

The concept of CZM is based on a presumption of a zone of softening ahead of a sharp crack tip in the material. Within this zone, the opening is resisted by cohesive tractions as illustrated in Figure 31. The important assumption in the bilinear CZM formulations is that, the material remains linear-elastic until it reaches its tensile strength (σ max). After the maximum limit is reached it degrades linearly to zero at finite displacement. This is the simplest and numerically most convenient traction-displacement curve. Other shapes could also be considered for the study; however due to time limitation, the current study was carried out only on the Exponential CZM model and Bilinear CZM model. Some of the observations and results related to Bilinear CZM model will be discussed in this subchapter [29].

Figure 31: Schematic showing crack bridging tractions in cohesive zone [29].

Bilinear CZM model is the simplest and numerically most convenient traction-displacement curve to implement because it is monotonic with no discontinuities. The area under the typical traction displacement curves is the absorbed energy which is given by

𝐺_𝑐 = ∫ 𝜎. 𝑑𝑢

𝛿𝑓

0

were σ is the interfacial stress, u is the crack opening displacement and δf is the displacement at failure. For an assumed shape of curve (see Figure 32), the stress at initiation σmax and displacement at failure can be set such that the energy absorbed per unit cracked area is equal to the material’s critical fracture energy, Gc thus preserving Griffith’s energy balance. The crack thus initiates once the maximum stress criteria are exceeded and has fully propagated when the stress is returned to zero. This gives CZM an advantage over other fracture

(45)

33

mechanics-based methods as mentioned in Section 2.6 since it can predict both initiation and propagation of a crack [29].

Figure 32: Various typical cohesive traction-displacement curves: (a) Triangular, (b) exponential, (c) trapezoidal, (d) perfectly plastic and (e) linear/Polynomial [29].

Implementations of interface elements have largely taken the form of planar two-dimensional elements, which can be either zero thickness with overlapping nodes or have a very small finite thickness. This is said to represent the thin resin-rich layer between plies. Interface elements can also be implemented in a discrete form as non-linear springs connecting adjacent nodes. Figure 33 shows schematically the implementation of such interface elements into a finite element mesh [29].

Figure 33: Schematic diagrams of different formulations of interface elements [29].

(46)

34

4.3.1.1 Implementing Bilinear CZM Model to existing model

Initial first step is taken to implement CZM material model. However, further work is needed to fully understand the post processing of CZM modeling and improve the module 3 to be more robust than what could be obtained within this thesis work.

The mode 1 dominated bilinear CZM model assumes that the separation of the material interfaces is dominated by the displacement jump normal to the interface as shown in Figure 34, while mode 2 or mode 3 dominated bilinear CZM models assume that the separation of the material interfaces is dominated by the displacement jump that is tangent to the interface as shown in Figure 34.

Figure 34: Mode 1 and mode 2 dominated Bilinear CZM law [24].

The interface element type Inter 205 [24] is implemented into the existing model as illustrated in Figure 35. Without failure analysis one could directly implement these elements to all ply interfaces. But that will then increase the computation time unnecessarily. To minimize the computation time, it is enough to implement the CZM mesh only at those interfaces which are closer to the critical ply interface location where there is higher possibility of IFF crack as observed in Module 2.

Interlaminar fracture energy has been tested for delamination growth in a 0°/0° ply interface for the material system mentioned in 3.2. This can be used to determine the input parameters for the bilinear material behavior with tractions and separation distances [24].

For realistic structural applications and loading, it is likely that there will be a component of mixed-mode loading [13]. Thus, mixed mode debonding which involves both normal separation and tangential slip is activated by inputting data items as given in Table 10 below [24].

(47)

35

Figure 35: Typical CZM mesh introduced into the existing model.

Table 10: Mixed mode bilinear material model input [24].

Constant Symbol Property

C1 σmax Maximum Normal Traction

C2 δnc Normal displacement jump at the completion of debonding C3 Tmax Maximum tangential traction

C4 δtc Tangential displacement jump at the completion of debonding C5 α Ratio of δn* to δnc or ratio of δt* to δtc

C6 β Non-dimensional weighting parameter

(48)

36

4.3.1.2 Post-processing for Progression Analysis

Post-processing of Module 3 is vast and not trivial like other analysis. The results are very sensitive to mesh and defined properties. Also, the solution converges only after several iterations and consideration of certain regularization. The first part in the module is to check whether the response of the model is the same as expected as per the input defined. One such illustration to check is shown in Figure 36.

Figure 36: Example illustrating bilinear behavior check.

The study of failure progression can be done either with help of cohesive interface stress or interface separation distance. One of the observations made here is that, when there is full separation the cohesive interface stress must approach zero while interface separation distance should further increase. For simplicity purposes, cohesive interface stress is studied in this thesis since the post processing of progression analysis involves several iterations and is time consuming. The simplest LC Tension case with the minimum number of sub step is illustrated further (see Figure 37 to Figure 43).

(49)

37

1. X-Component of interface stress

Figure 37: At δ = 0.9 mm.

Figure 38:At δ = 1.45 mm.

Figure 39: At δ = 2 mm.

(50)

38

2. XY-Component of interface stress

3. XZ – Component of interface stress

(51)

39

More detailed results for the critical sub-steps for both LC Tension and LC Compression can be found in Appendix 5 and Appendix 6.

(52)

40

5 Comparison with Experimental Results and Validation

Two tests were carried out for LC Tension and LC Compression respectively. The experimental rig used for testing is as shown in Figure 44. The test was carried out both experimentally and numerically as part of the NFFP5 Refact project, which involved composite damage study and process modeling.

Figure 44: Experimental test rig [7].

Average failure strains (measured by Digital Image Correlation) were compared and there was some correlation between experiments as illustrated in Figure 45 and Figure 46.

However, more tests need to be conducted to have better representative failure strain which might result in further tuning of all the 3 modules mentioned in the previous chapter.

(53)

41

Figure 45: Failure strains in Tensile tests.

Figure 46: Failure strains in Compressive tests.

It was observed in the experiment that the initiation of failure occurred at the centerline which was also found in the failure analysis discussed in Module 2. The model also predicts progression of failure comparable with the experiments. Expected failure displacements predicted with the help of Module 3 as mentioned in Appendix-5 and Appendix-6 are also comparable with experimental results as shown in Figure 47 and Figure 48.

0.58 0.592

1.491

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6

FE Results Test 1 Test 2

Failure Strain

Test Number

Tensile Test 1 and 2

0.58

0.311

0.498

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7

Failure Strain

Test Number

Compressive Test 1 and 2

(54)

42

Figure 47: Failure displacement in tensile tests

Figure 48: Failure displacement in compressive tests.

As mentioned in Section 4.1.1 there is a slight variation in geometry between the four test specimens, which will also affects the FE results mentioned above in the comparison study.

The FE results used in this study are based on the failure δ obtained after progression analysis. Initially, it was iterative process and computationally heavy. The average of the geometrical dimensions was used to run the simulation since Module 3 was not very robust.

0.8 0.8

1.2

0 0.2 0.4 0.6 0.8 1 1.2 1.4

Failure δin (mm)

Test Number

Tensile Test 1 and 2

0.6

0.7 0.7

0.54 0.56 0.58 0.6 0.62 0.64 0.66 0.68 0.7 0.72

Failure δin (mm)

Test Number

Failure Modeling of Curved Composite Beams