## Abstract

This study aimed to develop a method to build a ‘bridge’ between the macro fracture mechanics model and stochastic micromechanics-based properties so that the macro fracture mechanics model can be expanded to the fracture mechanics problem of functionally graded materials (FGMs) with stochastic mechanical properties. An analytical fracture mechanics model is developed to predict the stress intensity factors (SIFs) in FGMs with stochastic uncertainties in phase volume fractions. Considering the stochastic description of the phase volume fractions, a micromechanics-based method is developed to derive the explicit probabilistic characteristics of the effective properties of the FGMs so that the stochastic mechanical properties can be combined with the macro fracture mechanics model. A thought for choosing the samples efficiently is proposed so that the stable probabilistic characteristic of SIFs can be obtained with a very small sample size. The probability density function of SIFs can be determined by developing a histogram from the generated samples. The present method may provide a thought to establish an analytical model for the crack problems of FGMs with stochastic properties.

## 1. Introduction

Functionally graded materials (FGMs) are one class of non-homogeneous materials with a gradual transition between two or more phases, producing a variation in the physical characteristics of the composite at the macroscopic or continuum level to meet a desired functional performance. Compared with the traditional homogeneous materials, FGMs have many advantages. The bonding strength can be improved, and mismatch of material properties can be reduced by using FGMs. At the same time, there are still many challenging problems in the study of the mechanical behaviour of FGMs. Of them, facture problems have attracted great attention from the international research community.

Firstly, the effective properties must be determined before studying the fracture behaviour of FGMs. For actual FGMs, the distributions of their effective properties may be very complex. However, in order to obtain analytical solutions of the controlling differential equations, most of the published papers have assumed the effective properties to be very particular functions (e.g. exponential functions). Although this assumption is very particular, it may not be practical (Zuiker 1995). It was pointed out by Birman & Byrd (2007) that material properties evaluated according to theoretical models often disagree with the measured values of FGM constants, which indicates that a probabilistic approach to homogenization, accounting for uncertainty in the actual material distribution throughout the volume, may be justified. In classical micromechanics, there are some models for predicting effective properties, such as Eshelby's equivalent inclusion theory (Eshelby 1957), the rules of mixture (Mura 1991), the Hashin–Shtrikman model (Hashin & Strikman 1963), the self-consistent model (Hill 1965) and the Mori–Tanaka model (Mori & Tanaka 1973). For FGMs with locally pairwise interactions between particles, micromechanics thermoelastic models were developed by Yin *et al.* (2007). In the classical micromechanics models, the microstructural features (such as the phase volume fraction) must be defined locally with absolute certainty (Rahman & Chakraborty 2007). However, for actual FGMs, an issue of concern is the fact that manufacturing a desired gradient to exact specifications is often not possible, as considerable randomness can be observed from sample to sample (Ferrante & Graham-Brady 2005). A schematic of this can be seen in figure 1 (Ilschner 1996; Ferrante & Graham-Brady 2005). As shown in figure 1, Ferrante & Graham-Brady (2005) pointed out that the local volume fractions display random fluctuations around the average trend, and such fluctuations in the desired pattern also indicate that all the earlier-mentioned models assuming uniform gradation are only approximations of the actual microstructure. In recent years, the investigations of the effective properties for FGMs with stochastic uncertainties in phase volume fractions have attracted more and more attention (Buryachenko & Rammerstorfer 2001; Ferrante & Graham-Brady 2005; Rahman & Chakraborty 2007).

Once the effective properties are determined, the fracture behaviour of FGMs can be analysed. In past decades, there have been many research papers on the crack problems of FGMs with definite effective properties. Delale & Erdogan (1983) and Jin & Noda (1994) pointed out that the singularity of the crack-tip stress fields in FGMs with continuously varying properties has been proved to be the same as that in homogeneous materials. Some basic fracture mechanics concepts in FGMs have been studied by Jin & Batra (1996). In most of the previous analytical fracture models, the material properties were usually assumed to be very particular functions (mainly exponential functions) (Delale & Erdogan 1983; Chen & Erdogan 1996; Erdogan & Wu 1997; Guo *et al.* 2004, 2005; Guo & Noda 2008). Recently, Choi & Paulino (2008) investigated the thermoelastic contact problem of a FGM coating/substrate system, and the thermomechanical properties were assumed to be exponential functions so that the problem could be solved analytically. The problem of two collinear cracks in functionally graded piezoelectric materials under in-plane electromechanical loads was examined by . The dynamic behaviour of a Griffith crack situated at the interface of two bonded, dissimilar, functionally graded piezoelectric materials was considered by Rokne *et al.* (2009). Besides the above work, an analysis of a coupled plane elasticity problem of crack/contact mechanics for a coating/substrate system with functionally graded properties was performed by Choi & Paulino (2010). When the material properties are assumed to be exponential functions, the governing differential equations can usually be reduced to those with constant coefficients and then the exact solutions can be easily found. As for other arbitrary variations of elastic modulus, it is usually very difficult to obtain the exact solutions of the governing differential equations. In order to investigate the crack problem of FGMs with arbitrarily distributed properties, three types of multi-layered models have been proposed. The first model is the homogenous multi-layered model (Itou 2001; ; Wang *et al.* 2002), in which the FGMs are divided into sublayers with piecewise constant material properties. The second model is the linear multi-layered model (Huang *et al.* 2004; Zhong & Cheng 2008), in which the FGMs are divided into sublayers with elastic properties varying linearly in each sublayer. The third model is called the piecewise-exponential (PE) model (Guo & Noda 2007; Guo *et al.* 2008), in which the FGMs are divided into non-homogeneous layers, with each layer's properties varying exponentially; thus, the real material properties can be approached by a series of exponential functions. These methods are very significant for solving the crack problems of FGMs with general (but not stochastic) properties. There are some different characteristics between the above three types of models. In the first model, the material properties are still discontinuous between the sublayers (Guo & Noda 2007). Differently, the second and the third models have introduced the continuity of material properties between sublayers. Therefore, the influences of local material non-homogeneity on the crack-tip fields can be adequately allowed for. All the earlier-mentioned analytical fracture models assume the mechanical properties of the FGMs to be strictly definite. However, as stated already, the phase volume fractions in practical FGMs show random fluctuations around the average trend, which leads to stochastic features of the mechanical properties of FGMs. Thus, the corresponding fracture analyses of FGMs with stochastic mechanical properties will become very complex. As stated by Birman & Byrd (2007), when the probabilistic homogenization approach is implemented, a need will arise in the application of probabilistic mechanics to the analysis of the fracture characteristics of FGM structures. However, due to the complexity of the stochastic problems, there are very few published papers on the fracture problems of FGMs with stochastic properties. For example, by using the finite-element method, the probabilistic characteristics of the stress intensity factors (SIFs) in FGMs with stochastic mechanical properties were studied by Chakraborty & Rahman (2008). Nearly no analytical fracture model has been proposed to solve the crack problems of FGMs with stochastic properties. Actually, when the properties of FGMs are stochastic functions of position, it is very difficult to obtain the analytical solutions of the crack problem. Different from the numerical method of Chakraborty & Rahman (2008), the present study aimed to establish an analytical model—the extended piecewise-exponential (EPE) model for the crack problems of FGMs with stochastic mechanical properties.

The following analytical procedure includes:

— a stochastic description of phase volume fractions and stochastic micromechanics model for effective properties of FGMs;

— an EPE model for the crack problems of FGMs with stochastic mechanical properties;

— examples and discussions; and

— conclusions.

## 2. Stochastic micromechanics-based model for effective properties of the functionally graded materials

### (a) Description of the representative volume elements

A schematic of the microstructure of FGMs is shown in figure 2, in which the total zone is *A*=*A*_{1}∪*A*_{2}∪*A*_{3}. It includes two distinct material phases, phase 1 (white) and phase 2 (black). Either of them represents an isotropic and linear-elastic material.

According to Rahman & Chakraborty (2007) and Birman & Byrd (2007), a representative volume element (RVE) at one point of the FGMs characterizes material heterogeneity in the microscopic length scale. As shown in figure 2*a*, the particle–matrix region *A*_{1} comprises phase 1 (white) and phase 2 (black), and phase 1 is embedded in phase 2. Thus, the particles and matrix are phase 1 (white) and phase 2 (black) materials, respectively. Differently, as shown in figure 2*c*, the particle–matrix roles reverse in region *A*_{3} where the particles and matrix are phase 2 (black) and phase 1 (white) materials, respectively. In the transition region *A*_{2} shown in figure 2*b*, because phase 1 and phase 2 have similar volume fractions, the definition of the particle or matrix is ambiguous. This will be discussed in §2*d*.

Because the volume fractions are stochastic, they must be described by random fields. In the next part, the first objective is to establish a stochastic micromechanics model of the effective properties when the volume fractions *ϕ*_{p}(*x*) and *ϕ*_{m}(*x*) are stochastic.

### (b) Stochastic description of the volume fractions by random fields

The volume fractions of both component phases are denoted by *ϕ*_{p}(*x*) and *ϕ*_{m}(*x*). Either of them is bounded between 0 and 1. The constraint condition is
2.1
where the subscripts ‘p’ and ‘m’ refer to the particle and matrix, respectively.

When the volume fraction *ϕ*_{q}(*x*) (*q*=p or m) is a random variable, let *f*(*t*, *x*) denote the probability density function of *ϕ*_{q}(*x*). Thus, using the knowledge in Kaminski (2005), the mean *μ*_{ϕq}(*x*), the variance [*σ*_{ϕq}(*x*)]^{2} and the distribution function *F*_{ϕq}(*z*,*x*) of *ϕ*_{q}(*x*) can be given by
2.2
2.3
and
2.4
where the symbol *σ*_{ϕq}(*x*) denotes the standard deviation of *ϕ*_{q}(*x*).

### (c) Probabilistic characteristics of effective properties at particle–matrix regions

Considering an RVE at *x*∈*A*_{1}∪*A*_{3} (figure 2*a*,*c*), let *C*_{eff}(*x*) denote the effective elasticity tensor of FGMs. Then, *C*_{eff}(*x*) can be expressed as a function of the component materials' elasticity tensors and volume fractions, i.e.
2.5
where **C**_{m} and **C**_{p} are elasticity tensors of the particle and matrix, respectively. Substituting the constraint condition (2.1) into (2.5), we have
2.6

Let *C*_{ijkl}(*x*) denote the component of the effective elasticity tensor **C**_{eff}(*x*). It should be noted that *C*_{ijkl}(*x*) is only the function of *ϕ*_{p}(*x*) when the particles and matrix materials have been chosen. For convenience, we have rewritten equation (2.6) as
2.7
Considering the fact that the effective properties vary from sample to sample, and applying equation (2.7) to different samples, we have
2.8
Here, the superscript ‘(*J*)’ denotes the *J*th sample. For the *J*th sample of FGMs, and *ϕ*^{(J)}_{p}(*x*) denote the components of the effective elasticity tensor and the volume fraction of the particle, respectively. From equations (2.8), the effective properties of FGMs can be analysed using probability theory. If applying equations (2.8) directly, the sample size should be sufficiently large to obtain the probabilistic characteristics of the effective properties, which may cause very low efficiency in the following crack analysis. Therefore, a more efficient and practical method will be proposed in this part.

Because the phase volume fraction samples in equations (2.8) are subject to the same type of random field, by applying probability theory (Kaminski 2005), the distribution function of *C*_{ijkl}(*x*) can be written as
2.9
Here, the function *F*_{ijkl}(*z*,*x*)∈[0,1] denotes the probability when *C*_{ijkl}(*x*)<*z*. Substituting equation (2.7) into equation (2.9) gives
2.10
When *H*_{ijkl}[*ϕ*_{p}(*x*)] is a general continuous function, there must exist one domain *B*_{ijkl} in which *H*_{ijkl}[*ϕ*_{p}(*x*)]<*z*. Thus, *p*{*H*_{ijkl}[*ϕ*_{p}(*x*)]<*z*} can be given as
2.11

Here, *f*(*t*,*x*) denotes the probability density function of volume fraction *ϕ*_{p}(*x*). According to equations (2.10) and (2.11), the distribution function *F*_{ijkl}(*z*,*x*) of *C*_{ijkl}(*x*) in equation (2.11) can be written as
2.12
Then, the probability density function of *C*_{ijkl}(*x*) can be expressed as
2.13
Thus, similar to the definitions related to the volume fractions in equations (2.2) and (2.3), the mean *μ*_{ijkl}(*x*) and the variance [*σ*_{ijkl}(*x*)]^{2} of *C*_{ijkl}(*x*) can be given as
2.14
and
2.15
Here, the subscript ‘** ijkl**’ means the Einstein summation convention is not applied to ijkl, and the symbol

*σ*

_{ijkl}(

*x*) denotes the standard deviation of

*C*

_{ijkl}(

*x*). For practical FGMs,

*H*

_{ijkl}[

*ϕ*

_{p}(

*x*)] is usually a monotonically continuous function. The following two cases will be involved.

*Case 1.* *H*_{ijkl}[*ϕ*_{p}(*x*)] is a monotonically increasing function. Use of equation (2.11) gives
2.16
where denotes the inverse function of *H*_{ijkl}(*z*) and . Thus, equation (2.10) can be written as
2.17
If we define , the following relation can be obtained:
2.18
Substituting equation (2.18) into equation (2.17) gives
2.19
From equations (2.13) and (2.19), the probability density function of *C*_{ijkl}(*x*) can be obtained as
2.20
*Case 2.* *H*_{ijkl}[*ϕ*_{p}(*x*)] is a monotonically decreasing function. Using a procedure similar to that of case 1, we have
2.21
and
2.22
Thus, through the earlier-mentioned procedure, the probabilistic characteristics of *C*_{ijkl}(*x*) at the particle–matrix regions can be obtained.

### (d) Probabilistic characteristics of effective properties at the transition region

In the transition region *A*_{2}, it is not easy to identify the particle and matrix phases when the volume fraction of each phase is in the vicinity of 0.5. In this case, because the transition region *A*_{2} for a practical FGM lies between *A*_{1} and *A*_{3} and is much smaller than the regions *A*_{1} and *A*_{3}, the probability density function of the transition region can be approximately determined by using interpolation of the probability density functions of the particle–matrix regions. Thus, in the transition region *A*_{2}, we can choose two interpolation functions, *Q*_{1}(*x*) and *Q*_{2}(*x*). Then, the probability density function of *C*_{ijkl}(*x*) in the transition region *A*_{2} can be expressed as
2.23
where , *x*_{1} is the *x*-coordinate of the boundary between *A*_{1} and *A*_{2}, and *x*_{2} is the *x*-coordinate of the boundary between *A*_{2} and *A*_{3}. According to the probability theory, the following fundamental relations can be obtained:
2.24
Applying the above relations to equation (2.23) leads to
2.25
Using the following continuity conditions of the probability density functions:
2.26
and
2.27
and applying equations (2.23) and (2.25)–(2.27), the following relations can be determined:
2.28
2.29

Thus, both interpolation functions *Q*_{1}(*x*) and *Q*_{2}(*x*) can be defined as
2.30
and
2.31
It can be found that *Q*_{1}(*x*) and *Q*_{2}(*x*) satisfy the conditions (2.28) and (2.29). Thus, equation (2.23) can be rewritten as
2.32
Using equation (2.32), the probabilistic characteristics of *C*_{ijkl}(*x*) in the transition region *A*_{2} can be obtained.

## 3. The extended piecewise-exponential model for the crack problems of functionally graded materials with stochastic mechanical properties

### (a) Description of the crack problem

A cracked FGM plate with stochastic mechanical properties varying along the thickness direction will be considered. As shown in figure 3, an internal crack is located perpendicularly to the free surfaces with a crack length 2*c*. The thickness of the plate is *h*. Let *u*(*x*,*y*) and *v*(*x*,*y*) denote the displacement components in the *x*- and *y*-directions, respectively. The governing equations and constitutive equations are
3.1
3.2
Here, *σ*_{ij} (*i*,*j*=*x*,*y*) are stresses, *G*(*x*) is the shear modulus, *κ*(*x*)=3−4*υ*(*x*) (plane strain) and [3−4*υ*(*x*)]/[1+*υ*(*x*)] (plane stress). *υ*(*x*) is Poisson's ratio. The shear modulus and Poisson's ratio are stochastic functions that can be obtained by the method given in §2. The boundary conditions of the FGM plate are given by
3.3
3.4
3.5
3.6
and
3.7

### (b) Combination of the piecewise-exponential model and the stochastic mechanical properties

For each sample of the shear modulus *G*(*x*), the crack problem of the FGM plate in figure 3 can be studied by using the PE model (Guo & Noda 2007). Figure 4 shows the schematic of the layers in the PE model in which the FGM strip is divided into *M* layers, with each layer denoted by a subscript ‘*n*’ (*n*=1,2,…,*M*). For convenience, the layers containing the crack tips are marked by *n*_{1} and *n*_{2}, respectively. The nth layer lies between *x*=*h*_{n−1}and *x*=*h*_{n}. Particularly, *h*_{0}=0.

Now, the real shear modulus of the functionally graded strip with arbitrary properties is defined by
3.8
where ℜ(*x*) may be defined as a general function or tested from experiments. If the shear modulus of each layer in figure 4 is assumed to be
3.9
Applying the real material properties on both surfaces of each layer, we can obtain that
3.10
Solving equation (3.10), we have
3.11
Thus, when *M* is sufficiently large, the material properties of the FGM can be approached by a series of exponential functions.

For each layer, the constitutive equations and governing equations are
3.12
and
3.13
As mentioned already, the varying range of Poisson's ratio is small, and for simplicity, Poisson's ratio of each layer can be assumed to be constant. The boundary conditions are
3.14
3.15
3.16
3.17
and
3.18
The displacement continuity conditions can be written as
3.19
The stress continuity conditions can be expressed as
3.20
In order to obtain the singular integral equation, we introduce the following auxiliary function:
3.21
3.22
and
3.23
Applying the Fourier transform method to the equilibrium equations and using the constitutive equations, the expressions of displacement and stress components can be solved. Then, using the boundary and continuity conditions, the final singular integral equation regarding the unknown auxiliary function can be obtained as
3.24
where *k*_{n1}(*u*,*x*) and *h*_{n2}(*u*,*x*) are known expressions. More detailed information on the PE model can be found in Guo & Noda (2007). By using the method given by Erdogan & Gupta (1972) for solving the singular integral equations, equation (3.24) can be transformed into the following form:
3.25
where
3.26
3.27
3.28
and
3.29
More details can be founded in Erdogan & Gupta (1972). From equations (3.23) and (3.25), the unknowns *w*(*r*_{l}) can be solved. The SIF for an internal crack can be obtained by
3.30
and
3.31
From expressions (3.30) and (3.31), it can be found that the probabilistic characteristics of the SIFs are related to those of the shear modulus, *w*(−1) and *w*(1). It should be noted that *w*(−1) and *w*(1) can be determined from the singular integral equation (3.24) for each sample. Therefore, the probabilistic characteristics of *w*(−1) and *w* (1) are also related to those of the shear modulus.

However, if we calculate the samples one by one with the PE model, the probabilistic characteristics of the SIFs must be obtained by analysing a great number of samples. The aim of this study was to develop an analytical fracture mechanics model but not a simple numerical model for FGMs with stochastic mechanical properties. Therefore, it is not practical or efficient for the analytical model to calculate a great number of samples. Now, the key point lies in how to choose the samples as few as possible and obtain the stable probabilistic characteristic of the SIFs with high efficiency. We note that no matter how to choose the samples of the effective elastic properties, the probability distribution function drawn from the samples of the effective elastic properties should be consistent with equation (2.9). Moreover, it should be noted that the component *C*_{ijkl}(*x*) of the effective elasticity tensor is only related to the effective shear modulus and Poisson's ratio for isotropic FGMs studied in this paper. Because the varying range of Poisson's ratio is usually small, the stochastic characteristic of Poisson's ratio is omitted. Thus, only the probabilistic characteristic of the effective shear modulus needs to be considered. In order to make the probabilistic characteristic of the samples of the chosen effective shear modulus be consistent with equation (2.9) and at the same time, the samples can be chosen efficiently, the following method is proposed.

— From equation (2.9), we can obtain 3.32 Here,

*p*(*p*∈[0,1]) denotes the probability of*C*_{ijkl}(*x*)<*z*, and is the inverse function of the probability density function*F*_{ijkl}(*p*,*x*) of*C*_{ijkl}(*x*). It should be pointed out that the appropriate samples of the effective elastic component*C*_{ijkl}(*x*) can be chosen from appropriate samples*z*_{L}of*z*, while the appropriate samples*z*_{L}can be determined by choosing appropriate samples*p*_{L}of*p*. Thus, if appropriate samples*p*_{L}can be found, the appropriate samples of the effective elastic component*C*_{ijkl}(*x*) can be obtained.— On the basis of the earlier-mentioned analysis, the appropriate samples

*p*_{L}should be determined first. If we define 3.33 where*N*_{s}denotes the sample size (sample number) of the effective elastic component*C*_{ijkl}(*x*), it can be found that more samples can be concentrated in the zone near the peak probability density of*C*_{ijkl}(*x*). In order to state this more clearly, let us consider an example. If*N*_{s}=13, we have*p*_{L}=*L*/13 (*L*=1,2,…,13). The probability density function of the effective elastic component*C*_{ijkl}(*x*) at the position*x*can be obtained by the method in §2 as 3.34 and the schematic of its curve is shown in figure 5. According to equation (3.32), we have Because*p*_{L}is the probability of*z*<*z*_{L},*p*_{L}represents the area surrounded by the horizontal axis, the line of*z*=*z*_{L}and the curve of*g*_{ijkl}(*z*,*x*) in figure 5. It should be noted that the area between*z*_{L−1}and*z*_{L}(*L*=2,…,13) is a constant that can be obtained from equation (3.34) as*p*_{L}−*p*_{L−1}=*L*/*N*_{s}−(*L*−1)/*N*_{s}=1/*N*_{s}(=1/13 for the example in figure 5). This can lead us to the expected result that more samples*z*_{L}can be concentrated in the zone with a relatively larger probability density of the effective parameter*C*_{ijkl}(*x*) in figure 5.— The sample size

*N*_{s}can be determined so that the results of the SIFs exhibit stable probabilistic characteristics.

### (c) Probabilistic characteristics of the stress intensity factors

For each sample of the effective shear modulus *G*(*x*), we can obtain the SIFs using equations (3.30) and (3.31). For the ith sample, we have
3.35
and
3.36
According to probability theory (Kaminski 2005), the mean of the SIFs can be written as
3.37
The probability density function of the SIFs can be determined by developing histograms from the generated samples.

## 4. Examples and discussions

To check the earlier-mentioned method, a two-phase FGM with a crack shown in figure 3 will be considered. The microstructure of FGMs includes two distinct material phases. The macro properties vary along the thickness direction. As an example, *ϕ*_{p}(*x*) is defined as a one-dimensional beta random field (Ferrante & Graham-Brady 2005; Rahman & Chakraborty 2007; Chakraborty & Rahman 2008),
4.1
Here, *θ*_{1} and *θ*_{2} are distribution parameters, and is the gamma function. For practical cases, the mean *μ*_{ϕp}(*x*) and standard deviation *σ*_{ϕp}(*x*) of *ϕ*_{p}(*x*) at any point *x* can be obtained from actual tests. In this paper, in order to check the present method, *μ*_{ϕp}(*x*) and *σ*_{ϕp}(*x*) are directly assumed to be the following forms (Rahman & Chakraborty 2007):
4.2
4.3
Applying equations (2.2) and (2.3), the following expressions can be obtained:
4.4
and
4.5
Applying the Mori–Tanaka model (Mori & Tanaka 1973), the expression of the effective shear modulus of FGMs can be obtained as
4.6
and
4.7
where *G*_{p} and *G*_{m} denote the shear modulus of the particle and matrix constituents, respectively. The symbol *k*_{m} denotes the bulk modulus of the matrix. From equations (2.7) and (4.6), the function *H*_{ijkl}[*ϕ*_{p}(*x*)] for effective shear modulus can be expressed as
4.8
According to equation (2.10), the distribution function *F*(*z*,*x*) of the effective shear modulus gives
4.9
According to the theoretical formulas in §2, the distribution function, the probability density function, the mean *μ*_{G}(*x*) and the variance [*σ*_{G}(*x*)]^{2} of the effective shear modulus can be obtained as
4.10
4.11
4.12
and
4.13
where *σ*_{G}(*x*) denotes the standard deviation of the effective shear modulus.

Next, in order to check the influences of the ratio *G*_{p}/*G*_{m} on the probability density function of the effective shear modulus, the effective shear modulus of two-phase FGMs will be investigated for different modulus ratios *G*_{p}/*G*_{m}. For convenience, the effective shear modulus at the surface *x*=0 is kept constant. Figure 6 shows that the mean and the variance of the effective shear modulus *G*(*x*) increase gradually with the increase of *G*_{p}/*G*_{m}. It should be mentioned that the effective shear modulus *G*(*x*) at one point is still a random variable because the volume fraction *ϕ*_{p}(*x*) is a random variable. For example, when *x*/*h*=0.5, the probability density functions of *G*(*x*) for different values of *G*_{p}/*G*_{m} are shown in figure 7. Thus, the curve of the probability density function of *G*(*x*) at any point can be drawn similarly to those in figure 7 so that the influences of the ratio *G*_{p}/*G*_{m} on the probability density function of *G*(*x*) at any point can be found. The curves in figure 6*a*,*b* show that the ratio *G*_{p}/*G*_{m} may have obvious effects on the probabilistic characteristics of the effective modulus.

Next, the influences of the probabilistic characteristics of the stochastic volume fraction *ϕ*_{p}(*x*) on the probability density function of the effective shear modulus will be checked. For convenience, the ratio *G*_{p}/*G*_{m} is kept constant. Figure 8*a*,*b* displays the influences of the standard deviation *σ*_{ϕp}(*x*) of *ϕ*_{p}(*x*) on *G*(*x*). From figure 8*a*, it can be found that *σ*_{ϕp}(*x*) has very slight effects on the mean of *G*(*x*) at any point *x*. Differently, figure 8*b* shows that *σ*_{ϕp}(*x*) has obvious effects on the variance of *G*(*x*). As an example, the probability density functions of *G*(*x*) at *x*/*h*=0.5 are shown in figure 9 for different standard deviations of the volume fraction. Thus, from the probability density function of *G*(*x*) at one point, the influence of the standard deviation *σ*_{ϕp}(*x*) of *ϕ*_{p}(*x*) on the probabilistic characteristics of *G*(*x*) can also be found.

In the earlier-mentioned analyses, the probabilistic characteristics of the effective shear modulus *G*(*x*) with the variation of *G*_{p}/*G*_{m} and the stochastic distribution of the volume fraction *ϕ*_{p}(*x*) have been discussed. In the following part, the probabilistic characteristics of the SIFs will be analysed. Applying equations (3.32) and (3.33), we can obtain the distribution curves of the effective shear modulus *G*(*x*) from different samples. The schematics of these curves are shown in figure 10. For each curve of the shear modulus *G*(*x*) in figure 10, we can obtain the SIFs by using the earlier-mentioned PE model. As shown in figure 11*a*, for each crack length, it can be found that when the standard deviation of the volume fraction *ϕ*_{p}(*x*) and the ratio *G*_{p}/*G*_{m} is small (e.g. *G*_{p}/*G*_{m}=3), the SIFs for different samples vary very slightly. It is due to the fact that the standard deviation of *G*(*x*) is small when the standard deviation of the volume fraction *ϕ*_{p}(*x*) and the ratio *G*_{p}/*G*_{m} are small. Thus, the difference between the SIFs corresponding to different samples is very small so that the stochastic characteristics of the material properties nearly has no effect on the stress fields near the crack tips. Differently, the standard deviation of *G*(*x*) is large when the standard deviation of the volume fraction *ϕ*_{p}(*x*) and the ratio *G*_{p}/*G*_{m} are large. At this time, the stochastic characteristics of the material properties have an obvious effect on the stress fields near the crack tips. This case is shown in figure 11*b*.

Finally, the whole probabilistic characteristics of the SIFs will be discussed. In the following analysis, a central crack subjected to a constant compressing loading −*σ*_{0} at the crack face is considered. Poisson's ratio is assumed to be 0.3. The SIFs are normalized by In order to obtain the stable probabilistic characteristics of the SIFs, the sample size must be large enough to guarantee the convergence of the probability density function of the SIFs. In §3*b*, it has been stated how to choose samples to minimize the sample size. According to this method, a comparison of the probability density functions of the SIFs for different sample sizes (100, 125, 150 and 200) is shown in figure 12. It can be observed that the probability density functions of the SIFs for different sample sizes agree well. Therefore, when the sample size=100, it is large enough to obtain stable probabilistic characteristics of SIFs. Thus, with the sample size=100, the probability density function of the normalized SIFs for a crack with *c*/*h*=0.2 and *c*/*h*=0.25 are shown in figure 13*a*,*b*, respectively. It should be noted that with the earlier-mentioned method, a large number for sample size is not needed so that the present analytical model, the EPE model, can be applied to solve the crack problem of FGMs with stochastic mechanical properties efficiently.

## 5. Conclusion

A stochastic fracture mechanics model is proposed for predicting probabilistic characteristics of SIF for a crack problem of FGMs. The studied FGMs are subjected to statistical uncertainties in phase volume fractions. The original points of the present study include the following.

This study aimed to develop a method to build a ‘bridge’ between the macro fracture mechanics model and stochastic micromechanics-based properties so that the macro fracture mechanics model can be expanded to the fracture mechanics problem of FGMs with stochastic mechanical properties.

In the past decades, only a few analytical models have been developed to deal with the crack problems of FGMs with definite general mechanical properties (but not stochastic). For FGMs with stochastic properties stated in §2, nearly no work has been carried out on the analytical models of crack problems. This paper developed the work of Guo & Noda (2007) to deal with the crack problems of FGMs with stochastic mechanical properties. It should be emphasized that the present model is not a simple extension of the model of Guo & Noda (2007).

The new contributions of this paper include the following.

The original PE model cannot deal with the stochastic properties directly. In order to extend the PE model to the crack problems of FGMs with stochastic properties, one feasible form of the effective stochastic mechanical properties must be derived with consideration of the stochastic distributions of component materials. It should be emphasized that the present study aimed to develop an analytical model but not a numerical model (e.g. finite-element model) for the crack problems of FGMs with stochastic properties so that the probabilistic characteristic of the SIFs can be calculated efficiently with a small sample size. Therefore, different from Guo & Noda (2007), the present method should be able to not only deal with stochastic mechanical properties, but also obtain the stable probabilistic characteristics of the SIFs with high efficiency.

To realize the aim in paragraph (1) above, a micromechanics-based method is developed to obtain the explicit probabilistic characteristics of the effective properties of the FGMs so that the stochastic mechanical properties can be applied in the analytical model to obtain the SIFs. Although only the Mori–Tanaka model and beta distribution are used as examples in this study, other micromechanics models and stochastic distributions can be applied with this thought.

How to choose an appropriate sample size is the most common problem when studying the probability density function of SIFs. If the sample size is very large, the analytical model may lose its advantages. In order to apply the above analytical fracture model with high efficiency, a thought for choosing samples efficiently is proposed. With this thought, the stable probabilistic characteristic of SIFs can be obtained with very small sample sizes so that the present model can be applied with high efficiency.

## Acknowledgements

This work is sponsored by NSFC (11072067), NCET (08-0151), SRF for ROCS (SEM), Heilongjiang Science Fund for ROCS, Doctoral Programmes Foundation of Ministry of Education of China and Science Funds for Distinguished Young Scholar of Heilongjiang Province, China.

- Received March 11, 2012.
- Accepted April 16, 2012.

- This journal is © 2012 The Royal Society