Numerical Modelling of Damage Evolution and Failure Behavior of Continuous Fiber Reinforced Composites

© 2012 Wang and Zhang, licensee InTech. This is an open access chapter distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited. Numerical Modelling of Damage Evolution and Failure Behavior of Continuous Fiber Reinforced Composites


Introduction
Demands on high performance materials for use in a spectrum of structural and nonstructural applications are increasing to the point that monolithic materials cannot fully satisfy multifunctional requirements. One approach that has emerged is to develop advanced fiber composites whose properties are tailored for wide applications in engineering, such as turbine generators and ultrasonic aircrafts (Miracle, 2005). Composite materials are composed of reinforcements and matrices, in which mechanical behaviors can be seriously affected by their micro-structures. Therefore, it becomes necessary to deal with these materials from a micromechanical view point. Most importantly, an urgent task in the field of fiber-reinforced composites is to develop a relationship between the material structure and its performance under more severe loads and environments.
Damage and failure in fiber-reinforced composites evolves at several different length scales. At the smallest scale, pre-existing defects in fibers grow. Due to statistical distribution of such defect, the fiber strength exhibits large variability. Therefore, as increasing load is applied to the composite, the weakest fiber will first break. The loads carried by the broken fiber are redistributed among the remaining unbroken fibers and matrix as determined by the constitutive response of the fibers, matrix and interface. And this then causes other fibers near the failure site to fail and thus shed further load to intact fibers. Consequently, the failure of fibrous composites goes through a very complicated damage evolution, which is a combination of fiber fracture, matrix deformation and interfacial debonding and slipping around the fiber breaks, before it reaches ultimate failure. It is obvious that the connection between the microstructural scale and the macroscopic scale is nontrivial and involves mechanics, stochastics, and volume scaling (Curtin, 1999).
Local damage evolution at the scale of fiber diameter in composites under loading determines its fracture toughness, strength, and eventually lifetime. The mechanical behavior of such a composite depends on the evolution of multi-damages throughout the application of loading, and then modelling composites undergoing progressive damage becomes a complex procedure due to many mechanisms involved (Mishnaevsky & Brøndsted, 2008). The development of a micromechanical damage model has to take different aspects into account such as (i) the proposition of a micromechanics-based approach that describes the influence of the damage variable on material properties (Blassiau et al., 2008), (ii) the definition of a pertinent damage variable and its law of evolution (Kruch et al., 2006), and (iii) the use of an appropriate and efficient experimental technique for the evaluation of damage, such as neutron diffraction (Hanan et al., 2005) and acoustic emission (Bussiba et al., 2008), etc.
Reliability concerns in utilizing fibrous composites in structural application have motivated the development of many numerical and analytical failure models in the presence of multidamages. Different from the phenomenological approaches based upon the macroscopic level, the progressive model is needed to consider local damage mechanisms, such as fiber breakage, matrix deformation, interfacial debonding, etc (Kabir et al., 2006) and predict the dominant failure modes. This method seems to be more accurate but computationally complicated because it accounts for many failure mechanisms and is also related to damage accumulation correlated with material properties degradation. Recently computational micromechanics is also emerging as an accurate tool to study the mechanical behaviors of composites because of the sophistication of the modeling tools and the everincreasing power of digital computers (González, 2004). Within the framework, the macroscopic properties of a composite can be obtained by means of numerical simulation of the deformation and failure of the microstructure ).
In our previous work, the local cyclic shear plasticity of the interface around a broken fiber in ductile matrix composites under the in-phase and out-of-phase thermo-mechanical fatigue loads was analyzed by using the single-fiber shear-lag model (Zhang et al., 2002). A multifiber shear-lag model including matrix tensile modulus based upon an influence function superimposition technique was developed to simulate the nonlinear stress-strain response and the progressive failure of continuous fiber reinforced metal matrix composites under static tensile loading (Zhang & Wang, 2009) as well as thermomechanical fatigue loading (Zhang & Wang, 2010). This chapter will summarize our work and be organized as follows: an analytical model of the fibrous composite will be presented in Section II, in which fiber strength statistics, matrix behavior, and interfacial mechanics are explained in details. In Section III, we will introduce an influence function superimposition technique to derive stress profiles for any configuration of breaks, by considering local matrix plasticity and interface yield. In Section IV, numerical models combined with Monte-Carlo method will be developed to simulate progressive damage. In Section V, we will investigate failure behaviors of continuous fiber reinforced composites under cyclically thermomechanical loading. Finally, we will discuss limitations of the existing models, and aspects of the existing theories that require improvement.

Analytical model for multiple damage events
In the real multifiber composites of practical interest, the evolution of the fiber fragmentation during loading is, in principle, different because each individual fiber experiences a nonuniform stress due to the uniform applied stress plus stresses transferred from other broken fibers in the composite. The evolution of fiber damage thus depends crucially on the nature of the load transfer from broken or slipping fibers to unbroken fibers.
Consider an infinite, two-dimensional (2D) unidirectional composites reinforced with parallel, evenly-spaced fibers embedded in matrix material, in which the fibers and matrix regions are numbered in the serial number, shown in Fig. 1. The lamina is loaded in simple tension along the fiber direction. The width of fibers is D, and the width of each matrix region W can be related to the fiber volume fraction, Vf (Zhang et al., 2002). In many typical applications, fibrous composites are subjected to a cyclic mechanical loading with a superimposed variation in temperature. This type loading condition is referred to as thermomechanical fatigue (TMF) and can be regarded as one of the most severe types because the significant difference in the coefficient of the thermal expansion between the matrix and the fiber causes high thermal stresses and stress amplitudes raising irreversible deformation in the composites (Mall & Schubbe, 1994). In-phase and out-of-phase thermomechanical fatigue (TMF) loads are illustrated in Fig. 2. The tensile stress applied to composites varies between (σ∞) min and (σ∞) max while the temperature changes between Tmax and Tmin. For convenience, the state 'A' is used to represent the state of TMF loads with the maximum applied stress, and the state 'B' denotes the state of TMF load with the minimum mechanical load. It is clearly illustrated that in-phase conditions subject the composite to high stresses at hot temperatures and low stresses at low temperatures. Conversely, out-ofphase TMF conditions subject the composite to high stresses at low temperatures and vice versa since the peaks in the waveforms are 180 º apart (Williams & Pindera, 1995). Since the fiber strength exhibits large variability due to statistic distribution of defect, the weakest fiber break usually takes place at very early loading stage. After some loading level there are N fiber breaks occurred in the composite. The coordinates of the n-th fiber break, which locates in the kn-th fiber with vertical position n, are given by (kn, n ) (n = 1, 2,…, N).
If the interfacial strength is higher than the matrix yield stress in shear, the plastic deformation of the matrix will occur before fiber/matrix debonding (Beyerlein & Phoenix, 1996). This was observed in ductile resin and some metal matrix composites with a strong interface. Some experiments show that the microdamage and deformation modes include the fiber breakage, fiber pullout, debonding, and local plasticity around fiber breaks (Liu & He, 2001).
In order to analyze the complex combination of microdamage and deformation, we may propose a micromechanical model, shown in Fig. 3., in which a broken fiber, accompanying with its (reserve) tensile yielding matrix and its (reverse) shear yielding interface, and its debonding interface is called as a damage-plasticity event. Essentially the nonlinear cyclic behavior and thermomechanical fatigue failure of the composite is the result of the interactions among these damage-plasticity events under TMF loads.   In the loading phase the local interface shear yield and the matrix tensile yield around the fiber break could take place due to stress concentration as the load increases. The initial interface yield length caused during the loading phase is represented by, 2Ln, associated with the n-th fiber break. For fatigue case, when the load is released from state 'A' to state 'B' in the unloading phase, the yielded interface and the yielded matrix around the fiber breaks start the elastic unloading. The interface shear stress at broken sites relaxes faster than that at any other positions. Thus, the shear stress of broken sites will change sign. Subsequently, there are two possible cases: the complete elastic unloading and the reverse shear plastic yielding within a certain length, depending on the load range. If the load range is not such large that the shear stress of the broken site does not reach the reverse yield stress before state 'B', the unloading is completely elastic. In this case the unloadingreloading will be totally elastic, that is, shakedown occurs. If the load range is sufficiently large, the interface shear stress around the broken sites will change sign and will reach the reverse yield stress at some moment before the load reached state 'B'. The interface plasticity is modeled by the elastic, perfectly plastic shear stress-strain relation, ignoring the effect of the interface tensile stress, that is, where Gc and τs are the shear modulus of the interface and the yield shear stress, respectively; τ and γ denote the interfacial shear stress and shear strain, respectively. This interface constitutive law allows the interface to undergo the cyclic plasticity under TMF loading. The yield shear stress is related to the temperature by where c1, c2 and c3 are material constants.
With continuous mechanical unloading of composite, the reverse plastic yield of the interface will extend forward from the fiber break tip. Up to the loading moment 'B', the reverse plastic yield has occurred within a certain length, 2L'n. With the repeated loading and unloading processes, the cyclic interface shear plasticity takes place within the length, 2L'n, leading to debonding. In other words, the interface is debonded in the length of 2Ld because of cyclic plastic shear strain accumulation. It must be pointed out that material parameters in cyclic interface shear plasticity and debonding growth are difficult to measure directly. For the strong interface considered here, it is assumed that the interface shear yield is governed by the shear yield of matrix in a very thin layer. In other words, the material parameters of the interface are assumed to the same as those of the matrix.
Matrix tensile yield in the matrix regions neighboring with the fiber break sites could take place due to matrix tensile stress transferred from the tensile load released from the broken fiber. It is assumed that the M breaks of the N fiber breaks (M ≤ N) cause the associated matrix tensile yielding. Each matrix tensile yield zone takes the shape of a narrow strip spanning the matrix width, and is colinear with the associated fiber break. The constitutive behavior of matrix in the yield zones obeys the uniaxial cyclic elastic, perfectly plastic stressstrain relation. The yield tensile stress σs is related to the yield shear stress τs by von Mises yield criterion, that is,

Stress profiles for multiple damage events
Accurate prediction of the mechanical properties of the composite materials requires detailed knowledge of the micromechanical stress state in and around broken fibers as a function of the constituent material properties (Xia et al., 2002). Stress transfer between fiber and matrix is one of fundamental issues for a composite system (Hanan et al., 2003;Beyerlein & Landis, 1999). In order to obtain solutions for deformation in the composites containing multiple damage events, a theoretical model under the framework of 2D shearlag arguments will be used to derive stress profiles for any configuration of breaks in the presence of matrix tensile yield and interface shear yield. It is assumed that the tensile load in a fiber-reinforced composite loaded in the fiber direction is mainly carried by the fibers and matrix regions, where the load between fiber and matrix is transferred through the fiber/matrix interface shear stress. The transverse deformation in the fibers and matrix regions is neglected. The interface is assumed to be in as state of shear behavior as controlled by the axial displacements of the fiber and matrix materials. Shear-lag analysis of stresses and deformations, involving the cyclic thermoplasticity, must be carefully carried out by applying TMF loads to the composites incrementally.

Stress solution in loading phase
Let us consider the TMF loading phase. As the thermomechanical loads increase, N fiber breaks have appeared where each fiber break is accompanied by the interface yielding of length, 2Ln, and where M breaks (M ≤ N) cause the associated matrix tensile yielding simultaneously. The interfacial shear stress within the yield zone 0≤| -n|≤ Ln/lc associated with the n-th fiber-break has to be equal to the yield shear stress, it follows   n n c , , f o r / a n d 1 ,2,3, , , where τ (kn, ) is the interface shear stress at the location, (kn, ). Since the shear stress is constant (τs) within the interface yield region, the fiber stress can be obtained simply from the equilibrium condition of fiber, that is, , for 0 | | / and 1,2,3, , , where σ f (kn, ) is fiber stress at the location, (kn, ). Matrix tensile stress in the narrow matrix yield zone is In order to obtain solutions for fiber stress and stress in matrix outside the yield zones the superposition method (Beyerlein & Phoenix, 1996;Landis et al., 2000) has been further extended in our previous work (Zhang & Wang, 2009, 2010. Original superposition method of Beyerlein excludes the axial load-bearing capacity of matrix, which is suitable for polymer matrix composites. Our extended superposition method includes the axial load-bearing capacity of matrix, which is applicable to fiber reinforced metal matrix composites. The fundamental concepts behind influence superimposition techniques are to obtain analytical solutions for the influence of a 'unit' failure on the stress and displacement fields and then use a weighted superposition technique to obtain the solution for the full problem. Indeed, the units solutions for a single fiber-break and a single matrix-break, called the influence functions, are used as building blocks for the solutions for fiber stresses and matrix stresses in the composites. In the method, the narrow matrix yield zones is modeled as fictitious matrix breaks with surface tensile traction σs, and the interface yield length is considered as a fictitious debonding crack with surface shear traction τs. The M fictitious matrix breaks (matrix yield zones) can be decomposed into the M problems of single matrix-break with opening displacements, m t w (t = 1,2,3,…,M), that are yet to be determined by stress condition, eq.(5). For an interface yield length 0≤| -n|≤ Ln/lc associated with the n-th fiber break, the fiber segment with the known fiber stress given by eq. (4) can be modeled as a continuous distribution of fictitious fiber break with non-zero opening displacement   f n w  which is to be determined through the stress conditions, eqs. (3-4).
Therefore, the axial stresses in fibers and matrix regions and the interfacial shear stress, σ f (k, ), σ m (k, ) and τ (k, ) at any position (k, ) can be obtained by using the superposition method as follows, for 0 | | / and 1,2,3, , , ; Em are the Young's modulus of the matrix and K is the shear-lag parameter associated with fiber/matrix interfaces. By solving the integral eqs. It is worthy to emphasize that the stress analysis method mentioned above is exact under the framework of shear-lag argument. In that sense, moreover, the mechanical interactions between the multiple damage-plasticity events have been fully and exactly taken into account since the governing equations of shear lag model, the boundary conditions and the plastic yield conditions have been fully satisfied.

Stress solution in unloading phase
There are two possible cases. If the unloading is completely elastic, the solutions for the stresses and the displacements at the state 'B' can be attained by superposing the elastic solutions between A and B onto the solutions of the state 'A'. Otherwise, if the load range is sufficiently large, the local reverse yield of the interface and matrix will takes place before the thermomechanical load reaches state 'B'. We consider an unloading increment, during which the interface reverse yielding occurs in N ' (N ' ≤ N) lengths, 2L'n (n= 1, 2,…, N ' ), and matrix reverse yielding takes place in M ' (M ' ≤ M) narrow matrix zones. The interfacial shear stress within the reverse yield zone 0≤| -n|≤ L'n/lc has to be equal to the reverse yield shear stress, i.e.
  n n c , ,f o r0 | | ' / a n d 1 ,2 ,3 , , ' , Within the interface reverse yield region the fiber stress can be easily obtained from the equilibrium condition of fiber as follows     n n c 2 , , for 0 | | ' / and 1,2,3, , ', Matrix stress in the narrow matrix reverse yield zone is The incremental stress solutions for fibers and matrix regions in composite for the thermomechanical unloading increment can be obtained by using the superposition method mentioned in the previous section. The solutions for the total stresses after the thermomechanical unloading increment, then, can be obtained by superposing the incremental solutions onto the solutions before the unloading increment, that is, where the coefficients f j w (j = 1, 2,…, N') and m i w (i = 1, 2,…, M') are to be determined which belong to the incremental unloading. The subscript '0' represents the state before an incremental unloading. The undetermined interface reverse yield lengths 2L'n are governed by the interfacial shear stress at the end of reverse yield zone, n  = n + L'n /lc, and it follows from eqs. (12) and (17)

Local cyclic plasticity and debonding
The cyclic shear plasticity of the interface will take place only if L'n >0. Otherwise there is no alternating plastic shear and interface undergoes cyclic elastic deformation over the whole length after initial yielding of length Ln, that is, shakedown. Since the elastic, perfect-plastic constitutive relation has been assumed for the interface shear and matrix, the solutions for all cycles are the same before the new debonding takes place. The alternating plastic shear strain range can be derived from displacement solutions for loading and unloading condition. The total interface shear strain is evaluated by eq. (21a) where u f is the axial displacement of broken fiber and u m is the axial displacement of matrix adjacent to the broken fiber. The plastic shear strain can be easily computed by eq. (21b).
Debonding is an important mechanism that stimulates stress redistribution in composites.
Since plastic strain is dominant in the low cyclic fatigue regime (Llorca, 2002), the cyclic plastic shear strain range Δγp may be used to predict the debonding (Zhang et al., 2002) along with Coffin-Manson fatigue equations Δγp/2 = γf •(2Nf) c , where Nf is the cycles to interface debonding under the constant shear plastic strain range Δγp, material constants γf and c stand for the fatigue ductility coefficient and the fatigue ductility exponent, respectively. The debonding growth law based on Coffin-Mason equations represented by plastic shear strain range provides a means to account for the stable growth of debonding as the number of cycles increases, which may be relevant for the strong interface and/or the low applied load where the static debonding criterion fails to predict a further growth of debonding. The parameters in the debonding growth law are taken to the same as matrix since we assumed that the interface shear yield is governed by the shear yield of matrix in a very thin layer. This criterion for debonding based upon the shear-lag model ignores the effect of the interface tensile stress. This simplification is reasonable for the growth-dominated fatigue failure since the shear stress contributes to its growth and propagation along the fiber length. However, the interface tensile stress component could be important for the cause of the initiation of fiber/matrix debonding. Therefore, this simplification would predict a delayed debonding, leading to an overestimated fatigue lifetime for the initiation-dominated fatigue failure (high load). Although the composite stress range and the temperature range for both the in-phase and out-of-phase TMF loading conditions do not change with cycles, the local plastic shear strain ranges will change, due to new debonding, after some number of cycles. Therefore, one needs a damage accumulation fatigue rule to describe the debonding due to local varying amplitude fatigue. For simplicity, we assume that the fatigue life for interface debonding is governed by the linear damage accumulation rule ∑ Ni/(Nf)i = 1, where Ni is the number of the applied loading cycle leading to the constant plastic strain range (Δγp)i and (Nf)i is the cycles to interface debonding for the constant shear plastic strain range (Δγp)i. When the sum of the fractions from each step equals one, debonding is predicted (Liu, 2001) The stresses in composite with debonding can be obtained again by using the developed superposition method along with eqs. (22-25).

Statistical modelling
The catastrophic failure of fiber-reinforced composites is primarily dominated by the failure of fibers (Talrejia, 1995). The fibers typically exhibit variability in strength due to microflaws distributed randomly along the length. The variation of the fiber strength σf for length L can be characterized by the probability density function, F(σf), which is assumed to follow a two-parameter Weibull distribution, that is, where σ0 represents the scale parameter, or characteristic quantity of the material. L0 is the reference length when the characteristic fiber strength σ0 is measured. The Weibull modulus or shape parameter β controls the scatter of the fiber tensile strength in the distribution, experimentally found to describe a variety of materials. This scatter will become large with decreasing β. The parameters β and σ0 can be calculated by statistical method where E(σc) and D(σc) are the mean and variance of random variable, respectively.
In the approach we use here, we assume the total length of composite specimen to be modeled, 2LT, is divided into NL segments of equal length Δx such that Δx = 2LT/NL. The composite specimen to be modeled contains Nf fibers. For a given failure probability F the strength of the fiber segments with length Δx can be derived from the inversion of eq. (26) in which the failure probability F is a random number taken from the uniform distribution in the range [0,1], then the strength of each fiber segment can be obtained. By introducing the normalized fiber segment length L =Δx /D, eq. (29) can be written as where σ* acts as the scale parameter and σ* = σ0•(L0/D) 1/β .
The Monte-Carlo simulation technique coupled with the proposed analytical model is executed to simulate the mechanical failure process in fiber-reinforced composites under TMF loads. At the beginning of each simulation, all fiber segments are assumed to be intact. Their status changes from intact to fractured if the fiber stresses satisfy the failure criteria.
To employ the overall method described above, the algorithm is follows, a. Randomly assign each fiber element a strength according to a two-parameter Weibull distribution. b. From the displacement and the boundary condition, the stress for each segment of fiber, matrix, and interface are then obtained by using the stress analysis. c. Determine whether the fiber elements will break up or not, whether the matrix and interface around broken fibers will yield or not for each incremental load. If no new damage occurs at this loading step, we calculate the composite stresses. Otherwise, we recalculate the stress field by taking account of the new breakage, and repeat this step until no more damage occurs at the load level. d. Taking unloading. For the case except of complete elastic unloading, the cyclic plastic strain range along with Coffin-Manson equation is be used to predict the debonding e. Increase a new loading by a small increment and repeat step (b) and (c) and (d), until the stress-strain curve up to failure for the composite is obtained.
For each loading step at which the applied composite strain is given the overall composite stress can be evaluated by (1  ) , where Vf denotes the fiber volume fraction; the symbols f  and m  are used for the average fiber and matrix stress, respectively.
The average fiber stress is computed by and the average matrix stress is evaluated by The fatigue failure of the fiber reinforced ductile composites occurs as a result of accumulation of the large amounts of damage-plasticity events under cyclically thermomechanical loading.
As much more fiber breaks and the associated local thermoplasticity are accumulated, the composite as a whole will be unable to carry additional load and fail will ensue.

Predictions and concluding remarks
For illustrations, the Boron/Al continuous fiber reinforced composite is examined. The properties of the constituents are given below Ef = 400 GPa, Em = 70.2 GPa, αf = 6.3 µε/ º C, αm = 23.9 µε/ º C, D = 0.14 mm, Vf = 0.48, c = -0.65, γf = 0.42 Fig. 4. shows a simulated failure process of the fiber-reinforced composite under in-phase condition. In this figure, 4a-f indicate the damage configurations under some instantaneous stage, respectively. Because of large variability for fiber strength, the fiber element with the lowest strength is firstly broken at the early loading stage, see Fig. 4a. The high stress concentrations are generated in the matrix or the interface due to the fiber breakages, so that plastic yield appears for matrix or interface. As the applied load increases, more fiber breakages occur in the whole specimen, see Fig. 4b-c. Stress redistribution in the composite is caused by debonding after some cycles, shown in Figure 4d. Accumulation of severe damages can be observed until the composite completely fails, in Fig. 4e-f.   Fig. 4. Plastic yield onset in unloading and reloading is represented by the '+' symbols. Due to plastic deformation, the curve deviates from linear behavior beyond the elastic range at unloading and reloading state. It is distinctly seen that the composite strain accumulates with cycles under the constant stress amplitude. This attributes to the microdamage mechanism that a dominant 'critical cluster' of breaks is shaped, which leads to the failure of the rest of fibers.  In order to show accumulation of damage-plasticity in the composite, Figs. 6-8. illustrate the evolution of fiber break density, matrix yield density, and the total debonding length with TMF cycles, respectively. It is obvious that the deformation behavior is dependent on the load amplitude: the degree of damage exhibits a more aggravation with cycle times under higher amplitude. Fig. 9. shows the evolution of fiber cracking in the composites during reloading. Initial breaking takes place at about 750 MPa in the first cycles. The density of broken fiber increases as the number of load cycles increases. It should be emphasized that this is due to the stress redistribution induced by progressive debonding.  Fig. 10. is a plot of fatigue life by logarithm coordinate against maximum stress. The predicted S-N curves for the in-phase and out-of-phase TMF loads with a temperature range from 250 to 350 °C are plotted along with the experimental S-N curves tested by others (Nicholas, 1995;Liu & He, 2001).  In high stress level, the fatigue life for in-phase TMF conditions is considerably shorter than the life for the out-of-phase TMF conditions. On the other hand, in low stress level, the fatigue life for in-phase conditions is seen to be longer considerably than the life for out-ofphase conditions. This cross-over behavior is due to the difference of the dominated failure mode for this two TMF conditions. A micro-mechanism is presented in Figs. 11-12. Micro failure mechanism associated with the in-phase TMF is characterized by a fiber-dominated fracture with evidence of smaller matrix yielding. Evidence of larger matrix plasticity and less fiber fracture, suggestive of a matrix-dominated failure under the out-of-phase TMF, has been observed in the model simulation.

Conclusions
Modelling of the mechanical behavior of fibrous composite under progressive damages attracts much attention for a long time. The methodology developed in this paper provides an engineering tool to investigate damage evolution in fiber-reinforced composites, especially under cyclical thermomechanical loading. It has to be point out that we do not intend to make the quantitative comparison of the model with the experimental data. The predicted lifetime of fibrous composite is substantially higher than the measured, especially for the regime of high load level. There are several factors that may cause the quantitative discrepancy between the predicted and measured fatigue life. First, the model ignored the interface normal stress that dominates the initiation of debonding, leading to an estimate of delayed debonding and a longer lifetime of the composite. Second, the interface properties have been assumed to be the same as those of matrix; however, interfaces may contain inherent flaws or defects. This assumption will result in a longer lifetime of debonding growth.
As a final remark, it should be mentioned that the 2D modeling presented here does generalize the governing equations to include interactions with multiple damage events and develop the basic mechanics that are necessary to understanding failure mode of the composite produced by the failure of one or more of its components from a micromechanics perspective. The difference between predicted values and measured ones suggests, however, that consideration of the effect of variations in fiber strengths alone is not sufficient for predicting the variability of composite strength. A more precise calibration of the model, capable of explaining such effects as 3D fiber arrays and fiber/matrix interface sliding on composite strength, is the subject of our further study.