Multiscale analysis of fatigue crack initiation life for unidirectional composite laminates

Composite Structures 213 (2019) 271–283

Contents lists available at ScienceDirect

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

Multiscale analysis of fatigue crack initiation life for unidirectional composite laminates

T

⁎

M. Naderia, , J. Michopoulosb, N. Iyyera, K. Goelc, N. Phanc a

Technical Data Analysis, Inc. (TDA), 3190 Fairview Park Drive, Suite 650, Falls Church, VA 22042, USA Computational Multiphysics Systems Laboratory, US Naval Research Laboratory, Washington, DC 20375, USA c US Naval Air Systems Command, Patuxent River, MD 20670, USA b

A R T I C LE I N FO

A B S T R A C T

Keywords: Carbon fiber reinforced composite Multicontinuum theory Micromechanics Kinetic theory of fracture Finite element modeling

Three dimensional (3D) finite element analysis along with a damage model based on kinetic theory of fracture (KTF) are used to determine the onset of fatigue crack initiation in unidirectional carbon fiber reinforced composites. In addition, the Multicontinuum theory (MCT) and micromechanics are employed separately to enable the assessment of the stress state distribution at the length scale level of the constituents. While the MCT is implemented in the macroscale model formulation, micromechanical analysis of a representative unit cell has been also developed to capture the stress field in the micro length scale level. An ABAQUS™ user material (UMAT) subroutine has been developed and used for implementing the formulation. Numerical examples of various tests including tension, three point bending, four point bending, and short beam shear test of a unidirectional (UD) laminate are considered and the results are compared against actual test data. A detailed discussion, challenges and future plans on the applicability of the KTF are concluding this work.

1. Introduction The proliferation of composite materials has been growing intensively due to their superior weight to strength ratio, structural fatigue, and corrosion resistance properties. A major challenge in designing composite material structures against fatigue loading is the estimation of the onset of microcracks, which are responsible for a significant part of the fatigue life [1–5]. Consequently, a design framework and the associated methodology that can answer the question “when microcracks are initiated” is extremely crucial when the designer aims to predict the useful life the composite laminate especially in the high cycle fatigue region. Matrix microcracks are usually initiated from voids, defects, misalignment, or areas of stress concentrators at the early stage of fatigue life followed by the saturation stage. These matrix microcracks lead to degraded stiffness and strength that eventually result to a drop-off of the load carrying capacity. The load carrying capacity response and process of failure in composite laminates under fatigue loading conditions are both tightly connected to capturing the microcracks initiation onsets and associated locations. Therefore, capturing accurately initiation of microcracks is of paramount importance for robust design and modeling. Yet, despite considerable efforts related to progressive fatigue damage analysis of composites in the literature [1–19], to the best of

⁎

authors’ knowledge, few of these investigations have considered microcracks initiation based life predictions for composite laminates [1–5,9–11,13,14]. The majority of progressive fatigue damage models are focused on the propagation of cracks assuming either cracks already exist or manufacturing process defects/voids are always present [6–8,12,15]. The consequence of this assumption is the neglect of the possibility for crack initiation and its effect on life that can be significant [13]. To take into account fatigue crack initiation for life time prediction various approaches have been developed. These approaches have been mostly based on phenomenological models and rely on extensive test data. For example, for adhesively bonded joints, Quaresimin and Ricotta proposed a crack initiation law based on stress intensity factor combined with a virtual crack closure technique (VCCT)-based on a crack propagation method [13,14]. May and Hellett developed a combined crack initiation and propagation model using cohesive interface elements. Their damage initiation model is follows the classical Wohler (S-N) curve and the propagation approach is based on Paris law [1,10]. Nojavan et al. proposed a non Paris-law-based cohesive zone model for fatigue delamination. In this model the fatigue driven damage growth rate is directly calculated from calibrated parameters [11]. The above mentioned investigations rely mostly on S-N data. The drawbacks of using S-N curves are that a) a designer must recalibrate

Corresponding author. E-mail address: [email protected] (M. Naderi).

https://doi.org/10.1016/j.compstruct.2019.01.107 Received 6 November 2018; Received in revised form 8 January 2019; Accepted 29 January 2019 Available online 31 January 2019 0263-8223/ © 2019 Elsevier Ltd. All rights reserved.

Composite Structures 213 (2019) 271–283

M. Naderi et al.

the log-log fit each time the applied loading condition is altered; b) extensive experimental tests are required before designers attempt to model fatigue life. Jones et al. performed delamination analysis in polymer-matrix fiber composites under cyclic-fatigue loading in operational aircraft structures. Using fracture mechanics, they presented a methodology based on a variant of the Hartman-Schijve equation, to determine upper-bound fatigue crack growth rate curves. They mentioned that initiation value for energy release rate is strongly dependent on test protocol [20]. Recently, Sevenois et al. [15] developed a microscale based criterion to predict fatigue crack initiation for UD composite laminates subjected to multiaxial fatigue load. This damage initiation model is based on a power law equation along with three assumptions including linear elastic behavior between fiber and matrix, linear damage, and log-log type S-N curve data. However, this model requires a large amount of test data and had difficulty predicting tension-compression and compression-compression loading [5]. Most of the aforementioned models are macromechanical and they employ macroscale stress state and empirical parameters calibration for estimating damage and degradation while at the same time they ignore microscale effects and physic-based considerations. Alternatively, micromechanical-based models can elucidate underlying mechanisms contributing to fatigue failure and damage initiation at the constituent length scale level [21]. Various noteworthy studies have considered micromechanical modeling and simulation of fatigue failure for composite laminates [22–36]. Among them, Reifsnider and Gao have developed a fatigue failure criterion based at the micromechanical level to address failure of all constituent parts of composite laminate [22]. Akshantala and Talreja proposed a micromechanics based fatigue damage progression along with a semi-empirical fatigue failure criterion to predict the life of cross ply laminates under tensile loading [25]. Adibnazari et al. used a bridging micro-mechanical model to study stiffness and strength degradation of unidirectional laminates subjected to cyclic loading [28]. Krause developed a fatigue failure criterion to study the transverse crack initiation and evolution at the micromechanical level for fiber reinforced polymer under fatigue loading condition [29]. He used the hypothesis that fatigue response of composites is dominated by matrix. Using experimental tests on pure polymeric matrix, he derived a fatigue failure criterion for polymers under multiaxial loading conditions. The work is considering hysteresis energy as damage initiation quantity and is physically motivated. Fertig et al. [34–36] employed KTF (kinetic theory of fracture) as a physics based model to estimate matrix damage parameter and predict fatigue life of on-/off-axis laminates at both room and elevated temperature. Similar to Fertig [34–37], Hosseini Kordkheili et al. [33], Kapitzic et al. [38,39] and Dalgarno et al. [37] employed kinetic theory of fracture to calculate damage initiation in different composite laminates using minimal material data. Recently, some efforts have been invested in considering the presence of microcracking and sub-microcracking not explicitly but rather implicitly. This is accomplished in terms of homogenized equivalent scalar fields representing material properties degradation to capture damage as it has been demonstrated by Michopoulos et al. [40–42]. The above mentioned survey reveals that over the past decades, researchers have attempted to develop a physics-based fatigue models incorporating micromechanics in a manner that minimizes the number of required tests and associated experimental data every time the test conditions and load parameters are altered. While crack propagation at the macroscale has received significant attention, special attention needs to be focused on sub-microcrack and microcrack initiation as well as the relevant approaches that are considering a physics-based onset of sub-microcracks and microcracks. Yet, although the kinetic theory of fracture (KTF) is a physics-based failure criterion, it has not received much attention in the literature. The main advantage of the method is that it is physics-based and requires a minimal amount of material characterization data. Consequently, the focus of present work is to reexamine the potential of

kinetic theory of fracture (KTF) as a damage initiation model and identify its challenges and potential in the context of a multiscale framework approach. Hence, while using KTF, we plan to address microcracks initiation and answer questions such as “when is crack initiation onset?”, “what is the definition of a sub-microcrack and how is it connected to microcrack?” and “where and when to insert a sub-microcrack or microcrack in the model?”. Our plan is based utilizing a micromechanical approach for determining stress states at the constituent level as the input data for KTF-based model. KTF and micromechanics are implemented in finite element analysis (FEA) code to capture the onset of microcracks initiation during fatigue of unidirectional (UD) laminates. The underlying assumption is that microcracks initially occur in the matrix and the fiber/matrix interface. The body of the paper is organized as follows: First, the KTF and the relevant micromechanics are briefly-described. Then, a multiscale computational framework and the associated procedure are demonstrated in higher detail. Thirdly, numerical results compared against test data for a few specific testing conditions are followed by detailed discussion on the use of KTF and its challenges. The paper is finally closes with the overall conclusions and challenges for the future.

2. Background overview 2.1. Entropy/kinetic theory of fracture The kinetic theory of fracture developed by Coleman [43] and Zhurkov [44] states that the crack creation rate is a thermally activated process during cyclic loading. Damage initiates from small scale of a particular material constituent (atomic or molecular scale) and propagates to the macro/structural scale resulting final material rupture. In the current framework, we define the transition from the atomic or molecular length scale level to that of the microscale as the damage initiation level. Considering a unit cell as the smallest volume element representing the material ensemble of interest, the first sub-microcracks are anticipated to initiate from the matrix or the fiber matrix/interface rupture of atomic bonds or molecular chain. In order to associate the process of sub-microcracking with the concept of entropy one can assume the damage process as moving from an uncracked intact submicrostate to a cracked sub-microstate. It is assumed that the sub-microstates evolve at sub-microscale, which is a length scale that is some orders of magnitude smaller than fiber diameter (i.e. ∼0.1 of fiber diameter) where actual collections of molecular chains exist. Specifically, the term “sub-microscale” is referred to a length scale at a level greater than the molecular length scale and lower than the micro length scale. Once a group of broken molecular chains is saturated, and coalesces, sub-microcracks are formed in the unit cell. According to thermodynamics, entropy (S) can be expressed according to the following expression [45]:

S=−

∫ ρlnρdU ,

(1)

where, ρ is the probability density of sub-microstates transition from an uncracked sub-microstate to a cracked sub-microstate and U is the energy barrier which must be surpassed to initiate atomic bond rupture. Since thermal fluctuation energy of atoms is always stated in the distribution function form, the probability density of atoms or molecules to overcome the energy barrier U, is expressed as an exponential function as,

ρ=

1 U exp ⎛− ⎞ kT ⎝ kT ⎠

(2)

Hansen and Baker showed that the probability of microstates can be related to the characteristic frequency of atoms’ oscillation (ω0 ) and that the reaction rate (breakage or reformation) [45] can be expressed by 272

Composite Structures 213 (2019) 271–283

M. Naderi et al.

U K = ω0 exp ⎛− ⎞ ⎝ kT ⎠

constant [34,35]. It should be noted that this correction may not be enough as one compares two tests at different frequencies in which temperature rise is not significant in terms of changing polymer properties. Hence, one needs to decompose the cyclic load to take into account the effect of both stress range and mean stress. In what follows we describe how such an additive decomposition can be implemented.

(3)

where ω0 is oscillation frequency of atoms that is roughly kT/h with h as Plank’s constant (6.62606 × 10−34 J.s). This reaction rate K can therefore be written as follows [45]:

kT U exp ⎛− ⎞ h ⎝ kT ⎠

K=

(4)

2.3. Fatigue damage initiation

In the presence of applied stress, reaction rate can be intensified (called the reaction rate of bond breakage) or slowed down (called bond reformation rate) as expressed below for bond breakage rate (Kb)

Kb =

U − γσ kT ⎞, exp ⎛− h kT ⎠ ⎝

A critical question for fatigue of composites especially for high cycle fatigue is what is the proper definition of first initial crack/damage and what the length of the initial crack or damage may be. From a physical perspective, the onset of breakage of a group of molecular chains or/ and atomic bonds can be assumed to occur at the time when the first damage occurs. From an experimental perspective, the initiation may be considered as smallest detectable size of the crack/damage that can be identified via some experimental method. Depending on the measurement equipment, detectable crack size can be somewhere between a few tenth of a mm to a few mm [13]. However, from a multiscale analysis point of view, it depends on the length scales that are involved in the analysis. In the present work, the damage initiation size is in the order of magnitude of unit cell, which is some order of magnitude greater than the fiber diameter. Fig. 1 shows a schematic depicting cyclic loading for tension-tension (T-T), tension-compression (T-C), compression-tension (C-T), and compression-compression (C-C) loading. For simplicity, an in-phase saw-tooth loading is considered for analysis purposes. “R” is the load ratio, σmax, σmin, σm, σr are maximum stress, minimum stress, mean stress and stress range respectively. The total fatigue damage can be expressed as the sum of static damage Ds and fatigue damage Df. Static damage accounts for the damage caused by ramping the load to the mean stress, while fatigue damage is the damage produced by stress range in each cycle as follows;

(5)

where U is activation energy and γ is the activation volume which is approximately the volume of material participating in the rupture process. Also, the bond reformation rate (Kr) can be expressed as

Kr =

U + γσ kT ⎞. exp ⎛− h kT ⎠ ⎝

(6)

The thermal oscillation energy of atoms (kT) is competing against energy barrier (U-γα) to reach to a new equilibrium state. Once the thermal fluctuation energy exceeds the energy barrier, a large collection of atomic bonds and molecular chains fall apart. The coalescence of broken bonds in sub-microstate results in sub-microcrack. 2.2. Damage model In order to measure damage evolution at sub-microscale, a degradation law must be established as a function of bond breakage rate K. Hansen and Baker-Jarvis related the rate of sub-microcracks in matrix dN/dt to the concentration of weak regions NT and rupture rate of Kb according to [45]

D = Ds + Df .

dN = (NT − N ) Kb dt

(7)

Assuming zero initial damage, the static damage can be defined as

where N(t) is submicrocrack concentration with N(0) = 0 and NT represents the amount of weak regions at time zero. By normalizing Eq. (7) with respect to Nr that represents the amount of submicrocracks at rupture, it can be obtained that

dD = (DT − D) Kb; D (0) = 0 dt

kT U ⎡ Ds (t ) = DT ⎢1 − exp ⎜⎛− exp ⎛− ⎞ ⎝ kT ⎠ ⎝ h ⎣

1

(8)

dD =1 DT − D

∫0

t

γσ (τ ) ⎞ ⎞ ⎤ dτ , exp ⎛ ⎝ kT ⎠ ⎠ ⎥ ⎦

γσramp (τ ) ⎞

exp ⎛ ⎝ ⎜

σm τ t0

⎟

kT

⎠

⎤ dτ ⎟⎞ ⎥. ⎠⎦

(12)

with t0 as the duration of

(

)

(13)

The mean stress effect is not directly included in the damage evolution Eq. (10). One needs to distinguish damage caused by the stress range and mean stress. Therefore, in order to separate the effect of mean stress and stress range, Df can be decomposed to Df1 and Df2. Df1 is produced by an equivalent constant load of σ(τ) = σmax (for C-C) or σ(τ) = σmin (for T-T) during fatigue load, while Df2 is the damage produced by stress range σr, resulting to:

kT U Df 1 (t ) = DTf − (DTf − Ds ) exp ⎛− exp ⎛− ⎞ ⎝ kT ⎠ ⎝ h

Using Eq. (5) in the integration of Eq. (8), one can derive the following expression for damage evolution, ⎜

t

U ⎡ γσm γσm ⎛ kT ⎤⎞⎤ ⎛ ⎞ Ds (t ) = DT ⎡ ⎢1 − exp ⎜− h exp ⎝− kT ⎠ ⎢ kTt0 e kT − 1 ⎥ ⎟ ⎥. ⎦⎠⎦ ⎣ ⎝ ⎣

(9)

kT U D (t ) = DT ⎡1 − exp ⎛− exp ⎛− ⎞ ⎢ ⎝ kT ⎠ ⎝ h ⎣

∫0

For a ramped static load with σramp = ramp load application, we can obtain

where D = N/Nr, and DT = NT/Nr. The variable D can be considered as a damage measure parameter. It is assumed that initially there is no damage and at the time of sub-microcracks formation the limit of D approaches to unity. The parameter DT is obtained by enforcing the following condition [45] as supported with experimental observations based on the measurement of free radicals formed at materials rupture by Zhurkov and Kuksenko [46].

∫0

(11)

⎜

∫t

t

γσ (τ ) ⎞ ⎞ dτ exp ⎛ ⎝ kT ⎠ ⎠ ⎟

0

(14)

⎟

kT U Df 2 (t ) = DTf − (DTf − Ds ) exp ⎛− exp ⎛− ⎞ ⎝ kT ⎠ ⎝ h

(10)

⎜

where σ (τ ) is the time dependent magnitude of an effective stress measure that accounts for the multiaxial stress state. This equation can now be used to represent damage initiation evolution during fatigue loading. Nevertheless, it should be noted that Eq. (10) results in decreasing fatigue life with increasing mean stress, which contradicts testing data as reported in [34,35,38]. One remedy to address this issue, is to allow temperature rising during fatigue instead of maintaining it

∫t

t

0

γσ (τ ) ⎞ ⎞ dτ , exp ⎛ r ⎝ kT ⎠ ⎠ ⎟

(15) where damage DTf is obtained by enforcing Eq. (9) with Ds as lower limit of the integral of Eq. (9). For T-C or C-T load, D = Df2 without effect of Df1. The assumed sawtooth profile simplifies the stress function in Eqs. (14) and (15) and enables the included integrations. It should be mentioned here that for 273

Composite Structures 213 (2019) 271–283

M. Naderi et al.

Fig. 1. Schematic of static and fatigue damage for different load ratios: a) Tension-Tension (T-T), and Compression-Compression (C-C); b) Tension-Compression (TC), Compression-Tension (C-T).

Ψm = Cm ( ϕf ((Cc − Cf )(ϕm I + ϕf A))−1a + ηm − (ϕm I + ϕf A)−1ηc )

the sake simplicity no distinction is made here between damage caused by the tension and damage caused by compression load for the numerical examples presented here. However, in the Discussion section we will explain how one can distinguish two different tension and compression damages. The key element in the evaluation of damage in the KTF model is to the proper calculation of the equivalent stress of the constituents (fiber and matrix) in Eq. (9). Therefore, a micromechanical analysis is deemed necessary for assessing stress state in the constituents.

(19)

a = Cc ηc − ϕf C ηf − ϕm C m ηm f

(20)

where σM andσc are matrix and composite ply stress tensor respectively, ΔT is temperature difference between initial and current state and I is identity tensor. Cm , Cf , and Cc represent matrix, fiber, and composite ply stiffness matrices, respectively. Parameters ηm , ηf , and ηc are matrix, fiber, and composite ply thermal expansion vector, respectively. Parameters ϕf and ϕm are the fiber and matrix volume fractions, respectively. Ply level stresses at each integration point are obtained from a FEA of the ply and are introduced to MCT theory. In order to use constituent level stresses in kinetic theory, a scalar effective stress for matrix needs to be defined. Researches attempted to use different methods such as functions based on normalized fatigue strength, and modified von Mises stress [36,54–56]. Without loss of generality in this work, we simply adopted an effective stress measure that combines the stress perpendicular and parallel to the fiber similar to the work of Kapidzic et al. [38,39] as follows.

2.4. Micromechanical approach There are many approaches available in the literature to analytically evaluate the stress at each constituent using macro level stress [47–52]. In the present work, the multicontinuum theory (MCT) of Mayes and Hansen is adopted to extract stress at the micro level [50]. The aim of this work is not to judge the pros and cons of a micromechanics based analytical method. The accuracy of this method has been reported as reasonable in the World-Wide Failure Exercise (WWFE) [53]. Hence, we selected the MCT as an analytical based method to access stresses at micro level. A finite element based micromechanics analysis is considered here to assess the stress field of the constituents and to verify stresses obtained from analytical approach (MCT). It is noted that in a multiscale analysis framework, micromechanics based analysis has the advantage of lowering the computational cost while diminishing the accuracy.

σ=

σm2 ,22 + σm2 ,12

(21)

with σm,22 , σm,12 are the transverse normal and shear matrix component along transverse loading direction. The effective stress of Eq. (21) will be used in Eq. (9) for the evaluation of crack initiation damage.

2.5. Brief description of the multicontinuum theory

3. Numerical implementation

The idea of MCT is according to the constituent stress/strain field decomposition. MCT uses volume average of the microlevel stresses, calculated from constituent and composite constitutive [50,52]. The following equations are used in MCT.

The fatigue life prediction analysis requires the simultaneous solution of the constitutive equations and damage evolution law. A multiscale analysis scheme is employed in a manner that stress distribution of lower scale is used to approximate fatigue initiation life. The appropriate constitutive equations and damage evolution model are implemented in ABAQUS through the user subroutine UMAT. Both the lower scale FE and MCT based micromechanical analyses are connected to a macroscale model at each integration point through a UMAT routine. Since it is computationally expensive to simulate each cycle of loading, the cycle jump procedure is adopted assuming damage remains unchanged for a finite period of cycles increment (ΔN) [57]. To speed up the simulation, only the most critical elements are considered for

σM = Qm σc − Ψm (ΔT )

(16)

Qm = Cm (Cc (ϕm I + ϕf A))−1

(17)

A=−

ϕm (Cc − Cf )−1 (Cc − Cm) ϕf m

(18) 274

Composite Structures 213 (2019) 271–283

M. Naderi et al.

Fig. 2. Current fatigue numerical flow chart.

region of interest is considered for micro-macro analysis. The adopted frequency of the test was 10 Hz.

crack initiation region. The algorithmic implementation of the adopted multiscale numerical procedure is shown in Fig. 2. Macroscale stress at each integration point is initially calculated. Then it is checked whether any crack has or has not initiated in the model from previous increments. FE based micromechanics or MCT is subsequently employed to determine constituents’ stresses. In MCT calculated macrostress is the input for calculation of fiber and matrix contribution. While, in FE micromechanics, constituents’ stress is the outcome of small scale finite element simulation. In this case, macroscale strain distribution is considered as applied boundary condition for unit cell. Using equivalent stress, matrix damage evolution is calculated by Eqs. (11)–(14). It is noted that the focus of current work is initiation life in matrix failure. A small value for jump in the cycle is used to speed up the simulation process. Once microcracks are initiated in lower scale model, macroscale is immediately considered to contain a crack. We presumed that once an element of macroscale model (in the order of tenth of millimeter) has failed, a crack is initiated and the simulation is stopped (see Section 6 for more details). Throughout the entire simulation, it is assumed that initiation life coincides with final failure due to the expectation that the time for crack propagation is very short.

4.2. Three point bending (TPB) test Three point bending model set up as shown in Fig. 3b is performed on 90° laminate. The modeled domain was considered to have a length of 57.25 mm, thickness of 4.93 mm, and width of 6.35 mm with two support span of 50.8 mm. The interested region for multiscale analysis is in the middle and outermost layer of the beam. The beam is subjected to 10 Hz frequency of constant load ratio (R = 0.1) with different load severity levels ranging from 65% to 55%. Here the term “severity” is defined as the ratio the maximum applied fatigue load to the static load at failure. Test data are taken from the work of O’Brien et al. [58]. 4.3. Four point bending (FPB) test As shown in Fig. 3c, four point bending test is represented as a UD composite beam of length 57.25 mm, width of 6.35 mm and thickness of 3.3 mm. the span of roller supports is 50.8 mm while that of loading supports is 25.4 mm. The region of interest for multiscale analysis was defined to be in the middle and outermost layer of the beam. The beam is subjected to 10 Hz frequency of constant load ratio (R = 0.1) with different load severity levels ranging from 85% to 70%. Test data are taken from the work of O’Brien et al. [58].

4. Finite element models In order to analyze the potential of KTF method in predicting the onset of crack ignition life, a few numerical examples involving unidirectional laminates are chosen, including tension of 90°, 30°, and 45° unidirectional laminate, three point bending (TPB), four point bending (FPB) and short beam shear (SBS) tests. Fig. 3 depicts the respective models in terms of their geometry, dimensions and boundary conditions. The materials used were T800/2500 and IM7/8552 and the respective properties data are summarized in Table 1. Each macro model is meshed with C3D8 brick elements about 0.1 mm in size. Only one element is considered along the width assuming even distribution of damage along the width [10]. All analysis is based on isothermal material properties at room temperature.

4.4. Short beam shear (SBS) test The SBS is a composite beam of length 23 mm, thickness 2.5 mm and 5 mm width with as presented in Fig. 3d. To accurately capture deformation and stress distribution in SBS model, two rigid roller supports and one loading roller with a frictionless contact are modeled. Diameter of support roller and loading roller are 3 and 6 mm, respectively. The support rollers are spanned by 12.7 mm. Two interested regions for multiscale are one between 5th and 6th ply from the top 1.3 mm away from loading roller and the other region in the midplane and about 3.8 mm away between support and loading roller [10,59]. The beam is subjected to constant load ratio of R = 0.1 with different load severity levels ranging from 90% to 70%. Experimental test data and set up are reported in the work of May and Hallett [59]. The FE based micromechanics model is hexagon packed unit cell with periodic boundary conditions. At each integration point of macroscale interested region, the periodic boundary conditions are assigned

4.1. Tension test As shown in Fig. 3a, the plain on-/off-axis laminates (90°, 30°, and 45°) are fixed at one end while other end is subjected to fatigue loading. The plate is considered to have 10 mm height, 1 mm width and 0.1 mm thickness. In order to allow micro initiation in the model, only the 275

Composite Structures 213 (2019) 271–283

M. Naderi et al.

Fig. 3. Finite element models: a) UD laminate; b) three point bending (TPB); c) four point bending (FPB); d) short beam shear (SBS). Table 1 Material data for T800H/2500 and IM7/8552 [5,33,34,58,59,61,62]. T800H/2500

E11 (Gpa) E22 (Gpa) G12 (Gpa) G23 (Gpa) ν12 ν 23 U (KJ/mol) γ (KJ/Mpa/mol) α11 (10−6 °C−1) α 22 (10−6 °C−1)

IM7/8552

Laminate

Matrix

Fiber

Laminate

Matrix

Fiber

159 8.8 4.78 4 0.35 0.42 – – −1.06 26

4 4 1.5 1.5 0.35 0.35 125 1.4 26 38

245 12.8 15 5.5 0.23 0.25 – – −1.03 13.04

165 9 5.6 2.8 0.34 0.49 – –

4.08 4.08 1.5 1.5 0.38 0.38 160 1.0

276 19 27 7 0.2 0.2 – –

based on macroscale stress-strain field of full one cycle. The unit cell is modelled with wedge element to allow periodic mesh along the depth of the model. The unit cell has fiber volume fraction of 60% with radius of 3.5 μm. In what follows, PBC can be briefed as follows. More details are referred to [60]. f+ f− 0 ⎧ui − ui = 2a ∊ij ⎪ e+ e− 0 u − ui = 2a ∊ij ⎨ i ⎪ uiv + − uiv − = 2a ∊ij0 ⎩

(22)

Where ui f + and ui f − are displacement of opposite faces on unit cell, uie + and uie − are displacement of opposite edges on unit cell and uiv + and uiv − are displacement of opposite corner vertices on unit cell. 5. Results Fig. 4. a) Simulated stress-strain response of UD lamina under static loading compared with test data for T800/2500 [61], b) Simulated shear nonlinearity response against test data for IM7/8552 [66].

In this section the results of fatigue life simulation based on KTF are presented. First, the results of UD laminates response under static loading is verified against test data. Then, numerical verification examples for MCT and micromechanics are presented followed by 276

Composite Structures 213 (2019) 271–283

M. Naderi et al.

simulated fatigue initiation life results.

Table 2 Comparison of average stress obtained from MCT, FE unit cell and macro model.

5.1. Static test verification Fig. 4a and b shows the stress-strain response of UD 90°, 45°, and 30° laminate subjected to static loading for T800H/2500 and IM7/ 8552. Macroscale model of Fig. 3a is used to model static response. For T800H/2500, material non-linearity in shear is modeled using a third degree polynomial expression [63] (τyz = c1 γyz3 + c2 γyz2 + c3 γyz + c4 with c1 = 7.112e6; c2 = −6.892e5; c3 = 14.5e3; c4 = 0.33 Mpa). For IM7/8552, material non-linearity in shear is modeled using a Rambergτyz

τyz

( )

σ1 σ2 σ3

MCT (Mpa)

FE micromechanics (Mpa)

Macro stress at center element (Mpa)

143.3 21.85 61.8

150.2 27.5 66.7

100 ∼0 ∼0

temperature evolution due to heat condition is not considered. In addition, it should be noted that the simulations are restricted to mid stress level (mid stress is the stress resulting to a fatigue life approximately less than 105 cycles). Fatigue crack initiation life of T800H/ 2500 and IM7/8552 of four different test cases are modeled using KTF: on/off axis unidirectional laminate (UD), three point bending (TPB), four point bending (FPB), and short beam shear (SBS) test. The experimental fatigue data of T800H/2500 are obtained from the work of Kawai et al. [61] while test data of IM7/8552 are obtained from the work of Makeev et al. for UD laminate [70], O’Brien et al. for TPB, and FPB [58], and May and Hallett for SBS test [59]. Since the key parameters in KTF are of the activation energy U and activation volume γ they need to be calibrated based on creep or constant stress rate tests. In the absence of such tests, the simplest way to define these KTF parameters is to utilize available fatigue data of S-N curves for on/off axis UD laminates (90° or 45° laminate) and apply an optimization algorithm to approximate. Also, it is assumed that when microcracks are initiated, final failure occurs due to the short time assumed between initiation and final failure. However, this simplification demonstrated some drawbacks as they will be described later. In addition, since fatigue test data for composites shows a large scatter band, experimental data showing more repeated test data and larger scatter band would be beneficial to be considered for validation of KTF theory.

4.03

Osgood equation (γyz = 5 + 298 ) [64]. For both material, the simulation results are in good agreement with test data available in Refs. [61,65,66]. It should be noted that one might argue the rate dependency of shear stress-strain behavior for matrix material Fig. 4 [67,68]. Therefore, load-displacement curve to failure can be dependent on loading rate [68,69]. 5.2. MCT based micromechanics validation Before using an analytical micromechanics (here MCT), it is necessary to demonstrate that the implemented MCT microscale algorithm in UMAT correctly predicts the stresses at the constituents’ level. The FE model (90° ply) of Fig. 4a is subjected to 10 N tensile static load, which produces about 100 MPa stress along the loading direction and material is IM7/8552 of Table 1. As shown in Fig. 5, the center element in the model is considered for MCT and FE micromechanics analysis. FE micromechanics model is a hexagonally packed unit cell with periodic boundary conditions and with fiber volume fraction of 60%. The solution of macroscale model is used as the periodic boundary conditions on the unit cell. Multiscale stress analysis is performed and the results are summarized in Table 2. Fig. 5 shows the matrix stress (σ11) distribution in the unit cell and confirms that microscale stress distribution is different from macroscale stress. The average stress obtained from FE micromechanics is higher than the one obtained from MCT. Therefore, it should be noted that MCT can only capture the average unit cell stress. However, the FE micromechanics can provide both average and local stress distribution. MCT or other analytical based micromechanics are incapable of capturing local stress concentration in microscale especially where voids and defects are presented.

5.4. Fatigue simulation of T800H/2500 The UD laminates (90°, 45°, 30°) of Fig. 3a are subjected to cyclic loading with load ratios of R = −1, 0, 0.1, 0.5. We start with S-N fatigue data of 90° to calibrate U and γ. These parameters are used without any change to estimate fatigue life of different tests. As mentioned earlier, laminate properties and material data are summarized in Table 1. The associated FE model is shown in Fig. 3a. In order to initiate microcrack in the macro model, only one element (weakest) is considered to initiate a crack. It is noted that the crack propagation is not focus of the current work and out of our scope, the crack initiation onset is equivalent to final failure. Comparison between simulated and test data are presented in Fig. 6a for 90° UD laminate under fatigue load with different load ratios (R = −1, 0.1, 0.5). In these results, stress at the constituent level is captured using MCT. Fig. 6a depicts simulated

5.3. Fatigue initiation life prediction The main goal of this work is to predict fatigue crack initiation onset within the matrix of a composite. An important assumption here is that the time between crack initiation and final failure is very short. Also, it should be mentioned that fatigue properties are considered for isothermal conditions at room temperature (i.e. U and γ are constant) and

Fig. 5. MCT verification example along with matrix x direction stress distribution in FE micromechanics. Unit cell dimension and stress unit are based on μ unit. 277

Composite Structures 213 (2019) 271–283

M. Naderi et al.

of 45° angle ply (65, 60 and 55% of static strength) of R = 0.1 is plotted in Fig. 8c. As the fatigue process continues, damage accumulates until damage Di reaches the threshold value of one that is associated with the failed state of the element. 5.5. Fatigue simulation of IM7/8552 Fatigue crack initiation life of 90° IM7/8552 for four different tests including simple tension, TPB, FPB, and SBS are studied and the results are compared with test data. Predictions are performed according to MCT and FE micromechanics in which the results are based on the average stress of the unit cell. The evolution of the predicted number of cycles to initiation for UD 90 and 45° ply subjected to tension-tension fatigue load with the ratio of 0.1 is compared in Fig. 7a with relevant test data reported in [70]. The reported test data are reproduced in terms of load severity versus number of cycles to initiation. These experimental data are used for calibration of U and γ parameters of KTF. Some discrepancies are observed in the simulated data obtained from MCT and FE micromechanics. The reason is attributed to the effect of local stress concentration magnitude on the average stress of unit cell. As seen in Fig. 5 and Table 2, the average stress obtained from FE micromechanics is higher than that of obtained from MCT. However, the overall trend of simulated results compares well with test data. The results of TPB and FPB under tension-tension fatigue load with the ratio of 0.1 is presented in Fig. 7b and c. The simulation predictions are compared with experimental data taken from [58]. The error bars indicate the range of test data. Predictions are more conservative than the average values of the test data, yet they remain within the range of experimental data. Similar discrepancies between FE micromechanics and MCT predictions are observed as in the previous case. Fig. 8c depicts the comparison between simulation and test data of SBS test performed at constant fatigue load with the ratio of 0.1. As stated in [59], failure initiation takes place close to the contact of the loading roller with the specimen. Two possible damage initiation locations are selected for analysis according to the experimental observation [59]. One location is about 1.3 mm away from loading roller and around 5th and 6th layer from the top surface of the beam. The other region is about 3.8 mm away from loading roller and around midline of the beam. Rollers are modeled as rigid bodies to ensure a pronounced effect of their presence in terms of the stress distribution underneath them. As indicated in Fig. 8d, the predicted cycles to initiation for the first selected region (1.3 mm away from loading roller) appears to be more conservative than the test data. However, the predictions for the second selected region (3.8 mm away from loading roller) tend to over-predict the fatigue initiation life. This can be attributed to the fact that stress distribution close to roller tends to be more localized than away from roller due to the contact induced stress concentration. 6. Discussion

Fig. 6. Comparison of predicted fatigue life of T800H/2500 with test data [61] for a) 90° laminate at different load ratios of 0.5, 0.1, and −1; b) 90, 45, and 30° laminate at load ratio of 0.1; c) damage evolution at different load level for 45° laminate.

Overall, the current KTF based fatigue damage initiation life predictions are in a reasonable agreement within the range of test data with some slight under/over estimation. A crucial motivator for incorporating the KTF in our analysis is the fact that only a minimal set of material parameters are required while at the same time it can encapsulate the direct effect of temperature, creep and loading frequency. The very promising predictions of our approach still suggests that some further investigation and detailed discussions on KTF are still required. The presented results were based on some assumptions that in fact need further analysis for identifying the source(s) of their drawbacks. These assumptions involved the temperature effect, the effective stress selection, calibration of KTF parameters (U and γ), the definition of crack initiation onset, the low and high amplitude applied load, the compression load, and the KTF parameters for fiber/matrix interface.

fatigue life for 90° laminate at different load ratios of 0.5, 0.1, and −1. Within the range of mid stress level, the simulated data is in reasonable agreement with test data reported in [61]. Fatigue life predictions for different laminates (90, 45, and 30°) laminate under cyclic loading with ratio of R = 0.1 are presented in Fig. 6b. In the prediction simulations, the failure is considered to have occurred when the critical element with size of about one tenth of millimeter has failed. A comparison of simulation and test data shows reasonable trend considering assumptions mentioned previously. Also, fatigue damage initiation evolution Eq. (9) for three different load level 278

Composite Structures 213 (2019) 271–283

M. Naderi et al.

Fig. 7. Fatigue initiation life prediction for a) on/off axis UD laminate compared with experimental data [70]; b) three point bending compared with experimental data [58]; c) four point bending compared with experimental data [58]; d) short beam shear compared with experimental data [59].

level, temperature effect is not negligible. To qualitatively address temperature issue, a low dimensional and therefore simple temperature analysis based on viscoelastic formulation is presented. Composite materials under cyclic loading show temperature rise depending on the load level, test frequency and matrix or polymer type. Energy lost per cycle (area of hysteresis loop) is responsible for temperature rise during cycling. Following the work of Hanh and Kim [72], Xiao [73] and Daily and Broutman [74], and considering heat loss to the surrounding, temperature increase due to internal heat source for a plain specimen with no directional heat flux dependencies is given by [72,73]

q dT hA = − (T − T0) dt ρcp ρcp V

Fig. 8. Simulated temperature and damage evolution versus fatigue life for 90° UD T800h/2500 at load severity of 0.5, load ratio of 0.1, and different frequency of 10, 20, 50 Hz.

(23)

where T0 is ambient temperature, ρ is density, cp is specific heat, and h is convective heat transfer coefficient. A and V are surface area and volume of fatigue specimen. q is heat generation rate due to inelastic deformation of polymer and can be estimated using viscoelastic material theory. Assuming all input energy converts to heat, we can write

Some detailed discussions are referred to the work of Regal et al. [71]. In the following subsection we discuss some of above assumptions in more detail and highlight a plan for a future path.

q = wf

6.1. Temperature effect

(24)

with f as the frequency of the test and w as inelastic work. For current temperature analysis we assume that the inelastic work due to inelastic deformation is attributed to energy dissipation due to viscoelastic nature of polymeric materials during cyclic loading. The inelastic work is related to the damping property of the viscoelastic material, effective stress and strain through the following equation [75].

Temperature term in Eq. (6) conveys a physical feature of polymer during cyclic loading. For a fixed maximum applied stress, the higher the mean stress, the lower is the stress range and the smaller is temperature variation. Hence, lower molecular bond breakage or delay in microcracks initiation time is expected when compared with the low mean stress amplitude. To properly estimate fatigue life for any load history using KTF, temperature rise during cyclic loading should not be ignored. It is noted that current simulation was based on isothermal room temperature test data assuming temperature rise during fatigue is not significant on the matrix behavior. This may be true for low frequency and low stress level. However, for high frequency and high load

w = 0.5 ∗ φσe εe

(25)

where φ is damping factor, σe and εe are effective stress and strain in matrix. Having w and heat transfer coefficient h, transient temperature rise can be estimated as follows [72]. 279

Composite Structures 213 (2019) 271–283

M. Naderi et al.

ΔT =

q ⎛ ⎛ hA ⎞ ⎞ ⎜1 − exp ⎜ ρc V t ⎟ ⎟ ⎝ p ⎠⎠ ⎝

hA V

6.3. Calibration of KTF parameters (U and γ) (26)

Activation energy (U) and activation volume (γ) are the two crucial parameters in KTF. These parameters must be correctly estimated or calibrated through a series of constant stress or creep tests similar to those conducted by Zhurkov [44] and Hansen and Baker [45]. Throughout this work and most of work on fatigue life estimation using KTF, these parameters are calibrated from available S-N curves data of UD transverse or angle ply. However, one may argue the applicability and accuracy of calibrated parameters for the case of general applicability of KTF. Calibrating parameters from macroscale failure based SN data may be considered as an oversimplified method. S-N curves mostly represent the final break down of the material. Many sub-micro or micro cracks especially in matrix and fiber/matrix are already developed long before final rupture. Also, temperature history data pertained to same S-N curves are not available in the literature. Therefore, using macroscale S-N curve at failure without considering temperature effect can lead to inaccurate and oversimplified estimation of U, and γ that in turn may end up to generating unreasonable predictions for a more general multiscale fatigue load. Without loss of generality we advocate that it may be more reasonable to use small scale fatigue tests that can capture the onset of first initial damages developed in the matrix or matrix/fiber interface. Hence, the resulted S-N curves from these tests might be used to represent properly the onset of sub-micro or micro crack initiation depending on the accuracy of the device capturing physically distinct damage features. In addition, U, and γ are rate dependent parameters. A comprehensive understanding of behavior of U, and γ above room temperature is crucial. Most importantly, in composite materials, fiber/matrix interface is usually the weakest point in the absence of void/defect within just matrix. Interface initially fails and then sub-micro or micro cracks are developed in the matrix. Therefore, activation energy and volume at fiber/matrix interface might be different from U, and γ of pure matrix and one raises a question how to experimentally calibrate them at interface. We need to point out that as stated by Regel et al. [71], the physical meaning of term γ is related to activation volume in which elementary breakage takes place. The γ is the product of activation volume “Va” and overstress coefficient “q” as the ratio of local stress to average stress at the locations of local failure. The minimum activation volume is equal to volume per atoms in the material corresponding to the value of q = 1 without considering overstress. However, considering local overstress in molecular chain, the activation volume is some order of magnitude larger than volume per atoms in the materials.

Studies showed that activation energy U and activation volume γ of polymers are temperature dependent and will increase with temperature rise [76,77]. In fatigue life simulations based on KTF, Fertig [36] observed that around 80° temperature range, activation energy and activation volume increase up to 30% and 300%, respectively. For illustration purposes and in the absence of material data at temperature above the room temperature, we assume a linear increase in activation energy and activation volume up to 30% and 300%, respectively. However, this simplification may end up predicting a very conservative fatigue life. Therefore, additional experimental work is necessary to establish the relationship between temperature rise and activation energy/volume even with small range of temperature increase. Polymer properties beyond some temperature limit may be significantly different and therefore may have a drastic impact on the utilization of the KTF for fatigue initiation life predictions. In order to demonstrate the importance of temperature, simulations are performed for UD 90° ply under 57% load severity with the ratio of 0.1. Damage and temperature evolutions are simulated at different frequency of 10, 20, and 50 Hz. Temperature rise is calculated from Eq. (26) and the results are presented in Fig. 8. Density, specific heat, convective heat transfer and damping factor are set as ρ = 1600 kgm−3, Cp = 1100 J(kg°C)−1, h = 25 Wm−1 °C−1 [78] and φ = 2 (%). As anticipated, damage evolution is more intense once temperature rise becomes noticeable. The continuous temperature rise and shoot-up at the end is attributed to the coupling of the effective stress with the damage parameter Di that causes the inelastic energy to grow as the fatigue load application progresses in time. At higher frequency, temperature rise is evident and results in an increase of the activation volume and energy. Hence, this is consistent with the observation of a decrease in fatigue life initiation. Regardless of simplifications in temperature calculation and its effect of fatigue life, the importance of addressing temperature in KTF calculation is hardly disputable. In order to increase the accuracy of the calculation for any arbitrary loading and frequency, temperature must be considered, and in the future we plan to introduced a fully coupled temperature and stress analysis that satisfies energy and momentum conservation self-consistently. 6.2. Effective stress selection

6.4. Crack initiation onset

One of the most important parameters in fatigue life estimation using KTF is the definition of an equivalent stress in the matrix. Researches attempted to use different methods such as functions based on normalized fatigue strength, and modified von Mises stress [36,54–56]. For simplicity, in the current work we assumed microcracks in matrix are extended by stress perpendicular to fiber and stress along the fiber. However, this assumption might affect accuracy of time integral of bond rupture rate as well as the calibration of U and γ from transverse fatigue test data. For a general loading case where both tension and compression are present more accurate effective stress are necessary similar to the effective stress proposed by Fertig [36]. In the future, we plan to consider a volumetrically and distortionally decomposed strain energy density based stress measure. It is worthwhile mentioning here that the KTF presented here may break down at low and high stress level and the logarithm of time to break down is not linear anymore [45,79]. It is found that the kinetic nature of bond reformation becomes dominant at the low stress level and consequently time-to-rupture approaches to infinity. At very low stress levels molecular chain reformation or healing is active and cannot be neglected (Eq. (6)). In addition, for load histories containing compressive stress, bond reformation term of Eq. (6) cannot be neglected.

Fatigue induced damage can be considered as a multiscale problem that can greatly influence its macroscale manifestation that is very important for the designers community. Since any types of damages in composite start from the lower scale, a multiscale modeling technique appears to be a promising approach for understanding the underlying mechanisms of damage initiation and evolution. Also, it is extremely important to distinguish the definition of initiation onset and growth in the FE model. Fig. 9 presents a schematic of damage development from molecular scale to macro scale. Without referring to the molecular scale, one can adopt the appearance of sub-microcracks in the matrix or the matrix-fiber interface, once the density of broken molecular chain groups is saturated, at a scale corresponding to the tenth of fiber’s diameter. During fatigue loading, some of these sub-microcracks may arrest or block others while the rest may grow up to being fully developed as microcracks in the order of magnitude of a unit cell. As fatigue progresses, accumulation and growth of microcracks can be developed in the order of 100 μm which is about a representative volume structure (representative volume structure is a volume of material can mimic macroscale behavior). It is therefore fair to assume that after this stage, a tiny crack about tenth of mm is detectable and the macroscale 280

Composite Structures 213 (2019) 271–283

M. Naderi et al.

Fig. 9. Schematic of crack at different length scales and its order of magnitude.

macromolecular decomposition initiated by the active end primary free radicals in stressed polymers [44,71,80]. The electron paramagnetic resonance (EPR) tests showed that free radicals are formed in stressed polymers and microcracks are assumed to be initiated once generated free radicals approaches a critical value [80]. Zhurkov showed that the rate of free radicals exponentially grows under constant load rate for Nylon 6. However, time dependency of radical formation in constant load tests in which the rate of radical formation is negative cannot be justified. Other drawbacks in apply kinetic theory of fracture to polymer fracture is due fact that time to fracture equation in KTF theory is macroscopic based equation and how this equation is applicable for molecular chain rupture and radical formation [81]. Also, others argued that KTF theory cannot be hold for very short or very long experimental fracture tests [82–84]. In addition, atomic stress should be replaced with macroscopic stress in molecular chain rupture [81]. It is shown that KTF has very promising potential for being considered in a multiscale analysis framework for modeling microcrack initiation. Challenges and difficulties are discussed to elucidate the future research path for using KTF. Substantial effort is needed to enhance the understanding of KTF for damage initiation prediction. It is anticipated that, if used properly, the theory can provide stiffness and strength degradation during fatigue life without relying on S-N curves or any other phenomenological approaches which highly dependent on the loading conditions. The multiscale data-driven inversion methodology presented in [40–42] that was developed for deriving microscale material properties of healthy and damaged constituents can now be extended to identify the KTF parameters (U and γ) and for the case addressing the fatigue life of the relevant systems. Proper tests including constant stress and constant stress rate or creep tests will be conducted at macroscale in order to enable the characterization of the parameters U and γ. The appropriate data will be acquired by using advanced experimental techniques (i.e. full field strain measurement) and equipment (high resolution camera, thermography, etc.). Finally, considering above mentioned discussion on KTF, the authors’ attempts were to revisit KTF strength and drawbacks within a preliminary work. Further work still needs to be done to enhance the understanding of KTF and its applicability for modeling fatigue crack initiation life.

can start cracking. Stress contours resulting from FE micromechanics (Fig. 5) analysis revealed that stress concentrations appear to be high around fiber/ matrix interface and immediately under fatigue load density of broken groups of molecular chains is high. However, it takes time for the broken groups of molecular chains to coalesce and connect to each other before sub microcracks are initiated. For discussion purposes, we compare the initiation life obtained from MCT which uses average stress with FE micromechanics model which uses local stress distribution. Due to high stress concentration, elements around interface of fiber/matrix are failed in much shorter time than those corresponding to MCT-based initiation life time. The selected load severity and ratio of simulated results are 0.4 and 0.1 respectively. FE micromechanics predicts failure initiation around region with maximum stress concentration (see Fig. 6) less than about 50 cycles which is much less than fatigue life estimated by the MCT approach. The simulated results show noticeable difference between the two different definitions of crack initiation. Crack initiation definition in MCT is based on the macro model elements’ size and average stress of unit cell. While, the initiation definition in FE micromechanics is according to micro model elements’ size and localized stress distribution. Therefore, it is meaningful to observe that simulated fatigue initiation life of FE based is much shorter than that of MCT based. It is noted that the fatigue initiation life based on average stress obtained from MCT corresponds to the size of a crack in the order of magnitude of one tenth of millimeter (mm). On the other hand, fatigue initiation life based on localized stress obtained from FE micromechanics corresponds to the size of a crack in the order of magnitude of one tenth of micrometer (μm). In both cases, the input value for U and γ are obtained from S-N curves of macroscale fatigue tests. Depending on what the level of modeling and damage initiation definition are, the predictions can be different. Therefore, depending on the definition of damage initiation size and level of approximation, some information may be lost and this may result in over/under prediction. It is apparent that the damage initiation onset in KTF theory is where density of broken molecular chains is saturated and sub-microcrack in the order of tenth of fiber’s diameter are initiated. For the sake of computational cost, one can use the initiation onset as the failure of unit cell considering that material data U and γ are correctly calibrated using experiments accurately capturing moments of micro crack initiation life. It is noted that according to Zhurkov kinetic theory of fracture, submicrocracks formations are related to chain reactions of

7. Conclusions The KTF is used to calculate time to rupture of on/off axis UD 281

Composite Structures 213 (2019) 271–283

M. Naderi et al.

laminate subjected to cyclic loading. MCT and FE micromechanics are utilized to calculate average time history of stress in the matrix. Effective stress is based on the assumption that major failure type in matrix takes place perpendicular and parallel to the fiber. Various numerical examples including tension of UD laminate, three point bending, four point bending and short beam shear are chosen for verification of KTF to predict the initiation of damage onset in fatigue life. It is assumed that fatigue initiation life is almost equal to the time at final failure. Some issues and future plans are discussed in order to improve KTF usage in the context fatigue life estimation. Comparison of simulated results and predictions with experimental data at different frequencies reveals the importance of considering temperature in KTF theory. Although throughout the results section we assumed activation energy (U) and volume (γ) are constant, it seems in addition to temperature, they are rate dependent and further investigation about their calibration is necessary. The current work shows the potential and attractive features of KTF for analysis of fatigue crack initiation life prediction. Future works on the aforementioned issues are required for accurate predictions before the method is used for a general fatigue load history. The most important identified challenges are:

[14] Quaresimin M, Ricotta M. Stress intensity factors and strain energy release rates in single lap bonded joints in composite materials. Compos Sci Technol 2006;66:647–56. [15] Turon A, Costa J, Camanho P, Dávila C. Simulation of delamination in composites under high-cycle fatigue. Compos A Appl Sci Manuf 2007;38:2270–82. [16] El Kadi H, Ellyin F. Effect of stress ratio on the fatigue of unidirectional glass fibre/ epoxy composite laminae. Composites 1994;25:917–24. [17] Plumtree A, Cheng G. A fatigue damage parameter for off-axis unidirectional fibrereinforced composites. Int J Fatigue 1999;21:849–56. [18] Varvani-Farahani A, Haftchenari H, Panbechi M. An energy-based fatigue damage parameter for off-axis unidirectional FRP composites. Compos Struct 2007;79:381–9. [19] Kawai M, Itoh N. A failure-mode based anisomorphic constant life diagram for a unidirectional carbon/epoxy laminate under off-axis fatigue loading at room temperature. J Compos Mater 2014;48:571–92. [20] Jones R, Kinloch AJ, Michopoulos JG, Brunner AJ, Phan N. Delamination growth in polymer-matrix fibre composites and the use of fracture mechanics data for material characterisation and life prediction. Compos Struct 2017;180:316–33. [21] Talreja R. Damage and fatigue in composites–a personal account. Compos Sci Technol 2008;68:2585–91. [22] Reifsnider K, Gao Z. A micromechanics model for composites under fatigue loading. Int J Fatigue 1991;13:149–56. [23] Subramanian S, Reifsnider K, Stinchcomb W. A micromechanics model for the prediction of fatigue characteristics of off-axis unidirectional laminates. Appl Compos Mater 1994;1:69–74. [24] Abdelal G, Caceres A, Barbero E. A micro-mechanics damage approach for fatigue of composite materials. Compos Struct 2002;56:413–22. [25] Akshantala NV, Talreja R. A micromechanics based model for predicting fatigue life of composite laminates. Mater Sci Eng, A 2000;285:303–13. [26] Huang Z-M. Micromechanical modeling of fatigue strength of unidirectional fibrous composites. Int J Fatigue 2002;24:659–70. [27] Shokrieh M, Esmkhani M. Fatigue life prediction of nanoparticle/fibrous polymeric composites based on the micromechanical and normalized stiffness degradation approaches. J Mater Sci 2013;48:1027–34. [28] Adibnazari S, Farsadi M, Koochi A, Khorashadizadeh S. New approach for fatigue life prediction of composite plates using micromechanical bridging model. J Compos Mater 2015;49:309–19. [29] Krause D. A physically based micromechanical approach to model damage initiation and evolution of fiber reinforced polymers under fatigue loading conditions. Compos B Eng 2016;87:176–95. [30] Qian C, Westphal T, Nijssen R. Micro-mechanical fatigue modelling of unidirectional glass fibre reinforced polymer composites. Comput Mater Sci 2013;69:62–72. [31] Sayyidmousavi A, Bougherara H, Fawaz Z. A multiscale approach for fatigue life prediction of polymer matrix composite laminates. J Reinf Plast Compos 2015;34:1099–109. [32] Li W, Cai H, Li C, Wang K, Fang L. Micro-mechanics of failure for fatigue strength prediction of bolted joint structures of carbon fiber reinforced polymer composite. Compos Struct 2015;124:345–56. [33] Hosseini Kordkheili SA, Toozandehjani H, Soltani Z. A progressive multi-scale fatigue model for life prediction of laminated composites. J Compos Mater 2017;51:2949–60. [34] Fertig R, Kenik D. Physics-based fatigue life prediction of composite structures. NAFEMS World Congress 2011: NAFEMS. Boston 2011:12. [35] Fertig R, Kenik D. Predicting composite fatigue life using constituent-level physics. 52nd AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics and Materials Conference 19th AIAA/ASME/AHS Adaptive Structures Conference 13t2011. p. 1991. [36] Fertig R. Bridging the gap between physics and large-scale structural analysis: a novel method for fatigue life prediction of composites. Proceedings of the SAMPE Fall Technical Conference. 2009. [37] Dalgarno RW, Action JE, Robbins DH, Engelstad SP. Failure simulations of openhole IM7/977-3 coupons subjected to fatigue loading using Autodesk Helius PFA. J Compos Mater 2017;51:2119–29. [38] Kapidžić Z, Ansell H, Schön J, Simonsson K. Fatigue bearing failure of CFRP composite in biaxially loaded bolted joints at elevated temperature. Compos Struct 2015;127:298–307. [39] Kapidžić Z, Ansell H, Schön J, Simonsson K. Fatigue bearing failure of CFRP composite in bolted joints exposed to biaxial variable amplitude loading at elevated temperature. Compos Struct 2016;142:71–7. [40] Komninelli F, Iliopoulos A, Michopoulos JG. Towards identification of lower scale composite material properties from higher scale experimental data via inverse analysis of coupled multiscale models. ASME 2015 international design engineering technical conferences and computers and information in engineering conference. American Society of Mechanical Engineers; 2015. p. V01AT2A061–V01AT02A. [41] Michopolous J, Iliopoulos A, Hermanson J, Steuben J, Komnineli F. High performance data-driven multiscale inverse constitutive characterization of composites. In: Brancherie D, Feissel P, Bouvier S, Ibrahimbegovic A, editors. Mechanical engineering and solid mechanics series: from microstructure investigations to multiscale modeling: bridging the gap. London, UK: ISTE Ltd; 2017. p. 197–238 Hoboken, NJ: John Wiley & Sons Inc: 197–238 Chapter 8. [42] Michopoulos J, Iliopoulos A, Steuben J, Hermanson J. Multiaxial and multiscale implications of dissipative behavior of composites. IOP Conf Series: Mater Sci Eng 2018:012005. IOP Publishing. [43] Coleman BD. Time dependence of mechanical breakdown phenomena. J Appl Phys 1956;27:862–6. [44] Zhurkov S. Kinetic concept of the strength of solids. Int J Fract 1984;26:295–307.

a) Temperature and rate dependence consideration and their effects; b) Accurate calibration of U and γ for fiber, matrix and fiber/matrix interface. c) Molecular chain reformation behavior at low stress level and compressive stress. d) Definition of initiation size and onset. Acknowledgment M. Naderi and N. Iyyer of Technical Data Analysis Inc. acknowledge financial support from US Navy under Contract N68335-18-C-0076. J Michopoulos acknowledges support from Office of Naval Research through the Naval Research Laboratory’s core funding. References [1] May M, Hallett SR. An advanced model for initiation and propagation of damage under fatigue loading–part I: Model formulation. Compos Struct 2011;93:2340–9. [2] May M, Pullin R, Eaton M, Featherston C, Hallett SR. An advanced model for initiation and propagation of damage under fatigue loading–part II: Matrix cracking validation cases. Compos Struct 2011;93:2350–7. [3] Quaresimin M, Carraro P. Damage initiation and evolution in glass/epoxy tubes subjected to combined tension–torsion fatigue loading. Int J Fatigue 2014;63:25–35. [4] May M, Hallett SR. Damage initiation in polymer matrix composites under highcycle fatigue loading–a question of definition or a material property? Int J Fatigue 2016;87:59–62. [5] Sevenois R, Garoz D, Gilabert F, Spronk S, Van Paepegem W. Microscale based prediction of matrix crack initiation in UD composite plies subjected to multiaxial fatigue for all stress ratios and load levels. Compos Sci Technol 2017;142:124–38. [6] Bergmann H, Prinz R. Fatigue life estimation of graphite/epoxy laminates under consideration of delamination growth. Int J Numer Meth Eng 1989;27:323–41. [7] Blanco N, Gamstedt EK, Asp L, Costa J. Mixed-mode delamination growth in carbon–fibre composite laminates under cyclic loading. Int J Solids Struct 2004;41:4219–35. [8] Harper PW, Hallett SR. A fatigue degradation law for cohesive interface elements–development and application to composite materials. Int J Fatigue 2010;32:1774–87. [9] Henaff-Gardin C, Lafarie-Frenot M. The use of a characteristic damage variable in the study of transverse cracking development under fatigue loading in cross-ply laminates. Int J Fatigue 2002;24:389–95. [10] May M, Hallett SR. A combined model for initiation and propagation of damage under fatigue loading for cohesive interface elements. Compos A Appl Sci Manuf 2010;41:1787–96. [11] Nojavan S, Schesser D, Yang Q. An in situ fatigue-CZM for unified crack initiation and propagation in composites under cyclic loading. Compos Struct 2016;146:34–49. [12] Nojavan S, Schesser D, Yang Q. A two-dimensional in situ fatigue cohesive zone model for crack propagation in composites under cyclic loading. Int J Fatigue 2016;82:449–61. [13] Quaresimin M, Ricotta M. Fatigue behaviour and damage evolution of single lap bonded joints in composite material. Compos Sci Technol 2006;66:176–87.

282

Composite Structures 213 (2019) 271–283

M. Naderi et al.

[45] Hansen AC, Baker-Jarvis J. A rate dependent kinetic theory of fracture for polymers. Int J Fract 1990;44:221–31. [46] Zhurkov S, Kuksenko V. The micromechanics of polymer fracture. Int J Fract 1975;11:629–39. [47] Gotsis P, Chamis CC, Minnetyan L. Prediction of composite laminate fracture: micromechanics and progressive fracture. Compos Sci Technol 1998;58:1137–49. [48] Huang Z-M. A bridging model prediction of the ultimate strength of composite laminates subjected to biaxial loads. Compos Sci Technol 2004;64:395–448. [49] Jin K-K, Huang Y, Lee Y-H, Ha SK. Distribution of micro stresses and interfacial tractions in unidirectional composites. J Compos Mater 2008;42:1825–49. [50] Mayes JS, Hansen AC. Composite laminate failure analysis using multicontinuum theory. Compos Sci Technol 2004;64:379–94. [51] Huang Z. Tensile strength of fibrous composites at elevated temperature. Mater Sci Technol 2000;16:81–94. [52] Nelson EE, Hansen AC, Steven Mayes J. Failure analysis of composite laminates subjected to hydrostatic stresses: a multicontinuum approach. J Compos Mater 2012;46:2461–83. [53] Kaddour A, Hinton M, Soden P. A comparison of the predictive capabilities of current failure theories for composite laminates: additional contributions. Compos Sci Technol 2004;64:449–76. [54] Raghava R, Caddell RM, Yeh GS. The macroscopic yield behaviour of polymers. J Mater Sci 1973;8:225–32. [55] Bauwens J. Yield condition and propagation of Lüders’ lines in tension–torsion experiments on poly (vinyl chloride). J Polymer Sci Part A-2: Polymer Phys 1970;8:893–901. [56] Ha SK, Jin KK, Huang Y. Micro-mechanics of failure (MMF) for continuous fiber reinforced composites. J Compos Mater 2008;42:1873–95. [57] Lemaitre J. A course on damage mechanics. Springer Science & Business Media; 2012. [58] O'Brien TK, Chawan AD, Krueger R, Paris IL. Transverse tension fatigue life characterization through flexure testing of composite materials. Int J Fatigue 2002;24:127–45. [59] May M, Hallett SR. An assessment of through-thickness shear tests for initiation of fatigue failure. Compos A Appl Sci Manuf 2010;41:1570–8. [60] Barbero EJ. Finite element analysis of composite materials using AbaqusTM. CRC Press; 2013. [61] Kawai M, Suda H. Effects of non-negative mean stress on the off-axis fatigue behavior of unidirectional carbon/epoxy composites at room temperature. J Compos Mater 2004;38:833–54. [62] Chamis CC, Abdi F, Garg M, Minnetyan L, Baid H, Huang D, et al. Micromechanicsbased progressive failure analysis prediction for WWFE-III composite coupon test cases. J Compos Mater 2013;47:2695–712. [63] Falzon B, Apruzzese P. Numerical analysis of intralaminar failure mechanisms in composite structures. Part I: FE implementation. Compos Struct 2011;93:1039–46. [64] He Y, Makeev A. Nonlinear shear behavior and interlaminar shear strength of

[65] [66]

[67] [68] [69]

[70] [71] [72] [73] [74] [75] [76] [77] [78] [79] [80]

[81] [82] [83] [84]

283

unidirectional polymer matrix composites: a numerical study. Int J Solids Struct 2014;51:1263–73. Makeev A, Seon G, Lee E. Failure predictions for carbon/epoxy tape laminates with wavy plies. J Compos Mater 2010;44:95–112. He Y, Makeev A, Shonkwiler B. Characterization of nonlinear shear properties for composite materials using digital image correlation and finite element analysis. Compos Sci Technol 2012;73:64–71. Jones R, Hui D. Analysis, design and assessment of composite repairs to operational aircraft. Aircraft Sustainment and Repair. Elsevier; 2018. p. 325–462. Jones R, Sawyer J, Chiu W. Studies in the matrix dominated failures of composite joints. Compos Struct 1999;44:1–16. Chiu W, Galea S, Jones R, Pitt S. Load history dependence of graphite epoxy joints/ repairs: an experimental study. Proc 11th Int Conf on Composite Materials; 1997. p. 357–64. Makeev A, Nikishkov Y, Seon G, Lee E. Fatigue structural substantiation for thick composites. Proceedings of the 17th ICCM, Scotland, UK; 2009. Regel' VR, Slutsker AI, Tomashevskiĭ ÉE. The kinetic nature of the strength of solids. Phys Usp 1972;15:45–65. Hahn H, Kim RY. Fatigue behavior of composite laminate. J Compos Mater 1976;10:156–80. Xiao X. Modeling of load frequency effect on fatigue life of thermoplastic composites. J Compos Mater 1999;33:1141–58. Dally J, Broutman L. Frequency effects on the fatigue of glass reinforced plastics. J Compos Mater 1967;1:424–42. Findley WN, Davis FA. Creep and relaxation of nonlinear viscoelastic materials. Courier Corporation; 2013. Cook WD, Mayr AE, Edward GH. Yielding behaviour in model epoxy thermosets—II. Temperature dependence. Polymer 1998;39:3725–33. Su WF, Carr S, Brittain J. Thermally stimulated discharge current studies on lowtemperature relaxation in epoxy resin. J Appl Polym Sci 1980;25:1355–63. Amiri M, Naderi M, Khonsari M. An experimental approach to evaluate the critical damage. Int J Damage Mech 2011;20:89–112. Ettouney OM, Hsiao C. Time-dependent fracture strength of solid bodies. J Appl Phys 1988;64:4884–8. Zhurkov S, Zakrevskyi V, Korsukov V, Kuksenko V. Mechanism of submicrocrack generation in stressed polymers. J Polymer Sci Part A-2: Polymer Phys 1972;10:1509–20. Fava RA. Polymers physical properties. Academic Press; 1980. Henderson C, Graham P, Robinson C. A comparison of reaction rate models for the fracture of solids. Int J FractMech 1970;6:33–40. Curran DR, Seaman L, Shockey DA. Dynamic failure in solids. Phys Today 1977;30:46. Peterlin A. Fracture mechanism of drawn oriented crystalline polymers. J Macromol Sci Part B Phys 1973;7:705–27.