Predicting tensile behaviors of short flax fiber-reinforced polymer–matrix composites using a modified shear-lag model

Abstract

Natural fiber-reinforced composites are increasingly being used in the industry. The fiber–matrix interfacial properties of the composites are influenced by many factors, including chemical treatment of the natural fiber, type of polymer matrix, composites fabrication method, and process and the service environment of the composites. In this paper, a modified shear-lag model based on a cohesive fiber/matrix interface is proposed and applied to the analysis of the stress–transfer characteristics and the tensile properties of unidirectional short flax fiber-reinforced composites. The model takes into account of the interfacial shear stiffness, bonding strength between fiber end face and matrix, fiber aspect ratio and fiber volume fraction. 3D finite element models of the composites using a cohesive zone method are used to verify the accuracy of the modified shear-lag model. The fiber tensile strength and the composite tensile elastic modulus are significantly influenced by the interfacial shear stiffness, fiber aspect ratio, and fiber volume fraction. The bonding strength between the fiber end face and the matrix only has an effect when the interfacial shear stiffness is low. The predicted results from the modified shear-lag model show good agreement with the finite element analysis and experimental results in the literature. The modified cohesive shear-lag model provides a simple and effective method for analyzing fiber axial stress, shear stress in the fiber/matrix interface, and tensile elastic modulus of the final composite.

Keywords

Fiber-reinforced polymeric composites shear-lag model cohesive interface tensile elastic modulus fiber axial stress

Introduction

Natural fibers, such as flax, hemp, kenaf, and bamboo are being widely used as reinforcements for composite materials. Compared with traditional synthetic fibers, such as carbon and glass fibers, natural fibers are abundant, lightweight, biodegradable, recyclable, and available at a low price. In addition, the specific mechanical properties of natural fibers are comparable to those of glass fibers used as reinforcement materials.^1–3 This has led to increasing industrial application of natural fiber composites.^4–6 Short natural fiber-reinforced composites (SNFRCs) have been used in nonstructural parts in the automobile industry primarily because they can be easily fabricated by the rapid, low-cost injection-molding process.^7,8

There are a number of factors influencing the mechanical properties of SNFRCs, including the distributions of fiber length and orientation in the final composites, the intrinsic mechanical properties of the fiber and the matrix, and the fiber/matrix interfacial properties. The fiber/matrix interface plays an important role in the load transfer between the fiber and the matrix, which in turn determines the mechanical properties of fiber-reinforced composites, such as elastic modulus, tensile strength, and fracture toughness.^9–11 In comparison to glass fiber and carbon fiber composites, the interfacial shear properties between the natural fiber and matrix are generally poor and can easily deteriorate in humid environment.¹² Efficient and accurate models that predict the effect of fiber/matrix interface properties on the mechanical behavior of SNFRCs are useful to guide research aimed for improvement of natural fiber/matrix interface.

The classical shear-lag model (SLM) developed by Cox was the first micromechanics model for predicting the elastic properties of short-fiber composites based on computation of the stress distribution in the fiber, matrix, and fiber/matrix interface.¹³ Since then, modifications to the classical SLM have been made to improve prediction results for specific uses.^14–21 Nardone and Prewo¹⁹ took into account the load transferred from the matrix to the end faces of the fiber that had been ignored in Cox’s original SLM,¹³ in which fiber was considered to be either continuous or very long relative to diameter. A modified SLM developed by Ji and Zhao¹⁸ incorporated the influences of fiber shape on the tensile stress distribution in the fiber. Zhang and Qiu²¹ developed a modified SLM for fibers with irregular cross-sectional shapes to provide accurate calculations for interfacial shear stress through introducing an equivalent diameter of fiber and a dimensionless shear-lag constant in the classical SLM. Virk et al.²² introduced a fiber area correction factor (FACF) to transfer fiber properties from values based on apparent cross-sectional area (CSA) to values based on the statistically corrected true CSA for the fiber batch in the generalized rule of mixtures equation. They achieved accurate prediction of the modulus and strength of SNFRCs. McCartney¹⁶ introduced interfacial friction in SLM in the analysis of interfacial debonding during fiber pullout or pushout.

Finite element (FE) modeling is now widely used in the simulation of fiber-reinforced composites thanks to recent advancement in numerical simulation and computer software. A number of research groups have proposed a number of methods to characterize the fiber/matrix interface for use in FE modeling of fiber-reinforced composites, including the cohesive zone method (CZM),^23–25 the perfect bonding method,^26,27 and the frictional contact method.⁹ The FE models developed by Modniks and Andersons²⁶ and Modniks et al.²⁷ to investigate the influence of fiber volume fraction on the mechanical properties of flax fiber-reinforced polymer matrix composites assumed a perfect fiber/matrix adhesion, i.e., the interface will never fail when the composite is subjected to external loading. Guillebaud-Bonnafous et al.⁹ developed a 2D FE model of fragmentation test based on two different contact interactions in the interface, a perfect bonding and a frictional contact, and used the model to investigate the adhesion properties in hemp yarn–matrix composites. Modniks and Andersons²³ simulated the nonlinear deformation of a short flax fiber-reinforced composite accounting for the effect of the interfacial shear strength using the CZM. Shokrian et al.²⁵ established 3D FE models of a short glass fiber-reinforced composite using CZM to study the longitudinal and transverse tensile behavior of the composite material with different interfacial shear strength. Chen and Yan²⁴ introduced cohesive fiber/matrix interface into the classical SLM to analyze the stress distribution of fiber pullout and developed FE models based on CZM to simulate fiber pullout and to verify the accuracy of the refined SLM. Insofar, CZM-based FE modeling is the most effective and convenient approach to numerical simulation of the cohesive behavior of fiber/matrix interface.

As discussed above, most investigations on the micromechanical properties of short fiber–reinforced composites (SFRCs) using the SLM are based on the perfect bonding between the fiber and matrix.^14–21 However, the tensile modulus of the SNFRCs can be overestimated using SLM based on the assumption of perfect bonding for the poor adhesion between hydrophilic natural fibers and hydrophobic polymer matrix.²⁸ No study of the stress distribution of short fiber–reinforced composite systems using SLM with consideration of the cohesive behaviors of fiber/matrix interface has been identified, except for the work on fiber pullout from single fiber composites by Chen and Yan.²⁴ In this work, we used the experimental load–displacement curves microdroplet tests^12,29 to construct cohesion behaviors of the natural fiber/matrix interface and establish a modified SLM to analyze the elastic properties of the SNFRCs. Using this modified SLM, we investigated the effects of cohesive interfacial shear stiffness, elastic properties of the fiber and the matrix, fiber aspect ratio, and average fiber axial stress at the embedded end face of fiber on the stress–transfer characteristics and axial tensile properties of the composite. The axial stress in the fiber and the tensile elastic modulus of the composite are compared to published experimental data.

A modified SLM with a cohesive fiber/matrix interface

A cylindrical fiber was always adopted in the modeling of flax fiber–reinforced composites to simplify the analysis of stress transferring at the fiber–matirx interface and achieved remarkably accurate predictions of experimental data.^23,26,27 Figure 1 shows a general unit cell consisting of a cylindrical discontinuous fiber of length $2 l$ and diameter $2 r_{0}$ completely embedded in the cylindrical continuous matrix of length $2 L$ and diameter $2 R$ , and a cohesive interface. This unit cell is adopted in the modified SLM to analyze the load transfer from the matrix to the fiber. A uniform axial tensile stress $σ_{c}$ is applied at the ends of the unit cell along the z-direction. The embedded end face ( $z = 0$ and $z = 2 l$ ) between the fiber and matrix can be perfectly bonded or completely debonded. In the modified SLM, the fiber and matrix are assumed to be elastic and the interface ( $r = r_{0}$ ) between the fiber and matrix is modeled as a cohesive surface in pure shear mode. Besides, the outer surface of matrix ( $r = R$ ) is stress free.

Figure 1.

Schematic diagram of a unit cell with a cohesive fiber/matrix interface.

Experimental pullout curves^12,29 can be used as basis to construct a cohesive traction–separation function that describes the shear behavior of the cohesive interface, as shown in Figure 2. The cohesive traction–separation constitutive relation is expressed in equation (1).

\begin{matrix} τ = {\begin{matrix} K_{0} δ (0 \leq δ \leq δ_{0}) \\ τ_{1} (δ > δ_{0}) \end{matrix} \end{matrix}

(1)

where K₀ is the initial shear stiffness of the cohesive fiber/matrix interface, δ is the tangential displacement of the cohesive fiber/matrix interface.

Figure 2.

Cohesive traction–separation model: (a) microdroplet pullout test,²⁹ (b) typical pullout force–displacement curve,³⁰ and (c) cohesive traction–separation constitutive relation.

As the modified SLM with a cohesive fibre/matrix interface in this study is based on the classical SLM¹³ and the modified SLM,^14–17,24 the basic hypotheses inherent in a shear-lag analysis can also be invoked, except for the discontinuous displacement between the fiber and the matrix at the cohesive interface.

In our modified SLM, the material properties of the fiber and the matrix are assumed to be isotropic elastic and are described by Hooke’s laws:¹⁷

ɛ_{zz} = \frac{σ_{zz}}{E} - \frac{ν}{E} (σ_{rr} + σ_{θ θ}) and γ_{rz} = \frac{τ_{rz}}{G}

(2)

where E, G, and ν are the elastic modulus, shear modulus, and Poisson’s ratio of the isotropic elastic material, respectively.

The shear stress distribution in the fiber and matrix are assumed to be of the form:¹⁶

τ_{rz}^{f} (r, z) = τ (r_{0}, z) I (r) and τ_{rz}^{m} (r, z) = τ (r_{0}, z) O (r)

(3)

where

τ_{rz}^{f} (r, z)

and

τ_{rz}^{m} (r, z)

represent the shear stress in the fiber and the matrix, respectively.

τ (r_{0}, z)

is the axial shear stress at the cohesive interface. The functions

I (r)

and

O (r)

respectively define the shear stress in the fiber and the matrix along the radial r-direction.

According to the basic assumption of the SLM,^19,31 the axial tensile and shear strains of the fiber or the matrix can be given in terms of the axial displacement $w (r, z)$ in the fiber or the matrix, respectively, as below:

ɛ_{zz} = \frac{\partial w}{\partial z} and γ_{rz} = \frac{\partial w}{\partial r}

(4)

The displacement crossing the interface ( $r = r_{0}$ ) is discontinuous because a cohesive interface exists between the fiber and matrix. The relationship among the axial displacements of the fiber ( $w^{f} (r_{0}, z)$ ) and the matrix ( $w^{m} (r_{0}, z)$ ), and the tangential displacement of the cohesive fiber/matrix interface ( $δ (z)$ ) satisfies the following equation (5):²⁴

w^{f} (r_{0}, z) - w^{m} (r_{0}, z) = δ (z)

(5)

Combining equations (2) to (5), the relation among the average axial displacements of the fiber ( ${\bar{w}}^{f} (z)$ ) and the matrix ( ${\bar{w}}^{m} (z)$ ), and the tangential displacement of the cohesive fiber/matrix interface ( $δ (z)$ ) satisfies the following equation (6):²⁴

\begin{matrix} {\bar{w}}^{m} (z) - {\bar{w}}^{f} (z) + δ (z) \\ = {\frac{\bar{[rI (r)]}}{2 G_{f}} + \frac{\bar{[(R^{2} / r^{2} - 1) rO (r)]}}{2 G_{m}}} τ (r_{0}, z) \end{matrix}

(6)

where G_f and G_m are the shear modulus of the fiber and the matrix, respectively. Furthermore,

\bar{[rI (r)]}

and

\bar{[(R^{2} / r^{2} - 1) rO (r)]}

can be calculated from equation (7),^24,32 which defines the average of any variable

f (r)

over the circular cross section

(r_{i} \leq r \leq r_{0})

\begin{matrix} {\begin{matrix} \bar{f (r)} = \frac{1}{π (r_{0}^{2} - r_{i}^{2})} \int_{r_{i}}^{r_{0}} 2 π rf (r) d r \\ \bar{[rI (r)]} = \frac{1}{π r_{0}^{2}} \int_{0}^{r_{0}} 2 π r [rI (r)] d r \\ \bar{[(R^{2} / r^{2} - 1) rO (r)]} = \frac{1}{π (R^{2} - r_{0}^{2})} \int_{r_{0}}^{R} 2 π r [(R^{2} / r^{2} - 1) rO (r)] d r \end{matrix} \end{matrix}

(7)

Considering the load equilibrium of an infinitesimal element of the fiber and the matrix at a distance of z,³¹ the average axial stress of the fiber ( ${\bar{σ}}_{zz}^{f} (z)$ ), the axial force of the fiber ( $P^{f} (z)$ ), and the shear stress of the fiber ( $τ_{rz}^{f} (r_{0}, z)$ ) satisfy equation (8). Likewise, the average axial stress of the matrix ( ${\bar{σ}}_{zz}^{m} (z)$ ) and the shear stress of the matrix ( $τ_{rz}^{m} (r_{0}, z)$ ) satisfy equation (9).

\frac{d {\bar{σ}}_{zz}^{f} (z)}{dz} = \frac{{dP}^{f} (z)}{π r_{0}^{2} dz} = - \frac{2 π r_{0} τ_{rz}^{f} (r_{0}, z)}{π r_{0}^{2} dz} = - \frac{2 τ (r_{0}, z)}{r_{0}}

(8)

\frac{d {\bar{σ}}_{zz}^{m} (z)}{dz} = \frac{{dP}^{m} (z)}{π (R^{2} - r_{0}^{2}) dz} = - \frac{2 π r_{0} τ_{rz}^{f} (r_{0}, z)}{π (R^{2} - r_{0}^{2}) dz} = - \frac{2 τ (r_{0}, z)}{(R^{2} - r_{0}^{2}) / r_{0}}

(9)

Differentiate equation (6) with respect to z yields

{\bar{ɛ}}_{zz}^{m} (z) - {\bar{ɛ}}_{zz}^{f} (z) + \frac{d δ (z)}{dz} = {\frac{\bar{[rI (r)]}}{2 G_{f}} + \frac{\bar{[(b^{2} / r^{2} - 1) rO (r)]}}{2 G_{m}}} \frac{d τ (r_{0}, z)}{dz}

(10)

where

{\bar{ɛ}}_{zz}^{f} (z)

and

{\bar{ɛ}}_{zz}^{m} (z)

are the average axial strains of the fiber and the matrix on the cross section, respectively.

Averaging equation (2) over the fiber cross section and the matrix cross section, respectively, we obtain the average axial strain of the fiber ( ${\bar{ɛ}}_{zz}^{f} (z)$ ) and the matrix ( ${\bar{ɛ}}_{zz}^{m} (z)$ ) below:

\begin{matrix} {\begin{matrix} {\bar{ɛ}}_{zz}^{f} (z) = \frac{{\bar{σ}}_{zz}^{f} (z)}{E_{f}} - \frac{ν_{f}}{E_{f}} [{\bar{σ}}_{rr}^{f} (z) + {\bar{σ}}_{θ θ}^{f} (z)] \\ {\bar{ɛ}}_{zz}^{m} (z) = \frac{{\bar{σ}}_{zz}^{m} (z)}{E_{m}} - \frac{ν_{m}}{E_{m}} [{\bar{σ}}_{rr}^{m} (z) + {\bar{σ}}_{θ θ}^{m} (z)] \end{matrix} \end{matrix}

(11)

where the superscripts f and m denote the fiber and matrix, respectively. E and ν are the elastic modulus and Poisson’s ratio of the isotropic elastic material, respectively.

By adopting equation (12) proposed in Nairn,¹⁷ equation (11) can be simplified to equation (13).

\begin{matrix} | \frac{ν_{f}}{E_{f}} [{\bar{σ}}_{rr}^{f} (z) + {\bar{σ}}_{θ θ}^{f} (z)] | << | \frac{{\bar{σ}}_{zz}^{f} (z)}{E_{f}} | and \\ | \frac{ν_{m}}{E_{m}} [{\bar{σ}}_{rr}^{m} (z) + {\bar{σ}}_{θ θ}^{m} (z)] | << | \frac{{\bar{σ}}_{zz}^{m} (z)}{E_{m}} | \end{matrix}

(12)

{\bar{ɛ}}_{zz}^{f} (z) = \frac{{\bar{σ}}_{zz}^{f} (z)}{E_{f}} and {\bar{ɛ}}_{zz}^{m} (z) = \frac{{\bar{σ}}_{zz}^{m} (z)}{E_{m}}

(13)

As the total force on any cross section of the unit cell in the z-direction should be equal to the force acted on its end (Figure 1), we have

r_{0}^{2} {\bar{σ}}_{zz}^{f} (z) + (R^{2} - r_{0}^{2}) {\bar{σ}}_{zz}^{m} (z) = R^{2} σ_{c}

(14)

Differentiating equation (8) and substituting it together with equations (13) and (14) into equation (10) yield

\begin{matrix} \frac{d^{2} {\bar{σ}}_{zz}^{f} (z)}{{dz}^{2}} = \frac{[\frac{r_{0}^{2}}{R^{2} - r_{0}^{2}} + \frac{E_{m}}{E_{f}}]}{r_{0} E_{m} {\frac{\bar{[rI (r)]}}{4 G_{f}} + \frac{\bar{[(\frac{R^{2}}{r^{2}} - 1) rO (r)]}}{4 G_{m}}}} \\ \times {{\bar{σ}}_{zz}^{f} (z) - \frac{\frac{R^{2}}{r_{0}^{2}} σ_{c}}{[1 + (\frac{R^{2}}{r_{0}^{2}} - 1) \frac{E_{m}}{E_{f}}]} \\ - \frac{(\frac{R^{2}}{r_{0}^{2}} - 1) \frac{E_{m}}{E_{f}}}{[1 + (\frac{R^{2}}{r_{0}^{2}} - 1) \frac{E_{m}}{E_{f}}]} \frac{d δ (z)}{dz}} \end{matrix}

(15)

In conclusion, the governing equations for the SLM with a cohesive fiber/matrix interface are made up of equations (15), (8) and (1). The stress distribution in the fiber and matrix, and the axial tensile modulus of the composite can be computed from the above governing equations subjecting to appropriate boundary conditions. For clarity, the average of the axial tensile stress in the fiber ( ${\bar{σ}}_{zz}^{f} (z)$ ) and matrix ( ${\bar{σ}}_{zz}^{m} (z)$ ) will be respectively represented by $σ_{f} (z)$ and $σ_{m} (z)$ .

Analysis of the tensile modulus and stress of the composites

In this part, the modified SLM is applied to the analysis of the axial tensile properties of the natural fiber composites in the initial linear elastic stage. At this stage, the damage of the cohesive fiber/matrix interface does not appear and the shear stiffness of the cohesive interface is equal to K₀. By differentiating equation (1) and combining with equation (8), we find²⁴

\frac{d δ (z)}{dz} = \frac{1}{K_{0}} \frac{d τ (z)}{dz} = - \frac{1}{K_{0}} \frac{d τ (r_{0}, z)}{dz} = \frac{r_{0}}{2 K_{0}} \frac{d^{2} σ_{f} (z)}{{dz}^{2}}

(16)

where

τ (z)

is the shear stress in the cohesive fiber/matrix interface, and

τ (z) = - τ (r_{0}, z)

Substituting equation (16) into equation (15), it can be obtained that:

\begin{matrix} \frac{d^{2} σ_{f} (z)}{{dz}^{2}} = \frac{2 (\frac{r_{0}^{2}}{R^{2} - r_{0}^{2}} + \frac{E_{m}}{E_{f}})}{r_{0}^{2} [\frac{E_{m}}{2 G_{f}} \bar{[\frac{r}{r_{0}} I (r)]} + \frac{E_{m}}{2 G_{m}} \bar{[(\frac{R^{2}}{r^{2}} - 1) \frac{r}{r_{0}} O (r)]} + \frac{E_{m}}{r_{0} K_{0}}]} \\ \times {σ_{f} (z) - \frac{R^{2}}{r_{0}^{2}} \frac{1}{[1 + \frac{E_{m}}{E_{f}} (\frac{R^{2}}{r_{0}^{2}} - 1)]} σ_{c}} \end{matrix}

(17)

We also introduce the following scaling parameters:

η = z / r_{0} and ζ = \frac{l}{r_{0}}

where ζ is the aspect ratio of the fiber.

\begin{matrix} α = {\frac{2 (\frac{r_{0}^{2}}{R^{2} - r_{0}^{2}} + \frac{E_{m}}{E_{f}})}{[\frac{E_{m}}{2 G_{f}} \bar{[\frac{r}{r_{0}} I (r)]} + \frac{E_{m}}{2 G_{m}} \bar{[(\frac{R^{2}}{r^{2}} - 1) \frac{r}{r_{0}} O (r)]} + \frac{E_{m}}{r_{0} K_{0}}]}}^{1 / 2}, \\ and β = \frac{R^{2}}{r_{0}^{2}} \frac{1}{[1 + \frac{E_{m}}{E_{f}} (\frac{R^{2}}{r_{0}^{2}} - 1)]} \end{matrix}

Assuming $l \approx L$ , we obtain

V_{f} = \frac{r_{0}^{2}}{R^{2}}, β = \frac{1}{V_{f} + \frac{E_{m}}{E_{f}} (1 - V_{f})}

where V_f is the fiber volume fraction.

Then, equation (17) can be expressed as:

\frac{d^{2} σ_{f} (η)}{d η^{2}} = α^{2} [σ_{f} (η) - β σ_{c}]

(18)

and the interfacial shear stress is

τ (r_{0}, η) = - \frac{1}{2} \frac{d σ (η)}{d η}

(19)

Note that the unit cell in this paper is different from that by Chen and Yan,²⁴ which leads to a different load equilibrium condition in equation (14) and a different β value. On the other hand, the shear-lag parameter α is the same as that obtained by Chen. When the initial interfacial shear stiffness K₀ is significantly greater than $E_{m} / r_{0}$ , equation (18) is reduced to the classical SLM with a perfect bonding fiber/matrix interface.²⁴

The unit cell in Figure 1 is symmetrical about the middle point $z = l$ . Therefore, the tangential displacement of the cohesive fiber/matrix interface at the middle location δ is equal to zero. This defines the following boundary conditions:

σ_{f} (η = 0) = ϕ σ_{c}, τ (r_{0}, η = ζ) = 0

(20)

where the fiber end face coefficient ϕ defines the average fiber axial stress at the embedded end face (

η = 0

2 ζ

) in relation to the applied load

σ_{c}

. For the case of a completely free fiber end face (

η = 0

2 ζ

ϕ = 0

. For the case of a perfectly bonded fiber end face (

η = 0

2 ζ

ϕ > 0

.^19,24

Solving the second-order differential equation (18) subjecting to the boundary conditions in equation (20) yields

σ_{f} (η) = β σ_{c} {1 + (\frac{ϕ}{β} - 1) \frac{\cosh [α (η - ζ)]}{\cosh (α ζ)}}

(21)

σ_{m} (η) = \frac{R^{2}}{R^{2} - r_{0}^{2}} σ_{c} - \frac{r_{0}^{2}}{R^{2} - r_{0}^{2}} σ_{f} (η) = \frac{1}{1 - V_{f}} σ_{c} - \frac{V_{f}}{1 - V_{f}} σ_{f} (η)

(22)

τ (r_{0}, η) = - \frac{1}{2} α (ϕ - β) σ_{c} \frac{\sinh [α (η - ζ)]}{\cosh (α ζ)}

(23)

When the fiber volume fraction is very small ( $V_{f} \approx 0$ ), β can be simplified to $E_{f} / E_{m}$ . Comparing the $σ_{f} (η)$ obtained from equation (21) with the result of Zhao and Ji,³¹ we can see that when the fiber volume fraction is very small ( $V_{f} \approx 0$ ) and the fiber ends are well bonded ( $ϕ = 1$ ), there is only a slight difference in the expression of the shear-lag parameter α for introducing K₀ into α.

The average axial strain of the matrix $\bar{ɛ_{m}}$ can then be expressed as

\begin{matrix} \bar{ɛ_{m}} = \frac{\frac{r_{0}}{l} \int_{0}^{\frac{l}{r_{0}}} σ_{m} (η) d η}{E_{m}} \\ = \frac{β V_{f} σ_{c}}{E_{m} (1 - V_{f})} {(\frac{1}{β V_{f}} - 1) + (1 - \frac{ϕ}{β}) \frac{\tanh (α ζ)}{α ζ}} \end{matrix}

(24)

In order to predict the axial tensile Young’s modulus of the composite, it can be assumed that the average strains in the matrix and the composite are equal.³³

\bar{ɛ_{c}} = \bar{ɛ_{m}}

(25)

where

\bar{ɛ_{c}}

is the average axial strain of the composite.

The axial Young’s modulus of the composite E_c can thus be obtained:

E_{c} = \frac{σ_{c}}{\bar{ɛ_{c}}} = \frac{(1 - V_{f}) E_{m}}{β V_{f} {(\frac{1}{β V_{f}} - 1) + (1 - \frac{ϕ}{β}) \frac{\tanh (α ζ)}{α ζ}}}

(26)

Equation (26) describes the relationship between the axial Young’s modulus of the composite E_c and the Young’s modulus of matrix E_m. Equation (26) tells us that E_m, E_f, ϕ, K₀, and ζ can influence the axial tensile Young’s modulus of the composite.

FE model of the unit cell

A meshed 3D FE model of the unit cell containing a flax fiber embedded in the matrix and a thin interface layer of cohesive elements is presented in Figure 3. The two end faces of the fiber are free from the matrix. The axial displacement at the left end face of the matrix is constrained and a uniform displacement load is applied to the right end face of the matrix. The matrix and the fiber are meshed with linear hexahedral elements (C3D8R), as well as the interface layer meshed with cohesive elements (COH3D8).

Figure 3.

FE model of the unit cell.

The mechanical properties of the fiber and matrix are both considered as linear elastic, as shown in Table 1.²³ A traction–separation law, which is described in detail in equation (1), is applied to define the cohesive properties of the fiber/matrix interface with a interfacial shear strength of 4 MPa.³⁴ In order to investigate the interfacial shear stiffness on the tensile properties of the composite, a range of values for fiber/matrix interfacial shear stiffness have been used in computation.

Table 1.

The properties of the fiber and matrix.²³

Property	Young’s modulus E (GPa)	Poisson’s ratio v	Dimensions
PP matrix	E_m = 1.6	ν_m = 0.4
Flax fiber	E_f = 69	ν_f = 0.4	r₀ = 8 µm, ζ: 25, 75, 200

To investigate the effect of fiber volume fraction (V_f) on the tensile properties of the composite, the geometric model of the unit cell needs to be adjusted according to the V_f. As the radius of the flax fiber is constant, the length of the fiber can be calculated from equation (27). According to the assumption in Analysis of the tensile modulus and stress of the composites section, in the FE model of the unit cell, the length of the matrix is set to be slightly larger than the length of the fiber. Therefore, the radius of the unit cell (R) can be obtained from equation (28). In the FE models of the unit cell, the values of the V_f are assumed to be 0.1, 0.2, and 0.3.

L_{f} = 2 r_{0} / ζ

(27)

R = \sqrt{\frac{r_{0}}{V_{f}}}

(28)

where r₀, L_f, and ζ are the radius, length, and aspect ratio of the flax fiber, respectively.

Results and discussion

Axial tensile stress in the fiber

Different assumptions on the radial dependence of the shear stress in the fiber and matrix have been developed in SLMs.^13,14,16,17 In this paper, we adopt the assumption $I (r) = r / r_{0}$ and $O (r) = r_{0} R (R / r - r / R) / (R^{2} - r_{0}^{2})$ , which satisfies the interface continuity and boundary conditions for the shear stress.^16,17 The shear-lag parameter α as a function of $E_{m} / (r_{0} K_{0})$ for the different fiber volume fraction V_f is shown in Figure 4. The value of the shear-lag parameter α decreases with the increase of the ratio $E_{m} / (r_{0} K_{0})$ ; this is in agreement with the results of Chen and Yan.²⁴ Furthermore, the shear-lag parameter α for a low $V_{f} < 0.1$ decreases rapidly when the ratio $E_{m} / (r_{0} K_{0}) > 1$ . The shear-lag parameter α for a larger $V_{f} > 0.1$ begins decreasing rapidly when the ratio $E_{m} / (r_{0} K_{0}) > 0.1$ . However, the shear-lag parameter α remains nearly a constant when $E_{m} / (r_{0} K_{0}) < 0.1$ . It can be concluded that the initial shear stiffness of the cohesive fiber/matrix interface has a great influence on the shear-lag parameter α for a large V_f. That is to say that the interface properties have a more significant effect on the shear stress transmission between the fiber and matrix of the composite with a larger fiber volume fraction. In addition, when $E_{m} / (r_{0} K_{0}) < 0.1$ , the initial shear stiffness K₀ can be considered as large enough and the fiber/matrix interface can be treated as perfectly bonded.

Figure 4.

The shear-lag parameter as a function of $E_{m} / (r_{0} K_{0})$ for different V_f.

The longitudinal tensile stress in the flax fiber surrounded by the polypropylene (PP) matrix as functions of the fiber aspect ratio (ζ), fiber volume fraction (V_f), and the ratio $E_{m} / (r_{0} K_{0})$ is calculated from equation (21) and presented in Figure 5. As shown in Figure 5, the ratio of the fiber stress to the stress applied on the composite ( $σ_{f} (z) / σ_{c}$ ) increases from 1 at the end of the fiber to a maximum at the center of the fiber. This indicates that the fiber in the composite bears more tensile stress and cracks should initiate in the middle of the fiber. The maximum of the $σ_{f} (z) / σ_{c}$ ratio increases with the increase of the fiber aspect ratio, when $ζ \leq 5$ and remains constant when $ζ > 5$ (Figure 5(a)). The maximum of the $σ_{f} (z) / σ_{c}$ ratio increases with the decrease of the V_f(Figure 5(b)), in agreement with the results of Zhao and Ji.³¹ As shown in Figure 5(c), the maximum of the $σ_{f} (z) / σ_{c}$ ratio decreases with the increase of the $E_{m} / (r_{0} K_{0})$ when $E_{m} / (r_{0} K_{0}) > 10$ and reaches a constant when $E_{m} / (r_{0} K_{0}) < 1$ . When the $E_{m} / (r_{0} K_{0}) < 1$ , the part of fiber with a large relative tensile stress is more than 95% of the fiber length and thus fiber breakage can initiate at any position in this region, which is consistent with experimental results.³⁵

Figure 5.

Longitudinal distribution of relative tensile stress in the fiber as a function of (a)ζ, (b)V_f, and (c) $E_{m} / (r_{0} K_{0})$ .

Figure 6 compares the longitudinal distribution of relative tensile stress in the fiber calculated from equation (21) with finite element analysis (FEA) results. There is a strong agreement between the FEA results and our results based on the modified SLM with a cohesive fiber/matrix interface. The main difference is that the longitudinal tensile stress in the two ends of the fiber from FEA is higher than that from the modified SLM. This is caused by stress concentration, which has not been taken into account in equation (21). This supports the proposed modified SLM for the analysis of the fiber tensile stress.

Figure 6.

Comparison of the relative tensile stress predicted by the modified SLM with the FEA results.

According to the above discussion, it can be concluded that fiber aspect ratio ζ, fibre volume fraction V_f, and the $E_{m} / (r_{0} K_{0})$ have large influences on the longitudinal distribution of tensile stress in the fiber. The magnitude of the relative tensile stress in the fiber is large with a high aspect ratio and a strong initial shear stiffness of the interface. This means that fibers with high aspect ratio or strong initial shear stiffness can initiate a break at any position when the relative tensile stress and the interfacial shear strength are high.

Elastic modulus of composites

The longitudinal elastic modulus of unidirectional short flax fibre–reinforced PP matrix composite as functions of fiber aspect ratio (ζ), fiber volume fraction (V_f), and ratio $E_{m} / (r_{0} K_{0})$ is calculated from equation (26) and shown in Figure 7. The longitudinal elastic modulus of the composite (E_c) increases with the increase of the fiber aspect ratio and fiber volume fraction and the decrease of the ratio $E_{m} / (r_{0} K_{0})$ . Furthermore, E_c maintains a near constant when $E_{m} / (r_{0} K_{0}) < 0.1$ . That is to say, the initial shear stiffness of the fiber/matrix interface is large enough to effectively transmit the load from the matrix to the fiber when the ratio $E_{m} / (r_{0} K_{0}) < 0.1$ . It can be seen from Figure 7(a) that the E_c increases slowly with ζ for ζ greater than 100 when the ratio $E_{m} / (r_{0} K_{0}) < 1$ , whereas it increases almost linearly with ζ over the whole range when the ratio $E_{m} / (r_{0} K_{0}) > 10$ . As shown in Figure 7(b) and (c), the influence of the ratio $E_{m} / (r_{0} K_{0})$ on the E_c is greater at large values of V_f and ζ. This is because the fiber can bear more longitudinal tensile load with a higher interface shear stiffness (lower $E_{m} / (r_{0} K_{0})$ ). Therefore, enhancing the interface properties of the composite is a highly effective approach to improve the elastic modulus of short flax fiber–reinforced PP matrix composites, especially for composites containing a high volume fraction of long fibers (large fiber aspect ratio).

Figure 7.

Longitudinal tensile modulus of the composite as functions of (a) the fiber aspect ratio (ζ), (b) the fiber volume fraction (V_f), and (c) the ratio $E_{m} / (r_{0} K_{0})$ .

Figure 8 shows the influence of the bonding strength between the fiber end face and the matrix on the longitudinal tensile modulus of the composite. The top of the error bar represents the value of E_c at $ϕ = 1$ , which means high adhesion between the matrix and fiber end face, while the bottom of the error bar represents the value of E_c at $ϕ = 0$ , which means no adhesion between the matrix and fiber end face. The length of the error bar decreases with the increase of the fiber aspect ratio when $E_{m} / (r_{0} K_{0}) \leq 1$ , while it increases with the increase of the fiber aspect ratio when $E_{m} / (r_{0} K_{0}) = 100$ . That is to say, the bonding strength between the fiber end face and the matrix does not have a significant effect on the longitudinal tensile modulus for long fibers for composites with good fiber/matrix bonding. This result agrees with Nardone and Prewo.¹⁹ The underlining reason for this is that the fiber/matrix interface area for the composite with long fibers is far greater than the area of the fiber end face, and thus most of the load is transmitted to the fiber through the fiber side interface at high interfacial shear stiffness. At low interfacial shear stiffness, the load transmission through the fiber side interface is ineffective, and thus in relative terms, the load transmission between the fiber and matrix through the fiber end face becomes an important mechanism for load bearing.

Figure 8.

Influence of the bonding strength between the fiber end face and the matrix on the longitudinal tensile modulus of the composite.

Figure 9 presents the longitudinal tensile modulus of the composite as a function of the ratio $E_{m} / (r_{0} K_{0})$ in different fiber aspect ratios and fiber volume fractions calculated from equation (26) in comparison with FEA results. There is a good agreement between the FEA results and the results calculated from the modified SLM, although the longitudinal tensile modulus from the modified SLM is marginally higher than the FEA results. This is because the influence of the other directional principal stress on the z axial strain of the composite has been ignored in the modified SLM due to the simplifications made from equations (11) to (13).

Figure 9.

Comparison of the longitudinal tensile modulus of the composite predicted by the modified SLM with FEA results as a function of the ratio $E_{m} / (r_{0} K_{0})$ : (a) in different fiber aspect ratios and (b) in different fiber volume fractions.

The longitudinal tensile modulus of the composite with a high interfacial shear stiffness calculated from equation (26) is compared with that calculated from the Halpin–Tsai model,³⁶ which was used to predict the elastic modulus of natural fiber–reinforced thermoplastics in Facca et al.²⁸ As shown in Figure 10, the results from the modified SLM and Halpin–Tsai model match each other very well, indicating that the modified SLM has the same accuracy in predicting the tensile modulus of a composite with a perfectly bonded interface by assuming a large interfacial shear stiffness. The advantage of the proposed modified SLM is that it is also capable of predicting the elastic modulus of composites with any initial interfacial shear stiffness.

Figure 10.

Comparison of the longitudinal tensile modulus of the composite with a high interfacial shear stiffness predicted by the modified SLM and that by the Halpin–Tsai model.

Comparison to experiment results

Tensile stress and strain in fiber

A great deal of results have been published on the stress–transfer characteristics of the single fiber–reinforced polymeric composite in fragmentation experiments using Laser Raman spectroscopy.^35,37,38 Different volume fractions of short fiber in epoxy composites were loaded up to the levels of stress sufficient to cause interfacial failure and the fiber strain profiles in the fragmentation processes were determined at each level of applied stress.³⁵

Figure 11 compares the experimentally determined fiber tensile strain along the length of the fiber in different carbon fibre–reinforced epoxy composites at an applied strain of 1%³⁵ to the predicted strain using the modified SLM. In the fragmentation test, the applied strain was transferred to the fiber from the matrix through the fiber/matrix interface and both end faces of the fiber were free from any constraint by the matrix. Therefore, value of ϕ should be zero value in the calculation of the fiber tensile stress by the modified SLM for the fragmentation test. The fiber tensile strain $ɛ_{f} (η)$ along the length of the fiber therefore can be calculated using equation (29). Details of the modeling parameters used in the prediction, including the mechanical properties of the fiber and the matrix and the interfacial cohesive properties, are given in Figure 11.

ɛ_{f} (η) = β \frac{E_{c} ɛ_{c}}{E_{f}} {1 - \frac{\cosh [α (η - ζ)]}{\cosh (α ζ)}}

(29)

Figure 11.

Comparison between analytical predictions and experimental results of fiber strain distribution in fragmentation test of single fiber composites.

Figure 11 shows a good agreement between the experiment data and the results obtained from the modified SLM using appropriate shear stiffness of the cohesive interface. Combining the experiment data with the prediction results, we can find that the interfacial shear stiffness between the carbon fiber without electrolytic surface oxidative treatment and epoxy is much less than that of the carbon fiber with electrolytic surface oxidative treatment and epoxy. The low interfacial shear stiffness between the untreated carbon fiber and the matrix leads to the decrease of maximum strain in the fiber, as shown in Figure 11. Therefore, the modified SLM with a cohesive fiber/matrix interface can be used to analyze the stress–transfer characteristics of fiber-reinforced composite with different interfacial properties.

Elastic modulus of composite

The mechanical properties of short flax fiber–reinforced polymeric composites have been widely studied using FE model of a unit cell made of a single fiber embedded in matrix and equating the fiber volume fraction to that of the actual composites.^23,26 The tensile elastic modulus of the composites with different flax fiber volume fractions had been obtained from the FE models by Modniks and Andersons.²⁶ The Young’s modulus of the flax fiber and matrix used in their FE models were 69 GPa and 1.6 GPa, respectively. The length and the diameter of the flax fiber were 1.21 mm and 16 µm, respectively.

Table 2 compares the tensile elastic modulus of a short flax fiber–reinforced polymeric composite calculated from the modified SLM to that from the FE model by Chen and Yan.²⁴ In the FE model, the two end faces of the fiber were fully bonded with the matrix and the fiber–matrix interface was perfectly bonded. Therefore, a large interfacial shear stiffness corresponding to

E_{m} / (r_{0} K_{0}) = 0.01

and a large

ϕ = 2.7

are used in the modified analytical SLM. Table 2 shows that the tensile elastic modulus of the composite from our prediction is smaller than the FE results by Modniks at low fiber volume fraction (0.13). But difference between the two results narrows as the high fibre volume fraction increases. At the fiber volume ratio of 0.3, the difference between the two methods was only 0.5%. There are two main factors jointly influencing this difference between the SLM prediction and the FE results, the tensile strain in the simplifications and assumptions adopted in the modified SLM as shown in equation (13), and the stress transmission of the fiber end face with the matrix. The latter changes with the fiber volume fraction and fibre aspect ratio. At a lower fiber volume fraction, the stress transmission between the fiber end face and the matrix has a greater effect, which leads to an underestimate of the tensile stress in the composite and thus the FE results were greater than the SLM prediction. As the fiber volume fraction increases, the effect of the simplifications and assumptions in equation (13) gradually increases and eventually surpassed the effect of the fiber end face, resulting in an underestimate of the tensile strain of the composite. So the SLM prediction at high fiber fraction was slightly greater than the FE results. Overall, the modified SLM prediction improves with the fiber volume fraction. At a fiber volume fraction of 0.2, the SLM predictions were 3% lower than the corresponding FE results. At a fiber volume fraction of 0.3, the SLM predictions were 0.5% higher than the FE results.

Table 2.

Comparison of the tensile elastic modulus of a short flax fiber-reinforced polymeric composite.

Fiber volume fraction	Finite element results of Modniks (GPa)	Modified SLM Prediction (GPa)	Error
0.13	10.3	9.4	8.7%
0.21	15	14.6	2.7%
0.29	20.3	20.4	0.5%

SLA: shear-lag model.

Conclusions

In this paper, a modified SLM with a cohesive fiber/matrix interface is developed and applied to the analysis of the effect of interfacial properties on the elastic properties of unidirectional short flax fiber-reinforced polymeric composites. 3D FE models of the composite material accounting for the fiber/matrix interfacial properties are built using CZM to simulate the tensile properties of the composite. The tensile stress in the fiber and the tensile elastic modulus of the composite with different interfacial properties, fiber aspect ratios, and fiber volume fractions are computed using both the modified SLM and the FE models.

The analytical solutions show that the interfacial shear stiffness, fiber aspect ratio, and fiber volume fraction all have significant effects on the tensile stress of the fiber and the tensile elastic modulus of the composite. The maxima of the fiber relative tensile stress and the composite tensile elastic modulus increase with the increase of the interfacial shear stiffness and the fiber aspect ratio and then plateau after a critical interfacial shear stiffness and a critical fiber aspect ratio have been reached. The maximum of the fiber relative tensile stress decreases with the increase of the fiber volume fraction, while the tensile elastic modulus of the composite increases with the increase of the fiber volume fraction. The bonding strength between the fiber end face and the matrix has a great influence on the tensile elastic modulus of the composite at low fiber aspect ratio and low interfacial shear stiffness. At high fiber aspect ratio and high interfacial shear stiffness, this influence becomes quite small. The prediction of composite tensile elastic modulus using the modified SLM agrees well with that calculated from the Halpin–Tsai model when the product of the interfacial shear stiffness and the fiber radius ( ${rK}_{0}$ ) is far greater than the elastic modulus of the matrix (E_m).

There is a strong agreement between the predicted results using the modified SLM and the simulation results obtained from the FE models using CZM. The predicted fiber tensile stresses using the appropriate interfacial cohesive properties fitted published experimental results using Laser Raman Spectroscopy technique very well. The predicted tensile modulus of the composite using appropriate interfacial cohesive stiffness and fiber end face coefficient showed a good agreement with the simulation results using the FE models by assuming both the two fiber end faces are perfectly bonded with the matrix. The modified SLM provides a simple and effective method for the analysis of the stress–transfer characteristics and axial tensile properties of short fibre–reinforced composites accounting for the interfacial shear stiffness, fiber aspect ratio, fiber volume fraction and bonding strength between the fiber end face and the matrix.

Footnotes

Declaration of Conflicting Interests

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

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: The authors acknowledge the financial support from CSIRO Manufacturing, National Natural Science Foundation of China (No.51575417), Innovative Research Team Development Program of Ministry of Education of China (No. IRT_17R83), and 111 Project (No. B17034). Xiaoshuang Xiong also gratefully acknowledges the scholarship from China Scholarship Council (201606950038) which enabled the work at CSIRO in Australia.

References

Kong

Lee

Park

. Design and manufacturing of automobile hood using natural composite structure. Compos Part B Eng 2016; 91: 18–26.

Kabir

Wang

Lau

et al.

Chemical treatments on plant-based natural fibre reinforced polymer composites: an overview. Compos Part B Eng 2012; 43: 2883–2892.

Shah

Schubel

Clifford

. Can flax replace E-glass in structural composites? A small wind turbine blade case study. Compos Part B Eng 2013; 52: 172–181.

Ghasemnejad

Occhineri

Swift-Hook

. Post-buckling failure in multi-delaminated composite wind turbine blade materials. Mater Design 2011; 32: 5106–5112.

Scarponi

Messano

. Comparative evaluation between E-Glass and hemp fiber composites application in rotorcraft interiors. Compos Part B Eng 2015; 69: 542–549.

Scarponi

. Hemp fiber composites for the design of a Naca cowling for ultra-light aviation. Compos Part B Eng 2015; 81: 53–63.

Karus

Kaup

. Natural fibres in the European automotive industry. J Ind Hemp 2002; 7: 119–131.

Bledzki

Faruk

Sperber

. Cars from bio-fibres. Macromol Mater Eng 2006; 291: 449–457.

Guillebaud-Bonnafous

Vasconcellos

Touchard

et al.

Experimental and numerical investigation of the interface between epoxy matrix and hemp yarn. Compos Part A Appl Sci Manuf 2012; 43: 2046–2058.

10.

Brighenti

Carpinteri

Scorza

. Mechanics of interface debonding in fibre-reinforced materials. J Compos Mater 2015; 50: 2699–2718.

11.

Mokhtar

Alexandre

Alain

et al.

Damage analysis of composites reinforced with alfa fibers: viscoelastic behavior and debonding at the fiber/matrix interface. J Appl Polym Sci 2016; 133: 1–10.

12.

Zhang

Milanovic

Zhang

et al.

Effects of humidity conditions at fabrication on the interfacial shear strength of flax/unsaturated polyester composites. Compos Part B Eng 2014; 60: 186–192.

13.

Cox

. The elasticity and strength of paper and other fibrous materials. Br J Appl Phys 1952; 3: 72–79.

14.

Hsueh

. Elastic load transfer from partially embedded axially loaded fibre to matrix. J Mater Sci Lett 1988; 7: 497–500.

15.

Hsueh

. Interfacial debonding and fiber pull-out stresses of fiber-reinforced composites VII_improved analyses for bonded interfaces. Mater Sci Eng A-Struct 1992; 154: 125–132.

16.

Mccartney

. New theoretical model of stress transfer between fibre and matrix in a uniaxially fibre-reinforced composite. Proc R Soc of Lond 1989; 425: 215–244.

17.

Nairn

. On the use of shear-lag methods for analysis of stress transfer in unidirectional composites. Mech Mater 1997; 26: 63–80.

18.

Zhao

. Location of tensile fracture within rigid-brittle inclusions in a ductile flowing matrix. Tectonophysics 1993; 220: 23–31.

19.

Nardone

Prewo

. On the strength of discontinuous silicon carbide reinforced aluminum composite. Scr Metall 1986; 20: 43–48.

20.

Xia

Okabe

Curtin

. Shear-lag versus finite element models for stress transfer in fiber-reinforced composites. Compos Sci Technol 2002; 62: 1141–1149.

21.

Zhang

Qiu

. Modified shear lag model for fibers and fillers with irregular cross-sectional shapes. J Adhes Sci Technol 2003; 17: 397–408.

22.

Virk

Hall

Summerscales

. Modulus and strength prediction for natural fibre composites. Mater Sci Technol 2013; 28: 864–871.

23.

Modniks

Andersons

. Modeling the non-linear deformation of a short-flax-fiber-reinforced polymer composite by orientation averaging. Compos Part B Eng 2013; 54: 188–193.

24.

Chen

Yan

. A shear-lag model with a cohesive fibre–matrix interface for analysis of fibre pull-out. Mech Mater 2015; 91: 119–135.

25.

Shokrian

Shelesh-Nezhad

Soudmand

. 3D FE analysis of tensile behavior for co-PP/SGF composite by considering interfacial debonding using CZM. J Reinf Plast Comp 2015; 35: 365–374.

26.

Modniks

Andersons

. Modeling elastic properties of short flax fiber-reinforced composites by orientation averaging. Comp Mater Sci 2010; 50: 595–599.

27.

Modniks

Joffe

Andersons

. Model of the mechanical response of short flax fiber reinforced polymer matrix composites. Proced Eng 2011; 10: 2016–2021.

28.

Facca

Kortschot

Yan

. Predicting the elastic modulus of natural fibre reinforced thermoplastics. Compos Part A Appl Sci Manuf 2006; 37: 1660–1671.

29.

Zhi

Long

Miao

. Influence of microbond test parameters on interfacial shear strength of fiber reinforced polymer-matrix composites. Compos Part A Appl Sci Manuf 2017; 100: 55–63.

30.

Chen

Miao

Ding

. Influence of moisture absorption on the interfacial strength of bamboo/vinyl ester composites. Compos Part A Appl Sci Manuf 2009; 40: 2013–2019.

31.

Zhao

. Refinements of shear-lag model and its applications. Tectonophysics 1997; 279: 37–53.

32.

Nayfeh

. Thermomechanically induced interfacial stresses in fibrous composites. Fibre Sci Technol 1977; 10: 195–209.

33.

Nair

Kim

. Modification of the shear lag analysis for determination of elastic modulus of short-fiber (or whisker) reinforced metal matrix composites. J Appl Mech 1992; 59: 176–182.

34.

Awal

Cescutti

Ghosh

et al.

Interfacial studies of natural fibre/polypropylene composites using single fibre fragmentation test (SFFT). Compos Part A Appl Sci Manuf 2011; 42: 50–56.

35.

Galiotis

. A study of mechanisms of stress transfer in continuous- and discontinuous-fibre model composites by laser Raman spectroscopy. Compos Sci Technol 1993; 48: 15–28.

36.

Halpin

. Stiffness and expansion estimates for oriented short fiber composites. J Compos Mater 1969; 3: 732–734.

37.

Melanitis

Galiotis

. Interfacial micromechanics in model composites using laser Raman spectroscopy. P Roy Soc A-Math Phy 1993; 440: 379–398.

38.

Mehan

Schadler

. Micromechanical behavior of short-fiber polymer composites. Compos Sci Technol 2000; 60: 1013–1026.