main.dvi

Transcription

1 Copyright c 2005 Tech Science Press CMC, vol.2, no.3, pp , 2005 Micromechanical Analysis of Interphase Damage for Fiber Reinforced Composite Laminates Yunfa Zhang 1 and Zihui Xia 1,2 Abstract: In the present study, the initiation and evolution of the interphase damage and their influences on the global stress-strain relation of composite laminates are predicted by finite element analysis on a micromechanical unit cell model. A thin layer of interphase elements is introduced and its stress-strain relation is derived based on a cohesive law which describes both normal and tangential separations at the interface between the fiber and matrix. In addition, a viscous term is added to the cohesive law to overcome the convergence difficulty induced by the so-called snap-back instability in the numerical analysis. The matrix behavior is described by a recently developed nonlinear viscoelastic constitutive model. As application examples, glass fiber/epoxy unidirectional laminates under off-axis loadings are analyzed. One-quarter of the unit cell is used in the analysis accounting for the geometrical symmetry of the model, and the corresponding periodic boundary conditions for combined global shear and normal loading are derived. Results show that the initiation and evolution of the interphase damage can be well simulated and the predicted global stress-strain responses are in good agreement with the experimental results. keyword: Interphase damage, micromechanical modeling, composite laminates, off-axis loading, finite element analysis, viscoelastic analysis, periodic boundary conditions. 1 Introduction For a unidirectional laminate under tensile loading, the fiber/matrix interphase damage and matrix cracking are the two main damage modes [Smith (1987)]. In most cases, damage process in laminates initiates by the above two intralamina damage mechanisms. For multidirec- 1 Department of Mechanical Engineering, University of Alberta, Edmonton, Alberta, Canada, T6G 2G8 2 Corresponding author. Tel.: ; fax: ; zihui.xia@ualberta.ca tional laminates, intralamina damages may exist at the loading level much lower than the final failure load or even during the material curing process. The most common example is the loading of the 0 /90 cross-ply laminates along the 0 ply direction. In this case the damage occurs in the 90 plies at a rather low load, however, the laminates can still carry higher load in the 0 direction, although at a reduced stiffness, see, e.g., Hoover et al. (1997). Thus for an accurate analysis, it is imperative to consider the prevailing damage mechanisms, especially the intralamina ones. In most advanced composite materials, there is invariably an interphase between the two main constituents of the reinforced composites. The scale of the intrephase may be very small of the order of several microns [Daniel and Ishai (1994)]. Therefore, in many studies the interphase can be treated as an interface. Interphase/interface determines the stress transfer between the reinforcement and matrix, and thus it influences the damage process of laminates if the bonding is weak. Both macromechanical [e.g., Tian et al. (2004)] and micromechanical [e.g., Okada et al. (2004)] approaches can be used to model the mechanical response of the composite laminates, with or without an assumed damage occurrence. However, damages in the composite laminates initiate at a microscopic level and it would be desirable that a micromechanical approach be used, in which the fiber, matrix, and interphase/interface are explicitly accounted for in the model [Pagano and Yuan (2000), Zhang et al. (2005)]. Assuming a uniform distribution of fibers, a representative volume element (RVE) or a repeated unit cell (RUC) can be used to conduct micromechanical studies [Sun and Vaidya (1996)]. By using the micromechanical approach, not only the global properties of the composites can be obtained, but also various damages at microscopic level can be predicted provided proper damage initiation criteria and evolution rules are available. Some recent micromechanical studies of laminates involving matrix cracking or interphase damage are

2 214 Copyright c 2005 Tech Science Press CMC, vol.2, no.3, pp , 2005 reported e.g., in Zhang et al. (2005), Zhu et al. (1998), Yuan et al. (1997), Zhu and Achenbach (1991), among others. It is worthwhile to note that several researchers have introduced cohesive interface cracking models into the micromechanical RUC analysis. A significant advantage of the cohesive cracking model is that both the initiation and propagation of the damage can be properly simulated [Tvergaard (1990a, 1990b), Legarth, (2004), among others]. Most micromechanical analyses so far are confined to uniaxial loadings along the material principal directions, as reviewed by Pagano and Yuan (2000). For off-axis loadings, only limited attentions have been paid. Examples of recent studies of off-axis loadings are Zhu and Sun (2003) for thermoplastic matrix composites; Aghdam et al. (2001) for metal matrix composites, and Zhang et al. (2005) for thermoset matrix composites. In this paper, the initiation and evolution of interphase damage and their effect on the global stress-strain relation of composite laminates are predicted by finite element micromechanical analysis. A three-dimensional RUC model based on periodic fiber array is established, and appropriate periodic boundary conditions for combined shear and normal stress loading, with consideration of the geometrical symmetry of the RUC, are derived. To model both the normal and tangential separations of the fiber/matrix interface, a thin layer of interphase elements is introduced in the unit cell model. The stress-strain relation of the interphase element is derived based on a cohesive law [Tvergaard (1990a, 1990b)], which characterizes the dependence of the tractions on the displacement discontinuities across the interface. In addition, a viscous term is added to the cohesive law to overcome the convergence difficulty induced by the so-called snap-back instability in the numerical iterations. The matrix behavior is described by a recently developed nonlinear viscoelastic constitutive model [Xia et al. (2003)]. As application examples, glass fiber/epoxy unidirectional laminates under off-axis loading with angles of 90, 45, 30, and 20 are analyzed, respectively. Results show that the initiation and evolution of the interphase damage can be well simulated. And the predicted global stress-strain curves are also in good agreement with the experimental results. 2 Unit Cell Model and Periodic Boundary Conditions 2.1 Unit cell of a unidirectional laminate Assuming fiber distribution is periodic across the crosssection, a repeated unit cell (RUC) can be isolated from the composite laminates. The periodic fiber sequences commonly used are the square array and the hexagonal array, see Sun and Vaidya (1996), for example. In this paper, the square array of fiber distribution is assumed, resulting in a rectangular parallelepiped RUC containing one fiber, as shown in Fig.1. The cross section of the unit cell is a square with sides a, and thickness in the fiber direction is b (Fig.1b). The radius of the fiber, R, is determined by the fiber volume fraction V f of the composite such thatv f = πr 2/ a Periodic boundary conditions The analysis of a unit cell under general periodic boundary conditions can be found in Suquet (1987), Michel et al. (1999), Aboudi et al. (2001), among others, and the implementation in FEM analysis for multiaxial loading case can be seen in Xia et al. (2003), Xia et al. (2005). For a periodic array, the displacement field can be expressed as [Suquet (1987)]: u i (x 1,x 2,x 3 )=ε ij x j +u i (x 1,x 2,x 3 ) (1) where ε ij is the global strain applied to the periodic body, and u i (x 1,x 2,x 3 ) is the periodic part of the deformation. In addition, for a periodic RUC, the tractions on the opposite boundary surfaces should also meet the continuity condition, i.e. σ ij (P)n j (P)= σ ij (Q)n j (Q) (2) where P and Q are periodic points (with the same inplane coordinates) on the two opposite boundary surfaces, n is the unit outward normal vector to the surfaces, see Fig. 1(b). Generally, u i (x 1,x 2,x 3 ) is unknown prior to the solution, thus Eq. 1 cannot be directly applied as the displacement boundary conditions. However, since u i is periodic, i.e. its value is the same at the opposite surfaces of a RUC, therefore, u i (P) u i (Q)=[x j (P) x j (Q)]ε ij (3)

3 Micromechanical Analysis of Interphase Damage for Fiber Reinforced Composite Laminates 215 In a displacement FEM scheme, Eq. 3 can be implemented as displacement constraint equations. Furthermore, it has been shown in Xia et al. (2005) that in a displacement FEM scheme, a unique solution can be obtained by application of Eq.3 as the displacement boundary conditions and the traction continuity condition, Eq. 2, will be satisfied automatically. 2.3 Off-axis loading In this study, a unidirectional laminate under off-axis loading will be analyzed as an application example, Fig. 1(a). The off-axis tensile loading applied to a unidirectional laminate can be decomposed into a set of multiaxial loading in the principal material directions of a lamina, as shown in Fig.1(c). σ 11 = σcos 2 θ σ 22 = σsin 2 θ (4) σ 12 = σsinθcosθ Therefore the global stress and corresponding strain vectors can be written as: {σ} = {σ 11, σ 22, 0, σ 12, 0, 0} T (5) {ε} = {ε 11, ε 22, ε 33, ε 12, 0, 0} T (6) Note that Eq. 3 is the periodic boundary conditions for the entire RUC model. Referring to Fig.1, the geometry, material and loading of the present problem have three symmetries, viz.: (1) along x 1 axis (fiber direction), the geometry, material, stress/strain, and periodic part of displacement are invariant; (2) mirror symmetry about the plane x 3 = 0; and (3) π -rotation symmetry about the x 3 axis. In Appendix A, a detail derivation of the periodic boundary conditions for the one-quarter RUC, Fig. 1(c), under the combined shear and normal loading is provided considering the above three symmetries. In such a way, the computation time required can be significantly reduced, especially for the current nonlinear problem. The applied boundary conditions to the one-quarter model are summarized in the following: On planes x 1 = ±b/2 u 1 (b/2, x 2, x 3 ) u 1 ( b/2, x 2, x 3 )=ε 11 b u 2 (b/2, x 2, x 3 ) u 2 ( b/2, x 2, x 3 )=0 u 3 (b/2, x 2, x 3 ) u 3 ( b/2, x 2, x 3 )=0 (7) On plane x 2 = 0 u 1 (x 1, 0, x 3 )=ε 11 x 1 u 2 (x 1, 0, x 3 )=0 σ 23 (x 1, 0, x 3 )=0 On plane x 2 = a/2: u 1 (x 1, a/2, x 3 )=ε 11 x 1 +2ε 12 (a/2) u 2 (x 1, a/2, x 3 )=ε 22 (a/2) σ 23 (x 1, a/2, x 3 )=0 On plane x 3 = 0 u 3 (x 1, x 2, 0)=0 σ 23 (x 1, x 2, 0)=0 σ 31 (x 1, x 2, 0)=0 On plane x 3 = a/2: u 3 (x 1, x 2, a/2)=δ 33 σ 23 (x 1, x 2, a/2)=0 σ 31 (x 1, x 2, a/2)=0 (8) (9) (10) (11) Note that the constant δ 33 in the first of Eq.11 is not specified, therefore, it ensures that, on the plane x 3 = a/2, the total normal traction vanishes, i.e. Z σ 33 ds = 0 (12) Thus, the global stress component σ 33 = 0. In an off-axis loading of unidirectional laminates under strain control, the global strain in the loading direction (see Fig.1 (a) and (c)) is ε = ε 11 cos 2 θ +ε 22 sin 2 θ +2ε 12 cosθsinθ (13) For each time step, t, the strain increment is given by ε = ε t (14) where ε is the applied global strain rate. To simulate the off-axis loading, an iterative procedure is required to ensure that proper proportions of the increments of ε 11, ε 22, ε 12 are applied, so that Eq. 4 is satisfied at each step. The iteration procedure is as follows [Zhang et al. (2005)]: (i) For each time step t, we have the trial increments ofε 11, ε 22, ε 12, which satisfy Eq. 13.

4 216 Copyright c 2005 Tech Science Press CMC, vol.2, no.3, pp , 2005 (a) (b) (c) Figure 1 : Unidirectional laminates and unit cell model: (a) A unidirectional laminate under off-axis loading; (b) Unit cell for a square fiber array; (c) Quarter of the unit cell (ii) The solution gives the stress distribution in the unit cell, so the global stress components can be calculated from {σ} = 1 Z {σ(x 1, x 2, x 3 )}dv (15) V V where V is the volume of the unit cell. (iii) Equation 4 is checked and, if it is satisfied (within certain error limit), then one proceeds to the next step. If not, new increments of ε 11, ε 22, ε 12 are obtained and steps (i) to (iii) are repeated. For a small time step, it could be assumed that the increments of ε 11, ε 22, ε 12 are proportional to the corresponding increments of average stress components, then the new increments of ε 11, ε 22, ε 12 can be estimated from the average stresses, Eq.15. Numerical calculation indicated that the required increments of ε 11, ε 22, ε 12 could be obtained through a few iterations. 3 Cohesive Debonding Law and Cohesive Interphase Element 3.1 Cohesive law The fiber-matrix interphase/interface will be modeled by a cohesive zone model proposed by Tvergaard (1990a, b). The behavior to be captured by a cohesive law is that, as the cohesive surface separates the magnitude of the tractions at first increases, reaches a maximum and then it decreases with the increasing separation finally approaching zero. In the interface coordinate system shown in Fig. 1(b), let [u 1 ],[u 2 ] and [u 3 ] be the displacement differences across the interface and T 1, T 2 and T 3 are the corresponding tractions in the 1, 2, and 3 directions, respectively. To account for the combined effect of normal and shear decohesions, a non-dimensional parameter λ is defined as [Tvergaard (1990a,b)]: { ([u3 ) ] 2 ( ) [u2 ] 2 ( ) } [u1 ] 2 1 /2 λ = + + (16) δ n δ t1 where δ n, δ t1, δ t2 are displacement difference values in the three directions, corresponding to the complete separation. To determine the tractions a function F(λ) is chosen, F(λ)= 27 4 σ max(1 2λ +λ 2 ) for 0 λ 1 (17) And the interface tractions are given by the expressions T 3 = [u 3] δ n F(λ) T 2 = α [u 2] δ t2 F(λ) (18) T 1 = α [u 1] δ t1 F(λ) For the case of purely normal separation ([u 1 ]=[u 2 ] 0) the maximum traction is σ max, and the total separation occurs at [u 3 ]=δ n, and the work of separation per unit interface area is 9σ max δ n /16. Thus, to represent the interfacial debonding behavior of a given interface the values δ t2

5 Micromechanical Analysis of Interphase Damage for Fiber Reinforced Composite Laminates 217 of the five parameters,δ n, δ t1, δ t2, σ max and α have to be chosen such that the maximum traction and the work required for the separation for the combined modes are approximated as well as possible. 3.2 Cohesive interphase element The traction-displacement jump relation is relatively difficult to implement into a commercial finite element code, especially for three-dimensional cases. Instead, in this study the ideal interface (zero thickness) is represented by a very thin layer of interphase element with thickness of h.assuming uniform strains throughout the thickness, then, the strains can be expressed as a function of the displacement jumps as follows (in the local coordinate system shown in Fig. 1(b)): {ε 33 γ 23 γ 13 } = 1 h {[u 3][u 2 ][u 1 ]} (19) The critical strains are defined from the above equation as {ε c 33 γc 23 γc 13 } = 1 h {δ n δ t2 δ t1 } (20) Substituting Eqs. 19 and 20 into the cohesive law expressed by Eqs , and replacing the corresponding tractions T 3, T 2 and T 1 by the stresses σ 33, τ 23 and τ 13,a stress-strain constitutive law is obtained as follows: σ 33 = ε 33 ε c F(λ) 33 τ 23 = α γ 23 γ c F(λ) (21) 23 τ 13 = α γ 13 γ c F(λ) 13 where F λ = 27 2 σ max( 1 +λ) λ = 1 ( ε33 ε 33 λ ε c 33 ε c + γ γ c 23 γ 23 γ c 23 + γ 13 γ c 13 γ 13 γ c 13 ) (24) This cohesive constitutive law can be implemented into a FEM code as a stress-strain relationship. It should be noted that the constitutivelaw is defined in the coordinate system 1-2-3, in the FEM analysis, it should be transformed to the coordinate system (x 1, x 2, x 3 ), see Fig Snap-back instability during interphase damage In the simulation of damage, the snap-back of the stressstrain curve frequently results in convergence difficulties in the numerical calculations. Snap-back of the stressstrain curve indicates the strain decreasing with the decreasing of the stress [Gao and Bower (2004), Yang and Proverbs (2004)]. In this section, a simple interfacial model is considered to study the snap-back behavior and the numerical strategy to overcome the convergence difficulties. A one dimensional interfacial model is shown in Fig. 2. It consists of three plane strips representing the fiber, interphase and matrix, respectively. The fiber and the matrix are assumed to be elastic and the interphase model is described by the cohesive law. This one dimensional model can be solved either analytically or numerically. In the analytical solution, the interphase is idealized as an interface, and the cohesive law expressed by Eqs are used. While in the FEM solution, the interphase is represented by a thin layer of interphase element, and cohesive constitutive relations expressed by Eqs are used. First, let us examine the following analytical solution: λ = { (ε33 ε c 33 ) 2 + ( γ23 γ c 23 ) 2 + ( γ13 γ c 13 ) 2 } 1 /2 (22) The incremental expressions are obtained from Eq. 21 as σ 33 = ε 33 ε c F(λ)+ ε 33 F 33 ε c 33 λ λ τ 23 = γ 23 γ c F(λ)+ γ 23 F 23 γ c 23 λ λ (23) τ 13 = γ 13 γ c F(λ)+ γ 13 F 13 γ c 13 λ λ ε = λ(δ n /2R)+(σ/σ max )(σ max /2E m )(1 +E m /E f ) σ/σ max = 27 4 λ(1 λ)2 (25) Fig. 3 portrays the stress-strain curves for the above solution with E m =3450 MPa, E f =72400 MPa, σ max =60 MPa and different ratios of δ n /R. It can be seen that the stress-strain curves show the snap back instability when δ n /R = or By differentiating Eq. 25, the condition for this instability is, when (σ max /δ n )/(E /R) > 4/9 (26)

6 218 Copyright c 2005 Tech Science Press CMC, vol.2, no.3, pp , 2005 where E = E m E f /(E m +E f ). where η n, η t2, η t1 are viscous coefficients. In Fig. 4, the FEM numerical solution for δ n /R = with η n = 5.0 is compared with the analytical solution. It can be seen that the snap-back instability is avoided and the curve after the complete separation seems physically more realistic. In addition, the influence of the thickness of the cohesive interphase element, h, is also shown by comparing the results with two different thicknesses. It is seen that, the influence of the h is negligible as long as it is sufficiently small. Thus in the analyses of composite laminates to follow, the thickness of the interphase element will be taken to be 2% of the fiber radius. Figure 2 : One dimensional interphase model. Figure 3 : Global stress vs. global strain for a simple interphase model. In the FEM analyses, for the three ratios of δ n /R = 0.1, and 0.05the same global stress-strain curves as that in Fig. 3 are obtained. However, as for the ratios of δ n /R = and 0.01, the calculations can be continued only until the maximum global stress is reached and no converged solution can be obtained thereafter. To overcome this difficulty, when the snap-back occurs, an artificial viscous term is introduced, as suggested in Gao and Bower (2004). Equation 18 is modified as follows T 3 = [u 3] F(λ)+ η n d[u 3 ] δ n δ n dt T 2 = α [u 2] F(λ)+ η t2 d[u 2 ] δ t2 δ t2 dt T 1 = α [u 1] F(λ)+ η t1 d[u 1 ] δ t1 δ t1 dt (27) Figure 4 : Response with artificial viscousity and comparison with the analytical solution 4 Numerical Results of Laminates under Off-axis Loading The off-axis loading of an E-glass fiber/epoxy (Epon 828/Curing agent Z) unidirectional laminate has been analyzed as an application example. This example is chosen because a detailed experimental study [Ishai (1971)] revealed that the bonding between fiber and matrix in this material system was relatively weak and interfacial debonding was the dominant failure mechanisms of the unidirectional laminates under off-axis loading. Figure 5 shows the fractographs of the failed unidirectional laminate specimens under 90, 45 and 20 off-axis loadings, in which those fibers with relatively smooth surfaces indicate complete separations between the fibers and the epoxy matrix.

7 Micromechanical Analysis of Interphase Damage for Fiber Reinforced Composite Laminates 219 form of nonlinear viscoelastic model recently developed in Xia et al. (2003) is used to model the constitutive response of the epoxy matrix, Epon 828/Z. The material constants needed for the viscoelastic model are also listed in Tab. 1. The physical significances of the constants and the details of the viscoelastic model are documented in the above reference. The parameters for the cohesive model used in the calculations are:σ max = 72MPa, δ n /R = δ t1 /R = δ t2 /R=1.0 % and α = 0.9, where R is the radius of the fiber. The cohesive stress-strain relations for the interphase element and the viscoelastic constitutive model of the epoxy resin were implemented into the FEM code ADINA as user-defined material subroutines. The calculations were conducted on a SGI Origin 2000 computer system. To compare with the test results in Ishai (1971), uniaxial tensile loads at four off-axis angles (90, 45, 30, and 20 ) were applied to the unidirectional laminates at a constant strain rate of 10 3 s 1. Figure 5 : Fractographs of coupon specimens under offaxis tensile loadings [from Ishai (1971)]:(a) 90 off-axis loading; (b) 45 off-axis loading; (c) 20 off-axis loading In the current FEM analysis the glass fiber is assumed to be linearly elastic with elastic constants listed in Tab. 1. The epoxy matrix has a highly nonlinear viscoelastic behavior. Therefore, a viscoelastic constitutive relation for the matrix material is required for an accurate prediction of the response of the composite. A differential Figure 6 : Finite element mesh of the RUC model Figure 6 shows the finite element meshes of the RUC model with 8 node brick elements. A thin layer of interphase element around the fiber (h/r= 0.02, h and R are the thickness of the interphase element and the radius of the fiber, respectively) is shown in the inset. Note that in the fiber direction (x 1 direction), only one layer of element is sufficient, since all the stress and strain components are invariant along the fiber direction.

8 220 Copyright c 2005 Tech Science Press CMC, vol.2, no.3, pp , 2005 Table 1 : Constants of constitutive model for matrix Constituents E (MPa) ν α 1 τ 1 R E-glass fiber / / / Matrix E 1 (σ)= e σ MPa 4.1 Prediction of Global Stress-strain Curves Figure 7 shows the predicted global stress-strain curves and the comparison with the test results of Ishai (1971). All the calculations used the same set of material constants as shown in Table 1 and the off-axis angles consideredwere90,45,30,20,asinthetest. It can be seen that the stress-strain curves of the unidirectional laminates with four different angle orientations are well predicted by the present micromechanical analysis. In the calculations for all the four off-axis angles, η n /σ max = η t /σ max = 0.15 were used to avoid the convergence problems caused by the snap-back instabilities. Figure 7 : Global stress-strain curve of a unidirectional laminate under off-axis loadings. 4.2 Microscopic Deformation and Interphase Damage Initiation The microscopic deformation including the interphase damage behavior can be directly obtained from the current analysis. Figures 8(a) and (b) show the deformed RUC and the separation of the interphase for 90 and 45 off-axis loadings, respectively. The deformations of 30 and 20 loadings are similar to that of the 45 loading. As shown in Fig. 8(a), for the 90 (transverse) loading, the maximum separation of the interphase is normal to the interphase with [u 3 ] 0, [u 1 ]=[u 2 ]=0. For 45 off-axis loading, the interphase damage shows a mixed mode: opening normal to the interphase and sliding along the 1 direction, i.e. in this case,[u 3 ] 0, [u 1 ] 0, but [u 2 ]=0 as shown in Fig. 8(b). Note also that for the transverse (θ = 90 ) loading, the deformed surfaces of the RUC remains rectangular, while for 45 off-axis loading, the surfaces originally perpendicular to the fiber do not remain planes after the deformation. Figures 9 to 11 plot distribution of the tractions in the interphase for the four fiber angles at the instant prior to but close to the peak point of the stress-strain curves. For 90, 45, 30, and 20 off-axis angles, the corresponding global strains are ε =0.5%, 0.9%, 1.0%, and 1.0%, respectively. It can be seen from Figs. 9 to 11 that, along the circumferential direction of the interphase (see inset in Fig. 9 for angle φ ), the normal traction T 3 and sliding traction T 1 for all the four off-axis angles have maximum values at φ=0, thus indicates that the interphase debonding will initiate at φ = 0. However, the proportions of tractions T 3 and T 1 are different for different off-axis angles θ. When θ = 90 (transverse loading), the tractions T 2 and T 1 are very small at φ = 0, while T 3 /σ max is close to 1 as indicted by Fig. 9, thus the initiation of the damage is mainly caused by the traction T 3.Whenθ =45, T 3 /σ max 0.6, T 1 /σ max 0.7 (at φ=0 ), thus both tractions T 3 and T 1 contribute to the initiation of the interphase damage, resulting in the interphase separates with

9 Micromechanical Analysis of Interphase Damage for Fiber Reinforced Composite Laminates 221 (a) (b) Figure 8 : Deformed shape of the RUC and interphase damage for different off-axis angles: (a) 90 off-axis loading; (b) 45 off-axis loading Figure 9 : Distribution of T 3 Figure 10 : Distribution of T 2

10 222 Copyright c 2005 Tech Science Press CMC, vol.2, no.3, pp , 2005 (a) (b) Figure 12 : Traction distributions in interphase prior to and after damage: (a) Traction T 3 ; (b) Traction T 2 the combined normal opening and the sliding along the fiber circumference. 4.3 Evolution of interphase damage Figure 11 : Distribution of T 1 The propagation of the interphase damage can also be depicted by the interphase traction distributionsin the RUC. For transverse loading, for example, Fig. 12 shows the distribution of the interphase tractions at different loading levels. At ε= 0.5 %, the traction T 3 reaches the maximum value at φ = 0, thus the interphase starts separating at this point. The traction decreases with the increase of the separation, thus the normalized traction decreases to about 65% at global strain of 0.7 % and the total separation occurs at a global strain of 0.8 %, i.e. the traction drops to zero. With further increase of loading to 0.9 % and 1.1 %, the damaged zone propagates further along the circumferential direction. The interphase damage propagates to an angle of φ = 40 at the global strain of ε = 1.1%.

11 Micromechanical Analysis of Interphase Damage for Fiber Reinforced Composite Laminates Conclusions The interphase damage analysis has been successfully incorporated into the micromechanical finite element analysis for the fibrous composite materials. Both the macroscopic and microscopic responses of unidirectional laminates under off-axis loadings are well predicted based on the properties of the constituents and that of the interphase. In particular, the following conclusions can be drawn from the current investigation: The interphase element developed based on the cohesive law is capable of modeling both normal and tangential separations of the interface. Thus the entire response including damage initiation and evolution of the composite can be simulated. The interphase element is convenient for the implementation in a FEM code since an equivalent stress-strain relationship is used. Appropriate periodic boundary conditions should be applied to the repeated unit cell model of the composite. From the general periodic boundary conditions for an entire RUC model, the periodic boundary conditions for one-quarter of the unit cell accounting for the symmetry conditions are derived in a rigorous manner for the analysis of unidirectional laminates under off-axis loading. For unidirectional laminates under transverse loading, the interphase damage initiates as the normal opening while for the off-axis loadings combined mode of normal opening and in-plane sliding is predicted. Using the same set of material constants, the global stress-strain curves at different off-axis angles are predicted and they are in good agreement with the experimental results. Acknowledgement: The authors would like to express special thanks to Prof. Fernand Ellyin for helpful discussions and suggestions during the course of this investigation. The research is supported by the Natural Sciences and Engineering Research Council of Canada (NSERC) through grant to Zihui Xia. References Aboudi, J.; Pindera, M.-J.; Arnold, S. M. (2001): Linear thermoelastic higher-order theory for periodic multiphase materials. ASME Journal of Applied Mechanics, vol. 68, pp Aghdam, M. M.; Pavier, M. J.; Smith, D. J. (2001): Micro-mechanics of off-axis loading of metal matrix composites using finite element analysis. International Journal of Solids and Structures, vol. 38, pp Daniel,I.M.;Ishai,O.(1994): Engineering Mechanics of Composite Materials. Oxford University Press, Oxford. Gao, Y. F.; Bower, A. F. (2004): A simple technique for avoiding convergence problems in finite element simulations of crack nucleation and growth on cohesive interfaces. Modelling and Simulation in Materials Science and Engineering, vol. 12, pp Hoover, J. W.; Kujawski, D.; Ellyin, F. (1997): Transverse cracking of symmetric and unsymmetric glassfiber/epoxy-resin laminates. Composites Science and Technology, vol. 57, pp Ishai, O. (1971): Failure of unidirectional composites in tension. Journal of the Engineering Mechanics Division, Proceedings of the ASCE, vol. 97, pp Legarth, B. N. (2004): Unit cell debonding analyses for arbitrary orientation of plastic anisotropy. International Journal of Solids and Structures, vol. 41, pp Michel, J. C.; Moulinec, H.; Suquet, P. (1999): Effective properties of composite materials with periodic microstructure: a computational approach. Computer Methods in Applied Mechanics and Engineering, vol. 172, pp Okada, H.; Liu, C. T.; Ninomiya, T.; Fukui, Y.; Kumazawa, N. (2004): Analysis of particulate composite materials using an element overlay technique, CMES: Computer Modeling in Engineering & Sciences, vol. 6, pp Pagano, N. J.; Yuan, F. G. (2000): The significance of effective modulus theory (homogenization) in composite laminate mechanics. Composites Science and Technology, vol. 60, pp Smith, B. W. (1987): Fractography for continuous fiber composites. In: Engineered Materials Handbook, Volume 1, Composites. ASM International, Ohio, USA, pp Sun, C. T.; Vaidya, R. S. (1996): Prediction of composite properties from a representative volume element. Composites Science and Technology, vol. 56, pp Suquet, P. (1987): Elements of homogenization theory for inelastic solid mechanics. In: Sanchez-Palencia, E.;

12 224 Copyright c 2005 Tech Science Press CMC, vol.2, no.3, pp , 2005 Zaoui, A., (eds). Homogenization Techniques for Composite Media. Springer-Verlag, Berlin, pp Tian, J.; Gabbert, U.; Berger, H.; Su, X. (2004): Lamb Wave Interaction with Delaminations in CFRP Laminates, CMC: Computers, Materials, & Continua, vol. 1, pp Tvergaard, V. (1990a): Effect of fiber debonding in whisker-reinforced metal. Material Science and Engineering, vol. A125, pp Tvergaard, V. (1990b): Analysis of tensile properties for a whisker-reinforced metal matrix composite. Acta Metallurgica, vol. 38, pp Xia, Z.; Hu, Y.; Ellyin, F. (2003): Deformation behavior of an epoxy resin subject to multiaxial loadings, part II: constitutive modeling and predictions. Polymer Engineering and Science, vol. 43, pp Xia, Z.; Zhang, Y.; Ellyin, F. (2003): A unified periodical boundary conditions for representative volume elements of composites and applications. International Journal of Solids and Structures, vol. 40, pp Xia, Z.; Zhou, C.; Yong, Q.; Wang, X. (2005): On selection of repeated unit cell model and application of unified periodic boundary conditions in micro-mechanical analysis of composites. International Journal of Solids and Structures, (in press). Yang, Z. J.; Proverbs, D. (2004): A comparative study of numerical solutions to non-linear discrete crack modeling of concrete beams involving sharp snap-back. Engineering Fracture Mechanics, vol. 71, pp Yuan, F. G.; Pagano, N. J.; Cai, X. (1997): Elastic moduli of brittle matrix composites with interfacial debonding. International Journal of Solids and Structures, vol. 34, pp Zhang, Y.; Xia, Z.; Ellyin, F. (2005): Nonlinear viscoelastic micromechanical analysis of fibre reinforced polymer laminates with damage evolution. International Journal of Solids and Structures, vol. 42, pp Zhu, C.; Sun, C. T. (2003): Micromechanical modeling of fiber composites under off-axis loading. Journal of Thermoplastic Composite Materials, vol. 16, pp Zhu, H.; Sankar B.V.; Marney, R.V. (1998): Evaluation of failure criteria for fiber composites using single element micromechanics. Journal of Composite Materials, vol. 32, pp Zhu, H.; Achenbach, J. D. (1991): Radial matrix cracking and interphase failure in transversely loaded fiber composites. Mechanics of Materials, vol. 11, pp Appendix A: Periodic Boundary Conditions for a Quarter RUC under Off-axis Loading For the off-axis loading case shown in Fig. 13, the geometry, material, stress/strain, and periodic part of displacement are invariant along the x 1 direction. From Eq.1, the general form of the displacement field can be written as: u 1 (x 1, x 2, x 3 )=ε 11 x 1 +2ε 12 x 2 +u 1(x 2, x 3 ) u 2 (x 1, x 2, x 3 )=ε 22 x 2 +u 2(x 2, x 3 ) u 3 (x 1, x 2, x 3 )=ε 33 x 3 +u 3(x 2, x 3 ) (A1) (A2) (A3) Note that Eqs. A1 and A2 are obtained through adding a set of rigid body rotation about x 3 axis, u 1 = ε 12 x 2, u 2 = ε 12 x 1,toEq.A1. Appendix A:.1 Derivation of boundary conditions Eqs. 10 and 11 For off-axis loading (Eqs. 5 and 6), the full size RUC showninfig. A1(a)hasamirrorsymmetryofplane x 3 = 0, thus the displacement and stress components of two arbitrary symmetric points satisfy: u 3 (x 1, x 2, x 3 )= u 3 (x 1, x 2, x 3 ) σ 23 (x 1, x 2, x 3 )= σ 23 (x 1, x 2, x 3 ) σ 31 (x 1, x 2, x 3 )= σ 31 (x 1, x 2, x 3 ) And for the other displacement or stress components: u i (x 1, x 2, x 3 )=u i (x 1, x 2, x 3 ) σ ij (x 1, x 2, x 3 )=σ ij (x 1, x 2, x 3 ) (A4) (A5) (A6) Therefore, we have: 1) From Eqs. A4-A6, on plane x 3 = 0, u 3 = 0, σ 23 = σ 31 = 0, i.e. Eq.10. 2) From Eqs. A3 and A4, for an arbitrary point P on plane x 3 = a/2 and its symmetric point Q on x 3 = a/2, u 3 (x 1, x 2, a/2) u 3 (x 1, x 2, a/2)=ε 33 a u 3 (x 1, x 2, a/2)= u 3 (x 1, x 2, a/2) u 3 (x 1, x 2, a/2)=ε 33 (a/2) }

13 Micromechanical Analysis of Interphase Damage for Fiber Reinforced Composite Laminates 225 (a) (b) (c) Figure 13 : Reduce to one-quarter of the RUC under off-axis loading: (a) Unit cell; (b) Half of the unit cell; (c) Quarter of the unit cell i.e. the first equation of Eq.11. 3) From Eqs. 2, A5 and A6, we have, for σ 23 and σ 31, σ 23 (x 1, x 2, a/2)=σ 23 (x 1, x 2, a/2) σ 23 (x 1, x 2, a/2)= σ 23 (x 1, x 2, a/2) σ 23 (x 1, x 2, a/2)=0 σ 31 (x 1, x 2, a/2)=σ 31 (x 1, x 2, a/2) σ 31 (x 1, x 2, a/2)= σ 31 (x 1, x 2, a/2) σ 31 (x 1, x 2, a/2)=0 i.e. the second and third equations of Eq.11. Therefore the RUC reduced to one half of the original full size, Fig. A1 (b). Appendix A:.2 Derivation of boundary conditions Eqs. 8 and 9 TheproblemshowninFig.A1(a)hasalsoaπ -rotation symmetry about x 3 axis, thus the displacement and stress components of two arbitrary symmetric points meet: u 1 (x 1, x 2, x 3 )= u 1 ( x 1, x 2, x 3 ) u 2 (x 1, x 2, x 3 )= u 2 ( x 1, x 2, x 3 ) σ 23 (x 1, x 2, x 3 )= σ 23 ( x 1, x 2, x 3 ) σ 31 (x 1, x 2, x 3 )= σ 31 ( x 1, x 2, x 3 ) And for the other displacement or stress components: u i (x 1, x 2, x 3 )=u i ( x 1, x 2, x 3 ) } } (A7) (A8) (A9) (A10) σ ij (x 1, x 2, x 3 )=σ ij ( x 1, x 2, x 3 ) 1) For two arbitrary π -rotation symmetric points A and B on plane x 2 = 0: From Eq. A1, and note u 1 are independent of x 1, thus u 1 (x 1, 0, x 3 ) u 1 ( x 1, 0, x 3 ) = ε 11 x 1 ε 11 ( x 1 )+u 1(x 1, 0, x 3 ) u 1( x 1, 0, x 3 ) = 2ε 11 x 1 From Eq. A7 u 1 (x 1, 0, x 3 )= u 1 ( x 1, 0, x 3 ) Thus from the above two relations, we obtain u 1 (x 1, 0, x 3 )=ε 11 x 1 i.e. the first equation of Eq.8. Similarly from Eqs. A2 and A8 we can obtain: u 2 (x 1, 0, x 3 ) u 2 ( x 1, 0, x 3 )=0 u 2 (x 1, 0, x 3 )= u 2 ( x 1, 0, x 3 ) Thus u 2 (x 1, 0, x 3 )=0, i.e. the second equation of Eq. 8. From Eq. A9, σ 23 (x 1, 0, x 3 )= σ 23 ( x 1, 0, x 3 ), but σ 23 is also independent of x 1, thus σ 23 (x 1, 0, x 3 )=0

14 226 Copyright c 2005 Tech Science Press CMC, vol.2, no.3, pp , 2005 i.e. the last equation of Eq. 8. 2) For two arbitrary π -rotation symmetric points M and N on planes x 2 = ±a/2: From Eq. A1, note u 1 is independent of x 1 and u 1 is periodic, thus u 1 (x 1, a/2, x 3 ) u 1 ( x 1, a/2, x 3 ) = ε 11 x 1 ε 11 ( x 1 )+2ε 12 (a/2) 2ε 12 ( a/2) +u 1(x 1, a/2, x 3 ) u 1( x 1, a/2, x 3 ) = 2ε 11 x 1 +2ε 12 a +u 1(x 1, a/2, x 3 ) u 1( x 1, a/2, x 3 ) = 2ε 11 x 1 +2ε 12 a From Eq. A7, u 1 (x 1, a/2, x 3 )= u 1 ( x 1, a/2, x 3 ) Thus from the above two relations we have, u 1 (x 1, a/2, x 3 )=ε 11 x 1 +2ε 12 (a/2) i.e. the first equation of Eq.9. Similarly, from Eqs. A2 and A8, the second equation of Eq.9 can be obtained. Finally from Eq. A9 the last equation of Eq. 9 can be obtained.