- Email: [email protected]
The Mori–Tanaka stiffness tensor: diagonal symmetry, complex fibre orientations and non-dilute volume fractions
Mechanics of Materials 33 (2001) 531±544
www.elsevier.com/locate/mechmat
The Mori±Tanaka stiness tensor: diagonal symmetry, complex ®bre orientations and non-dilute volume fractions J. Schjùdt-Thomsen *, R. Pyrz Institute of Mechanical Engineering, Aalborg University, Pontoppidanstraede 101, DK 9220 Aalborg East, Denmark Received 3 April 2001; received in revised form 9 June 2001
Abstract This paper considers the diagonal symmetry of the stiness tensors predicted using the Mori-Tanaka (MT) approach. Since the MT approach may yield asymmetric stiness tensors this paper considers an alternative approach to ensure the diagonal symmetry. Furthermore, an extension of the MT approach to non-dilute volume fractions is considered. It is shown that the MT approach by Benveniste [Mech. Mater. 6 (1987) 147], may not yield diagonally symmetric stiness tensors in situations when a statistical ®bre orientation function is incorporated. Only when the inclusions are spherical and randomly distributed or fully aligned will the MT predicted stiness tensor be diagonally symmetric. The extension to non-dilute volume fractions is seen to lie within the Hashin±Shtrikman±Walpole (HSW) bounds and comparison to other approaches and experiments found in the literature shows very good agreement. Ó 2001 Elsevier Science Ltd. All rights reserved. Keywords: Mori±Tanaka stiness tensor; Diagonal symmetry; Complex ®bre orientations; Non-dilute volume fractions
1. Introduction Traditional continuum mechanics is a branch of the physical sciences concerned with the deformations and motions of continuous material media under the in¯uence of external eects. These media are in the form of a collection of material points interconnected by some internal forces. A basic assumption of continuum mechanics is that the mass, stress and strain are regarded as essentially uniform within an in®nitesimal material point and its surrounding neighbourhood. However, when dealing with non-homogeneous materials such as composite materials, cellular materials or polycrystals, the aforementioned assumption may no longer be valid, since the material point and its neighbourhood is generally not uniform, but consist of various constituents with dierent properties and shape. On the other hand the traditional continuum mechanic is a well established ®eld of science and thus applying continuum mechanics approaches to non-homogeneous media seems bene®cial. This requires the concept of a representative
*
Corresponding author. Tel.: +45-9635-9330; fax: +45-9815-1411. E-mail address: [email protected] (J. Schjùdt-Thomsen).
0167-6636/01/$ - see front matter Ó 2001 Elsevier Science Ltd. All rights reserved. PII: S 0 1 6 7 - 6 6 3 6 ( 0 1 ) 0 0 0 7 2 - 2
532
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
volume element (RVE), see e.g., Hill (1963) or Nemat-Nasser and Hori (1993), which is then used in the continuum theory as the material point. The RVE must include a large number of constituents and be statistically representative of the local continuum properties. This approach has been taken by several authors to calculate the overall elastic properties of various non-homogeneous materials, Hill (1963), Hashin (1963, 1965), Walpole (1969, 1985), Kwon and Dharan (1995), Ponte Castaneda and Willis (1995), Molinari and El Mouden (1996), Torquato (1997, 1998a,b) and Dvorak and Srinivas (1999). Non-elastic properties were considered by Hill (1965), Fotiu and NematNasser (1996), Bhattacharyya and Weng (1996), Ponte Castaneda (1996), Nebozhyn and Ponte Castaneda (1999), and Schjùdt-Thomsen and Pyrz (2000). Several other investigations have been carried out and references may be found in Nemat-Nasser and Hori (1993). One of the methods which has been used extensively is the Mori-Tanaka (MT) model, based on the work by Mori and Tanaka (1973) and Eshelby (1957). This model has proven to be quite accurate in predicting the eective properties of various materials with either random orientation or total alignment of the reinforcing phases. However, for composites containing inclusions of various shapes or orientations the MT predicted stiness tensor may violate the symmetry requirement that Cijkl Cklij (Ponte Castaneda and Willis, 1995; Li, 1999; Ferrari, 1991; Benveniste et al., 1991a). This is indeed a drawback of the MT model, but one should bear in mind that the MT model was originally proposed for composites containing inclusions of similar shape. Benveniste et al. (1991a) found that the MT stiness tensor is symmetric for twophase composites and for multi-phase composites if the inclusions are of similar shape and have the same orientation with respect to a ®xed reference frame. However, the same authors (Benveniste et al., 1991b), found that if the ``typical inclusion'' is a coated particle, then the implementation can be carried out within the framework of two-phase composites and thus gives diagonally symmetric stiness tensors (Benveniste et al., 1989; Chen et al., 1990). Li (1999) also concludes that it is the extension from two-phase composites to multi-phase composites with dierently shaped reinforcements, that causes the con¯icting results ± not the original MT assumption. In most of the work on composite materials using the MT approach the inclusions have been assumed to be either aligned or randomly dispersed in the matrix. However, in real short ®bre composites the ®bres are neither aligned nor dispersed at random directions. Fibre orientation is governed by the nature of the complicated ¯ow ®eld of the matrix during ®lling of the mould. The consequence of the ¯ow ®eld is that the ®bres are oriented at dierent directions, somewhere between the two extremes; random and fully aligned. Investigations considering various ®bre orientations, using dierent approaches have been considered by Pipes et al. (1982), Takao et al. (1982), Bozarth et al. (1987), Maekawa et al. (1989), Marzari and Ferrari (1992), Sayers (1992), Ngolle and Pera (1997) and Schjùdt-Thomsen and Pyrz (2000). Furthermore, Torquato (1997, 1998a,b), Molinari and El Mouden (1996) and Ponte Castaneda and Willis (1995) also considered spatial distribution and non-dilute volume fractions of the reinforcement phase. In light of the above observations the present paper focuses on a numerical study of the MT predicted stiness tensor for various ®bre orientations using a statistical ®bre orientation function. Furthermore, should the MT stiness tensor be non-symmetric an alternative approach is proposed which will yield a diagonally symmetric stiness tensor. Additionally, an approach is proposed which modi®es the MT approach to take non-dilute volume fractions into account. In the paper a boldface ordinary capital letter will denote a fourth-order tensor, A Aijkl , the inner product is AB Aijkl Bklmn and the inverse of a tensor is denoted as A 1 . 2. Eective stiness of a two-phase composite An estimate of the eective stiness of a two-phase composite can be obtained from Benveniste (1987) analysis of the MT model. Without going into detail the ®nal expression for the stiness is given as:
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
L MT L1 c2 f L2
L1 Tgc1 I c2 fTg
1
533
1
where indices ``1'' and ``2'' denote the matrix and reinforcing phase respectively and I is the fourth-order identity tensor. Curly brackets denote averaging over all possible directions. The tensor T relates the uniform strain in an ellipsoidal particle embedded in an all matrix medium to the average matrix strain. e 2 Te 1 ;
2
where T I SL1 1 L2
L1
1
3
and S is the Eshelby tensor, and thus the strain concentration tensor, A2 , is A2 T c1 I c2 T
1
4
In order to obtain the eective properties of a heterogeneous material the RVE is needed. Obviously the validity of this approach relies heavily on the de®nition of the RVE. The RVE should be the smallest possible unit which include all distinct elements of the microstructure that have a ®rst-order in¯uence on the properties and should represent the relative dispersion of the phases in a statistically representative manner (Nemat-Nasser and Hori, 1993). The present method is based on a transformation of fourth-order tensors followed by a spatial averaging over the RVE corresponding to the material point in the continuum, thus giving the overall averaged properties. The ®bre orientations in the RVE are now described through an orientation distribution function in the Euler angles, g h; /; b. The stiness Lijkl in the X1 X2 X3 coordinate system is related to the stiness Lpqrs in the local ®bre coordinate system through the transformation rule of fourth-order tensors as Lijkl aip ajq akr als Lpqrs ;
5
where aij is a function of the Euler angles, see Appendix A. In order to incorporate the stiness contribution from all orientations, aip ajq akr als Lpqrs is multiplied by the ®bre orientation distribution function and integrated over all possible orientations to obtain the weighted orientation average stiness as
Fig. 1. Euler angles de®ning the relation between the local coordinate system, x1 x2 x3 and the global coordinate system, X1 X2 X3 .
534
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
Z hLijkl i
2p
Z
0
2p
0
Z
p
g h; /; baip ajq akr als Lpqrs sin h dh d/ db;
0
6
where hi denotes averaging over all possible directions. g h; /; b is the orientation distribution function de®ned in the Euler coordinates h; /; b, Lpqrs is the stiness tensor of the RVE. The transformation matrix in terms of the Euler angles de®nes the orientation of the axes x1 x2 x3 ®xed in the RVE with respect to the axis X1 X2 X3 ®xed in the continuum, as shown in Fig.1. The concept of averaging over various orientations to estimate overall properties is not a new one. The approach is generally accepted and has been used by several authors in a variety of ®elds, for investigating the molecular orientation distribution of polymers (Bower, 1972), and to describe the electroelastic properties of piezoelectric materials (Huang and Kuo, 1996; Kuo and Huang, 1997). The major dierence between the present approach and the previous ones is that this work considers the full stiness tensor whereas the previous ones considered only the in-plane stiness properties of laminates or the Young's and shear moduli of the short ®bre composites. Furthermore, the present approach will yield a diagonally symmetric MT stiness tensor for non-trivial ®bre orientations, which has not been carried out previously. Reconsidering Eq. (1), the eective stiness of a two-phase material can be expressed in terms of an orientation dependent part and an orientation independent part as: 1 L MT pqrs Lpqrs fLpqrs g:
7
The stiness tensor Lpqrs in Eq. (6) is now taken to be the orientation dependent part of Eq. (7), for a composite with fully aligned ®bres. The averaging procedure will then provide the overall stiness of a composite with ®bre orientation described through the function g h; /; b. One of the objectives of the present work is to obtain a diagonally symmetric stiness tensor for arbitrary ®bre orientations. The proof that hLijkl i is indeed diagonally symmetric is quite straightforward. Since Lpqrs in Eq. (6) is calculated for a composite with fully aligned ®bres, Lpqrs Lrspq , i.e., it is symmetric, thus Lijkl aip ajq akr als Lpqrs ; Lklij akr als aip ajq Lrspq and it is easily seen that these two expressions are identical and the resulting eective stiness tensor, given by Eq. (6) is diagonally symmetric. 2.1. Orientation distribution function The orientation distribution function used in this work was used previously by Schjùdt-Thomsen and Pyrz (2000), and is given as (Maekawa et al., 1989) g h; /; b R hb ha
sin h
2P 1
sin h
cos h
2P 1
2Q 1
cos h
2Q 1
dh
;
8
where ha and hb are the upper and lower limits of the angle h being present in the distribution, 0 6 ha ; ha 6 p and P P 1=2, Q P 1=2. As can be seen from Eq. (8) the ®bre orientation function describes the orientation only in terms of h, meaning that the composite will be isotropic or transversely isotropic with respect to the X3 axis, depending on the values of P and Q.
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
535
3. Non-dilute volume fractions As was pointed out in Benveniste (1987), particle interaction is taken into account to a certain extent, since for a composite reinforced with spherical particles his reformulation of the MT-model coincides with the composites spheres assemblage of Hashin (1962). However, comparing the results from the MT approach with higher order approximations as the one proposed by Cherkaoui et al. (1995) the traditional MT approach gives smaller values as compared to the mentioned approach and experiments. The strain concentration tensors are usually calculated using Eq. (3), where the Eshelby tensor is calculated using the matrix material as the comparison medium. In order to extent the MT approach to non-dilute volume fractions the concentration tensor needs to be calculated dierently. This may be done in dierent ways. Molinari and El Mouden (1996) and Cherkaoui et al. (1995) formulated an integral equation of the unknown strain, leading to a set of linear equations from which the concentration tensors were found in a self consistent like manner. Another approach may be simply to replace L1 in Eq. (3) by another suitably chosen comparison medium, that is 1 T I SL0 1 L2 L0 ; 9 where L0 is the stiness tensor of the comparison medium, in which the Eshelby tensor must be calculated. An admissible choice of comparison medium is such that the overall stiness do not violate the bounds (Dvorak and Srinivas, 1999). In the present analysis the following approach is taken. Initially T is calculated as 1 T I SL1 1 L2 L1 and used to obtain the isotropic MT stiness, (Lisotropic ), with a volume fraction c02 less than the actual volume fraction c2 , that is 0 6 c02 6 c2 and used as the comparison medium stiness h i 1 1 T0 I SLisotropic L2 Lisotropic : 10 Now T0 is replacing T in Eq. (1). Using the actual volume fraction now gives the overall eective stiness. This approach is somewhat similar to the successive embedding procedure of Nemat-Nasser and Hori (1993, p. 330). Using Eqs. (1) and (4) gives the eective MT stiness in a shorter form as L MT L1 c2 L2
L1 A2 :
11
The eective stiness tensor determined from the non-dilute approach can be written more generally as L N 1 L1 c2 L2
N
L1 A2 ;
12
N A2
where is a function of c02 and N denotes the order of the approximation, i.e., for N 1 and c02 c2 the original MT stiness tensor is obtained and if N ! 1 and c02 c2 the self consistent method is obtained. For N 2 the stiness tensor from the non-dilute approach is obtained. 4. Results and discussion For illustrative purposes the matrix and ®bre properties were chosen to be those of polypropylene and glass ®bres, with the following properties: E 1 1:5 GPa, m 1 0:32, E 2 73 GPa and m 2 0:22 and the ®bre volume fraction, c2 0:35 in the subsequent calculations. However, it should be emphasized that the approach can ± of course ± be used for other material combinations.
536
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
Now two approaches are followed: (1) the overall stiness tensor, denoted as, LhMTi , is obtained by averaging the terms in curly brackets in Eq. (1) over all possible orientations, described the orientation function, Eq. (8) and (2) The stiness tensor predicted by the MT method is inserted into Eq. (6), thus giving the overall stiness for a composite with orientations described by Eq. (8), denoted as L . 4.1. Random spherical particles In this situation the particles ``orientation'' is described by taking the orientation parameters P and Q equal to Q P 0:5, Fig. 2, thus giving an isotropic material. In such cases the material eective properties are characterised through the eective bulk and shear moduli. For the sake of con®dence, however, the full 6 6 stiness matrices are shown. 2
LhMT i
4033:36 6 1674:13 6 6 1674:13 6 6 6 0 6 4 0 0 2
4033:36 6 1674:13 6 6 1674:13 6 L 6 6 0 6 4 0 0
1674:13 4033:36 1674:13 0 0 0
1674:13 4033:36 1674:13 0 0 0
1674:13 0 1674:13 0 4033:36 0 0 1179:61 0 0 0 0
1674:13 1674:13 4033:36 0 0 0
0 0 0 1179:61 0 0
0 0 0 0 1179:61 0
3 0 7 0 7 7 0 7 7; 7 0 7 5 0 1179:61
3 0 0 7 0 0 7 7 0 0 7 7: 7 0 0 7 5 1179:61 0 0 1179:61
From the two results it is seen that the stiness matrices are indeed identical and diagonally symmetric ± as they should be.
Fig. 2. Fibre orientation distribution function for random and partial aligned ®bres.
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
537
4.2. Partially aligned ellipsoidal particles The inclusions are now in the shape of ellipsoids elongated in the x3 direction and having aspect ratios, a3 =a1 a3 =a2 30. The ®bres are oriented according to Q 11 and P 0:5, which corresponds to partial alignment along the X3 axis and is shown in Fig. 2. The composite material is now transversely isotropic, with the X1 X2 plane as the isotropic plane. The eective stiness matrices are now 3 2 4965:42 2081:40 3362:94 0 0 0 7 6 2081:40 4965:42 3362:94 0 0 0 7 6 7 6 4481:31 4481:31 17417:51 0 0 0 hMTi 7; 6 L 6 7 0 0 0 3047:47 0 0 7 6 5 4 0 0 0 0 3047:47 0 0 0 0 0 0 1416:72 2
4573:82 6 1888:16 6 6 2948:46 L 6 6 0 6 4 0 0
1888:16 4573:82 2948:46 0 0 0
2948:46 2948:46 14532:12 0 0 0
0 0 0 2528:65 0 0
3 0 0 7 0 0 7 7 0 0 7: 7 0 0 7 5 2528:65 0 0 1342:82
Now the stiness matrices are dierent. The LhMTi stiness matrix is not diagonally symmetric and the asymmetry is clearly signi®cant as was also shown by Benveniste et al. (1991a) for circular discs and continuous ®bres. The main reason for the asymmetry may be due to the strain concentration tensors which are generally not diagonally symmetric (Dvorak and Bahei-El-Din, 1997). The reason is that the concentration tensors must satisfy the requirement that c1 A1 c2 A2 I (for a two-phase material). However, unless all the inclusions are aligned the concentration tensors may not satisfy this consistency requirement. To overcome this inconsistency the concentration tensors must be normalised (Nemat-Nasser and Hori, 1993). Walpole (1969) suggested this to be done by replacing the concentration tensors by r Ar hAr i 1 ; A
13
where hAr i indicates summation over all orientations, shapes and sizes. However, this normalisation has not been considered in this work, but it was used by Iwakuma and Nemat-Nasser (1984) and Nemat-Nasser and Obata (1986) in their self-consistent modelling of ®nite deformations of polycrystalline solids. Increasing the degree of alignment by taking P 0:5 and Q 100, as shown in Fig. 3 makes the asymmetry of MT predicted stiness tensor even more signi®cant as shown below 3 2 3199:57 1390:32 1636:86 0 0 0 6 1390:32 3199:57 1636:86 0 0 0 7 7 6 6 3366:90 3366:90 17873:38 0 0 0 7 7 LhMTi 6 6 0 0 0 1270:75 0 0 7 7 6 4 0 0 0 0 1270:75 0 5 0 0 0 0 0 825:88 and the L tensor is
538
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
Fig. 3. Fibre orientation distribution function for nearly aligned ®bres.
2
3721:75 6 1637:26 6 6 1870:93 L 6 6 0 6 4 0 0
1637:26 3721:75 1870:93 0 0 0
1870:93 1870:93 21048:15 0 0 0
0 0 0 1475:96 0 0
3 0 0 7 0 0 7 7 0 0 7: 7 0 0 7 5 1475:96 0 0 1042:24
The L tensor is again diagonally symmetric. Finally, as a last example the situation where the ®bres are practically fully aligned. A perfect alignment of the ®bres can be described through a Dirac delta function of the angle h, however, in this example the parameters of the orientation distribution function are P 0:5 and Q 250, and is shown in Fig. 3 3 2 3701:88 1633:87 1681:11 0 0 0 7 6 1633:87 3701:88 1681:11 0 0 0 7 6 7 6 1681:11 1681:11 21853:93 0 0 0 7 6 L 6 7 0 0 0 1288:58 0 0 7 6 5 4 0 0 0 0 1288:58 0 0 0 0 0 0 1034:00 and for comparison the results from the original MT approach for fully aligned inclusions 3 2 3696:27 1634:23 1539:30 0 0 0 7 6 1634:23 3696:27 1539:30 0 0 0 7 6 7 6 1539:30 1539:30 22431:67 0 0 0 MT; align 7: 6 L 6 7 0 0 0 1148:44 0 0 7 6 5 4 0 0 0 0 1148:44 0 0 0 0 0 0 1031:01 4.3. Non-dilute volume fractions Previously it has been shown that the averaging procedure yields diagonally symmetric stiness tensors, so this approach has also been used in the situation when non-dilute volume fractions are considered. In order to do some checking on the present approach the Hashin±Shtrikman±Walpole (HSW), bounds have
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
539
been calculated. The bounds are calculated using the connections between the MT model and the HSW bounds (Weng, 1990). For random spherical particles it was shown by Weng (1990) that the MT moduli will coincide with the HSW lower/upper bound, denoted as L and L , if the matrix is the softest/hardest of the phases, respectively. In this context softest/hardest means that Lr L0 is positive/negative semi-de®nite for all r and L0 is the stiness of the comparison medium. First, consider the situation when the inclusions are randomly dispersed spherical particles, thus the lower HSW bound coincides with the moduli of the MT model when the matrix stiness is taken to be that of polypropylene 3 2 4033:36 1674:13 1674:13 0 0 0 7 6 1674:13 4033:36 1674:13 0 0 0 7 6 7 6 1674:13 1674:13 4033:36 0 0 0 7 L 6 7 6 0 0 0 1179:61 0 0 7 6 5 4 0 0 0 0 1179:61 0 0 0 0 0 0 1179:61 and the upper bound is 2 38941:73 6 8458:58 6 6 8458:58 L 6 6 0 6 4 0 0
obtained when the matrix phase is taken to be that of the glass particles 3 8458:58 8458:58 0 0 0 7 38941:73 8458:58 0 0 0 7 7 8458:58 38941:73 0 0 0 7: 7 0 0 15241:57 0 0 7 5 0 0 0 15241:57 0 0 0 0 0 15241:57
The stiness prediction from the present approach is found to be 3 2 4707:56 1867:28 1867:28 0 0 0 7 6 1867:28 4707:56 1867:28 0 0 0 7 6 7 6 1867:28 1867:28 4707:56 0 0 0 7; 6 L 6 7 0 0 0 1420:13 0 0 7 6 5 4 0 0 0 0 1420:13 0 0 0 0 0 0 1420:13 where the stiness of the comparison medium is initially L0 L1 and in the next step L0 equals that of the MT approach for ®bre volume fractions equal to half the value of the actual composite. Taking c02 12 c2 is only done for illustrative purposes. Indeed, any choice of comparison medium leading to an eective stiness that does not violate the HSW bounds is admissible (Dvorak and Srinivas, 1999). 4.4. Randomly oriented ®bres Now the inclusions are taken to be ®bres with aspect ratios a3 =a1 a3 =a2 30 and the distribution parameters are P Q 0:5. The lower bound is now 3 2 7480:66 2801:32 2801:32 0 0 0 7 6 2801:32 7480:66 2801:32 0 0 0 7 6 7 6 2801:32 2801:32 7480:66 0 0 0 7 6 L 6 7 0 0 0 2339:67 0 0 7 6 5 4 0 0 0 0 2339:67 0 0 0 0 0 0 2339:67
540
and the upper bound 2 37234:66 6 7747:95 6 6 7747:95 L 6 6 0 6 4 0 0
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
7747:95 37234:66 7747:95 0 0 0
7747:95 0 7747:95 0 37234:66 0 0 14743:35 0 0 0 0
and the present approach yields 2 9524:65 3364:71 3364:71 6 3364:71 9524:65 3364:71 6 6 3364:71 3364:71 9524:65 L 6 6 0 0 0 6 4 0 0 0 0 0 0
0 0 0 3079:97 0 0
0 0 0 0 14743:35 0
3 0 7 0 7 7 0 7 7 0 7 5 0 14743:35
3 0 0 7 0 0 7 7 0 0 7: 7 0 0 7 5 3079:97 0 0 3079:97
From the two set of results above the eective bulk and shear moduli are calculated along with the upper and lower bounds on these properties. For the random spherical particles, the lower bound on the eective bulk modulus j 2460 MPa, the upper bound is j 18 620 MPa and the present approach yields j 2814 MPa. The shear moduli are: l 1180, l 15 241 and l 1420 MPa. For the random oriented ®bres the same properties are: the lower bound on the eective bulk modulus j 4361 MPa, the upper bound is j 17 577 MPa and the present approach yields j 5418 MPa and the shear moduli l 2340, l 14 743 and l 3080 MPa. It is clearly seen that in the situations considered the present approach does not violate the HSW bounds since j 6 j 6 j and l 6 l 6 l . The composite Young's modulus for the polypropylene/glass particle composite has been calculated and is shown in Fig. 4, for the full range of volume fractions. It is seen that the predicted Young's modulus does not violate the bounds. The curve of the Young's modulus from the original MT approach is not very smooth for larger volume fractions. This may be due to the fact that interaction is not accounted for suciently, whereas the curve from the non-dilute MT approach is much smoother, indicating that the use of a dierent comparison medium enhances the models ability to account for interaction.
Fig. 4. Upper and lower HSW bounds along with the predictions of the non-dilute MT method. (Glass particle/polypropylene composite.)
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
541
Fig. 5. Normalised Young's modulus. Theoretical data and experimental data are from Cherkaoui et al. (1995).
Finally, the results from this approach is compared to results from Cherkaoui et al. (1995), for an epoxy matrix reinforced with glass particles. The material properties are: E 1 2:9 GPa, m 1 0:4, E 2 73 GPa and m 2 0:2, c02 0:15 and c2 0:3. The results are shown in Fig. 5. Cherkaoui et al. (1995) found that the calculated composite Young's modulus is Ec 6:3 GPa and the c experimentally determined Young's modulus is Eexp 5:97 GPa. The present approach yields an eective Young's modulus, E 6:05 GPa, whereas the original MT approach gives EMT 5:3 GPa. Very good agreement is found for the present approach compared to other investigations. Finally, the Young's modulus for this composite is calculated and shown in Fig. 6. Comparing this ®gure with Fig. 4 and remembering that the MT approach is coincident with the lower bound it is believed that for larger mismatch between the constituents moduli (E 2 =E 1 48:7 in Fig. 4 and E 2 =E 1 25:2 in Fig. 6) the non-dilute MT approach yields more reasonable results.
Fig. 6. Upper and lower HSW bounds along with the predictions of the non-dilute MT method. (Glass particle/epoxy composite.)
542
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
5. Conclusions In the present paper the diagonal symmetry of the stiness tensors predicted by the MT model have been considered for various ®bre orientations. It was found that this method does not provide diagonally symmetric stiness tensors if the averaging is carried out only on the concentration tensor T as described in Benveniste (1987). However, averaging the stiness tensor for a composite with unidirectionally aligned ®bres according to a ®bre orientation function yields diagonally symmetric eective stiness tensors for non-trivial ®bre orientations and coincides with the predictions of the MT method for random spherical particles and fully aligned ®bres. The MT approach has been extended to incorporate the eect of non-dilute volume fractions in a manner similar to the successive embedding procedure and results have been compared to results in the literature and very good agreement is found. The non-dilute approach provides a ¯exible method for predicting the eective properties of various composite materials. Due to the complex nature of composite materials it is dicult to know the degree of interaction beforehand and thus which model may adequately describe the overall response. The advantage of this model is that it is not limited to any speci®c comparison medium but can be used as a data analysis tool. For instance, having measured the Young's modulus for dierent volume fractions of a certain composite it is possible to calculate c02 and use the non-dilute approach to predict the other elastic constants. This is in fact what has been done in the comparison between the present results and the results from Cherkaoui et al. (1995). Furthermore, the present approach is compared to the HSW bounds and the results are seen to lie within these bounds.
Appendix A The components of the transformation matrix aij in terms of the Euler angles, h; /; b, are given as: 8 a11 cos h cos / cos b sin / sin b; > > > > a12 cos h cos / sin b sin / cos b; > > > > > a13 sin h cos /; > > > < a21 cos h sin / cos b cos / sin b; aij a22 cos h sin / sin b cos / cos b; > > a23 sin h sin /; > > > > > a sin h cos b; 31 > > > > a sin h sin b; > : 32 a33 cos h:
References Benveniste, Y., 1987. On a new approach to the application of Mori±Tanaka's theory in composite materials. Mech. Mater. 6, 147. Benveniste, Y., Dvorak, G.J., Chen, T., 1989. Stress ®elds in composites with coated inclusions. Mech. Mater. 7, 305. Benveniste, Y., Dvorak, G.J., Chen, T., 1991a. On diagonal and elastic symmetry of the approximate eective stiness tensor of heterogeneous media. J. Mech. Phys. Solids. 39, 927. Benveniste, Y., Dvorak, G.J., Chen, T., 1991b. On eective properties of composites with coated cylindrically orthotropic ®bres. Mech. Mater. 12, 289. Bhattacharyya, A., Weng, G.J., 1996. Theoretical calculation of the stress±strain behaviour of dual-phase metals with randomly oriented spheroidal inclusions. Metall. Mater. Trans. A 27, 2359.
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
543
Bower, D.I., 1972. Investigation of molecular orientation distributions by polarised Raman scattering and polarised ¯ourescence. J. Polym. Sci. 10, 2135. Bozarth, J.M., Gillespie, J.W., McCullough, R.L., 1987. Fiber orientation and its eect upon thermoelastic properties of short carbon ®ber reinforced poly(etheretherketone) (PEEK). Polym. Comp. 8, 74. Chen, T., Dvorak, G.J., Benveniste, Y., 1990. Stress ®elds in composites reinforced with cylindrically orthotropic ®bres. Mech. Mater. 9, 17. Cherkaoui, M., Sabar, H., Berveiller, M., 1995. Elastic composites with coated reinforcements: a micromechanical approach for nonhomothetic topology. Int. J. Eng. Sci. 33, 829. Dvorak, G.J., Bahei-El-Din, A., 1997. Inelastic composite materials: transformation analysis and experiments. In: Suquet, P. (Ed.), Continuum Micromechanics. Dvorak, G.J., Srinivas, M.V., 1999. New estimates of overall properties of heterogeneous media. J. Mech. Phys. Solids 47, 899. Eshelby, J.D., 1957. The determination of the elastic ®eld of an ellipsoidal inclusion, and related problems. Proc. R. Soc. A 241, 376. Ferrari, M., 1991. Asymmetry and the high concentration limit of the Mori±Tanaka eective medium theory. Mech. Mater. 11, 251. Fotiu, P.A., Nemat-Nasser, S., 1996. Overall properties of elastic±viscoplastic periodic composites. Int. J. Plasticity 12, 163. Hashin, Z., 1962. The elastic moduli of heterogeneous materials. J. Appl. Mech. 29, 143. Hashin, Z., 1963. A variational approach to the theory of the elastic behaviour of multiphase materials. J. Mech. Phys. Solids 11, 140. Hashin, Z., 1965. On elastic behaviour of ®bre reinforced materials of arbitrary phase geometry. J. Mech. Phys. Solids 13, 119. Hill, R., 1963. Elastic properties of reinforced solids: some theoretical principles. J. Mech. Phys. Solids 11, 357. Hill, R., 1965. Continuum micromechanics of elasto-plastic polycrystals. J. Mech. Phys. Solids 13, 89. Huang, J.H., Kuo, W.S., 1996. Micromechanics determination of the eective properties of piezoelectric composites containing spatially oriented short ®bres. Acta Mater. 44 (44), 4889. Iwakuma, T., Nemat-Nasser, S., 1984. Finite elastic±plastic deformation of polycrystalline metals. Proc. R. Soc. Lond. A 394, 87. Kuo, W.S., Huang, J.H., 1997. On the eective electroelastic properties of piezoelectric composites containing spatially oriented inclusions. Int. J. Solids Struct. 34 (19), 2445. Kwon, P., Dharan, C.K.H., 1995. Eective moduli of high volume fraction particulate composites. Acta Metall. Mater. 43 (3), 1141. Li, J.Y., 1999. On micromechanics approximation for the eective thermoelastic moduli of multiphase composite materials. Mech. Mater. 31, 149. Maekawa, Z.I., Hamada, H., Yokoyama, A, 1989. Lamination theory of composite material with complex ®ber orientation distribution. In: Proceedings of the ICCS 5, pp. 701. Marzari, N., Ferrari, M., 1992. Textural and micromorphological eects on the overall elastic response of macroscopically anisotropic composites. J. Appl. Mech. 59, 269. Molinari, A., Mouden, Ess., 1996. The problem of elastic inclusions at ®nite concentration. Int. J. Solids Struct. 33 (20±22), 3131. Mori, T., Tanaka, K., 1973. Average stress in matrix and average elastic energy of materials with mis®tting inclusions. Acta Metall. 21, 571. Nebozhyn, M.V., Ponte Castaneda, P., 1999. The second order procedure: exact vs approximate results for isotropic, two-phase composites. J. Mech. Phys. Solids 47, 2171. Nemat-Nasser, S., Hori, M., 1993. Micromechanics: Overall Properties of Heterogeneous Materials. In: North-Holland Series in Applied Mathematics and Mechanics. North-Holland, Amsterdam. Nemat-Nasser, S., Obata, M., 1986. Rate dependent ®nite elasto±plastic deformation of polycrystals. Proc. R. Soc. Lond. A 407, 343. Ngolle, A., Pera, J., 1997. Microstructural based modelling of the elastic modulus of ®ber reinforced cement composites. Adv. Cem. Bas. Mat. 6, 130. Pipes, R.B., McCullough, R.L., Taggart, D.G., 1982. Behavior of discontinuous ®ber composites: ®ber orientation. Polym. Comp. 3, 34. Castaneda, P., Willis, J.R., 1995. The eect of spatial distribution on the eective behaviour of composite materials and cracked media. J. Mech. Phys. Solids 43 (12), 1919. Ponte Castaneda, P., 1996. Exact second-order estimates for eective mechanical properties of nonlinear composite materials. J. Mech. Phys. Solids 44, 827. Sayers, C.M., 1992. Elastic anisotropy of short ®bre reinforced composites. Int. J. Solids Struct. 29 (23), 2933. Schjùdt-Thomsen, J., Pyrz, R., 2000. Overall creep modelling of short ®bre reinforced composites with weakened interfaces and complex ®bre orientation distributions. Mech. Mater. 32, 349. Takao, Y., Chou, T.W., Taya, M., 1982. Eective longitudinal Young's modulus of misoriented short ®ber composites. J. Appl. Mech. 49, 536. Torquato, S., 1997. Eective stiness tensor of composite media- I. Exact series expansions. J. Mech. Phys. Solids 45 (9), 1421. Torquato, S., 1998a. Eective stiness tensor of composite media- II. Applications to isotropic dispersions. J. Mech. Phys. Solids 46 (8), 1411. Torquato, S., 1998b. Morphology and eective properties of disordered heterogeneous media. Int. J. Solids Struct. 35 (19), 2385. Walpole, L.J., 1969. On the overall elastic moduli of composite materials. J. Mech. Phys. Solids 17, 235.
544
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
Walpole, L.J., 1985. Evaluation of the elastic moduli of a transeversely isotropic aggregate of cubic crystals. J. Mech. Phys. Solids 33 (6), 623. Weng, G.J., 1990. The theoretical connection between Mori±Tanak's theory and the Hashin±Shtrikman±Walpole bounds. Int. J. Eng. Sci. 28, 1111.
www.elsevier.com/locate/mechmat
The Mori±Tanaka stiness tensor: diagonal symmetry, complex ®bre orientations and non-dilute volume fractions J. Schjùdt-Thomsen *, R. Pyrz Institute of Mechanical Engineering, Aalborg University, Pontoppidanstraede 101, DK 9220 Aalborg East, Denmark Received 3 April 2001; received in revised form 9 June 2001
Abstract This paper considers the diagonal symmetry of the stiness tensors predicted using the Mori-Tanaka (MT) approach. Since the MT approach may yield asymmetric stiness tensors this paper considers an alternative approach to ensure the diagonal symmetry. Furthermore, an extension of the MT approach to non-dilute volume fractions is considered. It is shown that the MT approach by Benveniste [Mech. Mater. 6 (1987) 147], may not yield diagonally symmetric stiness tensors in situations when a statistical ®bre orientation function is incorporated. Only when the inclusions are spherical and randomly distributed or fully aligned will the MT predicted stiness tensor be diagonally symmetric. The extension to non-dilute volume fractions is seen to lie within the Hashin±Shtrikman±Walpole (HSW) bounds and comparison to other approaches and experiments found in the literature shows very good agreement. Ó 2001 Elsevier Science Ltd. All rights reserved. Keywords: Mori±Tanaka stiness tensor; Diagonal symmetry; Complex ®bre orientations; Non-dilute volume fractions
1. Introduction Traditional continuum mechanics is a branch of the physical sciences concerned with the deformations and motions of continuous material media under the in¯uence of external eects. These media are in the form of a collection of material points interconnected by some internal forces. A basic assumption of continuum mechanics is that the mass, stress and strain are regarded as essentially uniform within an in®nitesimal material point and its surrounding neighbourhood. However, when dealing with non-homogeneous materials such as composite materials, cellular materials or polycrystals, the aforementioned assumption may no longer be valid, since the material point and its neighbourhood is generally not uniform, but consist of various constituents with dierent properties and shape. On the other hand the traditional continuum mechanic is a well established ®eld of science and thus applying continuum mechanics approaches to non-homogeneous media seems bene®cial. This requires the concept of a representative
*
Corresponding author. Tel.: +45-9635-9330; fax: +45-9815-1411. E-mail address: [email protected] (J. Schjùdt-Thomsen).
0167-6636/01/$ - see front matter Ó 2001 Elsevier Science Ltd. All rights reserved. PII: S 0 1 6 7 - 6 6 3 6 ( 0 1 ) 0 0 0 7 2 - 2
532
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
volume element (RVE), see e.g., Hill (1963) or Nemat-Nasser and Hori (1993), which is then used in the continuum theory as the material point. The RVE must include a large number of constituents and be statistically representative of the local continuum properties. This approach has been taken by several authors to calculate the overall elastic properties of various non-homogeneous materials, Hill (1963), Hashin (1963, 1965), Walpole (1969, 1985), Kwon and Dharan (1995), Ponte Castaneda and Willis (1995), Molinari and El Mouden (1996), Torquato (1997, 1998a,b) and Dvorak and Srinivas (1999). Non-elastic properties were considered by Hill (1965), Fotiu and NematNasser (1996), Bhattacharyya and Weng (1996), Ponte Castaneda (1996), Nebozhyn and Ponte Castaneda (1999), and Schjùdt-Thomsen and Pyrz (2000). Several other investigations have been carried out and references may be found in Nemat-Nasser and Hori (1993). One of the methods which has been used extensively is the Mori-Tanaka (MT) model, based on the work by Mori and Tanaka (1973) and Eshelby (1957). This model has proven to be quite accurate in predicting the eective properties of various materials with either random orientation or total alignment of the reinforcing phases. However, for composites containing inclusions of various shapes or orientations the MT predicted stiness tensor may violate the symmetry requirement that Cijkl Cklij (Ponte Castaneda and Willis, 1995; Li, 1999; Ferrari, 1991; Benveniste et al., 1991a). This is indeed a drawback of the MT model, but one should bear in mind that the MT model was originally proposed for composites containing inclusions of similar shape. Benveniste et al. (1991a) found that the MT stiness tensor is symmetric for twophase composites and for multi-phase composites if the inclusions are of similar shape and have the same orientation with respect to a ®xed reference frame. However, the same authors (Benveniste et al., 1991b), found that if the ``typical inclusion'' is a coated particle, then the implementation can be carried out within the framework of two-phase composites and thus gives diagonally symmetric stiness tensors (Benveniste et al., 1989; Chen et al., 1990). Li (1999) also concludes that it is the extension from two-phase composites to multi-phase composites with dierently shaped reinforcements, that causes the con¯icting results ± not the original MT assumption. In most of the work on composite materials using the MT approach the inclusions have been assumed to be either aligned or randomly dispersed in the matrix. However, in real short ®bre composites the ®bres are neither aligned nor dispersed at random directions. Fibre orientation is governed by the nature of the complicated ¯ow ®eld of the matrix during ®lling of the mould. The consequence of the ¯ow ®eld is that the ®bres are oriented at dierent directions, somewhere between the two extremes; random and fully aligned. Investigations considering various ®bre orientations, using dierent approaches have been considered by Pipes et al. (1982), Takao et al. (1982), Bozarth et al. (1987), Maekawa et al. (1989), Marzari and Ferrari (1992), Sayers (1992), Ngolle and Pera (1997) and Schjùdt-Thomsen and Pyrz (2000). Furthermore, Torquato (1997, 1998a,b), Molinari and El Mouden (1996) and Ponte Castaneda and Willis (1995) also considered spatial distribution and non-dilute volume fractions of the reinforcement phase. In light of the above observations the present paper focuses on a numerical study of the MT predicted stiness tensor for various ®bre orientations using a statistical ®bre orientation function. Furthermore, should the MT stiness tensor be non-symmetric an alternative approach is proposed which will yield a diagonally symmetric stiness tensor. Additionally, an approach is proposed which modi®es the MT approach to take non-dilute volume fractions into account. In the paper a boldface ordinary capital letter will denote a fourth-order tensor, A Aijkl , the inner product is AB Aijkl Bklmn and the inverse of a tensor is denoted as A 1 . 2. Eective stiness of a two-phase composite An estimate of the eective stiness of a two-phase composite can be obtained from Benveniste (1987) analysis of the MT model. Without going into detail the ®nal expression for the stiness is given as:
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
L MT L1 c2 f L2
L1 Tgc1 I c2 fTg
1
533
1
where indices ``1'' and ``2'' denote the matrix and reinforcing phase respectively and I is the fourth-order identity tensor. Curly brackets denote averaging over all possible directions. The tensor T relates the uniform strain in an ellipsoidal particle embedded in an all matrix medium to the average matrix strain. e 2 Te 1 ;
2
where T I SL1 1 L2
L1
1
3
and S is the Eshelby tensor, and thus the strain concentration tensor, A2 , is A2 T c1 I c2 T
1
4
In order to obtain the eective properties of a heterogeneous material the RVE is needed. Obviously the validity of this approach relies heavily on the de®nition of the RVE. The RVE should be the smallest possible unit which include all distinct elements of the microstructure that have a ®rst-order in¯uence on the properties and should represent the relative dispersion of the phases in a statistically representative manner (Nemat-Nasser and Hori, 1993). The present method is based on a transformation of fourth-order tensors followed by a spatial averaging over the RVE corresponding to the material point in the continuum, thus giving the overall averaged properties. The ®bre orientations in the RVE are now described through an orientation distribution function in the Euler angles, g h; /; b. The stiness Lijkl in the X1 X2 X3 coordinate system is related to the stiness Lpqrs in the local ®bre coordinate system through the transformation rule of fourth-order tensors as Lijkl aip ajq akr als Lpqrs ;
5
where aij is a function of the Euler angles, see Appendix A. In order to incorporate the stiness contribution from all orientations, aip ajq akr als Lpqrs is multiplied by the ®bre orientation distribution function and integrated over all possible orientations to obtain the weighted orientation average stiness as
Fig. 1. Euler angles de®ning the relation between the local coordinate system, x1 x2 x3 and the global coordinate system, X1 X2 X3 .
534
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
Z hLijkl i
2p
Z
0
2p
0
Z
p
g h; /; baip ajq akr als Lpqrs sin h dh d/ db;
0
6
where hi denotes averaging over all possible directions. g h; /; b is the orientation distribution function de®ned in the Euler coordinates h; /; b, Lpqrs is the stiness tensor of the RVE. The transformation matrix in terms of the Euler angles de®nes the orientation of the axes x1 x2 x3 ®xed in the RVE with respect to the axis X1 X2 X3 ®xed in the continuum, as shown in Fig.1. The concept of averaging over various orientations to estimate overall properties is not a new one. The approach is generally accepted and has been used by several authors in a variety of ®elds, for investigating the molecular orientation distribution of polymers (Bower, 1972), and to describe the electroelastic properties of piezoelectric materials (Huang and Kuo, 1996; Kuo and Huang, 1997). The major dierence between the present approach and the previous ones is that this work considers the full stiness tensor whereas the previous ones considered only the in-plane stiness properties of laminates or the Young's and shear moduli of the short ®bre composites. Furthermore, the present approach will yield a diagonally symmetric MT stiness tensor for non-trivial ®bre orientations, which has not been carried out previously. Reconsidering Eq. (1), the eective stiness of a two-phase material can be expressed in terms of an orientation dependent part and an orientation independent part as: 1 L MT pqrs Lpqrs fLpqrs g:
7
The stiness tensor Lpqrs in Eq. (6) is now taken to be the orientation dependent part of Eq. (7), for a composite with fully aligned ®bres. The averaging procedure will then provide the overall stiness of a composite with ®bre orientation described through the function g h; /; b. One of the objectives of the present work is to obtain a diagonally symmetric stiness tensor for arbitrary ®bre orientations. The proof that hLijkl i is indeed diagonally symmetric is quite straightforward. Since Lpqrs in Eq. (6) is calculated for a composite with fully aligned ®bres, Lpqrs Lrspq , i.e., it is symmetric, thus Lijkl aip ajq akr als Lpqrs ; Lklij akr als aip ajq Lrspq and it is easily seen that these two expressions are identical and the resulting eective stiness tensor, given by Eq. (6) is diagonally symmetric. 2.1. Orientation distribution function The orientation distribution function used in this work was used previously by Schjùdt-Thomsen and Pyrz (2000), and is given as (Maekawa et al., 1989) g h; /; b R hb ha
sin h
2P 1
sin h
cos h
2P 1
2Q 1
cos h
2Q 1
dh
;
8
where ha and hb are the upper and lower limits of the angle h being present in the distribution, 0 6 ha ; ha 6 p and P P 1=2, Q P 1=2. As can be seen from Eq. (8) the ®bre orientation function describes the orientation only in terms of h, meaning that the composite will be isotropic or transversely isotropic with respect to the X3 axis, depending on the values of P and Q.
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
535
3. Non-dilute volume fractions As was pointed out in Benveniste (1987), particle interaction is taken into account to a certain extent, since for a composite reinforced with spherical particles his reformulation of the MT-model coincides with the composites spheres assemblage of Hashin (1962). However, comparing the results from the MT approach with higher order approximations as the one proposed by Cherkaoui et al. (1995) the traditional MT approach gives smaller values as compared to the mentioned approach and experiments. The strain concentration tensors are usually calculated using Eq. (3), where the Eshelby tensor is calculated using the matrix material as the comparison medium. In order to extent the MT approach to non-dilute volume fractions the concentration tensor needs to be calculated dierently. This may be done in dierent ways. Molinari and El Mouden (1996) and Cherkaoui et al. (1995) formulated an integral equation of the unknown strain, leading to a set of linear equations from which the concentration tensors were found in a self consistent like manner. Another approach may be simply to replace L1 in Eq. (3) by another suitably chosen comparison medium, that is 1 T I SL0 1 L2 L0 ; 9 where L0 is the stiness tensor of the comparison medium, in which the Eshelby tensor must be calculated. An admissible choice of comparison medium is such that the overall stiness do not violate the bounds (Dvorak and Srinivas, 1999). In the present analysis the following approach is taken. Initially T is calculated as 1 T I SL1 1 L2 L1 and used to obtain the isotropic MT stiness, (Lisotropic ), with a volume fraction c02 less than the actual volume fraction c2 , that is 0 6 c02 6 c2 and used as the comparison medium stiness h i 1 1 T0 I SLisotropic L2 Lisotropic : 10 Now T0 is replacing T in Eq. (1). Using the actual volume fraction now gives the overall eective stiness. This approach is somewhat similar to the successive embedding procedure of Nemat-Nasser and Hori (1993, p. 330). Using Eqs. (1) and (4) gives the eective MT stiness in a shorter form as L MT L1 c2 L2
L1 A2 :
11
The eective stiness tensor determined from the non-dilute approach can be written more generally as L N 1 L1 c2 L2
N
L1 A2 ;
12
N A2
where is a function of c02 and N denotes the order of the approximation, i.e., for N 1 and c02 c2 the original MT stiness tensor is obtained and if N ! 1 and c02 c2 the self consistent method is obtained. For N 2 the stiness tensor from the non-dilute approach is obtained. 4. Results and discussion For illustrative purposes the matrix and ®bre properties were chosen to be those of polypropylene and glass ®bres, with the following properties: E 1 1:5 GPa, m 1 0:32, E 2 73 GPa and m 2 0:22 and the ®bre volume fraction, c2 0:35 in the subsequent calculations. However, it should be emphasized that the approach can ± of course ± be used for other material combinations.
536
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
Now two approaches are followed: (1) the overall stiness tensor, denoted as, LhMTi , is obtained by averaging the terms in curly brackets in Eq. (1) over all possible orientations, described the orientation function, Eq. (8) and (2) The stiness tensor predicted by the MT method is inserted into Eq. (6), thus giving the overall stiness for a composite with orientations described by Eq. (8), denoted as L . 4.1. Random spherical particles In this situation the particles ``orientation'' is described by taking the orientation parameters P and Q equal to Q P 0:5, Fig. 2, thus giving an isotropic material. In such cases the material eective properties are characterised through the eective bulk and shear moduli. For the sake of con®dence, however, the full 6 6 stiness matrices are shown. 2
LhMT i
4033:36 6 1674:13 6 6 1674:13 6 6 6 0 6 4 0 0 2
4033:36 6 1674:13 6 6 1674:13 6 L 6 6 0 6 4 0 0
1674:13 4033:36 1674:13 0 0 0
1674:13 4033:36 1674:13 0 0 0
1674:13 0 1674:13 0 4033:36 0 0 1179:61 0 0 0 0
1674:13 1674:13 4033:36 0 0 0
0 0 0 1179:61 0 0
0 0 0 0 1179:61 0
3 0 7 0 7 7 0 7 7; 7 0 7 5 0 1179:61
3 0 0 7 0 0 7 7 0 0 7 7: 7 0 0 7 5 1179:61 0 0 1179:61
From the two results it is seen that the stiness matrices are indeed identical and diagonally symmetric ± as they should be.
Fig. 2. Fibre orientation distribution function for random and partial aligned ®bres.
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
537
4.2. Partially aligned ellipsoidal particles The inclusions are now in the shape of ellipsoids elongated in the x3 direction and having aspect ratios, a3 =a1 a3 =a2 30. The ®bres are oriented according to Q 11 and P 0:5, which corresponds to partial alignment along the X3 axis and is shown in Fig. 2. The composite material is now transversely isotropic, with the X1 X2 plane as the isotropic plane. The eective stiness matrices are now 3 2 4965:42 2081:40 3362:94 0 0 0 7 6 2081:40 4965:42 3362:94 0 0 0 7 6 7 6 4481:31 4481:31 17417:51 0 0 0 hMTi 7; 6 L 6 7 0 0 0 3047:47 0 0 7 6 5 4 0 0 0 0 3047:47 0 0 0 0 0 0 1416:72 2
4573:82 6 1888:16 6 6 2948:46 L 6 6 0 6 4 0 0
1888:16 4573:82 2948:46 0 0 0
2948:46 2948:46 14532:12 0 0 0
0 0 0 2528:65 0 0
3 0 0 7 0 0 7 7 0 0 7: 7 0 0 7 5 2528:65 0 0 1342:82
Now the stiness matrices are dierent. The LhMTi stiness matrix is not diagonally symmetric and the asymmetry is clearly signi®cant as was also shown by Benveniste et al. (1991a) for circular discs and continuous ®bres. The main reason for the asymmetry may be due to the strain concentration tensors which are generally not diagonally symmetric (Dvorak and Bahei-El-Din, 1997). The reason is that the concentration tensors must satisfy the requirement that c1 A1 c2 A2 I (for a two-phase material). However, unless all the inclusions are aligned the concentration tensors may not satisfy this consistency requirement. To overcome this inconsistency the concentration tensors must be normalised (Nemat-Nasser and Hori, 1993). Walpole (1969) suggested this to be done by replacing the concentration tensors by r Ar hAr i 1 ; A
13
where hAr i indicates summation over all orientations, shapes and sizes. However, this normalisation has not been considered in this work, but it was used by Iwakuma and Nemat-Nasser (1984) and Nemat-Nasser and Obata (1986) in their self-consistent modelling of ®nite deformations of polycrystalline solids. Increasing the degree of alignment by taking P 0:5 and Q 100, as shown in Fig. 3 makes the asymmetry of MT predicted stiness tensor even more signi®cant as shown below 3 2 3199:57 1390:32 1636:86 0 0 0 6 1390:32 3199:57 1636:86 0 0 0 7 7 6 6 3366:90 3366:90 17873:38 0 0 0 7 7 LhMTi 6 6 0 0 0 1270:75 0 0 7 7 6 4 0 0 0 0 1270:75 0 5 0 0 0 0 0 825:88 and the L tensor is
538
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
Fig. 3. Fibre orientation distribution function for nearly aligned ®bres.
2
3721:75 6 1637:26 6 6 1870:93 L 6 6 0 6 4 0 0
1637:26 3721:75 1870:93 0 0 0
1870:93 1870:93 21048:15 0 0 0
0 0 0 1475:96 0 0
3 0 0 7 0 0 7 7 0 0 7: 7 0 0 7 5 1475:96 0 0 1042:24
The L tensor is again diagonally symmetric. Finally, as a last example the situation where the ®bres are practically fully aligned. A perfect alignment of the ®bres can be described through a Dirac delta function of the angle h, however, in this example the parameters of the orientation distribution function are P 0:5 and Q 250, and is shown in Fig. 3 3 2 3701:88 1633:87 1681:11 0 0 0 7 6 1633:87 3701:88 1681:11 0 0 0 7 6 7 6 1681:11 1681:11 21853:93 0 0 0 7 6 L 6 7 0 0 0 1288:58 0 0 7 6 5 4 0 0 0 0 1288:58 0 0 0 0 0 0 1034:00 and for comparison the results from the original MT approach for fully aligned inclusions 3 2 3696:27 1634:23 1539:30 0 0 0 7 6 1634:23 3696:27 1539:30 0 0 0 7 6 7 6 1539:30 1539:30 22431:67 0 0 0 MT; align 7: 6 L 6 7 0 0 0 1148:44 0 0 7 6 5 4 0 0 0 0 1148:44 0 0 0 0 0 0 1031:01 4.3. Non-dilute volume fractions Previously it has been shown that the averaging procedure yields diagonally symmetric stiness tensors, so this approach has also been used in the situation when non-dilute volume fractions are considered. In order to do some checking on the present approach the Hashin±Shtrikman±Walpole (HSW), bounds have
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
539
been calculated. The bounds are calculated using the connections between the MT model and the HSW bounds (Weng, 1990). For random spherical particles it was shown by Weng (1990) that the MT moduli will coincide with the HSW lower/upper bound, denoted as L and L , if the matrix is the softest/hardest of the phases, respectively. In this context softest/hardest means that Lr L0 is positive/negative semi-de®nite for all r and L0 is the stiness of the comparison medium. First, consider the situation when the inclusions are randomly dispersed spherical particles, thus the lower HSW bound coincides with the moduli of the MT model when the matrix stiness is taken to be that of polypropylene 3 2 4033:36 1674:13 1674:13 0 0 0 7 6 1674:13 4033:36 1674:13 0 0 0 7 6 7 6 1674:13 1674:13 4033:36 0 0 0 7 L 6 7 6 0 0 0 1179:61 0 0 7 6 5 4 0 0 0 0 1179:61 0 0 0 0 0 0 1179:61 and the upper bound is 2 38941:73 6 8458:58 6 6 8458:58 L 6 6 0 6 4 0 0
obtained when the matrix phase is taken to be that of the glass particles 3 8458:58 8458:58 0 0 0 7 38941:73 8458:58 0 0 0 7 7 8458:58 38941:73 0 0 0 7: 7 0 0 15241:57 0 0 7 5 0 0 0 15241:57 0 0 0 0 0 15241:57
The stiness prediction from the present approach is found to be 3 2 4707:56 1867:28 1867:28 0 0 0 7 6 1867:28 4707:56 1867:28 0 0 0 7 6 7 6 1867:28 1867:28 4707:56 0 0 0 7; 6 L 6 7 0 0 0 1420:13 0 0 7 6 5 4 0 0 0 0 1420:13 0 0 0 0 0 0 1420:13 where the stiness of the comparison medium is initially L0 L1 and in the next step L0 equals that of the MT approach for ®bre volume fractions equal to half the value of the actual composite. Taking c02 12 c2 is only done for illustrative purposes. Indeed, any choice of comparison medium leading to an eective stiness that does not violate the HSW bounds is admissible (Dvorak and Srinivas, 1999). 4.4. Randomly oriented ®bres Now the inclusions are taken to be ®bres with aspect ratios a3 =a1 a3 =a2 30 and the distribution parameters are P Q 0:5. The lower bound is now 3 2 7480:66 2801:32 2801:32 0 0 0 7 6 2801:32 7480:66 2801:32 0 0 0 7 6 7 6 2801:32 2801:32 7480:66 0 0 0 7 6 L 6 7 0 0 0 2339:67 0 0 7 6 5 4 0 0 0 0 2339:67 0 0 0 0 0 0 2339:67
540
and the upper bound 2 37234:66 6 7747:95 6 6 7747:95 L 6 6 0 6 4 0 0
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
7747:95 37234:66 7747:95 0 0 0
7747:95 0 7747:95 0 37234:66 0 0 14743:35 0 0 0 0
and the present approach yields 2 9524:65 3364:71 3364:71 6 3364:71 9524:65 3364:71 6 6 3364:71 3364:71 9524:65 L 6 6 0 0 0 6 4 0 0 0 0 0 0
0 0 0 3079:97 0 0
0 0 0 0 14743:35 0
3 0 7 0 7 7 0 7 7 0 7 5 0 14743:35
3 0 0 7 0 0 7 7 0 0 7: 7 0 0 7 5 3079:97 0 0 3079:97
From the two set of results above the eective bulk and shear moduli are calculated along with the upper and lower bounds on these properties. For the random spherical particles, the lower bound on the eective bulk modulus j 2460 MPa, the upper bound is j 18 620 MPa and the present approach yields j 2814 MPa. The shear moduli are: l 1180, l 15 241 and l 1420 MPa. For the random oriented ®bres the same properties are: the lower bound on the eective bulk modulus j 4361 MPa, the upper bound is j 17 577 MPa and the present approach yields j 5418 MPa and the shear moduli l 2340, l 14 743 and l 3080 MPa. It is clearly seen that in the situations considered the present approach does not violate the HSW bounds since j 6 j 6 j and l 6 l 6 l . The composite Young's modulus for the polypropylene/glass particle composite has been calculated and is shown in Fig. 4, for the full range of volume fractions. It is seen that the predicted Young's modulus does not violate the bounds. The curve of the Young's modulus from the original MT approach is not very smooth for larger volume fractions. This may be due to the fact that interaction is not accounted for suciently, whereas the curve from the non-dilute MT approach is much smoother, indicating that the use of a dierent comparison medium enhances the models ability to account for interaction.
Fig. 4. Upper and lower HSW bounds along with the predictions of the non-dilute MT method. (Glass particle/polypropylene composite.)
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
541
Fig. 5. Normalised Young's modulus. Theoretical data and experimental data are from Cherkaoui et al. (1995).
Finally, the results from this approach is compared to results from Cherkaoui et al. (1995), for an epoxy matrix reinforced with glass particles. The material properties are: E 1 2:9 GPa, m 1 0:4, E 2 73 GPa and m 2 0:2, c02 0:15 and c2 0:3. The results are shown in Fig. 5. Cherkaoui et al. (1995) found that the calculated composite Young's modulus is Ec 6:3 GPa and the c experimentally determined Young's modulus is Eexp 5:97 GPa. The present approach yields an eective Young's modulus, E 6:05 GPa, whereas the original MT approach gives EMT 5:3 GPa. Very good agreement is found for the present approach compared to other investigations. Finally, the Young's modulus for this composite is calculated and shown in Fig. 6. Comparing this ®gure with Fig. 4 and remembering that the MT approach is coincident with the lower bound it is believed that for larger mismatch between the constituents moduli (E 2 =E 1 48:7 in Fig. 4 and E 2 =E 1 25:2 in Fig. 6) the non-dilute MT approach yields more reasonable results.
Fig. 6. Upper and lower HSW bounds along with the predictions of the non-dilute MT method. (Glass particle/epoxy composite.)
542
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
5. Conclusions In the present paper the diagonal symmetry of the stiness tensors predicted by the MT model have been considered for various ®bre orientations. It was found that this method does not provide diagonally symmetric stiness tensors if the averaging is carried out only on the concentration tensor T as described in Benveniste (1987). However, averaging the stiness tensor for a composite with unidirectionally aligned ®bres according to a ®bre orientation function yields diagonally symmetric eective stiness tensors for non-trivial ®bre orientations and coincides with the predictions of the MT method for random spherical particles and fully aligned ®bres. The MT approach has been extended to incorporate the eect of non-dilute volume fractions in a manner similar to the successive embedding procedure and results have been compared to results in the literature and very good agreement is found. The non-dilute approach provides a ¯exible method for predicting the eective properties of various composite materials. Due to the complex nature of composite materials it is dicult to know the degree of interaction beforehand and thus which model may adequately describe the overall response. The advantage of this model is that it is not limited to any speci®c comparison medium but can be used as a data analysis tool. For instance, having measured the Young's modulus for dierent volume fractions of a certain composite it is possible to calculate c02 and use the non-dilute approach to predict the other elastic constants. This is in fact what has been done in the comparison between the present results and the results from Cherkaoui et al. (1995). Furthermore, the present approach is compared to the HSW bounds and the results are seen to lie within these bounds.
Appendix A The components of the transformation matrix aij in terms of the Euler angles, h; /; b, are given as: 8 a11 cos h cos / cos b sin / sin b; > > > > a12 cos h cos / sin b sin / cos b; > > > > > a13 sin h cos /; > > > < a21 cos h sin / cos b cos / sin b; aij a22 cos h sin / sin b cos / cos b; > > a23 sin h sin /; > > > > > a sin h cos b; 31 > > > > a sin h sin b; > : 32 a33 cos h:
References Benveniste, Y., 1987. On a new approach to the application of Mori±Tanaka's theory in composite materials. Mech. Mater. 6, 147. Benveniste, Y., Dvorak, G.J., Chen, T., 1989. Stress ®elds in composites with coated inclusions. Mech. Mater. 7, 305. Benveniste, Y., Dvorak, G.J., Chen, T., 1991a. On diagonal and elastic symmetry of the approximate eective stiness tensor of heterogeneous media. J. Mech. Phys. Solids. 39, 927. Benveniste, Y., Dvorak, G.J., Chen, T., 1991b. On eective properties of composites with coated cylindrically orthotropic ®bres. Mech. Mater. 12, 289. Bhattacharyya, A., Weng, G.J., 1996. Theoretical calculation of the stress±strain behaviour of dual-phase metals with randomly oriented spheroidal inclusions. Metall. Mater. Trans. A 27, 2359.
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
543
Bower, D.I., 1972. Investigation of molecular orientation distributions by polarised Raman scattering and polarised ¯ourescence. J. Polym. Sci. 10, 2135. Bozarth, J.M., Gillespie, J.W., McCullough, R.L., 1987. Fiber orientation and its eect upon thermoelastic properties of short carbon ®ber reinforced poly(etheretherketone) (PEEK). Polym. Comp. 8, 74. Chen, T., Dvorak, G.J., Benveniste, Y., 1990. Stress ®elds in composites reinforced with cylindrically orthotropic ®bres. Mech. Mater. 9, 17. Cherkaoui, M., Sabar, H., Berveiller, M., 1995. Elastic composites with coated reinforcements: a micromechanical approach for nonhomothetic topology. Int. J. Eng. Sci. 33, 829. Dvorak, G.J., Bahei-El-Din, A., 1997. Inelastic composite materials: transformation analysis and experiments. In: Suquet, P. (Ed.), Continuum Micromechanics. Dvorak, G.J., Srinivas, M.V., 1999. New estimates of overall properties of heterogeneous media. J. Mech. Phys. Solids 47, 899. Eshelby, J.D., 1957. The determination of the elastic ®eld of an ellipsoidal inclusion, and related problems. Proc. R. Soc. A 241, 376. Ferrari, M., 1991. Asymmetry and the high concentration limit of the Mori±Tanaka eective medium theory. Mech. Mater. 11, 251. Fotiu, P.A., Nemat-Nasser, S., 1996. Overall properties of elastic±viscoplastic periodic composites. Int. J. Plasticity 12, 163. Hashin, Z., 1962. The elastic moduli of heterogeneous materials. J. Appl. Mech. 29, 143. Hashin, Z., 1963. A variational approach to the theory of the elastic behaviour of multiphase materials. J. Mech. Phys. Solids 11, 140. Hashin, Z., 1965. On elastic behaviour of ®bre reinforced materials of arbitrary phase geometry. J. Mech. Phys. Solids 13, 119. Hill, R., 1963. Elastic properties of reinforced solids: some theoretical principles. J. Mech. Phys. Solids 11, 357. Hill, R., 1965. Continuum micromechanics of elasto-plastic polycrystals. J. Mech. Phys. Solids 13, 89. Huang, J.H., Kuo, W.S., 1996. Micromechanics determination of the eective properties of piezoelectric composites containing spatially oriented short ®bres. Acta Mater. 44 (44), 4889. Iwakuma, T., Nemat-Nasser, S., 1984. Finite elastic±plastic deformation of polycrystalline metals. Proc. R. Soc. Lond. A 394, 87. Kuo, W.S., Huang, J.H., 1997. On the eective electroelastic properties of piezoelectric composites containing spatially oriented inclusions. Int. J. Solids Struct. 34 (19), 2445. Kwon, P., Dharan, C.K.H., 1995. Eective moduli of high volume fraction particulate composites. Acta Metall. Mater. 43 (3), 1141. Li, J.Y., 1999. On micromechanics approximation for the eective thermoelastic moduli of multiphase composite materials. Mech. Mater. 31, 149. Maekawa, Z.I., Hamada, H., Yokoyama, A, 1989. Lamination theory of composite material with complex ®ber orientation distribution. In: Proceedings of the ICCS 5, pp. 701. Marzari, N., Ferrari, M., 1992. Textural and micromorphological eects on the overall elastic response of macroscopically anisotropic composites. J. Appl. Mech. 59, 269. Molinari, A., Mouden, Ess., 1996. The problem of elastic inclusions at ®nite concentration. Int. J. Solids Struct. 33 (20±22), 3131. Mori, T., Tanaka, K., 1973. Average stress in matrix and average elastic energy of materials with mis®tting inclusions. Acta Metall. 21, 571. Nebozhyn, M.V., Ponte Castaneda, P., 1999. The second order procedure: exact vs approximate results for isotropic, two-phase composites. J. Mech. Phys. Solids 47, 2171. Nemat-Nasser, S., Hori, M., 1993. Micromechanics: Overall Properties of Heterogeneous Materials. In: North-Holland Series in Applied Mathematics and Mechanics. North-Holland, Amsterdam. Nemat-Nasser, S., Obata, M., 1986. Rate dependent ®nite elasto±plastic deformation of polycrystals. Proc. R. Soc. Lond. A 407, 343. Ngolle, A., Pera, J., 1997. Microstructural based modelling of the elastic modulus of ®ber reinforced cement composites. Adv. Cem. Bas. Mat. 6, 130. Pipes, R.B., McCullough, R.L., Taggart, D.G., 1982. Behavior of discontinuous ®ber composites: ®ber orientation. Polym. Comp. 3, 34. Castaneda, P., Willis, J.R., 1995. The eect of spatial distribution on the eective behaviour of composite materials and cracked media. J. Mech. Phys. Solids 43 (12), 1919. Ponte Castaneda, P., 1996. Exact second-order estimates for eective mechanical properties of nonlinear composite materials. J. Mech. Phys. Solids 44, 827. Sayers, C.M., 1992. Elastic anisotropy of short ®bre reinforced composites. Int. J. Solids Struct. 29 (23), 2933. Schjùdt-Thomsen, J., Pyrz, R., 2000. Overall creep modelling of short ®bre reinforced composites with weakened interfaces and complex ®bre orientation distributions. Mech. Mater. 32, 349. Takao, Y., Chou, T.W., Taya, M., 1982. Eective longitudinal Young's modulus of misoriented short ®ber composites. J. Appl. Mech. 49, 536. Torquato, S., 1997. Eective stiness tensor of composite media- I. Exact series expansions. J. Mech. Phys. Solids 45 (9), 1421. Torquato, S., 1998a. Eective stiness tensor of composite media- II. Applications to isotropic dispersions. J. Mech. Phys. Solids 46 (8), 1411. Torquato, S., 1998b. Morphology and eective properties of disordered heterogeneous media. Int. J. Solids Struct. 35 (19), 2385. Walpole, L.J., 1969. On the overall elastic moduli of composite materials. J. Mech. Phys. Solids 17, 235.
544
J. Schjùdt-Thomsen, R. Pyrz / Mechanics of Materials 33 (2001) 531±544
Walpole, L.J., 1985. Evaluation of the elastic moduli of a transeversely isotropic aggregate of cubic crystals. J. Mech. Phys. Solids 33 (6), 623. Weng, G.J., 1990. The theoretical connection between Mori±Tanak's theory and the Hashin±Shtrikman±Walpole bounds. Int. J. Eng. Sci. 28, 1111.