Progressive failure mechanism of laminated composites under fatigue loading

Abstract

A cohesive zone model for delamination propagation in laminated composites under static and fatigue loading has been derived and validated with experimental data under different mode conditions. This study presents a new approach to quantify fatigue delamination degradation based on damage mechanics to evaluate the rate of fatigue damage ( $\partial D / \partial N$ ). The static damage evaluation and fatigue damage degradation are derived from damage surface concept. Both static and fatigue damage linked each other to establish fatigue crack growth formula in the laminated composites. A user-defined subroutine, UMAT, has been employed to develop and implement a damage model in ABAQUS. Two different specimens; a double cantilever beam and a single lap joint were used to investigate the effectiveness of the new method. The simulation results revealed that the developed model had good agreement with experimental data available in literature.

Keywords

Complex failure mechanisms crack damage properties Paris law energy dissipation

Introduction

Fatigue-driven delamination in laminated composites is one of the general failure modes in composite structure e.g. airplanes, automobiles and wind turbines, all of which undergo cyclic loading.¹ Fatigue is also a major contributor to failures in rotating machinery and has been addressed in many research work.^2,³ Life estimation of composite materials is considered a serious challenge for engineers because of the uncertainties associated with damage initiation sites and propagation direction in these composites. There is a real and urgent need for reliable analytical methods to reduce the time and cost of experimental testing, and to enable optimally designed structures.⁴

In composite structures, it is important to detect the damage at an early stage of failure and to know how the damage will grow during the service life of structures. A relationship between fatigue crack growth and stress intensity factor is generally written as a power law (i.e. the Paris law) for metals and composite materials.⁵

Different approaches have been employed to analyse fatigue delamination growth, such as linear elastic fracture mechanics e.g. the virtual crack closure technique (VCCT), and cohesive zone models. Although the VCCT technique is widely used to simulate fatigue crack propagation, it has some limitations regarding crack path and re-meshing requirements. An alternative to the VCCT to simulate interface crack propagation, such as delamination in composite materials, is the cohesive zone model.

Over the last twenty years many researchers have used finite element methods, including the cohesive zone model, rather than attempt analytical solutions to simulate crack growth

Maiti and Geubelle⁶ introduced the instantaneous interface stiffness degradation law under cyclical fatigue loading. The rate of change of cohesive stiffness was defined as a function of the number of cycles. They proposed a power law for the cohesive model that can be expressed in terms of discretized time steps. While the instantaneous cohesive stiffness was assumed by Serebrinsky and Ortiz⁷ to degrade exponentially with each unloading-reloading cycle.

Load envelope methods have been suggested by many researcher^8–14 to simulate fatigue damage growth based on experimental Paris law curves.

Gornet and Ijaz¹⁵ investigated high-cyclic elastic fatigue damage model for carbon fiber epoxy matrix laminates. A classical interface damage was used to predict the fatigue damage parameters for different mode-mixtures.

Nojavan et al.⁴ proposed a non-Paris law based fatigue cohesive zone model. The authors assumed a simple power law for fatigue delamination growth under pure Mode I and Mode II loading. Although their proposed model can predict mixed mode delamination growth from simple parameters, their research demonstrated that the suggested cohesive zone model can not be used for arbitrary damage onset and propagation. Therefore, a crack path must be pre-defined in order to arrange the elements along the crack path in the finite element model.

Fatigue loading

Composite structures can be subject to cyclic loads with different wave shapes such as sinusoidal, square, etc. The load fluctuates with a certain frequency (fn) with a constant amplitude and produces a periodical stress or strain. Figure 1 shows typical sinusoidal and square loading curves with relevant parameters i.e. average stress is ((σ_max+σ_min)/2), stress amplitude is (((σ_max−σ_min))/2), and R-ratio (σ_min/σ_max). High amplitudes of fatigue loading reduce the remaining life of the structure.¹⁶

Figure 1.

Sinusoidal and square stress functions.

In fracture mechanics, the crack growth rate depends on the stress intensity factor or the energy release rate, and propagates according to the Paris law. In the last few years many researchers have introduced algorithms into the cohesive zone model which include fatigue damage degradation. Belnoue et al.¹⁷ employed the Paris law to represent fatigue in an interface element based on change in strain energy release rate, ΔG. The crack growth rate, $\partial a / \partial N$ , was assumed to be a function of change in strain energy release rate, ΔG, in the crack tip within each cyclic loading and can be written as:

\partial a / \partial N = C (Δ G ⁄ Gc)^m

(1)

where G_C is the critical fracture energy, C and m are fitting parameters determined experimentally.

When the specimen is under mixed mode loading, i.e. Mode I and Mode II, the parameters C and m can be determined depending on linear or parabolic equations based on pure mode cases as in the following expressions¹²

Linear rule

C = \frac{G_{I}}{G_{I} + G_{II}} C_{I} + \frac{G_{II}}{G_{I} + G_{II}} C_{II}

(2)

m = \frac{G_{I}}{G_{I} + G_{II}} m_{I} + \frac{G_{II}}{G_{I} + G_{II}} m_{II}

(3)

Parabolic rule

\log C = \log C_{I} + (\frac{G_{II}}{G_{T}}) \log C_{m} + {(\frac{G_{II}}{G_{T}})}^{2} \log \frac{C_{II}}{C_{m} C_{I}}

(4)

m = m_{I} + m_{m} (\frac{G_{II}}{G_{I} + G_{II}}) + (m_{II} - m_{I} - m_{m}) {(\frac{G_{II}}{G_{I} + G_{II}})}^{2}

(5)

where the subscripts I and II indicates Mode I, Mode II, and mode ratio. m_m is a mixed mode parameter to be experimentally determined.

Fatigue delamination does not propagate if the strain energy rate does not reach the threshold value G_th. In this region no measurable fatigue crack growth occurs and it is called subcritical growth region, as shown in Figure 2 (Region I). The delamination grows linearly when the energy release rate during the cyclic loading exceeds the threshold value. This is known as the region of stable growth (Region II) when plotted on a log–log scale, as shown in Figure 2. The nonlinear or unstable growth (Region III) occurs if the strain energy release rate is near to the static fracture energy G_c, see Figure 2.¹³

Figure 2.

Typical fatigue delamination growth curve.¹³

The change in strain energy rate during fatigue loading can be calculated from the instantaneous maximum strain energy and the load ratio R as¹⁷

Δ G = (1 - R^{2}) G_{\max}

(6)

where R is the load ratio and is defined by the user input in the numerical model. G_max is calculated from the traction–separation curve of the cohesive zone model (CZM) for each fatigue cycle.

The interface element is easily simulated in a numerical model based on traction-separation responses, the fracture process zone ahead of a crack tip can be modelled using a CZM. Many different types of traction-separation relations have been investigated in the literature e.g. polynomial law,¹⁸ trapezoidal model,¹⁹ perfectly plastic relationship,²⁰ exponential law,²¹ and bilinear cohesive law.²² Among these the bilinear law gives good convergence of the finite element model and numerically is the simplest formulation.^13,23

It is worthwhile to mention that the static damage evolution should be included in modelling the delamination growth under fatigue conditions. The total damage parameter can be calculated during cyclic loading as:

D_{tot} = D + D_{f}

(7)

where D is static damage evolution and

D_{f}

is fatigue damage parameter.

The total damage $D_{tot}$ is used to update the cohesive traction at each element after each time-step. Complete failure of an element takes place once the total damage reaches unity and²⁴

τ = (1 - D) k δ static loading

(8a)

τ = (1 - D_{tot}) k δ fatigue loading

(8b)

where

τ

is cohesive traction, k is interface stiffness, and δ is relative displacement

Energy release rate at the crack tip can be calculated by integrating the traction-separation curve under cyclic loading, unless complete failure occurs. The fatigue loading causes stiffness degradation in the interface element, this leads to increase separation as total damage evolution accumulates. The highest strain energy release rate will be in those elements that lie in the process zone near to the numerical crack tip. The response of traction separation will follow a non-vertical path which matches the bi-linear law for static loading according to Turon et al.⁹ and Harper and Hallett.²⁴

Developed fatigue damage degradation

The quadratic nominal stress criterion to predict the onset of delamination under mixed mode loading is one of the most frequently adopted failure criteria.²⁵ This criterion has been successfully utilized by many researchers e.g. Hameed et al.,²⁶ Zou and Hameed,²⁷ Belnoue et al.¹⁷ and Kawashita and Hallett,¹² and is written as;

{F_{i} = (\frac{〈τ_{n}〉}{τ_{n c}})}^{2} + {(\frac{τ_{s}}{τ_{s c}})}^{2} + {(\frac{τ_{t}}{τ_{t c}})}^{2}

(9)

where

τ_{n}, τ_{s}, τ_{t}

are the normal, shear and tear stresses respectively and where the operator

〈x〉

is defined as;

〈x〉 = \{\begin{matrix} x, x > 0 \\ 0, x \leq 0 \end{matrix}

Damage growth under pure-mode loading is predicted using the energy release rate during loading and unloading, and fracture toughness:

F_{p} = \frac{G}{G_{c}}

(10)

The static damage evaluation is derived from damage surface concept and is written as:

Damage initiation criterion ( $F_{i}$ ) + Damage propagation criterion ( $F_{p}$ ) = 1

From this condition, the static damage evaluation can be expressed as

D = 1 - \sqrt{\frac{1 - F_{p}}{[\frac{k_{1}^{2} δ_{1}^{2}}{τ_{1 c}^{2}} + \frac{k_{2}^{2} δ_{2}^{2}}{τ_{2 c}^{2}} + \frac{k_{3}^{2} δ_{3}^{2}}{τ_{3 c}^{2}}]}}

(11)

The ratio between damaged ( $L_{d}$ ) and undamaged( $L_{u}$ ) parts was assumed by Turon et al.⁹ to equal the ratio of the dissipated energy over the representative interface part. Therefore, the damage growth criterion for pure mode can be written as:

F_{p} = \frac{G}{G_{c}} = \frac{L_{d}}{L_{u}}

(12)

Where $L_{u}$ is mesh size and $L_{d}$ is the length of damaged elements.

Substitute equation (12) into the static damage evaluation equation (11), to obtain the relation between damage evaluation and damaged length ratio;

D = 1 - \sqrt{\frac{1 - \frac{L_{d}}{L_{u}}}{[\frac{k_{1}^{2} δ_{1}^{2}}{τ_{1 c}^{2}} + \frac{k_{2}^{2} δ_{2}^{2}}{τ_{2 c}^{2}} + \frac{k_{3}^{2} δ_{3}^{2}}{τ_{3 c}^{2}}]}}

(13)

By differentiating the above equation with respect to the damaged length, we obtain a novel term $\frac{\partial D}{{\partial L}_{d}}$ , not found elsewhere in the literature and which can be written as:

\frac{\partial D}{{\partial L}_{d}} = \frac{1}{2} * \frac{1}{L_{u}} * \frac{1}{\sqrt{\frac{1 - \frac{L_{d}}{L_{u}}}{[\frac{k_{1}^{2} δ_{1}^{2}}{τ_{1 c}^{2}} + \frac{k_{2}^{2} δ_{2}^{2}}{τ_{2 c}^{2}} + \frac{k_{3}^{2} δ_{3}^{2}}{τ_{3 c}^{2}}]}}} * \frac{1}{[\frac{k_{1}^{2} δ_{1}^{2}}{τ_{1 c}^{2}} + \frac{k_{2}^{2} δ_{2}^{2}}{τ_{2 c}^{2}} + \frac{k_{3}^{2} δ_{3}^{2}}{τ_{3 c}^{2}}]}

(14)

The evaluation of the rate of fatigue damage for each cyclic loading ( $\partial D / \partial N$ ) on the traction–separation curve is a function of the novel term ( $\partial D / {\partial L}_{d}$ ) and the damage growth rate ( ${\partial L}_{d} / \partial N$ ) over the process zone. The expression for the required rate can be written as:

\frac{\partial D}{\partial N} = \frac{\partial D}{{\partial L}_{d}} * \frac{{\partial L}_{d}}{\partial N}

(15)

So,

\frac{\partial D}{\partial N} = \frac{1}{2} * \frac{1}{L_{u}} * \frac{1}{\sqrt{\frac{1 - \frac{L_{d}}{L_{u}}}{[\frac{k_{1}^{2} δ_{1}^{2}}{τ_{1 c}^{2}} + \frac{k_{2}^{2} δ_{2}^{2}}{τ_{2 c}^{2}} + \frac{k_{3}^{2} δ_{3}^{2}}{τ_{3 c}^{2}}]}}} * \frac{1}{[\frac{k_{1}^{2} δ_{1}^{2}}{τ_{1 c}^{2}} + \frac{k_{2}^{2} δ_{2}^{2}}{τ_{2 c}^{2}} + \frac{k_{3}^{2} δ_{3}^{2}}{τ_{3 c}^{2}}]} * \frac{{\partial L}_{d}}{\partial N}

(16)

The damage accumulates ahead of the crack tip in the fracture process zone to give additional delamination length. This new length of delamination is calculated based on the damage fatigue crack lengths, $L_{d}$ , of all elements which lie in the cohesive zone area. It is written as a summation function in the direction of delamination growth as:

a = \sum_{e \in L_{C Z}} L_{d}^{e}

(17a)

where

L_{C Z}

is the length of the cohesive zone and calculated as⁹

L_{C Z} = \frac{9 π}{32} \frac{G_{c}}{{(τ_{c})}^{2}} E

(17b)

E is the Young’s modulus of the material.

The equation (17a) is derived with respect to the number of cycles to obtain the crack growth rate, $\partial a / \partial N$ , as function of ( ${\partial L}_{d} / \partial N)$ .

\frac{\partial a}{\partial N} = \sum_{e \in L_{C Z}} \frac{{\partial L}_{d}^{e}}{\partial N}

(18)

The incremental rate of delamination is assumed to be equal for all elements within the cohesive zone area. Therefore, the summation sign can be rewritten in terms of the number of elements ( ${E N}_{C Z})$ in the process zone, $L_{C Z}$ .²⁸ The fatigue delamination growth rate can then be written as:

\frac{\partial a}{\partial N} = {E N}_{C Z} * \frac{{\partial L}_{d}}{\partial N}

(19)

where

{E N}_{C Z} = length o f cohesive zone (L_{C Z}) / element length (L_{u}) .

By substituting equation (19) into equation (16), the rate of fatigue damage ( $\partial D / \partial N$ ) is related to Paris’ law as:

\frac{\partial D}{\partial N} = \frac{1}{2} * \frac{1}{L_{u}} * \frac{1}{\sqrt{\frac{1 - \frac{L_{d}}{L_{u}}}{[\frac{k_{1}^{2} δ_{1}^{2}}{τ_{1 c}^{2}} + \frac{k_{2}^{2} δ_{2}^{2}}{τ_{2 c}^{2}} + \frac{k_{3}^{2} δ_{3}^{2}}{τ_{3 c}^{2}}]}}} * \frac{1}{[\frac{k_{1}^{2} δ_{1}^{2}}{τ_{1 c}^{2}} + \frac{k_{2}^{2} δ_{2}^{2}}{τ_{2 c}^{2}} + \frac{k_{3}^{2} δ_{3}^{2}}{τ_{3 c}^{2}}]} * \frac{1}{{E N}_{C Z}} * C {(\frac{Δ G}{G_{c}})}^{m}

(20)

In the UMAT subroutine, the fatigue damage degradation is updated in the model for each time-step as;

D_{f}^{new} = D_{f}^{old} + Δ N \frac{\partial D}{\partial N}

(21)

$Δ N$ is calculated as follows:

Δ N = frequency * step time

(22)

Double cantilever beam (DCB) model

Double cantilever beam (DCB) specimens as reported by Landry and LaPlante¹³ were adopted for this numerical simulation. The length and width of the specimen was 150 mm and 20 mm respectively with pre-crack length of 50 mm. The specimens have two unidirectional 16 ply-thick (3 mm thick) layers made from carbon/epoxy with material properties shown in Table 1. All simulations in this study were subjected to Mode I static and fatigue loading, and carried out using commercial finite element software ABAQUS. A two dimensional finite element model of DCB was implemented to reduce computational effort and save computer time. The layers were meshed using 4-node plane strain elements, whereas the cohesive elements were meshed with element type COH2D4 and placed between layers where the crack was expected to propagate.

Table 1.

Material properties of the carbon/epoxy layers used.¹³

E₁₁ (GPa)	E₂₂ (GPa)	E₃₃ (GPa)	$v_{12}$	$G_{12}$ (GPa)	$G_{I C}$ (N/m)
155	10.5	10.5	0.32	4.83	422

First, the static load was numerically applied on the DCB by specifying an opening displacement of each cantilever arm. The numerical results for energy as a function of displacement have been validated with experimental data reported by Landry and LaPlante.¹³ The response of strain energy with each initial displacement is in good agreement with experimental results, as shown in Figure 3.

Figure 3.

Comparison between developed model and experimental results of fracture toughness.

The experimental procedure was performed by Landry and LaPlante¹³ as follows; the specimens were cracked approximately 5 mm and then subsequently subjected to the quasi-static load. When opening displacement reached to 9.5 mm, the static load was stopped and then, sinusoidal displacement cycles with an amplitude of 8 mm were performed.

The novel approach to fatigue damage degradation was implement in ABAQUS via the UMAT subroutine. It was necessary to define the cohesive zone length in the model which was estimated based on traction values, fracture toughness and Young’s modulus of the cohesive material. Fatigue delamination growth is determined by relating the damage accumulation of cohesive elements to the Paris growth law, equation (20). The damage accumulates with each cyclic loading based on equation (21) to produce a new crack length. Figure 4 illustrates damage accumulation for the DCB under different numbers of cycles. It is clear that the delamination length increases when the number of cycles increases from 250 to 3500. The difference between the experimental results and developed model is determined as;

difference = \frac{experimental value - simulation value}{experimental value} %

(23)

Figure 4.

Delamination growth in DCB under different cyclic loading.

Figure 5.

Comparison between developed model predictions and experimental results for delamination growth.

For example, the difference at 250 cycle is calculated as;

difference = \frac{61 - 58}{61} % ≅ 5 %

In the same way, the difference at 3500 cycles can be determined and it is about 3% as shown in Figure 5.

The fatigue crack growth rate using developed can be represented by Paris plot as shown in Figure 6. The simulation results of mixed mode were validated against experimental data available in literature.²⁹

Figure 6.

Comparison of simulated fatigue crack growth rates and experimental data available in literature.²⁹

Single lap joint (SLJ) model

Adhesively bonded joints have many advantages e.g. lower weight, and better fatigue performance compared with classical mechanically fastened joints. Therefore, they are used in many structural applications in industries such as automotive and aerospace. The increase in demand for adhesively bonded joints has been accompanied by investigating and developing numerical models to predict the response of the structure, and to describe the behaviour of damage in the adhesive joints.³⁰ The cohesive zone model is a very effective way to model adhesively bonded joints using finite element analysis

Khoramishad et al.³¹ carried out static and fatigue loading of single lap joints (SLJ). The specimen consisted of two aluminium 2024-T3 layers bonded with FM 73 M OST toughened epoxy film adhesive. The substrates were pre-treated prior to bonding. This pre-treatment consisted of a chromic acid etch (CAE) and phosphoric acid anodising (PAA) followed by the application of BR 127 corrosion inhibiting primer, to maximise environmental resistance and bonding durability. The joints were cured at 120 °C, under 0.28 MPa pressure for 60 min. The dimensional details of the SLJ are shown in Figure 7. The overlap length, the width and the thickness of the bond line were 30, 12.5 and 0.2 mm, respectively.

Figure 7.

Single lap joint geometries.

The static tests were carried out by Khoramishad et al.³¹ on the single lap joints (SLJ) to find cohesive zone properties e.g. fracture energy. The displacement control was applied to the specimen at a rate of 0.1 mm/min. Average static strengths of 10.0 kN was observed for the SLJ. After that, the single lap joints (SLJ) specimen was tested under fatigue load condition at 5 Hz with a load ratio of 0.1.

The interface element length used over the cohesive zone depends on mesh size in the finite element model. It effects the sensitivity of the fatigue response according to Naghipour et al.³² Vassilopoulos,⁵ in previous studies investigated the effect of mesh element size on convergence of the numerical model and found that using a mesh size less than 1 mm for the interface element provided better solution convergence. Therefore, in the present simulation the element size adopted for the adhesive was 0.2 mm.

Fatigue modelling was performed at maximum fatigue load of 40% of the static final failure load which was measured experimentally at about 10 kN. When the specimen was loaded statically at 4 kN, as shown Figure 8, it can be seen that, initially, the damage parameter is less than unity for each of the three elements. When the fatigue loading was activated, the first and second elements needed 200 cycles for the rate of increase of damage to begin to accelerate rapidly whereas in the third element the rate of increase of damage began to accelerate after about 500 cycles. This was because the value of the damage parameter for the third element at the end of the static stage was less than for the first and second elements.

Figure 8.

Damage propagation of the first, second and third element.

The delamination growth was simulated based on the model developed here. It should be observed that the damage grew symmetrically on both side of the adhesively bonded joint as shown in Figure 9. The variation of the damage parameter under different numbers of cycles, as predicted by the present model is also illustrated in Figure 9. The elements having state variable (SDV3) equal to unity in the damage contours represent complete failure of the cohesive zone, signifying no fatigue resistance. Based on model proposed here the crack needed 2000 cycles to grow to 0.30 mm and propagated to 0.45 mm when the number of cycles reached 4000. The measured delamination lengths by Khoramishad et al.³¹ were 0.285 mm and 0.43 mm for 2000 and 4000 cycles respectively. Thus, the difference between numerical model and test was 5.2% at 2000 cycles and 4.6% at 4000 cycles.

Figure 9.

Delamination growth under different cyclic loading (a) developed model (b) experimental data by Khoramishad et al.³¹

Conclusions

The objective of this paper was to present a novel fatigue delamination accumulation model. The novel approach was tested numerically under two different load conditions; Mode I loading which was performed on a DCB specimen, and Mode II loading which was carried out on an adhesively bonded joint. In addition, the simulation was performed of the mixed-mode model and compared with experimental data available in the literature²⁹ as shown in Figure 6. The simulation findings were compared against experimental results available in literature, the maximum difference between them was 5.0% for the DCB specimen and 5.2% for the SLJ specimen. This is taken as confirmation that the proposed model converges with experimental data when tested using numerical simulation.

Footnotes

Declaration of Conflicting Interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) received no financial support for the research, authorship, and/or publication of this article.

ORCID iD

Ghalib R Ibrahim

References

Bak

BLV

Turon

Lindgaard

, et al. A simulation method for high-cycle fatigue-driven delamination using a cohesive zone model. Int J Numer Methods Eng 2016; 106: 163–191.

Ibrahim

Albarbar

Comparison between Wigner-Ville distribution- and empirical mode decomposition vibration-based techniques for helical gearbox monitoring. Proc Inst Mech Eng Part C J Mech Eng Sci 2011; 225: 1833–1846.

Ibrahim

Albarbar

Abouhnik

, et al. Adaptive filtering based system for extracting gearbox condition feature from the measured vibrations. Measurement 2013; 46: 2029–2034.

Nojavan

Schesser

Yang

QD.

A two-dimensional in situ fatigue cohesive zone model for crack propagation in composites under cyclic loading. Int J Fatigue 2016; 82: 449–461.

Vassilopoulos A P (ed). Fatigue and fracture of adhesively-bonded composite joints (No. 52). 2014. UK: Elsevier.

Maiti

Geubelle

PH.

A cohesive model for fatigue failure of polymers. Eng Fract Mech 2005; 72: 691–708.

Serebrinsky

Ortiz

A hysteretic cohesive-law model of fatigue-crack nucleation. Scr Mater 2005; 53: 1193–1196.

Robinson

Galvanetto

Tumino

, et al. Numerical simulation of fatigue-driven delamination using interface elements. Int J Numer Methods Eng 2005; 63: 1824–1848.

Turon

Costa

Camanho

, et al. Simulation of delamination in composites under high-cycle fatigue. Compos Part A Appl Sci Manuf 2007; 38: 2270–2282.

10.

Tumino

Cappello

Simulation of fatigue delamination growth in composites with different mode mixtures. J Compos Mater 2007; 41: 2415–2441.

11.

Moroni

Pirondi

A procedure for the simulation of fatigue crack growth in adhesively bonded joints based on the cohesive zone model and different mixed-mode propagation criteria. Eng Fract Mech 2011; 78: 1808–1816.

12.

Kawashita

Hallett

SR.

A crack tip tracking algorithm for cohesive interface element analysis of fatigue delamination propagation in composite materials. Int J Solids Struct 2012; 49: 2898–2913.

13.

Landry

LaPlante

Modeling delamination growth in composites under fatigue loadings of varying amplitudes. Compos Part B Eng 2012; 43: 533–541.

14.

Bak

BLV

Turon

Lindgaard

, et al. A benchmark study of simulation methods for high-cycle fatigue-driven delamination based on cohesive zone models. Compos Struct 2017; 164: 198–206.

15.

Gornet

Ijaz

A high-cyclic elastic fatigue damage model for carbon fibre epoxy matrix laminates with different mode mixtures. Compos Part B Eng 2011; 42: 1173–1180.

16.

Kiefer

2014. Simulation of high-cycle fatigue-driven delamination in composites using a cohesive zone model. London: Imperial College London.

17.

Belnoue

JPH

Giannis

Dawson

, et al. Cohesive/adhesive failure interaction in ductile adhesive joints part II: quasi-static and fatigue analysis of double lap-joint specimens subjected to through-thickness compressive loading. Int J Adhes Adhes 2016; 68: 369–378.

18.

Needleman

A continuum model for void nucleation by inclusion debonding. J Appl Mech Trans ASME 1987; 54: 525–531.

19.

Tvergaard

Hutchinson

JW.

The relation between Crack-Growth resistance and fracture process parameters in elastic plastic solids. J Mech Phys Solids 1992; 40: 1377–1397.

20.

Cui

Wisnom

MR.

A combined stress-based and fracture-mechanics-based model for predicting delamination in composites. Composites 1993; 24: 467–474.

21.

Xu X P and Needleman A. Numerical simulations of fast crack growth in brittle solids. Journal of the Mechanics and Physics of Solids 1994; 42: 1397–1434.

22.

Reedy

Mello

Guess

TR.

Modeling the initiation and growth of delaminations in composite structures. J Compos Mater 1997; 31: 812–831.

23.

Alfano

On the influence of the shape of the interface law on the application of cohesive-zone models. Compos Sci Technol 2006; 66: 723–730.

24.

Harper

Hallett

SR.

A fatigue degradation law for cohesive interface elements – eevelopment and application to composite materials. Int J Fatigue 2010; 32: 1774–1787.

25.

Ibrahim

Albarbar

A new approach to the cohesive zone model that includes thermal effects. Compos Part B Eng 2019; 167: 370–376.

26.

Hameed

Ibrahim

Albarbar

Effect of friction and shear strength enhancement on delamination prediction. J Compos Mater. Epub ahead of print 26 March 2020.

27.

Zou

Hameed

Combining interface damage and friction in cohesive interface models using an energy based approach. Compos Part A Appl Sci Manuf 2018; 112: 290–298.

28.

Skvortsov

Chernyakin

Glushkov

, et al. Simulation of fatigue delamination growth in composite laminates under mode I loading. Appl Math Modell 2016; 40: 7216–7224.

29.

Kenane

Benzeggagh

ML.

Mixed-mode delamination fracture toughness of unidirectional glass/epoxy composites under fatigue loading. Compos Sci Technol 1997; 57: 597–605.

30.

de Moura

MFSF

Goncalves

JPM.

Cohesive zone model for high-cycle fatigue of adhesively bonded joints under mode I loading. Int J Solids Struct 2014; 51: 1123–1131.

31.

Khoramishad

Crocombe

Katnam

, et al. Predicting fatigue damage in adhesively bonded joints using a cohesive zone model. Int J Fatigue 2010; 32: 1146–1158.

32.

Naghipour

Bartsch

Chernova

, et al. Effect of fiber angle orientation and stacking sequence on mixed mode fracture toughness of carbon fiber reinforced plastics: numerical and experimental investigations. Mater Sci Eng A Struct Mater Prop Microstruct Process 2010; 527: 509–517.