CFRP fibre metal laminates using a combined elastoplastic damage model and including delamination effects

CFRP fibre metal laminates using a combined elastoplastic damage model and including delamination effects

Composite Structures 114 (2014) 64–79 Contents lists available at ScienceDirect Composite Structures journal homepage: www.elsevier.com/locate/comps...

4MB Sizes 43 Downloads 103 Views

Composite Structures 114 (2014) 64–79

Contents lists available at ScienceDirect

Composite Structures journal homepage: www.elsevier.com/locate/compstruct

Progressive failure analysis of perforated aluminium/CFRP fibre metal laminates using a combined elastoplastic damage model and including delamination effects Jing-Fen Chen, Evgeny V. Morozov ⇑, Krishnakumar Shankar School of Engineering and Information Technology, University of New South Wales, Canberra, Australia

a r t i c l e

i n f o

Article history: Available online 5 April 2014 Keywords: Fibre–metal laminates Progressive failure analysis Combined elastoplastic damage model Delamination Layup configurations

a b s t r a c t In this work, a finite element (FE) model which includes in-ply and delamination damage effects is developed for progressive failure analyses of fibre–metal laminates (FMLs). The combined elastoplastic damage model recently developed by the authors is applied to represent the in-plane mechanical response of the composite layers. Delamination is taken into account by a cohesive zone method implemented using the cohesive elements available in Abaqus. The metal layers are simulated using the Abaqus built-in isotropic elastic–plastic model. The model is first employed to investigate the blunt-notched behaviour of a GLARE specimen, and then, applied to study the effects of layup configurations (i.e., unidirectional, cross-ply and quasi-isotropic) and the numbers of holes (i.e., 0, 1, 3, 6 and 9) on the mechanical response of perforated aluminium/CFRP FMLs. Also, the stress distribution and damage progression patterns together with the failure sequences for the quasi-isotropic FMLs with 1 hole and 9 holes under tensile loading are discussed. Ó 2014 Elsevier Ltd. All rights reserved.

1. Introduction Fibre–metal laminates (FMLs) are hybrid laminated structures composed of thin layers of metal alloy and composite materials. Incorporating composite layers into a metal alloy system greatly enhances the overall performance of the laminated systems which combine the beneficial features of both composites and metals. Experimental studies have shown that FMLs provide better tensile properties [1,2], fatigue and fracture resistance [2–4], and superior damage threshold energy [5] than conventional metal alloys. Also, these light-weight hybrid structures are capable of absorbing a significant amount of energy through localised fibre fracture and shear failure in the metal plies [2,6]. These features make FMLs materials of choice for aircraft structures, as evidenced by their successful applications in Airbus A380 [7] and Boeing 777 [8] aircraft. Although FMLs possess many attractive properties, their failure mechanisms are rather complex due to their inhomogeneous structures which are composed of constituents with significantly different properties that remain distinct in the final composition.

⇑ Corresponding author. Tel.: +61 2 6268 9542; fax: +61 2 6268 8276. E-mail address: [email protected] (E.V. Morozov). http://dx.doi.org/10.1016/j.compstruct.2014.03.046 0263-8223/Ó 2014 Elsevier Ltd. All rights reserved.

Several failure mechanisms, including plastic deformations of the aluminium plies, fibre fracture, matrix cracking and delamination between adjacent load-carrying plies, may contribute to the failure of a FML [9]. The influences of many factors, such as the number of notches [10] and their sizes [11–13], the stacking sequence of the composite layers [14,15], and delamination [11,16,17], on the ultimate strengths and failure modes of FMLs have been studied in the literature. Chang et al. [10] found that the lead fatigue crack growth rates in cross-ply GLARE 3-3/2 laminates containing multiple open holes increased because of the presence of multi-site damage. Wu et al. [11] found that the presence of a central hole in cross-ply GLARE 4-3/2 FMLs incurred a strength reduction of approximately 40%, even in a FML with the smallest hole diameter-to-width ratio of 0.1198. Yeh et al. [12] experimentally studied the effects of the notch sizes and constituents on the mechanical behaviour and ultimate notched strengths of hybrid boron/glass/aluminium FMLs. They also found that there was a decrease in the ultimate notched tensile strengths of specimens with larger notch sizes and the delamination area is larger for specimens with smaller hole diameters. Kawai and Arai [13] concluded from their off-axis experiments on cross-ply GLARE 3 that the notched strength of this material decreased with an increasing hole diameter, regardless of its loading direction with respect to its fibre orientation.

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

Also, the stacking sequences and layup configurations significantly affect the mechanical behaviour of FMLs as reported in publications [14,15]. Aktas and Dirikolu [14] found that pin-loaded carbon epoxy composite laminates with stacking sequences of     ½90 =45 =  45 =0 s had higher bearing strengths than those with     ½0 =45 =  45 =90 s . Seyed Yaghoubi et al. [15] studied the effects of layup configurations on the response of the GLARE 5-3/2 FMLs under low-velocity impact. They found that the one with the quasi-isotropic layup had the best impact resistance, followed by those with the angle-ply and cross-ply layups, while those with the unidirectional layup had the worst. In addition, investigations performed on notched and unnotched FMLs under tensile loading revealed that delamination is one of the dominant failure mechanisms of FMLs [11,16–18] and is likely to occur in the vicinity of the holes in notched FMLs. The presence of delamination causes stress redistribution and postpones fibre failure in the pre-preg layers of a FML and increases its notched strength [11,17,18]. From observations of failed notched GLARE specimens, the blunt-notched or residual strengths were lowest for those with only delamination, higher for those with only fibre failure and highest for those with both fibre failure and delamination. However, it should also be noted that as FMLs with built-in delamination due to poor impregnated fibres gave lower strengths, it was expected that a larger delamination area may cause a lower strength [18]. As delamination has a significant effect on the failure strengths of FMLs, modelling methods that do not take this into account may be not sufficiently accurate. Without considering delamination, the use of a three-dimensional (3D) fracture surface based on the ultimate strains of fibres and the aluminium alloy yields conservative predictions [17]. By disregarding delamination between the composite and metal layers, simulations performed by Van Rooijen et al. [19] resulted in overestimated bearing yield of the unidirectional GLARE 2-2/1 and cross-ply GLARE 3-3/2 with high edge distance-to-diameter ratios in which the influence of delamination was significant. Therefore, both in-plane failure mechanisms and out-of-plane delamination failure should be incorporated in any numerical model to realistically represent the behaviour of FMLs. As the mechanical response of FMLs is rather complex and affected by many factors, simulating their failure behaviour with greater accuracy remains a primary challenge. Currently, a common way of investigating the mechanical response and damage progression patterns of various failure modes, and estimating the notched strengths of FMLs, is to employ FE methods. The FE approach is one of the most flexible methods for studying the mechanical response of FMLs, irrespective of their geometries or layup configurations, under general loads and boundary conditions where analytical solutions are hardly applicable. The complete details of the stress distribution and damage progression patterns in the whole FE model can be obtained. It is necessary to incorporate progressive damage models in FE analyses to take into account the influence of various damage mechanisms on the residual strengths of notched FMLs. A few continuum damage mechanics (CDM)-based composite material constitutive models have been applied to progressive damage simulations of FMLs [20–22]. In these models, both damaging and linear elastic material behaviour have been assumed for modelling the composite layers of FMLs. Several investigations have adopted the Abaqus built-in fibre-reinforced material damage model to study the progressive failures of FMLs [12,20,23–25]. Although this model was able to provide reasonably accurate results in these studies, the values of fracture energies adopted for the damage evolution laws in this model were either those given by Lapczyk and Hurtado [20] which were not the actual values of the materials under consideration [20,23,24] or not reported at all [12]. Frizzell et al. [21,22] investigated the bearing

65

stress–strain response of pin-loaded glass FMLs and the effects of geometry on the bearing response of FML joints using a 3D continuum damage model. The above-mentioned two CDM-based models [20–22] are elastic-damage models that do not take into account the plasticity effects in composite materials. However, as mentioned in [26], some composite materials, such as AS4/PEEK, T300/1034-C and AS4/3501-6 composites, may exhibit nonlinearity or plasticity behaviour. The application of elasticdamage models (e.g. the Abaqus built-in fibre-reinforced material damage model) for predicting the mechanical response of such composite laminates may be insufficient when plastic deformations are present under loading [26]. In this work, a comprehensive FE model, which includes in-plane and delamination damage effects is developed for the progressive failure analyses of FMLs in Abaqus. The combined elastoplastic damage model [26], which is suitable for the progressive failure analyses of composite materials that exhibit plasticity behaviour, is applied to the numerical simulations through the user-defined material subroutine UMAT. Delamination occurring between different adjacent in-plane load-resisting layers is also taken into account by employing a cohesive zone model based on the COH3D8 cohesive elements available in Abaqus 6.10-3 [27]. The mechanical response of the metal layers is modelled using the isotropic elastic and isotropic hardening plasticity model available in Abaqus. The present FE modelling method is first employed to investigate the blunt-notched behaviour of a GLARE specimen, and then, applied to the progressive failure analysis of a 3/2 FML system which is made of three aluminium layers and two AS4/PEEK composite laminates which exhibit pronounced plasticity effects. To the best of the authors’ knowledge, little work has been carried out on the tensile failure behaviour and strengths of perforated FMLs containing multiple holes. Also, the effects of the layup configurations on the tensile failure behaviour of FMLs with different numbers of holes have seldom been reported in the literature. Thus, the objective of this study is to numerically investigate the tensile mechanical behaviour and damage growths of aluminium/CFRP FMLs with various layup configurations and different numbers of holes. The damage development patterns in the composite plies and the delamination behaviour in the FMLs are examined. This modelling approach is described in the following sections.

2. The finite element model The FE model developed in this work consists of a composite ply model, interface layer model and metal layer model as described in the following subsections. 2.1. Composite ply model – combined elastoplastic damage model The combined elastoplastic damage model for composite materials which exhibit plasticity and damage effects in their mechanical behaviour has been proposed, full details of which are given in [26]. However, for the integrity of this work, brief descriptions and relevant equations are provided. The combined elastoplastic damage model is presented for an elementary composite ply and consists of a plastic part which describes the plastic behaviour of composite materials under transverse and/or shear loading, the failure criteria that are used to predict the thresholds for damage initiation and growth, and damage evolution laws that account for the development of the failure process. 2.1.1. Stress–strain relationship For composite materials exhibiting plasticity response, the total strain tensor e is presented as a sum of the elastic and plastic strain parts ee and ep

66

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

e ¼ ee þ ep

ð1Þ

where the plastic strain ep represents all the irreversible deformations including those caused by microcracks. Damage accumulation and growth in the composite materials lead to the degradation of material properties which is taken into account by introducing damage variables in the stiffness matrix. The stress–strain relationships for damaged and undamaged composite materials are presented in the following forms:

r ¼ SðdÞ : ee ; r ¼ S0 : ee

ð2Þ

where bold-face symbols are used for variables of tensorial characteristics and symbol (:) denotes the inner product of two tensors with double contraction, e.g., ðSðdÞ : ee Þij ¼ SðdÞijkl eekl , where  are the summation convention is applied to the subscripts; r; r the Cauchy nominal stress tensor and the effective stress tensor (both are of the second order); S0 is the fourth-order constitutive tensor for linear-elastic undamaged unidirectional laminated composite; SðdÞ is the one for the associated damaged material. The explicit form of S0 is determined by elasticity theory for orthotropic materials; d represents the damage variable. In FE implementations, it is more convenient to use Voigt notation [28] to convert higher order symmetric tensor to lower order matrices. The Voigt form of the fourth order tensor SðdÞ adopted in this model is given as [26]

2 3 ð1  d1 ÞE01 ð1  d1 Þð1  d2 Þm021 E01 0 16 7 SðdÞ ¼ 4 ð1  d1 Þð1  d2 Þm012 E02 5 ð1  d2 ÞE02 0 D 0 0 0 Dð1  d3 ÞG12 ð3Þ where SðdÞ (printed in italic) is used to identify the Voigt form of SðdÞ; D ¼ 1  ð1  d1 Þð1  d2 Þm012 m021 ; parameters d1 ; d2 and d3 denote damage developed in the fibre and transverse directions, and under shear stress condition (they are scalar damage variables that remain constant throughout the ply thickness); E01 ; E02 ; G012 and m012 ; m021 are the elastic moduli and Poisson’s ratios of the undamaged unidirectional composite laminate. The damage variables are introduced as follows:

d1 ¼



d1t d1c

1 P 0 if r 1 < 0 if r

d2 ¼



d2t d2c

ð4Þ

where d1t ; d1c characterise the damage development caused by tension and compression in the fibre direction, and d2t ; d2c reflect the damage development caused by tension and compression in the transverse direction, respectively; d3 represents the damage effects on shear stiffness, it results from fibre and matrix cracking; d6 represents the damage effects on shear stiffness caused by matrix  1; r  2 are the effective stresses in the fibre and cracking only; r transverse directions. 2.1.2. Plastic model The irreversible deformations are allowed for by the plastic model which includes the yield criterion, plastic flow rule, hardening curve, and the hardening law. An equivalent form of the oneparameter plastic yield function for plane stress condition proposed in [29] is adopted to describe the irreversible strains exhibited by composites under transverse and/or shear loading:

 ; ~ep Þ ¼ Fðr

rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 2  þ 2ar  23 Þ  r ~ ð~ep Þ ¼ 0 ðr 2 2

r~ ð~ep Þ ¼ bð~ep Þn

ð5Þ

where a is a material parameter which describes the level of plastic deformation developed under shear loading compared with that  3 is the effective in-plane shear stress; under transverse loading; r

ð6Þ

where b and n are coefficients that fit the experimental hardening curve. These parameters together with the material parameter a are determined using an approach based on linear regression analyses of off-axis tensile tests performed on unidirectional composite specimens [29,30] The associated plastic flow rule and associated hardening law are assumed for composite materials. The corresponding evolution equations for the plastic strain ep and the equivalent plastic strain ~ep are expressed as:

e_ p ¼ k_ p @ r F; ~e_ p ¼ k_ p @ r~ F ¼ k_ p

ð7Þ

where k_ p P 0 is a nonnegative plastic consistency parameter (hereafter the following notations are used for the derivatives: @ x y ¼ @y=@x; y_ ¼ dy=dt); @ r F and @ r~ F represent the evolution directions of plastic strain ep and equivalent plastic strain ~ep , respectively. 2.1.3. Damage model In order to predict the damage initiation and propagation of each intralaminar failure of the material and evaluate the effective stress state, the damage initiation and propagation criteria fI are presented in the following form:

fI ð/I ; r I Þ ¼ /I  rI 6 0 I ¼ f1t; 1c; 2t; 2c; 6g

ð8Þ

where /I are the loading functions for different failure mechanisms adopted in the form of Hashin’s failure criteria [26,31]; r I is the damage threshold parameter which controls the size of the expanding damage surface and depends on the loading history. The development of damage in the material initiates when the value of /I exceeds the initial damage threshold r I;0 ¼ 1. Further damage growth occurs when the value of /I in the current stress state exceeds the value of r I in the previous loading history. The damage thresholds r I can be determined as:

rI ¼ maxf1; maxf/sI gg I ¼ f1t; 1c; 2t; 2c; 6g

2 P 0 if r 2 < 0 if r

d3 ¼ 1  ð1  d6 Þð1  d1t Þ

~ ð~ep Þ represents the isotropic hardening curve. An isotropic and r hardening curve proposed by Sun and Chen [29] is employed to represent the equivalent stress versus equivalent plastic strain hardening curve:

s 2 ½0; t

ð9Þ

Since damage is irreversible, the damage evolution rate should satisfy the following condition: d_ I P 0. The exponential damage evolution law is adopted for each damage variable and expressed in the following form [32]

dI ¼ 1 

1 expðAI ð1  r I ÞÞ I ¼ f1t; 1c; 2t; 2c; 6g rI

ð10Þ

where AI is the parameter that defines the exponential softening law. This parameter is determined by a numerical procedure presented in [26] aiming at alleviating the mesh dependency of the FE results which is a consequence of the use of any strain-softening material models. 2.1.4. Numerical implementation A strain-driven implicit integration procedure for the combined elastoplastic damage model has been developed using equations of continuum damage mechanics, plasticity theory and applying closest point projection return mapping algorithm. The simulation of the plasticity response is an incremental process that must be characterised by rate constitutive equations (e.g. Eq. (7)). To obtain the overall stress strain relationship, the numerical integration of the rate constitutive equations over a discrete sequence of time steps is required. In this way, the continuum problem is converted into a series of incremental problems. The solution to the nonlinear inelastic incremental constitutive problem under consideration is

67

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

achieved using the Newton–Raphson iterative method through which the nonlinear incremental problem is consistently linearised into a sequence of linear problems and subsequently solved [33]. The present incremental nonlinear inelastic problem is regarded as strain driven. The loading history is discretised into a sequence of time steps ½t n ; t nþ1 ; n 2 f0; 1; 2; 3; . . .g where each step is referred to as the ðn þ 1Þth increment. Driven by the strain increment Denþ1 , the discrete problem in the context of the backward Euler scheme for the combined elastoplastic damage model can be stated as ~ n ; rI;n ; dI;n g at  n ; rn ; r follows: for a given variable set fenþ1 ; epn ; ~epn ; r the beginning of the ðn þ 1Þth increment, find the updated ~ nþ1 ; r I;nþ1 ; dI;nþ1 g that fulfils all  nþ1 ; rnþ1 ; r variable set fepnþ1 ; ~epnþ1 ; r the governing constitutive equations that constitute the combined elastoplastic damage model at the end of the ðn þ 1Þth increment. To ensure the algorithmic efficiency of the Newton–Raphson method in the finite element analysis, a tangent operator consistent with the developed integration algorithm has been derived and implemented. The combined elastoplastic damage model is implemented in Abaqus using user-defined subroutine UMAT in which codes for the integration algorithm and tangent operator are written. Details concerning the numerical implementation and validation of the model can be found in [26]. Elements with a plane stress formulation, such as two-dimensional CPS4 plane stress elements, SC8R continuum shell elements, and S4, S4R conventional shell elements in which in-plane stresses are defined and out-of-plane transverse normal and shear stresses are equal to zero can be used with this model. In this study, the mechanical response of the composite layers is simulated with this model and SC8R continuum shell elements are used to discretise each composite layer. As previously discussed, the effect of delamination should be included in simulation approaches to obtain accurate predictions. A cohesive zone model is adopted for the modelling of delamination in the interface layers. 2.2. Interface layer model – cohesive zone model Many methods, such as the virtual crack closure technique (VCCT) [34,35] and cohesive zone approach [36], have been developed for the prediction of delamination. Compared with the VCCT, a pre-defined crack path and complex moving mesh technique are not required for a cohesive zone approach [37]. Thus, the latter has frequently been applied by researchers to study the onset and growth of delamination between composite layers and metal alloy layers [20–25,38]. In this work, a cohesive zone model implemented in COH3D8 cohesive elements available in Abaqus 6.10 is adopted for the simulation of delamination occurring in the adhesive layers. The use of cohesive elements involves the selections of a traction-relative displacement law, single- or mixed-mode delamination initiation criterion, single- or mixed-mode delamination propagation criterion, linear or exponential softening behaviour and element size. Details concerning the application of cohesive elements are discussed in this section. 2.2.1. Traction-relative displacement law Before the onset of delamination, an uncoupled linear-elastic traction-relative displacement law is specified to hold together the two surfaces of the composite layers above and below the adhesive interface. The uncoupled elastic traction-relative displacement law can be written as [39]

1

0

2

rz Kz C 6 B B szx C ¼ 6 0 A 4 @ szy

0

0 K zx 0

0

30

1

uz 7B C 7 B 0 5@ uzx C A K zy

uzy

ð11Þ

where x; y and z denote the three orthogonal directions of the local rectangular coordinate system assigned for the cohesive element’s mid-plane; z denotes the through-thickness direction and corresponds to the Mode I failure; the x and y directions span the midplane of the interface layer and correspond to the Mode II and III failures (the shear failures parallel and transverse to the fibre direction, respectively); rz is the out-of-plane normal stress; szx and szy are the out-of-plane transverse shear stresses; uz ; uzx ; uzy and K z ; K zx ; K zy are the corresponding relative displacements and penalty stiffnesses, respectively. The penalty stiffnesses for the traction-relative displacement law are defined as

8 0 > > < Ki K i ¼ ð1  di ÞK 0i > > : 0

ui 6 u0i u0i 6 ui 6 ufi ;

i ¼ z; zx; zy

ð12Þ

ui P ufi

where K 0i is the initial penalty stiffness; u0i and ufi correspond to the relative displacement of each single mode loading at delamination initiation and the point where delamination is completely formed; di is the damage variable that governs the delamination evolution law. Several methods for defining the initial penalty stiffnesses have been reported. A few constant values of interface stiffness have been applied to simulations of composite materials [40,41], e.g., 106 ; 107 and 108 N/mm3). Zou et al. [42] suggested that a reasonable value for the penalty stiffness could be 104 —107 mm1 times the inter-laminar shear and tensile strengths while other researchers assumed that it was a function of the interfacial material properties [43] or material properties at the laminate level [44]. In this study, the value of 106 N/mm3 is chosen for the three initial penalty stiffnesses required as it was determined by Davila et al. [40] based on a substantial number of simulations and considered to be suitable for modelling all the delamination modes. After the delamination initiation, the behaviour of the interface layer is controlled by a softening law for the delamination propagation. Traction-relative displacement laws with linear and exponential softening laws are available in Abaqus and, from the computational point of view, the one with the exponential softening law is deemed to be better because of its smoothness as it has only one vertex point whereas the bilinear law has two. In Abaqus, softening laws for delamination propagation can be defined based on effective displacement or fracture energy at complete delamination failure [39]. Ridha et al. [45] found that analysis applying the energy-based exponential softening law provided more accurate prediction of the failure strength of a repaired composite panel under tensile loading than that using the linear softening law. Therefore, the one with the energy-based exponential softening law is adopted in this work. For the energy-based softening laws (either linear or exponential softening laws), it is ensured in the main Abaqus computational procedure that the area under the traction-relative displacement curve is equal to the fracture energy defined. The evolution law of the damage parameter di for energy-based exponential softening law is expressed as [39]

dz ¼

Z

dzy ¼

uz

GIc 

u0z

Z

rz duz

uzy

u0zy

G0I

;

szy duzy GIIIc  G0III

dzx ¼

Z

uzx

u0zx

szx duzx GIIc  G0II

; ð13Þ

where GIc ; GIIc and GIIIc are the interfacial Mode I, II and III critical fracture energies, respectively; G0I ; G0II and G0III are the elastic energies at delamination initiations under single Mode I, II, and III loadings, respectively. It should be noted that in this energy-based

68

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

exponential softening law, the traction might not reduce immediately after delamination initiation as opposed to the displacement-based exponential softening law. The linear energy-based exponential traction-relative displacement law is illustrated in Fig. 1 for Mode I loading (the graphs reflecting Mode II and III loadings are similar). 2.2.2. Delamination initiation and propagation criteria In structural applications of composite materials, delamination growth often occurs under mixed-mode loading and delamination onset may happen before any of the traction components reach their corresponding maximum stress allowables. In this work, a quadratic nominal stress criterion is used to predict the onset of mixed-mode delamination. It allows for the effect of the interactions of the traction components on the delamination initiation. Considering that compressive normal traction does not normally contribute to the delamination onset, the quadratic nominal stress criterion is written in the following form as [36]

 2 hrz i

r0z



szx þ 0 szx

2

szy þ s0zy

!2 ¼1

ð14Þ

where hi ¼ ðjj þ Þ=2 is the Macaulay bracket; r0z ; s0zx , and s0zy denote the interface Mode I, II and III failure strengths, respectively. Once the mixed-mode delamination initiation criterion is met, delamination begins to propagate. The material stiffness of the cohesive layers starts to degrade for each mode in such a way that the total fracture energy absorbed per unit area under mixed-mode loading is governed by a delamination propagation criterion [46]. The mixed-mode propagation criterion defines the state of complete delamination for different mixed-mode ratios. A power law criterion and Benzeggagh–Kenane criterion are available in Abaqus. In this study, the power law criterion is adopted:



GI GIc

c

 þ

GII GIIc

c

 þ

GIII GIIIc

c

¼1

ð15Þ

where GIc ; GIIc and GIIIc are the interfacial Mode I, II and III critical fracture energies. For this criterion, the critical total fracture energy Gt;c at the completion of a mixed-mode delamination is the sum of the Mode I, II and III fracture energies, i.e. Gt;c ¼ GI þ GII þ GIII , where GI ; GII ; GIII reflect the works done by the normal and shear stresses in the cohesive layer and their corresponding relative displacements, respectively; The power law mixed-mode propagation criterion with c ¼ 1 and c ¼ 2 are frequently chosen in the literature. In this work, the one with c ¼ 1 is adopted as such selection has been found to be suited to predict the complete delamination failure and accurately capture the dependence of the mixed-mode fracture energy on the mode ratio of AS4/PEEK thermoplastic composite materials [41]. 2.2.3. Cohesive zone length Besides, it has been pointed out that the number of cohesive elements in the cohesive zone should be sufficient to accurately capture the stress distribution in the cohesive zone ahead of the

delamination tip and obtain accurate overall structural response [44,47]. For delamination simulations of thin orthotropic composite materials, Yang and Cox [48] suggested the effect of the specimen thickness should be included in the formulations for estimating the cohesive zone lengths of composite materials under different mode loadings. For cohesive zone models that specify non-zero tractions when the relative displacements are zero (i.e. rz ð0Þ > 0; szx ð0Þ > 0 and szy ð0Þ > 0), the formulations for estimating the cohesive zone lengths of a material under Mode I and II loadings are given as [48]

lcz;z ¼

E0I

!14

GIc 2

ðr0z Þ

3 4

h;

lcz;zx

vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi !ffi u u 0 GIIc t EII ¼ h 2 ðs0zx Þ

ð16Þ

where h is the half-thickness of the laminate; E0I and E0II are the equivalent elastic moduli under Mode I and II loadings, respectively. The procedure for determining the values of E0I and E0II can be found in [47]. According to the findings reported in [47,48], the cohesive zone length for Mode II loading is significantly longer than that for Mode I which indicates that the minimum number of cohesive elements in the cohesive zone is controlled by the Mode I cohesive zone length. The number of cohesive elements to be placed in the cohesive zone can be calculated as N e ¼ lcz =le . To accurately capture the stress distribution within the cohesive zone at the point of initial crack propagation, Harper and Hallett [47] suggested that at least two to three elements needed to be placed in this zone while Falk et al. [49] used two to five elements in their simulations and Davila et al. [40] used three elements. In this work, at least 9 cohesive elements are placed in the cohesive zones for the FMLs under consideration. For cohesive zone models which prescribe zero tractions at zero relative displacements (i.e. rz ð0Þ ¼ 0; szx ð0Þ ¼ 0 and szy ð0Þ ¼ 0), the cohesive zone lengths predicted using Eq. (16) will be smaller than the actual cohesive zone lengths [48]. Under such conditions, Yang and Cox [48] suggested that this equation can be used to determine a suitable element size, as it is only required that the element size be less than the calculated cohesive zone length. 2.3. Metal layer – isotropic elasto-plastic material model The metal layers in the FMLs are made of aluminium 2024-T3, the material response of which is assumed to be isotropic and can be modelled using the von Mises plastic model with isotropic hardening in which the isotropic elasto-plastic material model (also called the J 2 -material deformation theory) is available in Abaqus 6.10-3. The isotropic hardening behaviour is defined as the relationship between the yield stress and equivalent plastic strain in a tabular form. Considering the interactions of the aluminium layers with the AS4/PEEK carbon fibre-reinforced PEEK layers, the out-of-plane stresses may contribute to the plasticity response of the aluminium layers. To approximate this 3D stress state with greater accuracy, the aluminium layers are simulated using C3D8R solid elements.

Fig. 1. Linear energy-based exponential traction-relative displacement law for Mode I loading.

69

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

3. Numerical examples In order to validate the capability of the present modelling approach presented in the previous section, it is applied to investigate the blunt-notched behaviour of a GLARE specimen. Also, the FE modelling approach is employed to study the effects of layup configurations and numbers of holes on the mechanical response of perforated aluminium/CFRP 3/2 FMLs subjected to tensile loadings at their ends. 3.1. GLARE 4B-4/3-0.5 laminate The GLARE 4B-4/3-0.5 specimen is made of four aluminium 2024-T3 alloy layers alternating with three S2-glass/FM94-epoxy composite laminates. The composite layup between two alumin   ium layers of the specimen is ½90 =0 =90 . The GLARE 4B laminate has a length, width and thickness of 300, 50 and 3.125 mm, respectively. The thicknesses of each aluminium and composite layer are 0.5 and 0.125 mm. The diameter of central hole is 4.8 mm. The experimental investigation of this specimen was reported in [17]. The elastic material properties and strengths of the S2-glass/ FM94-epoxy composite are obtained from [50] and listed in Table 1. The cohesive material properties adopted for the simulation of the GLARE specimen is obtained from [51] and listed in Table 2. The values of fracture toughness parameters associated with fibre tensile and compressive failure G1t;c ¼ 81:5 N=mm and G1c;c ¼ 106:3 N=mm for IM7/8552 reported in [52] are adopted in this study. The value of the Mode II inter-laminar fracture toughness of S2-glass/ FM94-epoxy prepreg given as 9 N/mm in [53] is used as a substitute for the Mode II intra-laminar fracture toughness parameter G6;c . Using this value, the parameter G2c;c is calculated according to the  equation suggested in [54], i.e., G2c;c ¼ G6;c =cos53 ¼ 14:955 N=mm. The value of the Mode I intra-laminar fracture toughness parameter G2t;c is assumed to have the same value as that of parameter G2c;c . The parameters G1t;c ; G1c;c ; G2t;c ; G2c;c and G6;c are the material fracture toughnesses corresponding to each failure mechanism (see [26]). The plastic model parameters reported by [55] for S2-glass/8552-epoxy composite are adopted for the simulation of this laminate. All these parameters are listed in Table 1. For the

aluminium 2024-T3 alloy layers, the Young’s modulus is 73.8 GPa and Poisson’s ratio is 0.33. The isotropic hardening data are listed in Table 3. These material properties and input data for isotropic hardening were reported by Linde and de Boer [50]. In Abaqus, the yield stress at a given stress state is interpolated from the data listed in the table, and it remains constant when the equivalent plastic strain exceeds the last value provided in the tabular. Since necking usually occurs at 15% plastic strain [19], this type of behaviour is not defined in the FE model. Considering the symmetry of the geometry and loading of the laminate, only one-quarter of it is simulated. As previously discussed, an adequate number of cohesive elements should be placed in the cohesive zone. The element size in the region near the central hole in the simulation is 0.5  0.5 mm which should provide a sufficiently fine mesh to accurately capture the stresses in the cohesive zone ahead of the delamination tip. Since the laminate under consideration is subjected to in-plane tensile loading, delamination propagation is primarily caused by the interlaminar transverse shear stresses. The element size to be used in the simulation is governed by the cohesive zone length under Mode II loading. The cohesive zone length is calculated according to Eq. (16) and the procedure for determining the value of E0II presented in [47]. The other material properties required for calculating the cohesive zone length of this laminate are: E03 ¼ E02 ¼ 9:412 GPa, G013 ¼ G012 ¼ 5:548 GPa, G023 ¼ 3:0 GPa, m13 ¼ m12 ¼ 0:33; m23 ¼ 0:45. The half-thickness of the laminate is 1.56 mm. Therefore, the cohesive zone length for Mode II loading is determined as 12.3 mm. Correspondingly, there are at least 24 elements in the cohesive zone and an element size of 0.5 mm should provide a fine mesh for the laminate. In other region, a larger mesh size of 1  1 mm is used. The predicted and experimental stress versus strain curves of the laminate are shown in Fig. 2. As can be seen, the predicted curve agrees reasonably well with the experimental one at the initial elastic stage. At a later stage, as plasticity evolves in the aluminium layers, the mechanical behaviour of the GLARE 4B laminate is governed by the aluminium layers. The discrepancy between the experimental and predicted curves may result from the use of the hardening data which are obtained from other

Table 1 Material properties of S2-glass/FM94-epoxy composite and model parameters. E01

E02

G012

m012

s1t

s1c

s2t

s2c

53.98 GPa

9.412 GPa

5.548 GPa

0.33

2430 MPa

2000 MPa

50 MPa

160 MPa

50 MPa

G1t;c

G1c;c

G2t;c

G2c;c

G12;c

a

b

n

g

81.5 N/mm

106.3 N/mm

14.955 N/mm

14.955 N/mm

9 N/mm

1.4

540.396

0.25641

0.0002

ss

Where s1t ; s1c ; s2t ; s2c and ss are the material strength parameters; G1t;c ; G1c;c ; G2t;c ; G2c;c and G6;c are the material fracture toughness parameters corresponding to each failure mechanism; g is the viscosity coefficient [26].

Table 2 Material properties of GLARE 4B cohesive layers. K 0z ¼ K 0zx ¼ K 0zy 6

3

10 N/mm

r0z

s0zx

s0zy

GIc

GIIc

GIIIc

50 MPa

25 MPa

20 MPa

1.1 N/mm

1.1 N/mm

1.1 N/mm

Table 3 Isotropic hardening data for aluminium 2024-T3. Yield stress (MPa)

300

320

340

355

375

390

410

430

450

470

484

Equivalent plastic strain (%)

0.00

0.016

0.047

0.119

0.449

1.036

2.130

3.439

5.133

8.00

14.71

70

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

3.2.1. Laminate configurations and modelling aspects Each aluminium/CFRP 3/2 FML under consideration consists of three aluminium 2024-T3 alloy layers and two embedded AS4/ PEEK carbon fibre-reinforced thermoplastic laminates with one of    the following three layups: unidirectional ½04 , cross-ply ½0 =90 s    or quasi-isotropic ½0 =  45 =90 . The thicknesses of each aluminium alloy layer and AS4/PEEK layer are 0.305 and 0.127 mm, respectively, and that of the adhesive layer 0.001 mm. Schematic illustrations of the FMLs with different layups are presented in Fig. 4. The FMLs are 150 mm long (L), 40 mm wide (W) and 1.941 mm thick (t). The hole numbers contained in the perforated FMLs are 0, 1, 3, 6 and 9. The diameter D of these holes is 4 mm. The distances between neighbouring hole centres in the width and longitudinal directions are S1 ¼ 15 mm and S2 ¼ 50 mm, respectively. The configurations of these FMLs with various numbers of holes are shown in Fig. 5. In total, 15 numerical analyses were performed. The material properties and model parameters of the AS4/PEEK composites obtained from [26] and listed in Table 4, are used in the simulations. In this study, the value of 106 N=mm3 is chosen for the three initial penalty stiffnesses required as it has been successfully applied in delamination simulations of AS4/PEEK composites in [41]. The Mode I, II and III delamination failure strengths and the critical fracture energies for Mode I and II delamination failures of the cohesive layers are obtained from [41] and the Mode III critical fracture energy GIIIc from [56]. These cohesive material properties are listed in Table 5. The same material properties and hardening data as those used in the simulation of the aluminium layers in the GLARE 4B laminate (see Table 3) are applied to the present modellings. For laminates with a unidirectional or cross-ply layup, considering the symmetry of the laminate and its loading and to reduce simulation runtime, only one-eighth of it is modelled. The boundary conditions of the one-eighth single-hole quasi-isotropic FML FE model are illustrated in Fig. 6. Symmetric boundary conditions are

Fig. 2. Comparison of the experimental and predicted stress versus strain curves of the GLARE 4B laminate.

source (i.e. [50]) as they were not reported in [17]. The predicted   patterns of damage variables d1 and d2 in the first 0 and 90 layers of the GLARE 4B laminate at the end of the analysis are illustrated in Fig. 3. It follows from these graphs that the damage caused by  fibre breakage or matrix cracking in the 0 layer is localised in the vicinity of the hole, whereas the matrix damage spreads over  the whole 90 layer and fibre damage has not been initiated in this layer at all. 3.2. Perforated aluminium/CFRP 3/2 FMLs Also, the modelling approach considered in this study is used to investigate the effects of layup configurations and numbers of holes on the mechanical response of perforated aluminium/CFRP 3/2 FMLs subjected to tensile loadings.





Fig. 3. Predicted patterns of damage variables d1 and d2 in the first 0 and 90 layers of GLARE 4B laminate at the end of the analysis.

Fig. 4. Aluminium/CFRP laminates with various layups.

71

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

cohesive zone ahead of the delamination tip in the interface layers are well captured. 3.2.2. Results and discussion Comprehensive information related to the progressive failures of FMLs can be obtained from the appropriate FE analyses. In the following subsections, the effects of FMLs’ composite layup configurations and numbers of holes on the failure loads, notched tensile strengths, overall load versus displacement curves, stress distributions, damage progression patterns and failure sequences of the FMLs are discussed.

Fig. 5. Configurations of the aluminium/CFRP fibre–metal laminates.

prescribed for its three planes of symmetry. Each FML is subjected to tensile loading at its ends. A reference point RP is appointed to the right upper corner of its loaded surface (see Fig. 6) to which displacement-controlled tensile loading is assigned. A linear constraint equation is built between the reference point and loaded surface to constrain the displacements of all the nodes on the latter to have displacements equal to the applied displacement of the former. By doing this, the reaction force and corresponding displacement on the loaded surface can be conveniently output to this reference point. For laminates with a quasi-isotropic layup, full FE model of the whole laminate is built. The element size of the numerical simulations of the FMLs with unidirectional or crossply layups is 0.5 mm. As mentioned earlier, the number of elements placed in the cohesive zone ahead of the delamination tip should be sufficient to capture the stress distributions within this zone. Two to three elements should be sufficient to provide accurate predictions. Since the FMLs under consideration are subjected to in-plane tensile loading, as previously mentioned, delamination is governed by shear-type loading. The size of cohesive elements to be used in the simulation is governed by the cohesive zone length under Mode II loading. Similarly, the cohesive zone length under Mode II loading is calculated according to Eq. (16) and the procedure for determining the value of E0II presented in [47]. The other material properties required for calculating the cohesive zone length of the laminates are found in [41]: E03 ¼ E02 ¼ 10:3 GPa, G013 ¼ G012 ¼ 6:0 GPa, G023 ¼ 3:7 GPa, m13 ¼ m12 ¼ 0:32; m23 ¼ 0:45. The half-thickness of the laminates is 0.9705 mm. The cohesive zone length of these FMLs under Mode II loading is 4.63 mm, and at least 9 elements are placed in the cohesive zone under Mode II loading. For laminates with quasi-isotropic layup, an element size of 1 mm is used in the simulations. These two element sizes should ensure that accurate stress distributions within the

3.2.2.1. Overall load versus displacement curves and failure loads. The predicted load versus displacement curves of the aluminium/ CFRP FMLs (with various layups) containing a specific number of hole(s) are shown in Fig. 7. According to the numerical results, unidirectional FMLs provide the highest failure loads of all the FMLs with the three layups considered, irrespective of their numbers of hole(s). As observed in Fig. 7(a), the three un-notched laminated systems with different layups fail at approximately the same applied displacement. From this phenomenon, it is inferred that  the failure of fibres in the 0 composite plies leads to the immediate catastrophic collapses of these laminates. Examining damage developments in the three un-notched FMLs obtained from the simulations confirms that the occurrence of fibre damage in their  0 plies leads immediately to their catastrophic failures. As demonstrated in Fig. 7(b)–(e), for perforated FMLs containing specific numbers of holes, those with cross-ply and quasi-isotropic layups have close tensile loading capacities with the former providing higher tensile stiffnesses and the latter better ductilities. Also, FMLs with unidirectional and cross-ply layups fail at approximately the same applied displacements whereas those with quasi-isotropic layups do not reach their maximum failure loads at the same applied displacements. It follows from this phenome  non that the inclusion of þ45 and 45 composite layers contributes to additional ductility and the delay of failure in the notched quasi-isotropic FMLs. As many experiments have revealed, the  ½45 ns AS4/PEEK composite laminates experience significant plastic deformations under longitudinal loading [57,58]. The plastic deformations developed under shear loading are more pronounced than that developed under transverse loading. Due to the existence of holes, FMLs exhibit stress redistributions which  allow plastic deformations in the 45 composite layers to be well developed and postpone the failure of the notched laminates with quasi-isotropic configurations (compared with un-notched quasiisotropic FMLs in which the absence of notches results in poor stress relief). The predicted load versus displacement curves of the perforated AS4/PEEK-based FMLs with a specific layup and different numbers

Table 4 Material properties of AS4/PEEK and model parameters. E01

E02

G012

m012

s1t

s1c

s2t

s2c

ss

127.6 GPa

10.3 GPa

6.0 GPa

0.32

2023.0 MPa

1234.0 MPa

92.7 MPa

176.0 MPa

82.6 MPa

G1t;c

G1c;c

G2t;c

G2c;c

G12;c

a

b

n

g

128.0 N/mm

128.0 N/mm

5.6 N/mm

9.31 N/mm

4.93 N/mm

1.5

295.0274

0.142857

0.0002

Table 5 Material properties of the cohesive layers between the aluminium and AS4/PEEK composite layers. K 0z ¼ K 0zx ¼ K 0zy

r0z

s0zx ¼ s0zy

GIc

GIIc

GIIIc

106 N/mm3

80 MPa

100 MPa

0.969 N/mm

1.719 N/mm

2.01 N/mm

72

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

Fig. 6. Boundary conditions of the aluminium/CFRP laminates (the single-hole FML with the cross-ply layup is shown).

Fig. 7. Predicted load versus displacement curves of aluminium/CFRP FMLs with different numbers of holes.

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

of holes are shown in Fig. 8(a) (unidirectional), Fig. 8(b) (cross-ply) and Fig. 8(c) (quasi-isotropic). From these figures, it is apparent that, due to the presence of a single hole, the ultimate failure loads of unidirectional and quasi-isotropic FMLs decrease dramatically and the least of all three layups, respectively. These phenomena indicate that the unidirectional FMLs exhibit the highest notch

73

sensitivities and the quasi-isotropic FMLs the lowest. For multiple-hole FMLs with a specific layup, their failure loads are close. This is because the same number of holes (i.e. 3) are arranged in a row in the width direction and all rows and columns of holes are placed parallel to each other, the net tension areas of their cross sections normal to the loading direction are the same. Also, the distances between the parallel rows of holes in the longitudinal direction (i.e. S2 ) are sufficiently large to prevent adverse interactions of the stress distributions near the adjacent rows of holes as revealed in the stress distribution graphs of the composite and metal layers (e.g., in the quasi-isotropic FML with 9 holes illustrated in Figs. 17–19). Also, in Fig. 8(a)–(c), it is observed that the stiffnesses of multi-hole perforated FMLs with more holes are slightly less than those with fewer holes but their ductilities are slightly better because, in the former, delamination can occur in more damaged sites and lead to stress redistributions occurring in more areas which is beneficial for reducing the maximum stresses in their composite plies and as a consequence, their failures can be delayed. Also, it follows from Fig. 8(a)–(c) that increasing the numbers of holes in the same notched plane orthogonal to the loading direction of single-notched FMLs does not lead to as dramatic a decrease in their failure loads as in those of the un-notched FMLs. 3.2.2.2. Notched tensile strengths. In order to quantify the strength reductions, the gross and net notched tensile strengths as well as the gross and net notched tensile strength reduction ratios of the notched FMLs are calculated. The gross notched tensile strength rgross is defined as the value of the far field uniform tensile stress at failure and calculated using the ultimate failure load P u divided by the product of the specimen width W and thickness t as

rgross ¼

Pu Wt

ð17Þ

The net notched tensile strength is defined as the maximum value of the average net-section stress that a notched specimen is able to sustain and is calculated by the ultimate failure load P u divided by the net cross section length (the difference between the specimen width W and notched length nD, where n is the hole number and D is the hole diameter) and specimen thickness t as

rnet ¼

Fig. 8. Predicted load versus displacement curves of aluminium/CFRP FMLs with different layup configurations.

Pu tðW  nDÞ

ð18Þ

The gross and net notched strength reduction ratios are calculated as Rgross ¼ ðrult;0  rgross Þ=rult;0 and Rnet ¼ ðrult;0  rnet Þ=rult;0 , respectively, where rult;0 is the ultimate tensile strength of an unnotched laminate. The estimated gross and net notched tensile strengths of all the laminates as well as their corresponding strength reduction ratios are listed in Table 6. It is apparent from this table that the notch sensitivities of the aluminium/CFRP 3/2 FMLs are highest in those with unidirectional layups and lowest in those with quasi-isotropic layups. This finding is more obvious when the notched tensile strengths are expressed in terms of their net components. The gross notched tensile strength reduction ratio for the single-hole FML with unidirectional composite laminates is as high as 46.13% and its corresponding net counterpart is 40.14% whereas, for the single-hole FML with the quasi-isotropic layup, they are 17.0% and 7.77%, respectively. These data clearly demonstrate that the FMLs with quasi-isotropic layups exhibit the lowest notch sensitivities which is consistent with their highest ductilities due to the  inclusion of 45 composite plies. (see Fig. 7(b)–(e)). This phenomenon has also been observed from tests performed on cross-ply  GLARE 3 specimens loaded in the off-axis 45 direction [13]. For the multi-hole perforated FMLs with the same layups, as for their failure loads, their strength reduction ratios have similar values irrespective of their numbers of holes (3, 6 and 9). Also, it is

74

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

Table 6 Predicted failure loads and ultimate strengths of aluminium-based AS4/PEEK FMLs. Hole nos. Predicted failure loads P u (N) Gross tensile strengths (MPa)

No hole 1 hole 3 holes 6 holes 9 holes

rgross

Net tensile strengthsrnet (MPa)

Rgross (%)

Rnet (%)

Layup 1

Layup 2

Layup 3

Layup 1

Layup 2

Layup 3

Layup 1

Layup 2

Layup 3

Layup 1 Layup 2 Layup 3 Layup 1 Layup 2 Layup 3

95610.5 51509.3 43661.3 41550.3 39415.3

57691.2 35449.9 29426.9 28726.8 28756.8

43939.6 36471.2 30284.3 29478.5 29993.2

1231.46 663.44 562.36 535.17 507.67

743.06 456.59 379.02 370.0 370.39

565.94 469.75 390.06 379.68 386.31

1231.46 737.15 803.37 764.52 725.24

743.06 507.33 541.45 528.57 529.12

565.94 521.94 557.23 542.40 551.87

0.00 46.13 54.33 56.54 58.78

0.00 38.55 48.99 50.21 50.15

0.00 17.00 31.08 32.91 31.74

0.00 40.14 34.76 37.92 41.11

0.00 31.72 27.13 28.87 28.79

0.00 7.77 1.54 4.16 2.49

Layup 1: unidirectional; Layup 2: cross-ply; Layup 3: quasi-isotropic.

interesting to observe that their net notched tensile strengths are slightly higher than those of their single-hole counterparts. Accordingly, their net tensile strength reduction ratios are less. This phenomenon can be explained using the cohesive crack model applied by Camanho et al. [52] to explain the size effects in notched carbon/epoxy laminates. The cohesive crack model is described using a stress versus crack-opening displacement constitutive law. In the fracture process zone, the cohesive stress acting on the cracked surfaces is assumed to be a function of the crack-opening displacement. The area under the cohesive stress versus crack-opening displacement curve represents the critical fracture energy of the material Gc Rw (Gc ¼ 0 c rdw, where wc is the critical crack-opening displacement, i.e. the distance between the cracked surfaces where the cohesive stress vanishes). Several methods for determining the length of the fracture process zone LFPZ are reviewed by Turon et al. [44] and Carpinteri et al. [59] discussed the development of cohesive crack models. In this paper, for simplicity, a simplified cohesive constitutive model is assumed for the stress state in the cohesive zone. Considering a plate made of homogenised material and containing a central hole, the stress developed in the cohesive zone near the hole under loading is assumed to be constant and limited to its ultimate strength st . This simplified cohesive constitutive model is illustrated in Fig. 9(a). According to this simplified constitutive law, the stress distributions at failure on the assumed fracture planes of the specimens in the fracture process zones with one hole and three holes are illustrated in Fig. 9(b) and (c). It can be seen that, in the specimen with a single hole, the fracture process zones are only a small proportion of its width whereas, in the specimen with three holes, they cover a larger area. As the average stress in the fracture process zone is much higher than that in the remaining part of the specimen width, it is clear that the net notched tensile strengths of the multi-hole specimens are higher than those of the single-hole laminates.

According to the above discussion, calculating and measuring the net notched tensile strengths of notched FMLs is important. The values of this parameter for perforated FMLs with various hole numbers tend to be equal to a constant and are slightly higher in multi-hole perforated FMLs than in single-hole FMLs. When designing multi-hole perforated FML tension members, the net notched tensile strengths obtained from tests on those with a central hole could be adopted if there is a lack of data regarding those of multi-hole perforated FMLs. If the net notched tensile strength of a laminate is chosen, the design can be expected to be conservative but if the gross notched tensile strength is adopted, it tends to be unsafe. As FMLs with quasi-isotropic layups have the lowest notch sensitivities and highest ductilities, it is recommended that such laminates be chosen to avoid catastrophic failures in the case of unexpected damage occurring to them during their service lives. 3.2.2.3. Damage progression and stress distribution patterns. The damage progression and stress distribution patterns for the composite and metal layers in any given state can be obtained from the analysis results. In the interest of brevity, only those of quasi-isotropic FMLs with a single hole and 9 holes at failure loads are illustrated in Figs. 10–19. As composite plies behave elastically in the fibre direction, fibre  damage in the 0 composite plies of these two notched FMLs propagates along the net cross section (i.e. transverse to the loading direction) and is localised into the areas initiating from the notch tip, as illustrated in Fig. 10. In contrast, composite plies behave in a plastic and damaging way in the transverse direction, plastic deformations occurring in this direction are beneficial for reducing  their stress concentrations. Matrix damage in the 90 composite plies spreads over a larger area in the vicinity of the hole(s), as shown in Fig. 11. Fig. 12(a) illustrates the contour plots of the equivalent plastic strains in the middle metal layer of the single-hole laminate. Apparently, significant amounts of plastic deformation are

Fig. 9. Cohesive crack model and stress distributions along the fracture planes.

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

75



Fig. 10. Predicted fibre damage pattern in the first 0 composite layers of the two quasi-isotropic FMLs with 1 and 9 hole(s).



Fig. 11. Predicted matrix damage pattern in the first 90 composite layers of the two quasi-isotropic FMLs with 1 and 9 hole(s).

developed at the notch tips with the maximum equivalent plastic strain developed in the middle layer at failure as high as 13.20%. Fig. 12(b) shows the distribution patterns of the equivalent plastic strains in the middle metal layer of the laminate containing 9 holes. As can be seen, the plastic deformations developed in the aluminium layers localise in small areas in the notch tips throughout the loading process while the matrix damage developed in the  90 composite layers propagates quite far away from the hole edges (see Fig. 11). The maximum equivalent plastic strains in the aluminium layers corresponding to the maximum loads of all the FMLs are listed in Table 7 from which it is clear that substantial equivalent plastic strains develop in the vicinity of the holes in each notched FML. Only three related to the maximum loads (i.e. 15.19% in the single-hole unidirectional FML and 16.68% and 16.24% in the 6-hole and 9-hole cross-ply FMLs) slightly exceed 15% while the

others are in the range of the pre-defined isotropic hardening curve which indicates that, in most analyses, necking happens after the maximum loads are reached. Therefore, it is appropriate not to define necking behaviour in the aluminium layers in this study. In the un-notched FMLs, their maximum plastic equivalent plastic strains are rather small due to their poor stress redistribution capabilities. The failures of the FMLs under consideration are governed by the composite layers. The delamination patterns in the interface layers between the   90 and 45 plies of the quasi-isotropic FMLs with 1 and 9 holes at their failure loads are illustrated in Fig. 13(a) and (b). It follows from these graphs that a considerable amount of delamination has been developed in the vicinity of the notch tips in both FMLs.   The stress distributions in the 0 and 90 composite and middle aluminium plies of the two FMLs are illustrated in Figs. 14–19. It is evident that the presence of a hole or holes significantly affects the

76

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

Fig. 12. Predicted equivalent plastic strain in the middle aluminium layers of the two quasi-isotropic FMLs with 1 and 9 hole(s).





Fig. 13. Predicted delamination pattern in the interface layers between the 90 and 45 layers of the two quasi-isotropic FMLs with 1 and 9 hole(s).

Fig. 14. Predicted stress



r1 in the first 0 composite layer of the quasi-isotropic FML with a single hole.

stress distribution over a large proportion of a FML which emanates from the region of the hole(s). It follows from the stress distribution patterns of the quasi-isotropic FML with 9 holes

(Figs. 17–19) that the distance between two neighbouring hole centres in the longitudinal direction (i.e. S2 ) is sufficient to prevent the stresses in each layer in the vicinity of the holes from

77

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

Fig. 15. Predicted stress



r2 in the first 90 composite layer of the quasi-isotropic FML with a single hole.

Fig. 16. Predicted Mises stress in the middle aluminium layer of the quasi-isotropic FML with a single hole.



Fig. 17. Predicted stress

r1 in the first 0 composite layer of the quasi-isotropic FML with 9 holes.

Fig. 18. Predicted stress

r2 in the first 90 composite layer of the quasi-isotropic FML with 9 holes.



Fig. 19. Predicted Mises stress in the middle aluminium layer of the quasi-isotropic FML with 9 holes.

Table 7 Equivalent plastic strains in aluminium layers at maximum loads. Hole nos.

0 hole (%)

1 hole (%)

3 holes (%)

6 holes (%)

9 holes (%)

Unidirectional Cross-ply Quasi-isotropic

1.07 1.72 1.39

15.19 14.37 13.20

10.88 10.09 10.53

12.51 16.68 14.38

13.42 16.24 11.92

78

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

interacting with each other in contrast to those in the width direction which can result in adverse interactions. Corresponding to their high ductility and localised equivalent plastic strain distributions (as shown in Fig. 12(a) and (b)), the Mises stress distributions in the aluminium layers of both the quasi-isotropic FMLs are also localised into small areas in the vicinity of the holes, as demonstrated in Figs. 16 and 19. 3.2.2.4. Failure sequences. The failure sequences of various failure mechanisms in the two quasi-isotropic FMLs are examined and the points corresponding to the initiation of equivalent plastic strain in an aluminium layer, matrix cracking, delamination and fibre breakage are marked accordingly in Fig. 20(a) and (b). As can be seen, plasticity evolves in the aluminium layers prior to damage being initiated by other failure mechanisms, followed by  matrix damage in the 90 composite ply and fibre damage in the  0 composite ply after which delamination is initiated. The plasticity and damage accumulation processes happen at an early stage of the loading history, long before the ultimate failures of the notched quasi-isotropic FMLs, which indicates the significance of including progressive post-failure behaviour in a composite material model. Without including such effects, material models may greatly underestimate the ultimate failure loads of FMLs. The accumulation of plasticity in the aluminium layers of each FML causes a somewhat more pronounced stiffness loss at approximately the same applied displacement (see Fig. 7(a)–(e)).

Fig. 20. Predicted failure sequences of the two quasi-isotropic FMLs with 1 and 9 hole(s).



Although matrix damage development is observed in the 90 composite plies before the significant stiffness loss occurs in the two notched quasi-isotropic FMLs, such loss in the FMLs is not mostly attributed to matrix damage accumulation in these plies. This can be proven from the load versus displacement curve of the un-notched unidirectional FML in which no matrix damage occurs before its significant loss of stiffness. The fibre damage occurs in  the 0 plies and its initiation leads immediately to the catastrophic failure of the un-notched unidirectional FML.

4. Conclusions A FE model has been developed for progressive failure analyses of FMLs including in-ply and delamination damage effects. It consists of the combined elastoplastic damage constitutive model for simulating various in-plane progressive failure mechanisms and plasticity effects in composite laminates, a cohesive zone model implemented in cohesive elements for modelling delamination in the adhesive layers and an isotropic elastic–plastic model with isotropic hardening which represents the behaviour of the aluminium metal layers. Both the cohesive zone model and isotropic elastic– plastic model are available in the FE code Abaqus. Thus, interactions between the in-plane and out-of-plane failure mechanisms in the FMLs have been accounted for in the present FE model. This modelling method has been first applied to investigate the blunt-notched behaviour of a GLARE specimen. It has been shown that the present modelling approach provides reasonable prediction of the stress versus strain curve of the GLARE laminate. And then, it has been employed to study the effects of the layup configurations of embedded composite laminates and influences of the number of holes on the mechanical response of aluminium/CFRP FMLs. Concerning the effects of the layup configurations, it has been found that, of the FMLs with the same geometry, those having unidirectional layups provided the highest failure loads of all the three layups considered. The three un-notched FMLs (with three different layups) had poor stress redistribution capabilities and failed at approximately the same applied displacement because fibre failure dominated their failure mechanisms. Of these unnotched FMLs, the one having the quasi-isotropic layup has the lowest stiffness and offers the lowest failure load. Also, unidirectional and cross-ply FMLs (either notched or un-notched) with the same geometry failed at approximately the same applied dis placement. With the presence of holes and 45 layers, each quasi-isotropic FML provided better ductility than their notched unidirectional and cross-ply counterparts and offered similar failure loads with the latter. Regarding the influence of holes, they considerably reduced the tensile strengths of all the FMLs. Of the three layup configurations considered, the existence of one hole led to a dramatic decrease of the ultimate failure load in the unidirectional FML and the least in the quasi-isotropic FML. In comparison with the FMLs containing a single hole, the failure loads of all the FMLs with multiple holes were only slightly further reduced. The findings from the numerical simulations in this work show that the quasi-isotropic FMLs with holes exhibited the lowest notch sensitivities with net notched reduction ratio ranging from 1.54% to 7.77% which is consistent with their having the highest ductilities in comparison with FMLs with other layups. The stress distributions and damage progression patterns in the composite and metal layers of the quasi-isotropic FMLs with 1 hole and 9 holes, and their failure sequences were also discussed. Finally, it has been pointed out that it is crucial to include progressive failure models in the progressive failure analyses of FMLs as damage could accumulate long before the ultimate loads of FMLs are reached.

J.-F. Chen et al. / Composite Structures 114 (2014) 64–79

Acknowledgements One of the authors (J.F. Chen) gratefully acknowledges the financial support provided by the China Scholarship Council to enable this research to be undertaken. Also, the authors would like to acknowledge the support provided by the Intersect Australia Ltd. and the National Computational Infrastructure (NCI) in Canberra, Australia, which is supported by the Australian Commonwealth Government. References [1] Krishnakumar S. Fiber metal laminates – the synthesis of metals and composites. Mater Manuf Process 1994;9(2):295–354. [2] Cortes P, Cantwell WJ. The fracture properties of a fibre–metal laminate based on magnesium alloy. Composites Part B 2006;37(2–3):163–70. [3] Papakyriacou M, Schijve J, Stanzl-Tschegg SE. Fatigue crack growth behaviour of fibre–metal laminate GLARE-1 and metal laminate 7475 with different blunt notches. Fatigue Fract Eng Mater Struct 1997;20(11):1573–84. [4] Cortes P, Cantwell WJ. The tensile and fatigue properties of carbon fiberreinforced PEEK-Titanium fiber-metal laminates. J Reinf Plast Compos 2004;23(15):1615–23. [5] Vlot A, Kroon E, La Rocca G. Impact response of fiber metal laminates. Key Eng Mater 1998;141–143:235–76. [6] Cortes P, Cantwell WJ. The prediction of tensile failure in titanium-based thermoplastic fibre-metal laminates. Compos Sci Technol 2006;66(13):2306–16. [7] Botelho EC, Silva RA, Pardini LC, Rezende MC. A review on the development and properties of continuous fiber/epoxy/aluminum hybrid composites for aircraft structures. Mater Res 2006;9(3):247–56. [8] Vlot A, Gunnink JW, editors. Fibre metal laminates: an introduction. Dordrecht, London: Kluwer Academic Publishers; 2001. [9] Frizzell RM, McCarthy CT, McCarthy MA. An experimental investigation into the progression of damage in pin-loaded fibre metal laminates. Composites Part B 2008;39(6):907–25. [10] Chang PY, Yeh PC, Yang JM. Fatigue crack growth in fibre metal laminates with multiple open holes. Fatigue Fract Eng Mater Struct 2011;35(2):93–107. [11] Wu GC, Tan Y, Yang JM. Evaluation of residual strength of notched fiber metal laminates. Mater Sci Eng A-Struct 2007;457(1–2):338–49. [12] Yeh PC, Chang PY, Yang JM, Wu PH, Liu MC. Blunt notch strength of hybrid boron/glass/aluminum fiber metal laminates. Mater Sci Eng A-Struct 2011;528(4–5):2164–73. [13] Kawai M, Arai Y. Off-axis notched strength of fiber-metal laminates and a formula for predicting anisotropic size effect. Composites Part A 2009;40(12):1900–10. [14] Aktas A, Dirikolu MH. The effect of stacking sequence of carbon epoxy composite laminates on pinned-joint strength. Compos Struct 2003;62(1):107–11. [15] Seyed Yaghoubi A, Liu Y, Liaw B. Stacking sequence and geometrical effects on low-velocity impact behaviors of GLARE 5 (3/2) fiber-metal laminates. J Thermoplast Compos 2012;25(2):223–47. [16] Beumler T. Flying GLARE-a contribution to aircraft certification issues on strength properties in non-damaged and fatigue damaged GLARE structures. Ph.D. thesis; Delft University of Technology; 2004. [17] de Vries TJ. Blunt and sharp notch behaviour of Glare laminates. Ph.D. thesis; Delft University of Technology; 2001. [18] Vermeeren CAJR. The blunt notch behaviour of metal laminates: ARALL and GLARE. Tech Rep Report LR-617. Faculty of Aerospace Engineering, Delft University of Technology. The Netherlands; 1990. [19] Van Rooijen RGJ, Sinke J, De Vries TJ, Van Der Zwaag S. The bearing strength of fiber metal laminates. J Compos Mater 2006;40(1):5–19. [20] Lapczyk I, Hurtado JA. Progressive damage modeling in fiber-reinforced materials. Composites Part A 2007;38(11):2333–41. [21] Frizzell RM, McCarthy CT, McCarthy MA. Simulating damage and delamination in fibre metal laminate joints using a three-dimensional damage model with cohesive elements and damage regularisation. Compos Sci Technol 2011;71(9):1225–35. [22] Frizzell RM, McCarthy CT, McCarthy MA. Predicting the effects of geometry on the behaviour of fibre metal lamiante joints. Compos Struct 2011;93(7):1877–89. [23] Sugiman S, Crocombe AD, Katnam KB. Investigating the static response of hybrid fibre–metal laminate doublers loaded in tension. Composites Part B 2011;42(7):1867–84. [24] Sugiman S, Crocombe AD. The static and fatigue response of metal laminate and hybrid fibre–metal laminate doublers joints under tension loading. Compos Struct 2012;94(9):2937–51. [25] Nakatani H, Kosaka T, Osaka K, Sawada Y. Damage characterization of titanium/GFRP hybrid laminates subjected to low-velocity impact. Composites Part A 2011;42(7):772–81. [26] Chen JF, Morozov EV, Shankar K. A combined elastoplastic damage model for progressive failure analysis of composite materials and structures. Compos Struct 2012;94(12):3478–89.

79

[27] Abaqus. Dassault Systèmes Simulia Corp. Providence, RI, USA. 6.10-3 ed.; 2010. [28] Belytschko T, Liu WK, Moran B. Nonlinear finite elements for continua and structures. New York: John Wiley and Sons; 2000. [29] Sun CT, Chen JL. A simple flow rule for characterizing nonlinear behavior of fiber composites. J Compos Mater 1989;23(10):1009–20. [30] Winn VM, Sridharan S. An investigation into the accuracy of a one-parameter nonlinear model for unidirectional composites. J Compos Mater 2001;35(16):1491–507. [31] Hashin Z. Failure criteria for unidirectional fiber composites. J Appl Mech-T ASME 1980;47(2):329–34. [32] Faria R, Oliver J, Cervera MA. A strain-based plastic viscous-damage model for massive concrete structures. Int J Solids Struct 1998;35(14):1533–58. [33] Simo JC, Taylor RL. Consistent tangent operators for rate-independent elastoplasticity. Comput Method Appl M 1985;48:101–18. [34] Rybicki EF, Kanninen MF. A finite element calculation of stress intensity factors by a modified crack closure integral. Eng Fract Mech 1977;9(4):931–8. [35] Krueger R, Minguet PJ, O‘Brien TK. A method for calculating strain energy release rates in preliminary design of composite skin/stringer debonding under multi-axial loading. In: Composite structures: theory and practice. ASTM STP 1383; NASA Langley Research Center; 2000. p. 105–28. [36] Balzani C, Wagner W. An interface element for the simulation of delamination in unidirectional fiber-reinforced composite laminates. Eng Fract Mech 2008;75(9):2597–615. [37] Atas A, Mohamed GF, Soutis C. Modelling delamination onset and growth in pin loaded composite laminates. Compos Sci Technol 2012;72(10):1096–101. [38] Yamaguchi T, Okabe T, Yashiro S. Fatigue simulation for titanium/CFRP hybrid laminates using cohesive elements. Compos Sci Technol 2009;69(11– 12):1968–73. [39] 29.5.6 Defining the constitutive response of cohesive elements using a traction-separation description in Abaqus Analysis User’s Manual. Dassault Systèmes Simulia Corp, Providence, RI, USA. 6.10 ed.; 2010. [40] Davila CG, Camanho PP, de Moura MF. Mixed-mode decohesion elements for analyses of progressive delamination. In: Proceedings of the 42nd AIAA/ASME/ ASCE/AHS/ASC structures. No. AIAA-01-1486 in structural dynamics and materials conference. Seattle, Washington; 2001. p. 1–12. [41] Camanho PP, Davila CG, de Moura MF. Numerical simulation of mixed-mode progressive delamination in composite materials. J Compos Mater 2003;37(16):1415–38. [42] Zou Z, Reid SR, Li S, Soden PD. Modelling interlaminar and intralaminar damage in filament-wound pipes under quasi-static indentation. J Compos Mater 2002;36(4):477–99. [43] Daudeville L, Allix O, Ladeveze P. Delamination analysis by damage mechanics: some applications. Compos Eng 1995;5(1):17–24. [44] Turon A, Davila CG, Camanho PP, Costa J. An engineering solution for mesh size effects in the simulation of delamination using cohesive zone models. Eng Fract Mech 2007;74(10):1665–82. [45] Ridha M, Tan VBC, Tay TE. Traction-separation laws for progressive failure of bonded scarf repair of composite panel. Compos Struct 2011;93(4):1239–45. [46] Pinho ST, Iannucci L, Robinson P. Formulation and implementation of decohesion elements in an explicit finite element code. Composites Part A 2006;37(5):778–89. [47] Harper PW, Hallett SR. Cohesive zone length in numerical simulations of composite delamination. Eng Fract Mech 2008;75(16):4774–92. [48] Yang Q, Cox B. Cohesive models for damage evolution in laminated composites. Int J Fract 2005;133(2):107–37. [49] Falk ML, Needleman A, Rice JR. A critical evaluation of cohesive zone models of dynamic fracture. In: van der Giessen E, Forest S, Kubin L, editors. Proceedings of the 5th European mechanics of materials conference on scale transitions from atomistics to continuum plasticity, J Phys IV France, Vol. 11. Delft, the Netherlands: European Mechanics of Materials Conference; 2001. p. 43–50. [50] Linde P, de Boer H. Modelling of inter-rivet buckling of hybrid composites. Compos Struct 2006;73(2):221–8. [51] Hashagen F, de Borst R. Numerical assessment of delamination in fibre metal laminates. Comput Methods Appl Mech Eng 2000;185(2–4):141–59. [52] Camanho PP, Maimí P, Dávila CG. Prediction of size effects in notched laminates using continuum damage mechanics. Compos Sci Technol 2007;67(13):2715–27. [53] Moriniere FD, Alderliesten RC, Sadighi M, Benedictus R. An integrated study on the low-velocity impact response of the GLARE fibre–metal laminate. Compos Struct 2013;100:89–103. [54] Maimí P, Camanho PP, Mayugo JA, Dávila CG. A continuum damage model for composite laminates: Part II – Computational implementation and validation. Mech Mater 2007;39(10):909–19. [55] Tsai JL, Wang H. Modeling nonlinear rate dependent behaviors of composite laminates. J Chin Inst Eng 2007;30(1):141–8. [56] Chai H. Interlaminar shear fracture of laminated composites. Int J Fract 1990;43(2):117–31. [57] Lafarie-Frenot MC, Touchard F. Comparative in-plane shear behaviour of longcarbon-fibre composites with thermoset or thermoplastic matrix. Compos Sci Technol 1994;52(3):417–25. [58] Sun CT, Yoon KJ. Elastic–plastic analysis of AS4/PEEK composite laminate using a one-parameter plasticity model. J Compos Mater 1992;26(2):293–308. [59] Carpinteri A, Cornetti P, Barpi F, Valente S. Cohesive crack model description of ductile to brittle size-scale transition: dimensional analysis vs. renormalization group theory. Eng Fract Mech 2003;70(14):1809–39.