An analytical model for stress analysis of short fiber composites in power law creep matrix

International Journal of Non-Linear Mechanics 57 (2013) 39–49

Contents lists available at ScienceDirect

International Journal of Non-Linear Mechanics journal homepage: www.elsevier.com/locate/nlm

An analytical model for stress analysis of short fiber composites in power law creep matrix M. Mondali a,n, A. Abedian b a b

Department of Mechanical and Aerospace Engineering, Islamic Azad University, Science and Research Branch, Tehran, Iran Department of Aerospace Engineering, Sharif University of Technology, P.O. Box 11365-8639, Azadi Avenue, Tehran, Iran

art ic l e i nf o

a b s t r a c t

Article history: Received 23 July 2012 Received in revised form 27 May 2013 Accepted 8 June 2013 Available online 18 June 2013

The creep deformation behavior of short fiber composites has been studied by an approximate analytical model. A perfect fiber/matrix interfacial bond is assumed and a power law function is considered for describing the steady state creep behavior of the matrix material. The results obtained from the proposed analytical solution satisfy the equilibrium and constitutive creep equations. Also, a parametric study was undertaken to define the effects of geometric parameters on the steady state creep strain rate of short fiber composites. The present model is then validated using the results of finite element method. The predicted strain rate and stress components by the proposed analytical approach exhibit good agreement with the finite element results. & 2013 Elsevier Ltd. All rights reserved.

Keywords: Short fiber composites Steady state creep Stress analysis Finite element modeling Geometric factors

1. Introduction In recent years, the high-temperature creep behavior of short fiber composites has become a topic of considerable interest, primarily because these materials have a high potential for use in structural applications at elevated temperatures [1–11]. Kelly and Tyson [1,2] have used a frictional energy dissipation approach for strength behavior in metal matrix composites. They showed that the continuous fiber reinforced metal matrix composites display significantly better creep resistance and stress rupture life than the corresponding unreinforced matrix material [1–3]. Tension creep tests of 6061Al/20% SiC whisker composites have been conducted in the temperature range 505–644 K by Nieh [4] in 1984. He reported that the stress exponent for the steady state creep rate of the composite is approximately 20.5 and is not essentially a function of temperature while the stress exponent of the matrix material is 4. Nardone and Strife [5] examined the effects of stress and temperature on the creep behavior of SiC whisker reinforced 2124Al composite. They showed that the stress exponent changes from 8.4 at 450 K to 21 at 561 K. Dragon et al. [6] have investigated the creep behavior of other kind of metal matrix composites such as lead reinforced with nickel fibers. A review of creep data of discontinuous SiC–Al composites has been conducted by Mohamed et al. [7] which shows that the value

n

Corresponding author. Tel./fax: +98 21 44865239. E-mail address: [email protected] (M. Mondali).

0020-7462/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ijnonlinmec.2013.06.007

of the stress exponent is high and also variable. The conventional high temperature creep rupture tests on a commercially produced metal matrix composite and a wrought unreinforced equivalent alloy were performed by Greasly [8]. He reported that at typical creep service loads the performance of the reinforced alloy is superior. Also, some researchers have studied the mechanisms of longitudinal creep deformation and damage in laminated composites [9]. Ma and Tjong [10] have investigated the tensile creep behavior of SiC-whisker and particle reinforced 2124Al matrix composites and unreinforced 2124Al alloy at 573–673 K. Their results exhibited that the unreinforced 2124Al alloy has the stress exponents of 4.7–8.5 while the particulate and whisker-reinforced composites have the stress exponents of 9.4–17.1 and 11.5–18.2, respectively. Mondali et al. [11] highlighted the shortcomings of the previous FEM studies and combines the analytical findings with the power of FEM modeling. More specifically, they incorporated the fiber/matrix debonding parameter into the modeling and showed that in contrary to the available analytical results, the value of debonding parameter could only change in the range of 0–0.5. The increasing number of applications of short fiber composites makes it more important to understand and predict their creep characteristics and deformation mechanisms. It is well known that short fiber composites display significantly better creep resistance than the corresponding unreinforced matrix alloy. As it will be fully discussed in following sections, the unit cell of a short fiber composite can be divided into two separate regions which are considered as two continuous fiber composite models. Therefore,

40

M. Mondali, A. Abedian / International Journal of Non-Linear Mechanics 57 (2013) 39–49

a proper formulation for continuous or long fiber composite models should be derived, first. Note that in comparison with short fiber composites, creep modeling of continuous or long fiber composites is simpler because of the end stress of fibers in short fiber composites can hardly be determined. In recent years, some researchers have tried to analyze the short fiber composites using a new developed shear lag model namely imaginary fiber technique which is only presented for elastic problems [12–20]. A variety of analytical and numerical models have been developed to predict the second stage creep strain rate of short fiber composites [11,21–31]. These models mainly include the shear lag models [21–30], and the numerical models based on the finite element method (FEM) [11,31]. Compared with the numerical models, the shear lag model is simplest mathematically. However, because of simple assumption of the matrix uniform deformation rate, it cannot describe the variation of the stress distributions in the matrix and the fiber during the steady state creep deformation. Except for stress transfer from the matrix to the fiber, the shear lag model cannot well predict the composite creep strain rate. Therefore, the application of the shear lag model to short fiber composites has become limited over the time. In general, solution of the analytical model involves with two equilibrium equations and four constitutive equations [32]. In the shear lag models, only one of the equilibrium equations is used to solve the problem. Thus, the average axial stress in the fiber and matrix and the interfacial shear stress (i.e. sfz ; sm z , and τi ) as a function of fiber axial direction (z) could be obtained. In fact, instead of the constitutive equations, the simple velocity relations of the flowing matrix are used for defining the creep strain rate of the composite. Therefore, the shear-lag model cannot predict the stress components and also the displacement rates of the composite as functions of radial (r) and axial (z) coordinates. The above discussions show that a relatively simple analytical solution is still highly desired for describing the steady state creep behavior of short fiber composites. For such an approximate analytical solution the stress field components in the matrix and the fiber should satisfy all the constitutive and equilibrium equations considering the imposed boundary conditions. In general, with the axisymmetric assumption for problem geometry, the equilibrium and constitutive equations cannot be exactly satisfied even for the composites with elastic behavior of the matrix. However, some researchers have solved the elastic problem considering various approximations and simplifications made in the shear lag model [12–18]. Therefore, developing an analytical solution for creep behavior in short fiber composites is very complicated and should be obtained by FEM or other numerical methods. In 2009, Mondali et al. [33] developed a new analytical model based on the shear-lag theory for stress analysis and prediction of the steady state creep deformation in short fiber composites. In this new analytical model, a perfect bonding at the fiber/matrix interface was assumed and the steady state creep behavior of the matrix was described by an exponential law. Note that generally the creep behavior of the most metals and alloys at high temperatures is described by a power law. Consequently, some experiments data obtained for the metal matrix composites with matrix creeping are analyzed on the basis of the power law [25]. Hence, it is very essential to develop a new analytical solution for defining the second stage creep strain rate of short fiber composites when the steady state creep behavior of the matrix has been described by a power law. The main objective of the present study is thus to develop such an analytical model for stress analysis of short fiber composites in power law creep matrix. This model should satisfy the constitutive, and equilibrium equations and the overall boundary conditions. Since the geometrical shape and arrangement of fibers have a strong effect on the creep strain rate

of short fiber composites, a parametric study is undertaken to study the effects of these factors on the steady state creep strain rate of short fiber composites. Here, an axisymmetric unit cell representing a fiber with its surrounding matrix as two coaxial cylinders is assumed. For verification of the solution method, the Al6061/20%SiC composite [31] is selected as a case study and the results will be compared with the FEM results.

2. Problem definition 2.1. Composite model The cylindrical unit cell shown in Fig. 1 has been used by many researchers [12–18,25–27,29] to model a short fiber composite. In this model a cylindrical fiber with a radius a and a length 2l is surrounded by a coaxial cylindrical matrix with an outer radius b and a length 2l′. The volume fraction and aspect ratio of the fiber 2 0 are defined as f ¼ ða2 lÞ=ðb l Þ and s ¼ l=a, respectively. Also, 0 k ¼ l a=lb is a parameter related to the geometry of the unit cell. Axial stress s0 is uniformly applied on the end faces of the unit cell, i.e. z ¼ 7 l′. The cylindrical coordinate system (r; θ; z) is used with the origin shown in Fig.1. Due to symmetry in geometry, loading, and boundary conditions, the analysis is performed only on half of the unit cell, i.e. 0 ≤ z ≤l′. In this analysis, the following assumptions are made: (i) Steady state condition of stress is assumed. (ii) Elastic deformations are small and negligible as compared to creep deformations. (iii) The fibers behave elastically during the analysis and the steady state creep behavior of the matrix is described by a power law as given in Eq. (1), ε_ e ¼ Asne

ð1Þ

(iv) Perfect bonded at fiber/matrix interface is considered.

2.2. Equilibrium and constitutive equations The governing equilibrium equations for the axisymmetric problem considering the cylindrical coordinates (r; θ; z) are obtained as ∂sz ∂τrz τrz þ þ ¼0 ∂z ∂r r

ð2aÞ

∂sr ∂τrz sr −sθ þ þ ¼0 ∂r ∂z r

ð2bÞ

The generalized constitutive equations and the strain ratedisplacement rate relations for creep deformation of the matrix material in r; θ and z directions [32,33] are ε_ r ¼

ε_ e ∂u_ ½2sr −sθ −sz  ¼ ∂r 2se

Fig. 1. Representation of the unit cell.

ð3aÞ

M. Mondali, A. Abedian / International Journal of Non-Linear Mechanics 57 (2013) 39–49

ε_ θ ¼

ε_ e u_ ½2sθ −sr −sz  ¼ r 2se

ð3bÞ

ε_ z ¼

_ ε_ e ∂w ½2sz −sr −sθ  ¼ ∂z 2se

ð3cÞ

γ_ rz ¼

_ ∂u_ ∂w 3_εe þ ¼ τrz ∂z ∂r se

where the equivalent stress, se , is given by i1 1 h se ¼ pffiffiffi ðsr −sθ Þ2 þ ðsθ −sz Þ2 þ ðsz −sr Þ2 þ 6τrz 2 2 2

41

ð3dÞ

ð4Þ

and sr ; sθ ; sz and τrz are the stress components in the directions indicated by subscripts. Furthermore, ε_ e is the equivalent strain rate and is defined by Eq. (5), in terms of the strain rate components, pffiffiffi 2 ½ð_εr −_εθ Þ2 þ ð_εθ −_εz Þ2 þ ð_εz −_εr Þ2 þ 6_γ rz 2 1=2 ε_ e ¼ ð5Þ 3 where ε_ r ; ε_ θ ; ε_ z and γ_ rz are the strain rate components in the directions indicated by subscripts. Note that only the equilibrium equations (i.e. Eq.(2a, 2b)) are applicable to both the effective fiber and the matrix. Also, based on the available literature, solutions for the stress transfer problem between the fiber and the matrix require extensive numerical analysis to satisfy strictly the equilibrium equations, the stress– displacement rate relations, and the essential boundary conditions. However, simple analytical solutions are attainable when some of the above equations are satisfied approximately, and the rest of equations are satisfied exactly. Note that most important novelty of this research work is to developing such a simple analytical solution for short fiber composites in power law creep matrix. 2.3. Boundary conditions The basis of the proposed solution method here is somehow similar to the Hsueh's elastic solution known as the imaginary fiber technique [12–15]. In this approach, the composite is divided into two separate regions (i.e. regions I and II), see Fig. 2. For elastic problems of short fiber composites, the matrix connecting the fiber end to the surface of the unit cell is considered as an imaginary fiber, which has the elastic matrix properties. Hence, the stress transfer problem can be solved using two similar shearlag models in regions I and II. However, in the present problem shown in Fig. 2, region I contains a real elastic fiber with matrix creeping while in region II both the imaginary fiber and matrix are under creep. Therefore, a new proposed shear-lag model in region I must be developed which is not applicable in region II. The details of this approach will be discussed in next section. Here, since the unit cell should remain a right circular cylinder during the deformation, the iso-displacement condition is used for the present problem, see Fig. 2. Note that another type of the boundary conditions for stress analysis problems in fiber composites is the stress-free condition [17,18] which is not appropriate here. According to the iso-displacement condition, the overall boundary conditions including the surface conditions, the interface continuity conditions and the axial force equilibrium conditions, can be proposed as the following The surface conditions, i.e. r ¼ 0 ; r ¼ b ; z ¼ 0 and z ¼ l′, are given by Eq. (6a–6e). Note that here radial iso-displacement rate condition is used instead of the radial iso-displacement condition [18,33]. _ zÞ ¼ 0 uð0; _ zÞ ¼ u_ b uðb;

l ≤z ≤l′ 0 ≤z ≤l′

Fig. 2. Schematic illustration of two separate regions and the unit cell edges in the 2nd stage creep and elastic states.

_ 0Þ ¼ 0 wðr;

a ≤r ≤b

ð6cÞ

τrz ð0; zÞ ¼ τrz ðb; zÞ ¼ 0

0 ≤z ≤l′

ð6dÞ

τrz ðr; 0Þ ¼ τrz ðr; l′Þ ¼ 0

0 ≤r ≤b

ð6eÞ

As it is shown in Fig. 2, u_ b in Eq. (6b) is the prescribed radial displacement rate on the outer surface of the unit cell. Also, the boundary conditions on the fiber–matrix interface (i.e. at r ¼ a; 0 ≤z ≤l) are assumed as [18,33] _ zÞ ¼ 0 uða;

ð7aÞ

_ wða; zÞ ¼ 0

ð7bÞ

f τm rz ða; zÞ ¼ τrz ða; zÞ ¼ τi

ð7cÞ

f sm r ða; zÞ ¼ sr ða; zÞ ¼ sp

ð7dÞ

In reality, due to the elastic behavior of the fiber during the creep deformation of the composite, any matrix slipping on the fiber at the interface is prohibited. Therefore, both displacement rate components in the radial and axial directions at the interface are zero. Moreover, based on the shear-lag theory, the interfacial shear stress of the fiber and the matrix are assumed to be equal. The axial force equilibrium conditions on the cross-section of the unit cell in an average manner is presented by 2

2

b s0 ¼ a2 sfz þ ðb −a2 Þsm z

ð8Þ

where the superscripts m and f denote the respective matrix and fiber and the bar sign on the stress symbol (e.g. sfz ) denotes the average value over the cross-sections of the fiber or the matrix. Since the matrix only creeps, no index of form is given to _ w _ and ε_ e for avoiding complications. u;

3. Problem solution steps

ð6aÞ ð6bÞ

The exact analytical solution of the boundary-value problem defined by Eqs. (2) and (3), can hardly be obtained for the

42

M. Mondali, A. Abedian / International Journal of Non-Linear Mechanics 57 (2013) 39–49

composite cylinder model shown in Fig. 1. For solving such a complex problem, the creep behavior of region I is first determined without considering region II. Then using the appropriate boundary conditions, the creep strain rate of the full composite model is calculated. Accordingly, using the first equilibrium equation (i.e. Eq. (2a)) the shear stress in region I as a function of r and the interfacial shear stress (i.e. τi ðzÞ) are calculated for both the fiber and the matrix. Integrating Eq. (2a) over the cross-section of the fiber and then dividing by its cross sectional area π a2 yield the shear-lag equation as [19]: dsfz 2 ¼ − τi a dz

ð9Þ

where sf is the average value of sfz over the cross-section of the fiber, and τi is the interfacial shear stress, which is given in Eq. (7c). Since the axial stress gradient in fiber and matrix is assumed to be only a function of z, radial dependence of the shear stresses for the fiber and the matrix in region I can be obtained as [19] τfrz ¼ τm rz ¼

r τ a i a 2

b −a2

! 2 b −r τi r

Integrating Eq. (19) with respect to r from a to r and using Eq. (7b) will yield

ð13Þ

ð15Þ

Considering the above equations and also Eq. (3d), the shear strain rate (_γ rz ) at the fiber–matrix interface and at the outer surface can be obtained as given by Eq. (16a,16b), respectively. pffiffiffinþ1 γ_ rz jr ¼ a ¼ A 3 τi n ð16aÞ γ_ rz jr ¼ b ¼ 0

ð18Þ

ð11Þ

Substituting the above equation into Eq. (4) and using (7c) will give the equivalent stress at the interface as pffiffiffi se ða; zÞ ¼ 3τi ð14Þ Now, substituting Eq. (14) into (1) yields pffiffiffi ε_ e ða; zÞ ¼ Að 3 τi Þ n

∂u_ ¼0 ∂z

ð10Þ

Now, incorporating the boundary conditions on the fiber–matrix interface and substituting Eq. ((7a)–(7d)) into (3a)–(3c) and (4a)– (4c) give sm z ða; zÞ ¼ sp

The second approximation is introduced by assuming u_ to be independent of z coordinate in the matrix of region I. This assumption has been used by the researchers who have successfully solved the elastic problem by imaginary fiber technique [12–16,19,33]. This assumption is more accurate when the matrix creeps in the neighborhood of the elastic fibers. Therefore, u_ must satisfy the following condition,

Now substituting Eq. (3d) into (17) and then using Eq. (18) gives ! pffiffiffi nþ1 2 _ 3 ∂w aA b ¼ −r τi n ð19Þ 2 ∂r r ðb −a2 Þ

Eqs. (10 and 11) satisfy continuity of the shear stress at the interface (i.e. Eq. (7c)) as well as the free surface condition (i.e. Eq. (6d)). However, the interfacial shear stress (i.e. τi ðzÞ) still remains unknown. In the following, using the equation for matrix shear strain rate (i.e. Eq. (3d)) and applying the appropriate boundary conditions, _ is first calculated. Subsequently, the axial displacement rate (i.e. w) the matrix incompressibility concept during creep deformation is _ Then, incorporated to calculate the radial displacement rate (i.e. u). by the use of appropriate assumptions, the interfacial shear stress function can be found. Here, substituting Eq. (11) into Eq. (3d) gives ! 2 3_εe b aτi γ_ rz ¼ −r ð12Þ 2 se r b −a2

sm θ ða; zÞ ¼

the mentioned assumption, substituting Eqs. (14) and (15) into (12) will result ! pffiffiffi nþ1 2 3 aA b γ_ rz ¼ −r τi n ð17Þ 2 r ðb −a2 Þ

ð16bÞ

The first assumption of the solution method is introduced by substituting ε_ e ða; zÞ and se ða; zÞ instead of ε_ e and se into Eq. (12) for calculating γ_ rz in the matrix of region I as a function of τi and ‘r'. According to the mentioned approximation, γ_ rz at the interface (r ¼ a) and the outer surface (r ¼ b), i.e. Eq. (16a,b), is completely _ in the satisfied and will yield acceptable functions for u_ and w matrix of region I (i.e. a ≤ r ≤b). Also, as it will be shown later, the introduced errors by the assumption are negligible. Now, using

_ zÞ ¼ wðr;

pffiffiffi nþ1   aA 2 r 2 −a2 3 τi n b lnðr=aÞ− 2 2 ðb −a2 Þ

ð20Þ

Next, adding up the Eq. (3a–3c) gives the incompressibility condition as _ ∂u_ u_ ∂w þ þ ¼0 ∂r r ∂z

ð21Þ

Now, Eq. (20) is substituted into Eq. (21) and the obtained equation is then solved by the use of Eq. (7a). Therefore, u_ can be given by " !# pffiffiffi nþ1 2 anA 2 r 2 −a2 2b −a2 dτ 3 _ zÞ ¼ − uðr; þ r τi n−1 i b rlnðr=aÞ− 2 4 r dz 2ðb −a2 Þ ð22Þ Then substituting Eq. (22) into (18) gives τi n−1

dτi ¼C dz

ð23Þ

where parameter C could be found by substituting Eq. (6b) into (22) as u_ C ¼ − pffiffiffin−1 b 3 nA D1

ð24Þ

where D1 is a geometrical factor related to the fiber and the unit cell radii and can be written as " # 4 3 a 4b ln b=a 2 2 D1 ¼ −ð3b −a Þ ð25Þ 2 8b b −a2 Z

Also, substituting Eq. (6e) into (23) will give Z z τi n−1 dτi ¼ Cdz

τi 0

0

ð26Þ

Therefore, the shear stress at the interface is obtained as τi ðzÞ ¼ ðnCzÞ1=n

ð27Þ

Then, substituting Eq. (27) into (20) will result _ zÞ ¼ wðr;

pffiffiffi nþ1   anAC 2 r 2 −a2 3 z b lnðr=aÞ− 2 2 ðb −a2 Þ

ð28Þ

M. Mondali, A. Abedian / International Journal of Non-Linear Mechanics 57 (2013) 39–49

Next, with the use of Eq. (24), the axial displacement rate can be rewritten as   3au_ b r 2 −a2 2 _ zÞ ¼ − z ð29Þ b lnðr=aÞ− wðr; 2 2 D1 ðb −a2 Þ Also, by substituting Eqs. (23) and (24) into (22), the radial displacement rate can be written as " !# 2 3au_ b r 2 −a2 2b −a2 2 _ ¼ þ r ð30Þ b rlnðr=aÞ− uðrÞ 2 4 r 2D1 ðb −a2 Þ Now, using the conditions of constant radial displacement rate at the outer surface of the unit cell, i.e. u_ b , the constant volume of the model during creep deformation, and the rule of mixture equation, u_ b and the average axial stresses in the fiber and matrix (i.e. sfz and sm z ) will be obtained. Finally, using the rest of equations and applying the appropriate boundary conditions, the stress m m components in the matrix (i.e. sm r ; sθ ; sz ) and sp ðzÞ, which is the interfacial radial pressure, will be found. Considering Eqs. (29) and (30), it is clear that u_ b is still unknown. Here, the rule of mixture equation (i.e. Eq. (31)) is used to calculate the mentioned unknown. f f sz

þ

m ð1−f Þsz

¼ s0

ð31Þ f sz

m sz ,

and indicates where duple bar sign on the stress symbol, i.e. the average of the axial stress along both r and z directions. For f example, sz can be expressed as Z l 1 f sz ¼ sf dz ð32Þ l 0 z sfz ,

first. Substituting For this purpose, it is necessary to obtain Eq. (27) into (9) and then integrating the obtained equation with respect to z from l to z, give i 2C 1=n h ðnzÞ1þ1=n −ðnlÞ1þ1=n þ sfz jz ¼ l ¼− ð33Þ aðn þ 1Þ   is the end stress of the fiber. where sfz  z¼l In the previous shear-lag models, the end stress of fiber was usually assumed [20,23,34–36]. However, in this study the end stress can be obtained using the shear-lag theory in region II and applying the appropriate boundary conditions between two regions. Hence, the interfacial shear stress in region II must be determined, first. The behavior of this function is assumed to be linear in region (r ¼ a and l ≤z ≤ l′). Here, Eq. (27) will give the shear stress value at ðz ¼ lÞ while it is zero at (z ¼ l′), i.e. Eq. (6e). sfz

43

Therefore, the shear stress will be determined as τ′i ðzÞ ¼

ðl′−zÞ τ ðl′−lÞ l

ð34Þ

where the prime sign indicates that the equation is written for region II, which is so called the imaginary fiber region, see Fig. 3. Also, τl can be given by τl ¼ ðnClÞ1=n

ð35Þ

Therefore, Eq. (9) for the mentioned region can be written as ds′fz 2 ¼ − τ′i a dz

ð36Þ

where s′fz is the average axial stress in the imaginary fiber region. Substituting Eq. (34) into (36) and then integrating the result with respect to z from l to l' and applying the appropriate boundary condition (i.e. s′fz jz ¼ l′ ¼ s0 ) will yield   sfz 

z¼l

¼

ðl′−lÞ ðnClÞ1=n þ s0 a

ð37Þ

Now, by substituting Eq. (37) into (33), the average axial stress in the fiber can be obtained as sfz ¼ −

 2C 1=n ðnClÞ1=n  ðnzÞ1þ1=n þ ðn þ 1Þl′ þ ðn−1Þl þ s0 aðn þ 1Þ aðn þ 1Þ

ð38Þ

Also, using Eq. (38) in (8) gives the average axial stress in the matrix for region I, i.e. a ≤r ≤b and 0 ≤ z ≤l, as sm z ¼

2aC 1=n 2

ðn þ 1Þðb −a2 Þ

ðnzÞ1þ1=n −

aðnClÞ1=n 2

ðn þ 1Þðb −a2 Þ

  ðn þ 1Þl′ þ ðn−1Þl þ s0

ð39Þ Substituting Eq. (38) into (32) and then integrating with respect to z from 0 to l will give  ðnClÞ1=n l f þ s0 l′− sz ¼ ð40Þ 2n þ 1 a Since, sm z given by Eq. (39) only presents the axial matrix stress in region I, this unknown could not be calculated with the similar f m way as used for sz . This is because sz belongs to the total matrix of the composite model (i.e. regions I and II). Due to incompressibility condition of the unit cell, in the steady state creep stage the volume of the unit cell remains constant. Therefore, the relationship between the unit cell displacement rates at r ¼ b (i.e. u_ b ) and

Fig. 3. The analytical and FEM curves of the matrix shear and equivalent stresses at the fiber/matrix interface and the outer surface.

44

M. Mondali, A. Abedian / International Journal of Non-Linear Mechanics 57 (2013) 39–49

_ l′ ) could be found from Eq. (41) [33]. z ¼ l′(i.e. w 2l′u_ b _ l′ ¼ − w b

ð41Þ

Consequently, the composite strain rate can be obtained as ε_ c ¼

_ l′ 2u_ w ¼− b l′ b

ð42Þ m

Now, the average axial stress in the matrix (i.e. sz ) can be found by substituting Eq. (42) into (1) as  2u_ 1=n m sz ¼ − b ð43Þ Ab Then, by substituting Eqs. (40) and (43) into (31) and using Eq. (24), the only unknown parameter (i.e. u_ b ) can be obtained as 2 !1=n 3−n  Ab 4 f l bl 5 sn 1þ l′− ð44Þ u_ b ¼ − pffiffiffin−1 0 2 að1−f Þ 2n þ 1 2 3 D 1

Finally, substituting Eq. (44) into (42) gives 2 !1=n 3−n  f l bl 5 sn l′− ε_ c ¼ A 41 þ pffiffiffin−1 0 að1−f Þ 2n þ 1 2 3 D1

ð45Þ

The radial, circumferential and axial stresses in the matrix can be defined by the use of Eq. (3a,3b) as  2se ∂u_ u_ sm þ sm þ 2 ð46aÞ r ¼ z ∂r r 3_εe sm θ ¼

 2se ∂u_ u_ þ2 þ sm z r 3_εe ∂r

ð46bÞ

Now, using the above equations and also Eqs. (11), (2b), (21) result   _ ∂sm 4se ∂2 w ∂u_ u_ ∂ 2se z þ 2 þ − ∂r r ∂r 3_εe ∂r 3_εe ∂r∂z ! 2 a b dτi −r ¼0 ð47Þ þ 2 dz b −a2 r Finally, by substituting se as a function of ε_ e in Eq. (47) and incorporating Eqs. (1) and (3a–3c) and then integrating with respect to r from a to r and using Eqs. (13, 27), sm z ðr; zÞ is obtained as sm z ðr; zÞ ¼ sp þ hðr; zÞ

ð48Þ

where hðr; zÞ ¼

h i =AÞ1=n dr ð2_εr þ ε_ θ Þ ∂r∂ 2ð_εe3_ εe   aC r 2 −a2 2 − 2 ðnCzÞ−1þ1=n b lnðr=aÞ− 2 2 b −a R

r 4ð_εe =AÞ1=n ∂_εz a 3_εe ∂r

dr−

R

r a

ð49Þ

Also, the equations for ε_ r ; ε_ θ ; ε_ z and ε_ rz as functions of both r and zdirections can be easily obtained. Now, the average axial stress over the cross-section of the matrix is defined as [19] Z b 1 sm sm ð50Þ z ðzÞ ¼ z ðr; zÞð2π rÞdr 2 πðb −a2 Þ a Substituting Eq. (48) into (50) and then integrating with respect to r from a to b gives the interfacial radial pressure (i.e. sp ) as Z b 2 rhðr; zÞdr ð51Þ sp ¼ sm z − 2 b −a2 a where sm z was obtained as a function of z coordinate using Eq. (39).

Consequently, substituting the above equation into Eq. (48) gives the axial stress in the matrix of region I as Z b 2 m sm rhðr; zÞdr ð52Þ z ðr; zÞ ¼ sz þ hðr; zÞ− 2 b −a2 a Finally, substituting Eq. (52) into (46a, b) gives the radial and circumferential stresses in matrix of region I. Note that the integration of Eq. (52) can be done numerically for the selected test case.

4. Results and discussions To validate the present analytical model, the SiC/Al6061 composite is chosen as a test case. The elastic properties of material constituents are considered as Em ¼ 68:3 GPa, Ef ¼ 470 GPa, νm ¼ 0:345, and νf ¼ 0:17 where the subscripts m and f denote the matrix and fiber, respectively [31]. For the purposes of comparison, the finite element simulation using the commercial FEM code of ANSYS is carried out. The model geometry is chosen as shown in Fig. 1 and the surface conditions are applied as presented in Eq. (6a–6e). The axisymmetry approach with nonlinear quadratic element of PLANE 183 is used for FEM analysis. This element is a higher order eight-node element and has the capability of creep modeling. Results are initially presented for the base case simulation of the short fiber composite with a fiber volume fraction of 20%. The fibers have an aspect ratio of 5, as denoted by Nieh in TEM studies [4], and the unit cell aspect ratio is also set to 5 as suggested by the scanning electron microscopy measurements performed by Morimoto et al. [27]. A simulated stress of 80 MPa is applied at a temperature of 561 K. Under these conditions, the steady state creep constants of the matrix, A and n, in Eq. (1) are considered as A ¼ 2:0  10−13 and n ¼ 4:0, which are also in accordance with the suggestions made in [4,31]. Here, Norton implicit creep equation can be employed for creep modeling of the matrix material in ANSYS software. 4.1. Verification of the stress field components To illustrate the results, it is important to verify the most critical point of the solution method, first. This is the relationship between the equivalent stress at the interface ðse ða; zÞÞ and the interfacial shear stress ðτi Þ along the fiber length as given by Eq. (14). The results of the analytical and FEM calculations of the shear and equivalent stresses along the fiber length at the interface and also outer surface of the unit cell are shown in Fig. 3. Interestingly, an excellent correlation between the results is observed. Also, the relationship between se ða; zÞ and τi described by Eq. (14), is clearly noticeable. Note that the mentioned relationship was normally assumed in the studies performed by other researchers, while it is analytically proved here. For example, the following assumption, i.e. se ða; zÞ ¼ 2τi , was made by Morimoto et al. [27]. The ability of the proposed method is then confirmed by comparison of the analytically calculated se and τrz at the outer surface of the unit cell, i.e. at r ¼ b, with the FEM results, see Fig. 3. As expected, at the outer surface, τrz ¼ 0 and se are constant along the fiber length. The analytical and FEM calculations of the axial, circumferential, and radial stress components of the matrix along the fiber length at r ¼ a are shown in Fig. 4a-c. As it is clear from the figures, the values of these stress components are exactly the same. This is expected from Eq. (15), which is proved by the FEM results, as well. Fig. 4 also present the mentioned stress components at the outer surface of the unit cell in region I, i.e. at r ¼ b. According to the figures, the analytical results are in good agreements with the FEM results. However, the values of the stress components at r ¼ b are not the same.

M. Mondali, A. Abedian / International Journal of Non-Linear Mechanics 57 (2013) 39–49

45

550 6061 Al/20%SiC 561 K 80 MPa

50

0

-50 Analytical @ r = a Analytical @ r = b FEM @ r = a FEM @ r = b

-100

0.2

0.4

0.6

0.8

350 250

Analytical, fiber Average FEM, fiber @ r =0 FEM, fiber @ r =a/4 FEM, fiber @ r =a/2 FEM, fiber @ r =3a/4 FEM, fiber @ r =a Analytical, matrix Average FEM, matrix Average

150 50 -50 -150

-150 0

6061 Al/20%SiC 561 K 80 MPa

450

Axial stress in fiber (MPa)

Axial stress in matrix (MPa)

100

1

Normalized axial position, (z / l)

0

0.2

0.4

0.6

0.8

1

Normalized axial position, (z / l) Fig. 5. The analytical and FEM results of the fiber and matrix axial stresses.

Hoop stress in matrix (MPa)

100 6061 Al/20%SiC 561 K 80 MPa

50

not much far from the reality. Based on the figure, a good correlation between the FEM and the mentioned assumption exists. 4.2. Load transfer mechanism in short fiber composites

0

-50 Analytical @ r = a Analytical @ r = b FEM @ r = a FEM @ r = b

-100

-150 0

0.2

0.4

0.6

0.8

1

Normalized axial position, (z / l)

Radial stress in matrix (MPa)

100 6061 Al/20%SiC 561 K 80 MPa

50

0

-50 Analytical @ r = a Analytical @ r = b FEM @ r = a FEM @ r = b

-100

-150 0

0.2

0.4

0.6

0.8

1

Normalized axial position, (z / l) Fig. 4. The analytical and FEM curves of the matrix stresses at the fiber/matrix interface and the outer surface: (a) axial stresses; (b) circumferential stresses; and (c) radial stresses.

To show the validity of the method for calculating the axial stress in the fiber, the FEM axial stress at five different locations in the domain 0 ≤r ≤a along the fiber length and the average analytical axial stress at the fiber cross section are presented in Fig. 5. As it is seen, the calculated average analytical values exactly match the FEM results. Moreover, the analytically calculated average axial stress in the matrix matches quite well the FEM results. Also, the average axial stress at the fiber end, i.e. 0 ≤ r ≤ a and z ¼ l, is in good agreement with the FEM results. This is due to the approximation made for the shear stress in region II at r ¼ a. Fig. 3 shows that the linear assumption made for τ′i in Eq. (34) is

For better understanding of the load transfer mechanism in short fiber composites, a comparison of the average axial stress in elastic and second stage creep states is performed using the finite element method. The development of axial stress within an aligned short fiber composite in both cases is illustrated in Fig. 6. It shows the magnitude of axial stress as a function of radial and axial position within the unit cell upon initial loading and in steady state creep stage. The numbers printed around the surface plots indicate the magnitude of the stresses in certain important areas of the model. Upon initial loading, the fiber bears a large fraction of the load, with the stress reaching a value of 197 MPa or over twice the applied stress of 80 MPa at ðz ¼ 0Þ. Also, the axial stress in the matrix at ðz ¼ 0Þ is about 19 MPa. Note that there is a large stress concentration of 249 MPa at the sharp corner of the fiber. While the presence of this stress concentration is expected, its value may not be very accurate because the mesh is not refined enough in this area to capture the rapid increase of stress with position. However, as it was shown in Figs. 5 and 6, this stress concentration in the fiber should not affect the stress components in other regions of the unit cell. As the matrix creeps around the fiber, it forces the extension of the fiber because of the perfect bond at the fiber/matrix interface. This results in an increase by nearly a factor of two in the axial stress, i.e. 395 MPa, at the center of the fiber in the steady state creep stage. Consequently, in order to satisfy the equilibrium equation, i.e. Eq. (8), the axial stresses in the matrix around the center of the fiber become compressive. Therefore, a maximum compressive stress of 88 MPa appears in the matrix. While this compressive stress is counter-intuitive for the applied loading conditions, it is a result of the constraint imposed by the other fibers through the boundary conditions on the unit cell. The increasing of the absolute values of average axial stress in the fiber during second stage creep can be explained by comparing the elastic interfacial shear stress with the steady state interfacial shear stress, see Fig. 6. As the figure shows, the elastic shear stress is almost zero along half of the fiber length. However, due to the matrix creeping, more loads are transferred to the fiber and as a result the interfacial shear stress in the first half of the fiber length starts to increase. This change is such that the exponential nature of the stress in the elastic case [12–18] is replaced by a function as

46

M. Mondali, A. Abedian / International Journal of Non-Linear Mechanics 57 (2013) 39–49

in Eq. (27). The average axial stress in fiber, which is expressed by the area under the shear stress curve, can be obtained by integrating Eq. (9) with respect to z from l to z. Consequently, higher average axial stress could be found in fiber in the second stage creep compared to the elastic case. 4.3. Verification of the displacement rates In the present study, according to Eq. (44b), the composite creep strain rate could be found using radial displacement rate on the outer surface of the unit cell, i.e. u_ b . Note that finding all the stress and displacement rate components in region II, which can be obtained by the use of imaginary fiber technique, requires a creep–creep solution of the problem. However, solution of such a problem with creeping both fiber and matrix is not presently in hand. Therefore, as a check point for validity of ε_ c , the distributions of radial and axial displacement rates in region I are verified by the FEM results. Fig. 7 illustrates the u_ distribution in the domain a ≤r ≤b calculated by both FEM and the analytical method. Based on the assumption made by Eq. (18), u_ is independent of z and hence the radial displacement rate is the same for any axial location along the fiber length. This is also evident in Eq. (30), which is expressed for region I (a ≤ r ≤ b). As it is seen, an excellent correlation exists _ is a function of between the results. Also, according to Eq. (29), w both r and z coordinates. Therefore, the analytical and FEM calculations of axial displacement rate in region I (a ≤r ≤b) are presented at z ¼ l=2, see Fig. 7. As it is seen, the results obtained by both methods for γ_ rz , i.e. Eq. (3d), coincide exactly at r ¼ a and _ values at r ¼ a þ ðb−aÞ=2 along r ¼ b. This is also backed by the w the fiber length, which is an arbitrary radial location in the domain (a ≤r ≤b), see Fig. 8. However, the relatively small discrepancy seen _ along the radial and axial directions is attributed to the in w assumption made for γ_ rz distribution in Eq. (17). Note that the

Normalized radial & axial displacement rates

Fig. 6. FEM evolution of the axial and shear stresses within the aligned fiber composite during the elastic and the second stage creep states.

14 6061 Al/20%SiC 561 K 80 MPa

12 10

Analytical, radial disp. rate @ 0 < z < l Analytical, axial disp. rate @ z = l / 2 FEM, radial disp. rate @ 0 < z < l FEM, axial disp. rate @ z = l / 2

8 6 4 2 0 -2 -4

0

0.2

0.4

0.6

0.8

1

Normalized radial position, (r - a) / (b - a) Fig. 7. Analytical and FEM results of the matrix radial and axial displacement rates normalized by the absolute value of u_ b vs. the normalized radial position.

mentioned small difference has no effect on the prediction of global behavior of the composite. Also, the calculated radial displacement rate on the outer surface coincides well with the FEM results, see Figs. 7 and 8.

4.4. Verification of the composite creep strain rate Now, the predictions of the composite creep strain rate made by the proposed analytical method and FEM are compared with each other and also with the available experimental data [4,31]. As Fig. 9 shows, the results obtained by the analytical method under different applied loads (s0 ) are in good agreement with the FEM results. Eqs. (53) and (54) show the creep relationships

M. Mondali, A. Abedian / International Journal of Non-Linear Mechanics 57 (2013) 39–49

Normalized radial & axial displacement rates @ r = a + (b - a) / 2.

16 Analytical, radial disp. rate Analytical, axial disp. rate FEM, radial disp. rate FEM, axial disp. rate

6061 Al/20%SiC 561 K 80 Mpa

14 12 10 8 6 4 2 0 -2 0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

Normalized axial position, (z / l) Fig. 8. Analytical and FEM results of the matrix axial and radial displacement rates at r ¼ ða þ bÞ=2 normalized by the absolute value of u_ b vs. the normalized axial position.

47

discrepancy. First, the assumption of aligned short fibers is not representative of the actual case of SiC-Al composites where fibers are usually misoriented. Second, the perfect bond at the fiber/ matrix interface is not consistent with the microstructure of discontinuous SiC–Al composites. Third, the accumulation of damage, i.e. the formation of voids or the occurrence of debonding, at the interface would accelerate the creep strain rate of the composite and hence it would introduce a different stress exponent into the creep law. In most cases, interfacial debonding is an important mechanism to increase the creep strain rate of short fiber composites [7,11,27,31]. It has been shown both analytically and numerically [11,18] that under the condition of progressive interfacial debonding with increasing stress, the creep prediction by the modeling methods will highly improve. However, since the primary goal here is the development of the analytical approach, this effect is not discussed in the present research work. 4.5. Effect of geometric factors

1.E-04

1.E-05

Creep strain rate (1/sec)

n = 21 1.E-06

n=4

1.E-07

1.E-08 n=4

1.E-09

Unreinforced 6061Al Matrix [4,22] 6061Al/20%SiC composite [4,22] Analytical prediction FEM prediction

1.E-10

1.E-11 10

100

1000

Creep stress (MPa) Fig. 9. Comparison of the analytical and FEM predictions with composite creep data for a 6061Al/20%SiC at 561 K, [4,31].

obtained by the analytical method and FEM modeling, respectively. −15

4:0

ð53Þ

ε_ c ¼ 1:12  10−15 s0 4:0

ð54Þ

ε_ c ¼ 1:47  10

s0

However, based on Fig. 9, the available experimental results are far from the calculations made here. This difference is better understood by comparing Eq. (55), which represents the experimentally obtained data for 6061Al/20%SiC composite [4,31], with Eqs. (53) and (54). ε_ c ¼ 6:82  10−46 s0 20:5

ð55Þ

As Eqs. (53) and (54) show, in the analytical and FEM models, the stress exponent for the composite creep law in the absence of creep in the fibers is found to be equal to that for the matrix material. This characteristic contrasts with available experimental evidence, i.e. Eq. (55), which shows that the stress dependence of creep rate in the composite is much stronger than that in the matrix. Several explanations are proposed to account for this

Since creep deformation is a major limiting factor for application of composites in high temperature environments, it is useful to have a detailed analytical prediction of the creep strain rate in these materials as a function not only of the applied load, but also as a function of all important geometric parameters. That is because the geometrical shape and arrangement of fibers in composites have a strong effect on the creep strain rate of these materials. In this section, a parametric study was undertaken to define the dependence of the geometric parameters on the steady state creep rate of short fiber composites. For all subsequent analyses, an average axial stress of 80 MPa is applied at 561 o K to the composite. Under these conditions, the predicted creep strain rate is 6:02  10−8 sec−1 which is in good agreement with the experimental results, see Fig. 9. Here, three of important geometric factors in composites, i.e. the fiber aspect ratio s ¼ l=a, the factor k ¼ ððl′=bÞ=ðl=aÞÞ that somehow is related to the unit cell aspect ratio, and the fiber volume fraction 2 0 f ¼ ðða2 lÞ=ðb l ÞÞ, are considered. For illustrative purpose, substituting the mentioned geometric factors into Eq. (45) and performing some simplifications give " pffiffiffi  f ðs= 3Þ1þ1=n 1 2 ε_ c ¼ A 1 þ ðk =f Þ1=3 − 1 þ 2n 1−f !1=n 3 −n 8 ð1−ðf kÞ2=3 Þ 5 sn  ð56Þ 0 4lnðf kÞ−1=3 −ðf kÞ4=3 þ 4ðf kÞ2=3 −3 Eq. (56) is used to perform the intended parametric study. Note that a similar study has been performed in [31] where reinforcement spacing factor, i.e. b=l, considered instead of k parameter. However, the mentioned study, which is based on finite element method, could not easily address the parameters effects. In the present work, each of these important geometric variables is varied in turn over a specific range while the other two parameters are considered to be constant. Mathematically, this procedure separates the effect of a particular geometric variable on the creep rate, assuming that its effect is independent of the effect of the other variables, i.e. variable separation, ε_ c ¼ f ðf ; s; kÞ ¼ f 1 ðf Þ f 2 ðsÞ f 3 ðkÞ

ð57Þ

The range of values investigated and the typical base case values are listed in Table 1. As stated in Table 1, the aspect ratio parameter s can have any value from zero to infinity. In fact, no restrictions are imposed on this parameter by the geometry and so each unit cell corresponding to parameter s is meaningful. However, this is not the case with the other two parameters

48

M. Mondali, A. Abedian / International Journal of Non-Linear Mechanics 57 (2013) 39–49

Table 1 Range of geometrical parameters studied. Parameter

f

Base case value [4,22] Range

0.2 2

0-k f or k o 1 0-1=k f or k 41

k

s

1 pffiffiffi f -1=f

5 0-∞

1

I

0.9

k=√f l = l'

Fiber volume fraction,f

0.8

or III

0.7

Base case value; Design Point [4,22]

0.6 0.5

k=1/f or a=b

II

0.4 0.3 0.2 0.1

Path B

Path A

0 0

1

2

3

4

5

6

7

8

9

10

k

with increasing f . Therefore, the distance between the fibers will decrease, which produces a greater constraint on the matrix material between the fibers and consequently ε_ c will reduce. As Fig. 11 shows, the FEM results match quite well the analytical predictions. Note that for the FEM cases with very small or very high values of fiber volume fraction, since it is not possible to have a close control over the size and the number of elements, the discrepancy between the analytical and FEM calculations is more noticed. For path “B”, as Fig. 10 shows, parameter k can only change in the range of 0:447 ok o5. Fig. 12 shows the changes in ε_ c as a function of parameter k. As it is seen, with increasing k the creep strain rate will decrease. This indicates that the composite resistance to creep deformation is increased with increasing the k parameter. The effect of increasing the fiber aspect ratio on the creep strain rate of the composite is illustrated in Fig. 13. As expected, the creep strain rate of the composite decreases with increasing the fiber aspect ratio. In fact increasing either the fiber aspect ratio while holding the volume fraction constant reduces the spacing between fibers. Therefore, a higher degree of constraint is imposed on the matrix material which could cause a high reduction in the steady state creep rate of the composite. Interestingly, the FEM results are in good agreement with the analytical calculations performed in this research work. Note that for the finite element analysis it is not easily possible to perform the calculations for any value of

Fig. 10. Feasible and unfeasible regions of the geometric factors.

1.E+00

FEM Analytical

1.E-01 1.E-02

Composite creep rate Normalized by matrix creep rate

Composite creep rate Normalized by matrix creep rate

1.E+00

1.E-03 1.E-04 1.E-05 1.E-06 1.E-07 1.E-08 1.E-09 1.E-10

6061 Al / SiC Path A ; k = 1, s = 5 80 MPa, 561 K

1.E-11 1.E-12 1.E-13 1.E-14 0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

1.E-02

FEM

1.E-03 1.E-04 1.E-05 1.E-06 1.E-07 1.E-08 6061 Al / SiC Path B ; f = 0.2, s = 5 80 MPa, 561 K

1.E-09 1.E-10

1.447

2.447

3.447

4.447

5

k Fig. 12. Effect of k parameter on the analytical and FEM predictions of the steady state creep strain rate of composite.

1.E+00

Composite creep rate Normalized by matrix creep rate

namely k and f . As it is known, the fiber volume fraction f changes from zero to one. Considering this fact and thep relationships for k ffiffiffi and f , the limiting values for k are found to be f o k o 1=f . Since always b 4 a and l′ 4 l, it is easily provedp that ffiffiffi for any value of fiber volume fraction, the fact l′ 4 l dictates f ok and b 4 a dictates k o 1=f . Therefore, as shown in Fig. 10, curves l ¼ l′ and a ¼ b divide the graph f vs. k to three regions of I, II, and III, where only region II is geometrically feasible. In fact, regions I and III provide unit cells with l 4l′ and a 4 b which are geometrically unacceptable. In this study, the base case value in the feasible region II is given in Table 1 and the results are presented for two paths of “A” and “B”, see Fig. 10. Note that for both paths the value of fiber aspect ratio s is equal to 5. Also, the results are presented for the design point with constant values of k and f and varying parameter s from zero to infinity. Fig. 11 illustrates the effect of changing volume fraction of fiber on the creep strain rate of the composite. As it is seen, ε_ c reduces to zero by increasing fiber volume fraction. Here, considering the relationships between k, s and f results that both l′ and b decrease

Analytical

1.E-11 0.447

Fiber volume fraction, f Fig. 11. Effect of fiber volume fraction on the analytical and FEM predictions of the steady state creep strain rate of composite.

1.E-01

FEM Analytical

1.E-01 1.E-02 1.E-03 1.E-04 6061 Al / 20%SiC k = 1 , 80 MPa 561 K

1.E-05 1.E-06 0

5

10

15

20

25

30

Fiber aspect ration, s Fig. 13. Effect of fiber aspect ratio on the analytical and FEM predictions of the steady state creep strain rate of composite.

M. Mondali, A. Abedian / International Journal of Non-Linear Mechanics 57 (2013) 39–49

fiber aspect ratio, while the analytical method provides an accurate tool for such calculations.

5. Conclusions In the present work, a closed form solution for the steady state creep deformation of short fiber composites is proposed. As Eq. (45) shows, a simple expression for creep strain rate of the composite was obtained which is only a function of creep constants of the matrix material and also geometrical parameters of the unit cell. Also, this analytical model is perfectly able to calculate the stress field through the composite. Considering the obtained results the following conclusions could be made. First, the matrix stress and displacement rate components provided by the proposed analytical solution completely satisfy the equilibrium and constitutive creep equations. Moreover, the results match well the overall boundary conditions including the surface conditions, the interface continuity conditions, and the axial force equilibrium conditions. Second, the main advantage of the present method compared to the previous available analytical methods is in its capability of predicting the stress and displacement rate fields in the composite and also compatibility of the obtained results with the FEM predictions. Third, the analytical model can simulate a reduction in the creep rate of the composite but cannot adequately predict the high stress exponent of the material. This dictates the necessity for including of damage at the fiber/matrix interface and also the material imperfections into the analytical modeling. Forth, the simulated creep rate is strongly influenced by the geometric arrangement of the reinforcing phase of the composite. In particular the volume fraction and the aspect ratio of the reinforcement and the k parameter have a strong effect on determining the creep rate of the composite. Fifth, one of the advantages of the present analytical method compared to the FEM modeling is in its simplicity to have a parametric study which defines the effects of geometric parameters on the steady state creep behavior of short fiber composites. Note that for the limiting values of k parameter, as shown in Fig. 10, FEM modeling preparation is very complicated and also the obtained solution consists of high degree of approximation. Also, in creep FEM modeling convergence is a major issue and the process is very time consuming due to high nonlinearity involved. Finally, the most important application of the present work is to determine the creep strain rate of some MMC composites which their stress exponent in power law creep formulation is similar to stress exponent of the matrix material. The mentioned expression was also emphasized by Dragon and Nix for Ag/W composite which its stress exponent for both the matrix and the composite is 4.3 [31]. References [1] A. Kelly, W.R. Tyson, Tensile properties of fibre-reinforced metals: copper/ tungsten and copper/molybdenum, Journal of the Mechanics and Physics of Solids 13 (1965) 329–338. [2] A. Kelly, W.R. Tyson, Tensile properties of fibre reinforced metals—I. Creep of silver-tungsten, Journal of the Mechanics and Physics of Solids 14 (1966) 177–184. [3] M. Doruk, A.S. Yue, Creep behavior of fiber reinforced metal matrix composites, Metallurgical and Materials Transactions A7 (1976) 1465–1468. [4] T.G. Nieh, Creep rupture of a silicon-carbide reinforced aluminum composite, Metallurgical and Materials Transactions A15 (1984) 139–146. [5] V.C. Nardone, J.R. Strife, Analysis of the creep behavior of silicon carbide whisker reinforced 2124 Al (T4), Metallurgical and Materials Transactions A18 (1987) 109–114.

49

[6] T.L. Dragon, J.J. Schlautman, W.D. Nix, Processing and creep characterization of a model metal matrix composite: lead reinforced with nickel fibers, Metallurgical and Materials Transactions A22 (1991) 1029–1036. [7] F.A. Mohamed, K.T. Park, E.J. Lavernia, Creep behavior of discontinuous SiC–Al composites, Materials Science and Technology A150 (1992) 21–35. [8] A. Greasly, Creep of dispersion reinforced aluminum based metal matrix composite, Materials Science and Technology 11 (1995) 163–166. [9] S.W. Schwenker, D. Eylon, Creep deformation and damage in a continuous fiber-reinforced Ti–6Al–4V composite, Materials Science and Technology A27 (1996) 4193–4204. [10] Z.Y. Ma, S.C. Tjong, The high temperature creep behavior of 2124 aluminum alloys with and without particulate and SiC-whisker reinforcement, Composites Science and Technology 59 (1999) 737–747. [11] M. Mondali, A. Abedian, S. Adibnazari, FEM study of the second stage creep behavior of Al6061/SiC metal matrix composite, Computational Materials Science 34 (2005) 140–150. [12] C.H. Hsueh, Interfacial debonding and fiber pull-out stresses of fiber- reinforced composites, Materials Science and Engineering 123 (1) (1990) 1–11. [13] C.H. Hsueh, Interfacial debonding and fiber pullout stresses of fiber-reinforced composites part VII: improved analysis for bonded interfaces, Materials Science and Engineering 154 (1992) 125–132. [14] C.H. Hsueh, A modified analysis for stress transfer in fiber-reinforced composites with bonded fiber ends, Journal of Materials Science 30 (1995) 219–224. [15] C.H. Hsueh, R.J. Young, X. Yang, P.F. Becher, Stress transfer in a model composite containing a single embedded fiber, Acta Materialia 45 (4) (1997) 1469–1476. [16] A.H. Nayfeh, W.G. Abdelrahman, Micromechanical modeling of load transfer in fibrous composites, Mechanics of Materials 30 (1998) 307–324. [17] Z. Jiang, X. Liu, G. Li, J. Lian, A new analytical model for three-dimensional elastic stress field distribution in short fibre composite, Materials Science and Engineering A 366 (2004) 381–396. [18] A. Abedian, M. Mondali, M. Pahlavanpour, Basic modifications in 3D micromechanical modeling of short fiber composites with bonded and debonded fiber end, Computational Materials Science 40 (2007) 421–433. [19] X.L. Gao, K. Li, A shear-lag for carbon nanotube-reinforced polymer composites, International Journal of Solids and Structures 42 (2005) 1649–1667. [20] Z. Jiang, J. Lian, D. Yang, S. Dong, An analytical study of the influence of thermal residual stresses on the elastic and yield behaviors of short fiberreinforced metal matrix composites, Materials Science and Engineering A 248 (1998) 256–275. [21] S.T. Mileiko, Steady state creep of a composite with short fibres, Journal of Materials Sciences 5 (1970) 254–261. [22] D. McLean, Viscous flow of aligned composites, Journal of Materials Science 7 (1972) 98–104. [23] A. Kelly, K.N. Street, Creep of discontinuous fibre composites. II. Theory for the steady-state, Proceedings of the Royal Society of London A328 (1972) 283–293. [24] H. Fukuda, T.W. Chou, An advanced shear-lag model applicable to discontinuous fiber composites, Journal of Composite Materials 1 (15) (1981) 79–91. [25] H. Lilholt, Creep of fibrous composite materials, Composites Science and Technology 22 (1985) 277–294. [26] M. McLean, Creep deformation of metal–matrix composites, Composites Science and Technology 23 (1985) 37–52. [27] T. Morimoto, T. Yamaoka, H. Lilholt, M. Taya, Second stage creep of silicon carbide whisker/6061 aluminum composite at 573K, Journal of Engineering Materials and Technology 110 (1988) 70–76. [28] J.R. Pachalis, J. Kim, T.W. Chou, Modeling of creep of aligned short-fiber reinforced ceramic composites, Composites Science and Technology 37 (1990) 329–346. [29] Y.S. Lee, T.J. Batt, P.K. Liaw, Stress analysis of a composite material with short elastic fibre in power law creep matrix, International Journal of Mechanical Sciences 32 (10) (1990) 801–815. [30] Y.R. Wang, T.W. Chou, Analytical modeling of creep of short fiber reinforced ceramic matrix composite, Journal of Composite Materials 26 (9) (1992) 1269–1286. [31] T.L. Dragon, W.D. Nix, Geometric factor affecting the internal stress distribution and high temperature creep rate of discontinuous fiber reinforced metals, Acta Metallurgica et Materialia 38 (10) (1990) 1941–1953. [32] J.T. Boyle, J. Spence, Stress Analysis for Creep, Butterworths Press, Southampton58. [33] M. Mondali, A. Abedian, A. Ghavami, A new analytical shear-lag based model for prediction of the steady state creep deformations of some short fiber composites, Materials and Design 30 (4) (2009) 1075–1084. [34] V.C. Nardone, K.M. Prewo, On the strength of discontinuous silicon carbide reinforced aluminum composites, Scripta Metallurgica 20 (1986) 43–48. [35] T.W. Clyne, A simple development of the shear lag theory appropriate for composites with a relatively small modulus mismatch, Materials Science and Engineering A 122 (1989) 183–192. [36] V.M. Karbhari, D.J. Wilkins, An engineering modification to the shear-lag model as applied to whisker and particulate reinforced composites, Scripta Metallurgica 25 (1991) 707–712.