Brett A. Bednarcyk

Glenn Research Center, Cleveland, Ohio

Jacob Aboudi

Tel Aviv University, Tel Aviv, Israel 

Steven M. Arnold

Glenn Research Center, Cleveland, Ohio

The Effect of General Statistical Fiber Misalignment 

on Predicted Damage Initiation in Composites

NASA/TM—2014-218313

July 2014
	
	
	








	



	

!

	
	
"#"	




		


	
	

		


"$
	

%	



	&	

		


'	
	
	



	
(

	)
	

			


	




	
*
	

)	

+

")
		*



*	
	",

. 5$6$79:7$%)
	




		
	
;
	 


				


	
	
<
	
	

"






	







		


		
	+	*
	


	

		





	<


		


. 5$6$7=5=%)(9= 

		
	"
	 


&	>#	 

	
	*
	#	
	






(




<"

. $%)$%))5%)

"+

	 


	
		

. $%?5)5$59:7$%$


		




	"


	
	
	




	
	


	"

. 5$79:7$%


	
	
	

	

 

	
	
	
;





	*;

	

. 5$6$7)7%5+

	


	


		
'



&	
	



	

&
	&

			

?
	

	
	




	
	


*,

. 	
	




http://www.sti.nasa.gov

. 5+
"
	>

help@sti.nasa.gov

. ?<"
	>




	

(#@@BHKLKHLOQB

. 

	

(#

 @@BHKLKHLOQR

. U	
,


	

(#

 $	
		

	



KVVL	(	

6
	=(RVQKWHVBRQ
Brett A. Bednarcyk

Glenn Research Center, Cleveland, Ohio

Jacob Aboudi

Tel Aviv University, Tel Aviv, Israel 

Steven M. Arnold

Glenn Research Center, Cleveland, Ohio

The Effect of General Statistical Fiber Misalignment 

on Predicted Damage Initiation in Composites

NASA/TM—2014-218313

July 2014


	



	


X)	$	 

$%
@@VBL
	



$	
		

	



KVVL	(	

6
	=(RVQKWHVBRQ



	

	

LBQV*)


<	YRRBVR

	
",ZZ***


*
	#*

	"?
	
	
	


X)	$	

Level of Review,
	"	*"


NASA/TM—2014-218313 1 

The Effect of General Statistical Fiber Misalignment on 

Predicted Damage Initiation in Composites 

Brett A. Bednarcyk 

National Aeronautics and Space Administration 

Glenn Research Center 

Cleveland, Ohio 44135 

Jacob Aboudi 

Tel Aviv University 

Tel Aviv, Israel 69978 

Steven M. Arnold 

National Aeronautics and Space Administration 

Glenn Research Center 

Cleveland, Ohio 44135 

Abstract 

A micromechanical method is employed for the prediction of unidirectional composites in which the 

fiber orientation can possess various statistical misalignment distributions. The method relies on the 

probability-weighted averaging of the appropriate concentration tensor, which is established by the 

micromechanical procedure. This approach provides access to the local field quantities throughout the 

constituents, from which initiation of damage in the composite can be predicted. In contrast, a typical 

macromechanical procedure can determine the effective composite elastic properties in the presence of 

statistical fiber misalignment, but cannot provide the local fields. Fully random fiber distribution is 

presented as a special case using the proposed micromechanical method. Results are given that illustrate 

the effects of various amounts of fiber misalignment in terms of the standard deviations of in-plane and 

out-of-plane misalignment angles, where normal distributions have been employed. Damage initiation 

envelopes, local fields, effective moduli, and strengths are predicted for polymer and ceramic matrix 

composites with given normal distributions of misalignment angles, as well as fully random fiber 

orientation. 

1.0 Introduction 

It is well-known that a micromechanical analysis can provide the effective properties of composite 

materials from the knowledge of the constituent properties, their geometric arrangement, and their 

detailed interactions. However, the impact of fiber misalignment on the composite response, which often 

can occur in modern composites (Barwick and Papathanasiou, 2003), is typically not accounted for in 

design and analysis. A notable exception is the effect of longitudinal fiber misalignment on the 

compressive failure response of composites, which has been studied extensively in the literature (c.f., 

Budiansky, 1983; Wisnom (1990); Budianksy and Fleck, 1993; Schultheisz and Waas, 1996; Waas and 

Schultheisz, 1996; Liu et al., 2004; Basu et al., 2006; Numayr and Al Rjoub, 2013; Zhou et al., 2013; 

Ataabadia et al, 2014). Such misalignment can occur due to manufacturing deficiency and defects, curved 

surfaces, and fiber waviness, which was investigated by Kugler and Moon (2001). A standard approach to 

incorporating the effect of fiber misalignment in a unidirectional composite is to consider the anisotropic 

ply properties and perform a standard transformation to a desired fiber orientation (for instance, an inplane rotation, c.f. Jones (1975), Herakovich (1998)). While such a macromechanical approach will 

provide the effective properties of the ply, the associated local field distributions within the constituent 

materials are not available (as is always the case with the macromechanical approach). It is precisely 
NASA/TM—2014-218313 2 

these local fields that dictate and drive the damage initiation and progression, yielding, failure, and other 

nonlinearities, which impact the macroscopic composite response and the performance of structures 

composed of composite materials. Alternatively, if one utilizes an appropriate micromechanical analysis, 

the local field variations within the composite constituents can be predicted based on the knowledge of 

the strain (or stress) concentration tensors, which are naturally given by the micromechanics theory. 

Consequently, the effect of fiber misalignment on the local fields in the composite can be captured 

through the use of these concentration tensors. 

In the present investigation, the micromechanical analysis known as the High-Fidelity Generalized 

Method of Cells (HFGMC) (Aboudi et al., 2013) is enhanced to capture the effects of possible fiber 

misalignment in unidirectional composites. In the latter reference, HFGMC has been shown to provide a 

reliable and robust micromechanical method for predicting the properties and nonlinear response of a 

wide array of composites. Further, it was shown that the strain (or stress) concentration tensors are 

fundamental to predicting both the local and global composite response and thus they are used herein to 

predict the effect of fiber misalignment on the composite behavior. In particular, the concentration tensors 

are employed in an averaging procedure to determine the effect of a known distribution of misalignment 

angles (i.e., orientation distribution function), or completely randomly fiber orientation, on the local fields 

in the composite. Consequently, the initial yield and damage surfaces, field distributions in the 

constituents, as well as the standard effective properties, can be predicted for the composite. 

Averaging based on orientation has been employed extensively within the microsphere model (Miehe 

et al., 2004), which was originally developed for analyzing rubber (Miehe and Göktepe, 2005; Göktepe 

and Miehe, 2005). Since then, the microsphere model has been used extensively in the analysis of 

biological materials. Menzel and Waffenschmidt (2009) used the model, with an evolving orientation 

distribution function, to simulate remodeling in soft tissues. Murtada et al. (2010) used the microsphere 

approach to analyze smooth muscles, wherein a specialized distribution of the muscle contractile fiber 

orientation, as a function of stretch, was employed. Alstrué et al. (200!	-periodic von Mises 

orientation distribution function, and then a Bingham orientation distribution function (Alstrué et al., 

2010) in applying the microsphere model to blood vessels. Waffenschmidt et al. (2012) used the 

microsphere approach, with evolving orientation density functions, to model bone remodeling. Li et al. 

(2013) have recently applied the microsphere approach to structural composite materials, wherein a von 

Mises fiber orientation distribution has been employed. Orientational averaging was also recently 

employed by Modniks and Andersons (2013) to model the nonlinear mechanical response of short-flaxfiber-reinforced composites using an analytical approach. While the approach and equations were general 

in terms of fiber misalignment distribution, results focused on a uniform distribution. 

The remainder of this paper is organized as follows. A brief description of the HFGMC methodology 

is provided, followed by details on the incorporation of the statistical fiber misalignment distributions. 

Results are then given exhibiting the impact of fiber misalignment on the properties, damage initiation 

envelopes, and critical values of the local field variables. Two classes of composite materials are 

examined: a polymer matrix composite (PMC) system, namely, graphite/epoxy, and a ceramic matrix 

composite (CMC) system. While the PMC system exhibits extreme mismatch in properties between the 

constituents, the CMC system has low mismatch between the fiber and matrix, but includes a very 

compliant fiber interfacial layer. 

2.0 High-Fidelity Generalized Method of Cells Micromechanical Model 

The HFGMC micromechanical model is employed herein to predict the effective behavior of 

unidirectional composites with known fiber misalignment. This theory has been fully described by Aboudi 

et al. (2013). The continuously reinforced (i.e., doubly periodic) version of HFGMC is briefly outlined in 

the following. The doubly periodic microstructure considered is shown in Figure 1(a) in terms of the global 

coordinates (x2,x3). The repeating unit cell, Figure 1(b), defined with respect to local coordinates (y2,y3), of 

such a composite is divided into Nand N subcells in the y2 and y3 directions, respectively. Each subcell is 

labeled by the indices () with = 1,,Nand  = 1,,N, and may contain a distinct homogeneous 
NASA/TM—2014-218313 3 

Figure 1.—(a) A multiphase composite with doubly-periodic microstructures defined with respect to global coordinates 

(x2, x3). (b) The repeating unit cell is represented with respect to local coordinates (y2, y3). It is divided into Nand 

N subcells, in the y2 and y3 directions, respectively. (c) A characteristic subcell () with local coordinates 

( , ) ( )

3

( )

2

 y y whose origin is located at its center. 

material. The dimensions of subcell () in the y2 and y3 directions are denoted by hand l, respectively. A 

local coordinate system ( , ) ( )

3

( )

2

 y y is introduced in each subcell whose origin is located at its center. The 

local (subcell) constitutive equation of the elastic, anisotropic material is given by, 

( ) ( ) ( )    ij Cijkl kl (1) 

where ( ) ij , ( ) kl , and ( ) Cijkl are the components of the stress, strain, elastic stiffness tensors, 

respectively. 

The basic assumption in HFGMC is that the displacement vector ( ) ui in each subcell is expanded 

into quadratic forms in terms of its local coordinates ( , ) ( )

3

( )

2

 y y , as follows,

( )

(02)

2

( )2

3

( )

(20)

2

( )2

2

( )

(01)

( )

3

( )

(10)

( )

2

( )

(00)

( )

4

3

2

1

4

3

2

     1     



	













  



	













i  ij j  i  i  i   i Wi

l

W y

h

u x W y W y W y (2) 

where ij is the applied (external) average strains, and the unknown terms ( )

( )

 Wi lm must be determined 

from the fulfillment of the equilibrium conditions, the periodic boundary conditions, and the interfacial 

continuity conditions of displacements and tractions between subcells. The periodic boundary conditions 

ensure that the displacements and tractions at opposite surfaces of the repeating unit cell are identical. A 

principal ingredient in the HFGMC micromechanical analysis is that all these conditions are imposed in 

the average (integral) sense. 
NASA/TM—2014-218313 4 

As a result of the imposition of these conditions, a linear system of algebraic equations is obtained, 

which can be represented in the following form: 

KU  f (3) 

where the matrix K contains information on the geometry and properties of the materials within the 

individual subcells (), and the displacement vector U contains the unknown displacement coefficients 

( )

( )

 Wi lm , which appear on the right-hand side of Equation (2). The vector f contains information on the 

applied average strains ij . The solution of Equation (3) enables the establishment of the following 

localization relation which expresses the average strains    ij in the subcell () to the externally applied 

average strains ij in the form, 

ij  Aijkl kl

( ) ( ) (4) 

where ( ) Aijkl is the strain concentration tensor components, of the subcell (). The final form of the 

effective constitutive law of the multi-phase thermo-inelastic composite, which relates the average 

stresses ij and strains kl , is established as follows: 

ij  Cijklkl  (5) 

In this equation  Cijkl are components of the effective stiffness tensor, which are given by, 

( ) ( )

1 1

1  



 

 



 ijpq pqkl

N N

ijkl h l C A

HL

C (6) 

Next, the components of stress concentration tensor, ( ) Bijkl , which relate the average stresses in the 

subcell, ( ) ij , to the average (global) stresses, ij , are determined. By combining Equations (1), and (4), 

the subcell stresses are given by, 

ij  Cijpq Apqklkl

( ) ( ) ( ) (7) 

Then, using Equation (5), one obtains, 

ij  Bijkl kl

( ) ( ) (8) 

where, 

    Bijkl  Cijpq ApqrsSrskl

( ) ( ) ( ) (9) 

and  Srskl are components of the effective compliance tensor. 

Consequently, with the present information provided by the HFGMC, the local elastic fields 

throughout the composite constituents are known, as are the composite effective (homogenized) 

properties. 
NASA/TM—2014-218313 5 

3.0 Incorporation of Fiber Misalignment Distribution 

3.1 Transformation to Misalignment Angles 

Suppose that the fibers of the unidirectional composite are oriented in the x1-direction. In order to 

express a rotation of the fiber direction by angle  in the x1-x2 plane about the x3-axis (see Figure 2), the 

following transformation matrix must be employed, 





















 

  



0 0 1

sin cos 0

cos sin 0

(10) 

Subsequently, a rotation about the new x2-axis, namely /

2 x , by an angle  is performed (see Figure 2) 

which can be expressed by the following transformation matrix, 





















  

 



sin 0 cos

0 1 0

cos 0 sin

 (11) 

The full transformation is then given by, 











































3

2

1

/

3

/

2

/

1

x

x

x

x

x

x

T (12) 

where T   . 

Figure 2.—Transformation from the original (un-prime) to rotated (prime) 

coordinate system through the rotation angles  and . int x1 is an 

intermediate x1 coordinate direction. 
NASA/TM—2014-218313 6 

Consequently, the standard approach for transforming the fourth-order effective stiffness tensor  Cijkl , 

given the misalignment angles  and , is as follows, 

  Cijkl  TipTjqTkrTlsCpqrs/ (13) 

Referring to Equation (6), it can be observed that, in the framework of micromechanics analysis, 

transforming the effective stiffness tensor  Cijkl implies the transformation of, 

( ) ( ) ( ) Zijkl  Cijpq Apqkl (14) 

where [ ] ( ) ( )   Z Zijkl can be referred to as a mixed concentration tensor, as it relates the subcell stresses 

to the global strains, see Equation (7). Because this is a fourth-order tensor, its transformation is also 

given in the form of Equation (13). 

Equations (13) and (6) provide, 

( )/

1 1

( )

1 1

/ 1 1 



 





 

  

 

  ijkl

N N

ip jq kr ls pqrs

N N

ijkl h l Z

HL

h l T T T T Z

HL

C (15) 

In order to determine the local stresses in the subcells in the global coordinates including the effects 

of fiber misalignment, through the use of Equation (7), one obtains, 

ij  Zijkl kl

( ) ( )/ (16) 

The global stress-strain relations for the composite, including a fiber misalignment, can be expressed as, 

ij ijkl kl   S  / (17) 

where / / 1 [ ]    S  C . Substituting (17) into (16) leads to the final form of the local stresses in the subcells 

in the global coordinates, including the effects of fiber misalignment. 

ij  Zijpq Spqklkl  Bijkl kl

   ( )/ / ( ) ( )/ (18) 

3.2 Probability-Weighted Averaging 

Now consider a known statistical distribution of fiber misalignment, for instance, over a composite 

part, coupon, or structure. Herein a normal fiber misalignment distribution has been adopted with 

probability density function given by, 

  













 

	


 




    

2

2

1

exp 2

1

s

x

s

p x (19) 

where x is the random variable and  and s are the mean and standard deviation, respectively. Yurgartis 

(1987) and Jelf and Fleck (1992) have shown that local fiber misalignment angles in composites can be 

reasonably approximated by a normal distribution. Figure 3 shows an example of a normal distribution 

probability density function versus x/s for the case where  = 0. 
NASA/TM—2014-218313 7 

Figure 3.—Normal distribution probability density function, where the random 

variable has been normalized by the standard deviation and with the mean 

equal to zero. 

In order to account for the known fiber misalignment variability with respect to angles  and , see 

Equations (10) and (11), the following operation has been employed, 

             









  Z Z p p d d ijkl ijkl , ( ) ( )/ (20) 

providing the probability-weighted average of the mixed concentration tensor. It should be noted that, in 

Equation (20), the transformed mixed concentration tensor is a function of both transformation angles, 

and . Therefore, Equation (20) is valid for any choice of probability density functions for either angle. 

Observing Equation (15), it is clear that the probability-weighted average of the effective stiffness 

tensor of the composite is given by, 

( )

1 1

1 



 

 



 ijkl

N N

ijkl h l Z

HL

C (21) 

Finally, from Equation (16), the subcell stresses in the presence of probability-weighted averaging, 

are given by, 

ij  Zijkl kl

( ) ( ) (22) 

The global stress-strain relations for the composite, including probability-weighted averaging, can be 

expressed using Equation (17) as, 

ij ijkl kl   S   (23) 

-5 -4 -3 -2 -1 0 1 2 3 4 5

p(x)

x/s
NASA/TM—2014-218313 8 

where 1 * *  S  C . Substituting Equation (23) into Equation (22) gives the final form of the subcell 

stresses, as follows, 

ij ijpq pqkl kl Bijkl kl   Z S    ( ) ( )  ( ) (24) 

3.3 Fully Random Averaging 

In the case where the fibers within the composite are fully randomly oriented, the overall behavior of 

the composite is isotropic. The transformation of the mixed concentration tensor, Z()

, can be carried out 

using the three Eulerian angles, E, E, and E, Aboudi et al. (2013). A fourth-order tensor transformation 

is applied to relate the original and transformed systems. In tensorial notation, this transformation is given 

by (Arridge, 1985), 

( ) ( ) 

i j k  i i j j k k ZijkZ a a a a (25) 

where ai i are given by 

E

E E

E E

E E

E E E E E

E E E E E

E E

E E E E E

E E E E E

a

a

a

a

a

a

a

a

a

 

  

   

  

       

      

  

       

      

 

 

 

 

 

 

 

 

 

cos

sin sin

cos sin

sin sin

sin cos sin cos cos

cos cos sin sin cos

sin cos

cos sin sin cos cos

cos cos cos sin sin

3 3

3 2

3 1

2 3

2 2

2 1

1 3

1 2

1 1

(26) 

The mean value of ( ) Zi j k for random orientation is given by, 

      

 

 





    sin d ( )

0

2

0

2

0 8

( ) 1

i j k 2 i i j j k k ZijkZ d d a a a a (27) 

Therefore, the effective isotropic stiffness tensor of the composite with fully random fiber orientations is 

given by Equation (21), and the subcell stresses are given by Equation (24). Note that Christensen and 

Waals (1972) used a fully random averaging approach, along with the concentric cylinder model, to 

predict the elastic properties of composites with random oriented fibers. Also, Luo and Daniel (2003) 

used a fully random averaging approach in conjunction with the Mori-Tanaka (1973) micromechanical 

method to predict the effective properties of nanocomposites. 

4.0 Practical Application 

The equations presented above are useful for assessing the impact of fiber misalign on, not only the 

effective properties, but also the local fields within composite materials. Given the local stress and strain 

fields, initiation of damage or other detrimental nonlinear mechanisms can be predicted for use as 

estimates of the composite allowables. This section provides a brief step-by-step example of how the 

presented methodology can be applied, prior to the results, which are presented in the next section. 
NASA/TM—2014-218313 9 

Given a composite material with known fiber volume fraction and constituent material properties, the 

first step in applying the presented methodology is to discretize the composite repeating unit cell into a 

number of subcells. For example, consider a unidirectional graphite/epoxy composite with a fiber volume 

fraction of 0.6. The transversely isotropic graphite fiber properties are as follows: axial Young's modulus 

= 276 GPa, transverse Young's modulus = 15 GPa, axial shear modulus = 15 GPa, axial Poisson's 

ratio = 0.2, and transverse Poisson's ratio = 0.2. The isotropic epoxy matrix properties are: Young's 

modulus = 3.42 GPa and Poisson's ratio = 0.34. A 40 by 40 subcell discretization of the repeating unit 

cell can be considered, as shown in Figure 4. 

Next, the fiber misalignment present in the composite must be specified. Fully random fiber 

alignment is a special case, as described in Section 3.3. Otherwise, the methodology presented above 

considers a normal distribution of fiber misalignment, see Equation (19). The mean, , and standard 

deviation, s, of the in plane () and out-of-plane () fiber angles within the composite fully characterizes 

the fiber misalignment. For example, with the mean fiber angles representing perfectly aligned fibers, 

 =  = 0, the standard deviations might be s = 3° and s = 2°. Then, based on Equation (21), the 

effective elastic properties of the this misaligned composite can be predicted, as can the effective stresses 

throughout the composite, based on Equation (24), given specified global composite level stresses. The 

results for the specified example are that the predicted effective axial Young's modulus of the composite 

is reduced from 164.7 GPa (in the perfectly aligned case) to 163.2 GPa (reduction of 0.9 percent), 

whereas, in response to a global uniaxial applied stress of 1 MPa in the 0° fiber direction, the maximum 

von Mises stress in the matrix increases from 21.3 to 23.3 Pa (increase of 9.4 percent). Hence, a slight 

reduction in axial modulus and a moderate reduction in the axial stress allowable might be warranted for 

use in design of such a composite with this fiber misalignment distribution. 

Figure 4.—Graphite/epoxy repeating unit cell 

divided into 40 by 40 subcells. The fiber 

volume fraction shown is 0.6. 
NASA/TM—2014-218313 10 

5.0 Results and Discussion 

Results are shown for a polymer matrix composite (PMC) and a ceramic matrix composite (CMC) 

that exhibit the effects of fiber misalignment on the composite effective properties, local field 

distributions within the constituents, and the initial damage surfaces (envelopes). Both probabilityweighted averaging of the misalignment and the fully random fiber orientation case are presented. 

In order to generate the initial damage surfaces, two criteria were used, given the externally applied 

stresses, ij , or strains, ij , on the composite. The first is the traditional von Mises criteria, given by, 

eq  ij ij  Y ( ) ( ) ( ) ˆ ˆ 2

3 (28) 

where ˆ 3 ( ) ( ) ( ) ij  ij  kk !ij

   are the subcell deviatoric stresses, !ij is the Kronecker delta, and Y is the 

yield stress in simple tension. The second damage initiation criterion, suggested by Lemaitre and 

Chaboche (1990), is used to express the critical strain energy release rate associated with loss of stiffness 

in brittle materials, which are typically highly dependent on the state of triaxial stress, 3 h  kk . The 

		<"
7
	RQQV!!" 

           

 

2

1 3 1 2

3

2





	















   "   " 



  

eq

h Rv (29) 

where ( ) " is the Poisson's ratio of the isotropic subcell and ( ) h is the subcell average hydrostatic 

pressure. The damage initiation criterion is then given by, 

    eq Rv  cr   2

1

(30) 

where cr is a material parameter that specifies the stress at which damage initiates for a uniaxial state 

test. Thus, if the triaxiality function is 1 and cr = Y, this damage initiation criterion is equivalent to the 

von Mises criterion, Equation (28). 

To generate a damage initiation surface in a given global stress plane, ij – kl (e.g., 11 – 22 ), 

the stress components are expressed as, 

ij  Rcos# kl  Rsin # (31) 

where R is the radial distance from the origin at a point located on the initial damage surface in the 

ij – kl stress plane and # is the corresponding polar angle. For a given #, the local (subcell) stress 

components are readily determined (in terms of R) by substituting Equations (31) into Equation (8). These 

expression for the local stresses can then be substituted into either initiation criterion, Equation (28) or 

Equation (30). This provides the value of R/Y or R/cr corresponding to damage initiation at the specified 

polar angle #. 

As stated earlier, normal distributions for the fiber misalignment angles has been presently adopted, 

see Equation (19). In order to perform the integration including the normal distribution in Equation (20), 

finite limits of integration must be determined. Assuming that the average fiber misalignment angle (for 

both  and ) is  = 0, Figure 5 shows the probability density functions for various value of the standard 

deviation, s. It is clear that it is sufficient to utilize $4s as the limits of the infinite integrations in Equation 

(20) since wider limits will have negligible contributions (see Figure 3). 
NASA/TM—2014-218313 11 

Figure 5.—Normal probability density functions versus a misalignment angle are 

shown for various values of the standard deviation, s, where the mean is taken to 

be zero. 

5.1 Polymer Matrix Composite 

The first application considers the 0.6 unidirectional graphite/epoxy composite described previously 

in Section 4.0. The same constituent properties and same 40 by 40 subcell discretization of the repeating 

unit cell, as shown in Figure 4, have been considered. 

The HFGMC theory has been extensively validated versus experiment data, and extensively verified 

versus detailed finite element models, in both the linear and nonlinear regimes for polymer, metal, and 

ceramic matrix composites (c.f., Aboudi et al., 2013; Bednarcyk et al., 2004; Pineda et al., 2013; 

Matzenmiller and Kurnatowski, 2010). In order to verify the fully random fiber orientation averaging 

procedure, described in Section 3.3, comparison has been made to results for the effective Young's 

modulus and effective Poisson's ratio based on the equations presented by Christensen and Waals (1972). 

These equations are as follows, 

3[2 (8 12 7) 2( )]

[ (4 8 4) ][ (4 4 1) 6( )]

12 23 12 23

2

11 12

12 23 12 23

2

12 23 11 12

2

11 12

3      

         



 "  "   

 "  "   "  "     E K G GE K E K G G E D (32) 

4 (16 24 14) 4( )

(4 16 6) 4( )

12 23 12 23

2

11 12

12 23 12 23

2

11 12

3      

     



 "  "   

 "  "   

"  E K G G

E K G G

D (33) 

where  E11 ,  "12 ,  K23 ,  G12 , and  G23 are the effective composite axial Young's modulus, axial Poisson's 

ratio, transverse plane strain bulk modulus, axial shear modulus, and transverse shear modulus, 

respectively. Note that these equations assume transverse isotropy, with five independent material 

constants. Christensen and Waals (1972) used the concentric cylinder assemblage (CCA) model to predict 

the effective composite properties for use in Equations (32) and (33). Figure 6 compares the predicted 

fully random (isotropic) effective graphite/epoxy composite properties predicted by the present fully 

0

2

4

6

8

10

12

-75 -60 -45 -30 -15 0 15 30 45 60 75

p(

)

 (deg.)

s = 2

s = 4

s = 6

s = 8

s = 10

s = 15

s = 20

s = 25
NASA/TM—2014-218313 12 

random averaging approach implemented within HFGMC with the Christensen and Waals (C&W) (1972) 

equations (32) and (33), as a function of fiber volume fraction, normalized by the isotropic matrix 

properties. In addition to the CCA model, predictions are shown for the Self-Consistent Scheme (Hill, 

1965), the Mori-Tanaka (1973) Method, the Method of Cells (MOC) (Aboudi, 1989), and HFGMC, 

wherein the C&W averaging equations have been used. Note that, in the case of the CCA model, as done 

by Christensen and Waals (1972), the lower bound estimates were used for the composite transverse 

properties. As shown, the present averaging approach matches closely with these other predictions. It 

should be noted, however, that these other, simpler micromechanics theories cannot provide a reliable 

estimate of the local field distributions, and they all lack shear coupling (see Aboudi et al., 2013). Slight 

(<1 percent) differences exist between the present fully random averaging approach and the C&W 

averaging equations when using HFGMC with the square fiber packing arrangement represented by the 

employed repeating unit cell (see Figure 4) because, as is well-known, such a packing arrangement leads 

to six rather than five independent material constants associated with transverse isotropy. As mentioned 

previously, the C&W averaging equations assume transverse isotropy, whereas the present averaging 

approach includes no such assumption. While the present averaging approach, involving concentration 

tensors, produces accurate composite effective properties, the strength of the method involves its ability to 

predict local fields. 

The remainder of the PMC results focus on the 0.6 volume fraction composite shown in Figure 4. 

Figure 7 exhibits the predicted initial damage surfaces for the von Mises criterion, Equation (28), and the 

damage-based criterion, Equation (30), for various values of the standard deviation of the misalignment 

angles, s and s under normal loading in the fiber (x1) and transverse (x2) directions. Also shown in the 

figure are the envelopes for the fully random case. In the absence of misalignment (s = 0), the envelopes 

exhibit the expected shape in that the strength in the fiber direction (x1) is more than 40 times that in the 

transverse (x2) direction. As can be seen in Figure 7(a) and (b), the effect of an in-plane misalignment 

angle, , distribution, rotates the initial damage envelope such that the composite becomes weaker in 

the axial direction and stronger in the transverse direction. The fully random envelope, which represents 

an isotropic material, is nearly elliptical. The effect of the hydrostatic stress, which is included in the 

damage-based criterion through the trixiality function, Rv, can be clearly observed by comparing 

Figure 7(a) and (b). 

Figure 6.—Comparison of the present fully random averaging approach within HFGMC with the Christensen and 

Waals (C&W) (1972) Equations (32) and (33) using various micromechanics theories for the predicted normalized 

fully random isotropic graphite/epoxy composite properties as a function for fiber volume fraction. (a) Effective 

Young's modulus, E3-D. (b) Effective Poisson's ratio, "3-D. 

0

2

4

6

8

10

12

14

16

18

20

0.0 0.2 0.4 0.6 0.8 1.0

E3-D/Ematrix

Fiber Volume Fraction

Self-Consistent Scheme (C&W Eq.)

Mori-Tanaka (C&W Eq.)

MOC (C&W Eq.)

CCA (C&W Eq.)

HFGMC (Present Averaging)

HFGMC (C&W Eq.)

0.55

0.60

0.65

0.70

0.75

0.80

0.85

0.90

0.95

1.00

0.0 0.2 0.4 0.6 0.8 1.0

"3-D/"matrix

Fiber Volume Fraction

Self-Consistent Scheme (C&W Eq.)

Mori-Tanaka (C&W Eq.)

MOC (C&W Eq.)

CCA (C&W Eq.)

HFGMC (Present Averaging)

HFGMC (C&W Eq.)

(a) (b) 
NASA/TM—2014-218313 13 

Figure 7.—Initial damage surfaces for the graphite/epoxy composite. (a) von Mises criterion, Equation (28), in the 

11–22 plane, where an in-plane fiber misalignment angle, , has been considered. (b) Damage-based criterion, 

Equation (30), in the 11–22 plane, where an in-plane fiber misalignment angle, , has been considered. (c) von 

Mises criterion, Equation (28), in the 11–22 plane, where both an in-plane fiber misalignment angle, , and an outof-plane misalignment angle, , have been considered. (d) Damage-based criterion, Equation (30), in the 11–22 

plane, where both an in-plane fiber misalignment angle, , and an out-of-plane misalignment angle, , have been 

considered. 

Figure 7(c) and (d) show the initial damage surfaces in the presence of equal in-plane, , and out-ofplane, , misalignment angle distributions. It can be seen that the effect of incorporating misalignment out 

of the x1-x2 plane is to reduce the size of the failure envelopes significantly as the reinforcing influence of 

the fiber is lost for both the axial and transverse directions. 

Figure 8 shows the initial damage surfaces for the damage based criterion in the 11–12 plane for 

various values of the standard deviation of the misalignment angles. Also included is the fully random 

case. Note that the surfaces based on the von Mises criterion are not shown as they are quite similar to 

those shown in Figure 8 since the hydrostatic stress effect is relatively small. In fact, the greatest value of 

the triaxiality function, Rv (see Equation (29)), throughout the matrix material for all points on all damage 

-60

-40

-20

0

20

40

60



-6 -4 -2 0 2 4 6

11/Ym

22/Ym

Note: Scales are not equal s = 0%

fully 

random

s = 10%

s = 25%

-60

-40

-20

0

20

40

60

-6 -4 -2 0 2 4 6 11/cr-m

22/cr-m

Note: Scales are not equal s = 0%

fully 

random

s = 10%

s = 25%

-60

-40

-20

0

20

40

60



-6 -4 -2 0 2 4 6

11/Ym

22/Ym

Note: Scales are not equal s = s = 0%

fully 

random

s = s = 5%

s = s = 25%

s = s = 10%

-60

-40

-20

0

20

40

60

-6 -4 -2 0 2 4 6 11/cr-m

22/cr-m

Note: Scales are not equal s = s = 0%

fully 

random

s = s = 10%

s = s = 25%

s = s = 5%

(a) (b) 

(c) (d) 
NASA/TM—2014-218313 14 

surfaces shown in Figure 8 is 1.16, whereas, for the damage surfaces shown in Figure 7, this value is 5.15. 

It can be observed in Figure 8 that the initial damage envelopes both rotate and reduce size as the fiber 

misalignment angle distribution standard deviations are increased. It should be noted that, contrary to the 

cases shown in Figure 7, where additional fiber misalignment causes only minor strengthening in the 

transverse direction, in Figure 8, it can be observed that significant longitudinal shear strengthening 

occurs with additional in-plane misalignment (as one might expect). 

It is of interest to examine the variation of various composite properties with respect to the in-plane 

misalignment angle . Figure 9 shows a number of such properties, normalized with respect to their 

values in the absence of misalignment. This figure shows that, while the effective axial Young's modulus,  E1 , decreases with additional misalignment, the transverse Young's modulus,  E2 , and the axial shear 

modulus,  G12 , increase. The axial and transverse Young's moduli are relatively insensitive to small 

amounts of misalignment, whereas the axial shear modulus is very sensitive. It is worth mentioning that, 

by applying an in-plane transformation, the initially orthotropic composite material becomes monoclinic 

and the effective stiffness tensor component /

1212 C , see Equation (13), is greatly influenced by this 

transformation. As a result of the employed averaging procedure, see Equation (20), the additional 

monoclinic terms in the effective stiffness matrix vanish (assuming the mean misalignment angle,  = 0), 

but the strong effect remains on  C1212 , which is the effective axial shear modulus,  G12 . Note, however, 

that, as mentioned in the Introduction, the dependence of the effective elastic properties on the 

misalignment could also be determined using a macromechanical analysis (provided the presented 

probability-weighted averaging procedure were performed). The influence on damage initiation surfaces, 

however, could not be determined. 

Also shown in Figure 9 are the variations of the areas of the damage initiation envelopes predicted 

using the von Mises, Avm, and damage-based, AD, criteria in the 11–22 plane. The decreasing trends in 

these values indicate the overall loss of strength of the material in this plane due to fiber misalignment. 

However, for very large standard deviations of the misalignment angle, this trend appears to reverse. 

Figure 8.—Initial damage surfaces for the graphite/epoxy composite. (a) Damage-based criterion, Equation (30), in 

the 11–12 plane, where an in-plane fiber misalignment angle, , has been considered. (b) Damage-based 

criterion, Equation (30), in the 11–12 plane, where both an in-plane fiber misalignment angle, , and an out-ofplane misalignment angle, , have been considered. 

-60

-40

-20

0

20

40

60

-4 -2 0 2 4 11/cr-m

12/cr-m

Note: Scales are not equal s = 0%

fully 

random

s = 10%

s = 25%

-60

-40

-20

0

20

40

60

-4 -2 0 2 4 11/cr-m

12/cr-m

Note: Scales are not equal s = s = 0%

fully 

random

s = s = 10%

s = s = 25%

s = s = 5%

(a) (b) 
NASA/TM—2014-218313 15 

Figure 9.—Effect of in-plane fiber misalignment on the graphite/epoxy 

composite effective axial shear modulus  G12 , transverse Young's modulus 

 E2 , axial Young's modulus  E1 , area of the von Mises criterion damage 

initiation envelope Avm, area of the damage-based criterion damage 

initiation envelope AD, and strength, all normalized with respect to their 

pristine values. 

A final curve plotted in Figure 9 shows the impact of fiber misalignment on the predicted strength 

of the composite loaded in axial tension, as dictate by the maximum von Mises stress, ( ) eq , see 

Equation (28), in the matrix constituent material. That is, this ‘strength’ corresponds to the global uniaxial 

tensile stress for damage initiation in the matrix, predicted based on the maximum local von Mises stress 

throughout the matrix. It can be readily observed that the impact of fiber misalignment on this axial 

strength, which can only be predicted by a micromechanical analysis, is substantial (for example, at 

s = 4%, 16 percent reduction), and much greater than the impact on the axial Young's modulus,  E1 (for 

example, at s = 4%, only 1.2 percent reduction). Initially, as small amounts of misalignment are 

introduced, this curve decreases rapidly with increasing s due to increased stress (particularly shear) in 

the matrix. At greater values of s, the normal distribution flattens (see Figure 5), and the strength curve 

in Figure 9 levels off, approaching the value for fully random fiber in the x1-x2 plane (0.124). 

In order to further investigate effect of the fiber misalignment on the composite, Figure 10 shows 

surface plots of AD, Avm, the aforementioned von Mises stress based strength, and the effective axial Young's 

modulus (  E1 ), all normalized by their pristine values as in Figure 9. By adding the effect of out-of-plane 

misalignment, , it can be observed that the effect of the latter misalignment is greater than that of in-plane 

misalignment on AD and Avm (Figure 10(a) and (b)). As for the strength and effective axial Young's modulus, 

it appears that the effects of both in-plane and out-of-plane misalignments are similar. It is also clear that, 

across all values of both misalignment angle standard deviations, the effect on the effective axial Young's 

modulus, which could be calculated using the macro approach, is significantly less than the other measures 

of composite performance. The fall off in the effective Young's modulus is also much less steep. 

0.0

0.5

1.0

1.5

2.0

2.5

3.0

0 5 10 15 20 25

Normalized Property

s

G12*

E2*

E1* Avm

Strength

AD
NASA/TM—2014-218313 16 

 

Figure 10.—Surface plots as a function of the standard deviations s and s. (a) The area, AD, of the damage initiation 

envelopes predicted using the damage-based criterion in the 11–22 plane. (b) The area, Avm, of the damage 

initiation envelopes predicted using the von Mises criterion in the 11–22 plane. (c) The strength of the composite 

loaded axially, as dictate by the maximum von Mises stress, ( ) eq . (d) The effective axial Young's modulus,  E1 . 

Figure 11 shows the normalized von Mises equivalent stress distributions in the composite with no 

fiber misalignment (pristine) subjected to global uniaxial stress loading in the fiber and transverse 

directions. Also shown is the corresponding distribution for the case where the fibers are fully randomly 

oriented, whose response, as mentioned previously, is isotropic, which yields identical equivalent stress 

distributions irrespective of the orientation of the applied uniaxial loading. Comparing Figure 11(a) to 

Figure 11(c), it is clear that the stresses in the fiber increase and the stresses in the matrix decrease by 

introducing random fiber orientations. The opposite result can be observed when comparing Figure 11(b) 

to Figure 11(c). Thus the effect of random fiber orientation on the composite behavior is significant. 

0

2 4 6 8

10

15

20

25

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

0 2 4 6 8

10 15 20 25

s (%)

AD

s (%)

0.9-1

0.8-0.9

0.7-0.8

0.6-0.7

0.5-0.6

0.4-0.5

0.3-0.4

0.2-0.3

0.1-0.2

0-0.1

0

2 4 6 8

10

15

20

25

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

0 2 4 6 8

10 15 20 25

s (%)

Avm

s (%)

0.9-1

0.8-0.9

0.7-0.8

0.6-0.7

0.5-0.6

0.4-0.5

0.3-0.4

0.2-0.3

0.1-0.2

0-0.1

0

2 4 6 8

10

15

20

25

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

0 2 4 6 8

10 15 20 25

s

Strength

s

0.9-1

0.8-0.9

0.7-0.8

0.6-0.7

0.5-0.6

0.4-0.5

0.3-0.4

0.2-0.3

0.1-0.2

0-0.1

0

2 4 6 8

10

15

20

25

0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

0 2 4 6 8

10 15 20 25

s (%)

E1*

s (%)

0.9-1

0.8-0.9

0.7-0.8

0.6-0.7

0.5-0.6

0.4-0.5

0.3-0.4

0.2-0.3

0.1-0.2

0-0.1

(a) (b) 

(c) (d) 
NASA/TM—2014-218313 17 

 

Figure 11.—von Mises equivalent stress distributions normalized by the applied global uniaxial stress for the 

graphite/epoxy composite. (a) Pristine composite subjected to axial stress 11 . (b) Pristine composite subjected to 

transverse stress 22 . (c) Fully random fiber distribution subjected to either axial or transverse stress. 

5.2 Ceramic Matrix Composite 

A second application considers a 0.26 fiber volume fraction unidirectional SiC/SiC composite. The 

isotropic SiC fiber properties are as follows: Young's modulus = 400 GPa and Poisson's ratio = 0.17. The 

isotropic SiC matrix properties are: Young's modulus = 335 GPa and Poisson's ratio = 0.18. In addition, 

there is a BN-based compliant interfacial layer whose volume fraction is 0.1. Its isotopic properties are as 

follows: Young's modulus = 5 GPa and Poisson's ratio = 0.22. A 28 by 28 subcell discretization of the 

repeating unit cell, which now includes an explicit BN interfacial phase, has been considered, as shown in 

Figure 12. 

(a) (b) 

(c) 
NASA/TM—2014-218313 18 

Figure 12.—SiC/SiC CMC repeating unit cell 

divided into 28 by 28 subcells. The fiber 

volume fraction is 0.26 and the interfacial 

layer volume fraction is 0.1. 

Figure 13.—Initial damage surfaces for the SiC/SiC CMC composite. Damage-based criterion, Equation (30), 

in the 11–22 plane, where (a) an in-plane fiber misalignment angle, , has been considered, and (b) both 

an in-plane fiber misalignment angle, , and an out-of-plane misalignment angle, , have been considered. 

Figure 13 shows the predicted initial damage surfaces for the CMC for the damage-based criterion, 

Equation (30), for various values of the standard deviation of the misalignment angles, s and s, under 

normal loading in the fiber (x1) and transverse (x2) directions. Also shown in the figure are the envelopes 

for the fully random case. The CMC behaves quite differently than the PMC previously considered. In the 

longitudinal (fiber) direction, the CMC exhibits much less mismatch in properties between the fiber and 

matrix, and the interfacial compliant layer, which is in parallel with the fiber and matrix, has only a minor 

effect. In contrast, in the transverse direction, the interfacial compliant layer (now in series with the fiber 

and matrix) introduces a great deal of property mismatch. This is exactly opposite to the PMC, wherein 

the longitudinal direction property mismatch it much larger (80:1) compared to the transverse mismatch 

-1

-0.5

0

0.5

1

-1 -0.5 0 0.5 1 11/cr-m

22/cr-m

s = 0%

fully 

random

s = 10% s = 25%

-1

-0.5

0

0.5

1

-1 -0.5 0 0.5 1 11/cr-m

22/cr-m

s = s = 0%

fully 

random

s = s = 10%

s = s = 25%

(a) (b) 
NASA/TM—2014-218313 19 

(4:1). As a result, the effect of relatively small amounts of fiber misalignment on the CMC is considerably 

less compared to the PMC. Furthermore, also in contrast to the previously considered PMC, the addition 

of out-of-plane fiber misalignment, see Figure 13(b), does not exhibit a significant effect. 

Figure 14 shows the variation of various composite properties with respect to the in-plane 

misalignment angle , normalized with respect to their values in the absence of misalignment. Note that 

the scale of the ordinate axis in Figure 14 is significantly smaller than that in Figure 9 for the PMC. As 

such, it is clear that the effect of misalignment on the various properties for the CMC is considerably less 

pronounced. The strength remains the property most affected by the fiber misalignment, but the trend 

reversals in the Avm and AD, observed in Figure 9 are now absent. 

Finally, the normalized local von Mises equivalent stress distributions in the CMC are shown in 

Figure 15. Figure 15(a) and (b) show the stress distributions in the pristine composite loaded in the axial 

and transverse directions, respectively. The compliant layer can be clearly observed. The induced stress 

concentrations are significantly higher in the transversely loaded case, Figure 15(b), due to the compliant 

layer's influence, which isolates the fiber. Figure 15(c) shows the fully random CMC case, which, 

compared to the axially loaded pristine case, Figure 15(a), exhibits greater stress in the matrix and lower 

stress in the fiber. 

Figure 14.—Effect of in-plane fiber misalignment on the SiC/SiC CMC composite 

effective axial shear modulus  G12 , transverse Young's modulus  E2 , axial 

Young's modulus  E1 , area of the von Mises criterion damage initiation envelope 

Avm, area of the damage-based criterion damage initiation envelope AD, and 

strength, all normalized with respect to their pristine values. 

0.70

0.75

0.80

0.85

0.90

0.95

1.00

1.05

1.10

1.15

1.20

0 5 10 15 20 25

Normalized Property

s (%)

G12*

E2*

E1*

Avm

Strength

AD
NASA/TM—2014-218313 20 

 

Figure 15.—von Mises equivalent stress distributions normalized by the applied global uniaxial stress for the SiC/SiC 

CMC. (a) Pristine composite subjected to axial stress 11 . (b) Pristine composite subjected to transverse stress 

22 . (c) Fully random fiber distribution subjected to either axial or transverse stress. 

6.0 Conclusions 

A micromechanical formulation has been presented that is capable of predicting the effects of statistical 

distributions of fiber misalignment on the behavior of composite materials. The key unique feature is the 

ability of the method to predict not only effective composite properties, but also local fields. This is possible 

because this approach is based on the establishment of the composite concentration tensors, which include 

the effects of the fiber misalignment. It is thus possible (although not done herein) to predict the complete 

nonlinear response of the composite, including interfacial debonding and progressive damage. This 

approach was applied to the case of linear elasticity and perfect bonding, excluding progressive damage. 

However, in the absence of fiber misalignment, imperfect debonding and progressive damage have been 

addressed by Aboudi et al. (2013). In all cases presented, normal distributions of in-plane and out-of-plane 

fiber misalignment angles have been considered. It has been shown that the effect of these fiber 

misalignments are more significant in PMCs as compared to CMCs due to the fact that PMCs have 

significantly greater property mismatch. In the case of the PMC, the substantial variations of the initial 

damage surfaces with the standard deviation of the misalignment angles have been illustrated. Fully random 

(a) (b) 

(c) 
NASA/TM—2014-218313 21 

fiber orientation has also been considered and shown to have a major impact on the initial damage surfaces 

and local field distributions when compared to the pristine composite. 

It should be noted that the presented approach to capturing the effects of fiber misalignment is 

applicable to any micromechanical approach that can accurately predict the strain concentration tensor of 

a composite, including the finite element method. Due to the proven ability of the HFGMC model to 

analyze various nonlinear effects, the proposed method could also be extended to consider material 

nonlinearity such as damage and inelasticity within the constituents and to consider the thermalmechanical response of composites. Additional potential applications include analysis of composites with 

given fiber paths, short/nano fiber composites with known fiber alignment/waviness (c.f., Yu et. al, 2011) 

and multiscale analysis within laminated plate theories and finite element analyses, which could be used 

to capture misalignment effects on the response of woven composites. Application to biological materials, 

as has been done with the microsphere model, is also possible. 

References 

Aboudi, J. (1989) “Micromechanical Analysis of Composites by the Method of Cells” Applied Mechanics 

Reviews, 42, 193-221. 

Aboudi, J., Arnold, S.M., and Bednarcyk, B.A. (2013) Micromechanics of Composite Materials A 

Generalized Multiscale Approach, Elsevier, New York. 

Alastrué, V., Martinez, M.A., Doblaré, M., and Menzel, A. (2009) “Anisotropic micro-sphere-based finite 

elasticity applied to blood vessel modelling,” Journal of the Mechanics and Physics of Solids, 57, 178-

203. 

Alastrué, V., Sáeza, P., Martinez, M.A., Doblaré, M. (2010) “On the use of the Bingham statistical 

distribution in microsphere-based constitutive models for arterial tissue,” Mechanics Research 

Communications, 37, 700–706. 

Arridge, R.G.C. (1985) An Introduction to Polymer Mechanics, Taylor and Francis, London. 

Ataabadia, A.K., Hosseini-Toudeshkyb, H., and Rada, S.Z. (2014) “Experimental and analytical study on 

fiber-kinking failure mode of laminated composites,” Composites: Part B, 61, 84-93. 

Barwick, S.C. and Papathanasiou, T.D. (2003) “Identification of Fiber Misalignment in Continuous Fiber 

Composites,” Polymer Composites, 24, 475-486. 

Basu, S., Waas, A.M., and Ambur, D.R. (2006) “Compressive failure of fiber composites under multi-axial 

loading,” Journal of the Mechanics and Physics of Solids, 54, 611-634. 

Bednarcyk, B.A., Arnold, S.M., Aboudi, J., and Pindera, M.-J. (2004) “Local field effects in titanium matrix 

composites subject to fiber-matrix debonding,” International Journal of Plasticity, 20, 1707-1737. 

Budiansky, B. (1983) “Micromechanics,” Computers and Structures, 16, 3-12. 

Budianksky, B. and Fleck, N.A. (1993) “Compressive Failure of Fibre Composites,” Journal of the 

Mechanics and Physics of Solids, 41, 183-211. 

Christensen, R.M. and Waals, F.M. (1972) “Effective Stiffness of Randomly Oriented Fibre Composites,” 

Journal of Composite Materials, 6, 518-532. 

Göktepe, S. and Miehe, C. (2005) “A micro-macro approach to rubber-like materials–Part III: the microsphere model of anisotropic Mullins-type damage,” Journal of the Mechanics and Physics of Solids, 53, 

2259-2283. 

Herakovich, C.T. (1998) Mechanics of Fibrous Composites, John Wiley and Sons, New York. 

Hill, R. (1965). “A Self-Consistent Mechanics of Composite Materials,” Journal of the Mechanics and 

Physics of Solids, 13, 213-222. 

Jelf, P.M. and Fleck, N.A. (1994) “The Failure of Composite Tubes Due to Combined Compression and 

Torsion,” Journal of Materials Science, 29, 3080-3084. 

Jones, R.M. (1975) Mechanics of Composite Materials, Hemisphere Publishing Corp., New York. 

Kugler, D. and Moon, T.J. (2001) “Identification of the Most Significant Processing Parameters on the 

Development of Fiber Waviness in Thin Laminates,” Journal of Composite Materials, 36, 1451-1479. 

Lemaitre, J. (2001) Handbook of Materials Behavior Models, Elsevier, New York. 
NASA/TM—2014-218313 22 

Lemaitre, J. and Chaboche, J.L. (1990) Mechanics of Solid Materials, Cambridge University Press, 

Cambridge. 

Li, Y, Simon, J.-W., Reese, S. (2013) “Microsphere-based homogenization for unidirectional fiber 

reinforced composites,” 5th GACM Colloquium on Computational Mechanics, Hamburg, Germany. 

Liu, D., Fleck, N.A., and Sutcliffe, M.P.F. (2004) “Compressive strength of fibre composites with random 

fibre waviness,” Journal of the Mechanics and Physics of Solids, 52, 1481-1505. 

Luo, J.J. and Daniel, I.M. (2003) “Characterization and Modeling of Mechanical Behavior of Polymer/Clay 

Nanocomposites,” Composites Science and Technology, 63, 1607-1616. 

Matzenmiller, A. and Kurnratowski, B. (2010) “A Comparison of Micromechanical Models for the 

Homogenization of Microheterogeneous Elastic Composites,” Advances in Mathematical Modeling and 

Experimental Methods for Materials and Structures The Jacob Aboudi Volume, R. Gilat and L. BanksSills (Eds.), Springer, New York. 

Menzel, A. and Waffenschmidt, T. (2009) “A microsphere-based remodelling formulation for anisotropic 

biological tissues,” Philosophical Transactions of the Royal Society A, 367, 3499-3523. 

Miehe, C., Göktepe, S., and Lulei, F. (2004) “A micro-macro approach to rubber-like materials–Part I: the 

non-affine micro-sphere model of rubber elasticity,” Journal of the Mechanics and Physics of Solids, 52, 

2617-2660. 

Miehe, C. and Göktepe, S. (2005) “A micro-macro approach to rubber-like materials–Part II: the microsphere model of finite rubber viscoelasticity,” Journal of the Mechanics and Physics of Solids, 53, 2231-

2258. 

Modniks, J. and Andersons, J. (2013) “Modeling the non-linear deformation of a short-flax-fiber-reinforced 

polymer composite by orientation averaging,” Composites: Part B, 54, 188-193. 

Mori, T. and Tanaka, K. (1973) “Average Stress in Matrix and Average Energy of Materials With Misfitting 

Inclusions,” Acta Metallurgica, 21, 571-574. 

Murtada, S.I., Kroon, M., and Holzapfel, G.A. (2010) “Modeling the dispersion effects of contractile fibers 

in smooth muscles,” Journal of the Mechanics and Physics of Solids, 58, 2065-2082. 

Numayr, K.S. and Al Rjoub, Y.S. (2013) “Two analogous methods for estimating the compressive strength 

of fibrous composites,” Composites: Part B, 30, 290-296. 

Pineda, E.J., Bednarcyk, B.A., Waas, A.M., and Arnold, S.M. (2013) “Progressive failure of a unidirectional 

fiber-reinforced composite using the method of cells: Discretization objective computational results,” 

International Journal of Solids and Structures, 20, 1203-1216. 

Schultheisz, C.R. and Waas, A.M. (1996) “Compressive Failure of Composites, Part I: Testing and 

Micromechanical Theories,” Progress in Aerospace Sciences, 32, 1-42. 

Shuart, M.J. (1989) “Failure of Compression-Loaded Multidirectional Composite Laminates,” AIAA 

Journal, 27, 1274-1289. 

Waas, A.M. and Schultheisz, C.R. (1996) “Compressive Failure of Composites, Part II: Experimental 

Studies,” Progress in Aerospace Sciences, 32, 43-78. 

Waffenschmidt, T., Menzel, A., and Kuhl, E. (2012) “Anisotropic density growth of bone—A 

computational micro-sphere approach,” International Journal of Solids and Structures, 49, 1928–1946. 

Wisnom, M.R. (1990) “The Effects of Fibre Misalignment on the Compressive Strength of Unidirectional 

Carbon Fibre/Epoxy,” Composites, 21, 403-407. 

Yu, J., Lacy, T.E., Toghiani, H., Pittman Jr., C.U., and Hwang, Y. (2011) “Classical micromechanics 

modeling of nanocomposites with carbon nanofibers and interphase,” Journal of Composite Materials, 

45, 2401-2413. 

Yurgartis, S.W. (1987) “Measurement of Small Angle Fiber Misalignments in Continuous Fiber 

Composites,” Composites Science and Technology, 30, 279-293. 

Zhou, H.W., Yi, H.Y., Gui, L.L., Dai, G.M., Peng, R.D., Wang, H.W., and Mishnaevsky Jr., L. (2013) 

“Compressive damage mechanism of GFRP composites under off-axis loading: Experimental and 

numerical investigations,” Composites: Part B, 55, 119-127.