The interface crack problem for a functionally graded coating-substrate structure with general coating properties

The interface crack problem for a functionally graded coating-substrate structure with general coating properties

Accepted Manuscript The interface crack problem for a functionally graded coating-substrate structure with general coating properties Licheng Guo , S...

3MB Sizes 0 Downloads 58 Views

Accepted Manuscript

The interface crack problem for a functionally graded coating-substrate structure with general coating properties Licheng Guo , Shanqiao Huang , Li Zhang , Pengfei Jia PII: DOI: Reference:

S0020-7683(18)30133-1 10.1016/j.ijsolstr.2018.03.025 SAS 9950

To appear in:

International Journal of Solids and Structures

Received date: Revised date: Accepted date:

3 January 2018 23 March 2018 25 March 2018

Please cite this article as: Licheng Guo , Shanqiao Huang , Li Zhang , Pengfei Jia , The interface crack problem for a functionally graded coating-substrate structure with general coating properties, International Journal of Solids and Structures (2018), doi: 10.1016/j.ijsolstr.2018.03.025

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

ACCEPTED MANUSCRIPT

The interface crack problem for a functionally graded coating-substrate structure with general coating properties Licheng Guoa*, Shanqiao Huanga, Li Zhanga, Pengfei Jiaa a

CR IP T

Department of Astronautic Science and Mechanics, Harbin Institute of Technology, Harbin 150001, China ∗ Corresponding author. Tel.: +86 451 86403725; fax: +86 451 86403725. E-mail address: [email protected] (L.C. Guo) Abstract A piecewise-exponential model (PE model) is developed for the interface crack problem of a functionally graded coating-substrate structure in which the coating

AN US

mechanical properties are generally continuous while the substrate is homogeneous. By this way, the coating properties can be approached by a series of exponential functions without losing the continuity. Through a series of mathematical manipulations, the problem is reduced into a group of singular integral equations that

M

can be solved by numerical methods. After analyzing the mix-mode stress intensity factors (SIFs) and the strain energy release rates (SERRs) of some typical examples

ED

under plane loadings, the influences of the geometric parameters and coating properties variations on the fracture behaviors of the interface crack are presented.

PT

Those analyses can greatly benefit the designs and the manufactures of FGM coating-substrate structures.

CE

Keyword: Functionally graded coating-substrate structure; General mechanical properties; Piecewise-exponential model; SIFs and SERRs

AC

1. Introduction

Due to the continuous variation of material properties, the functionally graded

materials (FGMs) avoid the mismatch between different materials shown in traditional composite materials successfully and improve the strength and toughness. Nowadays the application of FGMs has been greatly expanded and the research on the fracture behaviors of FGMs has become an important aspect of FGMs’ performance assessing and designation. As FGMs are usually used in the form of coating covering on the homogeneous substrate, the cracks are likely to initiate at the interface between 1

ACCEPTED MANUSCRIPT

the coating and the substrate due to the manufacture processing or external loadings. So the fracture behaviors of this type of cracks possess significant research value. The crack problems of FGMs have attracted the attentions of many researchers. Erdogan and his coauthors made great contribution to this field in the past years using analytical methods. Delale and Erdogan (1983) investigated the internal crack problem for an infinite FGM plate with Young’s modulus that varied exponentially

CR IP T

alone the crack direction. This paper revealed that the influence of the variation of Poisson’s ratio on stress intensity factors (SIFs) is rather insignificant. Erdogan (1985) analyzed the anti-plane fracture problem in bonded nonhomogeneous half planes where the derivative of shear modulus was discontinuous. He found that the

AN US

discontinuity has little effect on the singularity of stress field at the crack tip. The edge crack problem of a FGM half plane under general loading was studied by Dag and Erdogan (2002). It is worth noting that the fracture behaviors of surface crack under sliding contact loading, conducted by a rigid circular stamp, were analyzed in

M

this work. Furthermore, Shbeeb et al. (1999a, b) conducted investigation on the multiple cracks problem in an infinite nonhomogeneous plate analytically. Since

ED

FGMs were used as thermal insulation materials initially, the investigations on thermal fracture problems are also very important, on which a series of analytical

PT

analyses have been finished by Noda and Jin (1993; 1994a; 1994b). The functionally graded piezoelectric material (FGPM) is also an important field of the applications of

CE

FGMs. Yan and Jiang (2009a; 2009b) conducted analytical investigations on the crack propagations in different FGPMs structures under electromechanical loadings. On the

AC

other hand, numerical methods are also important approaches for the crack problems of FGMs. Kim and Paulino (2002) used three finite element methods to calculate the mix-mode SIFs of two-dimensional fracture problems of FGMs. Gao et al. (2008) developed a boundary element method to investigate mode-I edge crack problem of FGMs. Yildirim et al. (2005) conducted three-dimensional analyses on the semi-elliptical crack problems in a FGM coating bonded to a homogeneous substrate under mode-I mechanical and transient thermal loadings. In order to simplify the calculation of the governing partially differential 2

ACCEPTED MANUSCRIPT

equations, the mechanical properties of FGMs are usually defined by particular functions, mainly exponential functions, in the analytical analyses of fracture problems. Obviously, such a practice is unable to deal with the general cases in engineering applications. So several representative models have been proposed to approximate the mechanical properties of FGMs, such as the homogeneous multi-layered model (HM model) (Wang et al., 2002) and the linear multi-layered

CR IP T

model (LM model) (Huang et al., 2004). However, the HM model loses the continuity of mechanical properties, and the LM model makes the governing partially differential equations difficult to solve for arbitrarily oriented crack problems. Under such cases, the piecewise-exponential model (PE model), proposed initially by Guo

AN US

and Noda (2007), possesses high computing efficiency and wider scope of application. Approaching the material properties by a series of exponential functions, both the continuity of material properties and the convenience of solving the governing partially differential equations are taken into consideration. Therefore, the PE model

M

was adopted to solve many fracture problems of FGMs. It was used further to investigate the fracture behaviors of a vertical crack in a FGM plate with arbitrary

ED

thermomechanical properties under steady thermal loading (Guo et al., 2008) and thermal shock (Guo and Noda, 2014), respectively. Bai et al. (2013) involved Laplace

PT

integral transform into the PE model and expand the PE model to study the transient fracture problem of a FGM plate. Guo et al. (2012) came up with an analytical model

CE

for the mode-I crack problems of FGMs with stochastic mechanical properties, which were described by the stochastic of the phase volume fractions, on the base of the PE

AC

model. In their work, a new method of choosing samples was proposed, which can help obtain the stable probabilistic characteristic of the SIFs efficiently. It is obvious that the PE model is a powerful approach to solve the fracture problems of FGMs under various cases. However, the PE model has not been expanded to the interface crack problems of FGM coating-substrate structures. Many scholars were also engaged in the studies on the interface crack problems of graded coating-homogeneous substrate structure, since this kind of structure is a main form of the applications of FGMs. Chen and Erdogan (1996) studied the 3

ACCEPTED MANUSCRIPT

interface crack problem of a FGM coating bonded to a homogeneous substrate, where the shear modulus of the coating varied exponentially. The influences of several parameters, such as Poisson’s ratio, the nonhomogeneity constants and the geometric parameters, on the fracture behaviors of interface crack were investigated. Then the interface crack problem of the same structure under concentrated loading was investigated (Guo et al., 2004). The studies on the anti-plane interface crack problem

CR IP T

under static and dynamic loadings were performed by Li et al. (2016). Furthermore, the interface crack problem of another coating-substrate structure was investigated (Shbeeb and Binienda, 1999), where the top and bottom surfaces of a FGM interphase were bonded to a homogeneous coating and substrate respectively and the crack was

AN US

along the interface between the FGM layer and substrate. Choi and Paulino (2010) analyzed a coupled problem of fracture/contact mechanics of the same structure. In their work, the fracture behaviors of the interface crack under the sliding contact loading imposed by a rigid flat punch were investigated. However, in most analytical

M

researches on interface crack problems of nonhomogeneous coating-substrate structure the material properties of FGM layers were assumed to be very particular

ED

functions, mainly exponential functions. This status requires a new model able to deal with the FGM coating with general mechanical properties so that the studies on the

PT

interface crack problems of FGM coating-substrate structures can be greatly closer to real applications.

CE

In this paper, a PE model is developed for the interface crack problem of a functionally graded coating-substrate structure in which the coating mechanical

AC

properties are generally continuous while the substrate is homogeneous. After a series of mathematical manipulations, the fracture problem is transformed into a group of singular integral equations, which can be solved by the Lobatto-Chebyshev method (Binienda and Arnold, 1995; Theocaris and Ioakimidis, 1977). Then the mix-mode stress intensity factors (SIFs) and strain energy release rates (SERRs) are obtained, and the influences of variation of FGM mechanical properties and geometric parameters of the structure on the fracture behaviors are presented. 4

ACCEPTED MANUSCRIPT

CR IP T

2. Problem formulation

AN US

Fig. 1 2D Geometry of the crack between a FGM coating with general properties and a homogeneous substrate.

The interface crack problem for the coating-substrate structure is illustrated graphically in Fig. 1. The FGM coating with general mechanical properties and the homogeneous substrate are bonded together while a crack is along the interface

M

between the coating and substrate. The interface crack problem is analyzed as a plane

ED

elastic fracture problem. The coating-substrate structure is considered as an isotropic strip with limited thickness along y-axis and infinite length along x-axis. The thicknesses of the coating and the substrate are 𝑕𝑀 and 𝐻0 respectively. As the

PT

displacement components in x-direction and y-direction are expressed as 𝑢(𝑥, 𝑦) and

CE

𝑣(𝑥, 𝑦) respectively, the constitutive relations of the strip can be written as 𝜇(𝑥,𝑦)

𝜕𝑢(𝑥,𝑦)

𝜇(𝑥,𝑦)

𝜕𝑥 𝜕𝑣(𝑥,𝑦)

AC

𝜎𝑥𝑥 (𝑥, 𝑦) = 𝑘(𝑥,𝑦)−1 *,𝑘(𝑥, 𝑦) + 1-

𝜎𝑦𝑦 (𝑥, 𝑦) = 𝑘(𝑥,𝑦)−1 *,𝑘(𝑥, 𝑦) + 1-

{

𝜕𝑦 𝜕𝑢(𝑥,𝑦)

𝜎𝑥𝑦 (𝑥, 𝑦) = 𝜇(𝑥, 𝑦),

𝜕𝑦

+ ,3 − 𝑘(𝑥, 𝑦)+ ,3 − 𝑘(𝑥, 𝑦)+

𝜕𝑣(𝑥,𝑦)

+

𝜕𝑦 𝜕𝑢(𝑥,𝑦) 𝜕𝑥

+

(1)

𝜕𝑣(𝑥,𝑦) 𝜕𝑥

-

where 𝜇(𝑥, 𝑦) is the shear modulus, 𝑘(𝑥, 𝑦) = ,3 − 𝜈(𝑥, 𝑦)-/,1 + 𝜈(𝑥, 𝑦)- for plane stress case , 𝑘(𝑥, 𝑦) = 3 − 4𝜈(𝑥, 𝑦) for plane strain case and 𝜈(𝑥, 𝑦) is the Poisson’s ratio. The equilibrium equations are 𝜕𝜎𝑥𝑥 𝜕𝑥 {𝜕𝜎 𝑦𝑦 𝜕𝑦

+ +

5

𝜕𝜎𝑥𝑦 𝜕𝑦 𝜕𝜎𝑥𝑦 𝜕𝑥

=0 =0

(2)

ACCEPTED MANUSCRIPT For convenience,𝜎 𝐶 and 𝜎 𝑆 are set to the stress fields of the coating and the substrate respectively. Then the boundary conditions of the coating-substrate structure can be expressed as (3)

{

𝑆 (𝑥, 𝜎𝑦𝑦 −𝐻0 ) = 0 𝑆 (𝑥, 𝜎𝑥𝑦 −𝐻0 ) = 0

(4)

and the continuity conditions on the interface are

𝐶 (𝑥, + ) 𝑆 (𝑥, − ) 𝜎𝑦𝑦 0 = 𝜎𝑦𝑦 0 { 𝐶 𝑆 (𝑥, − ) − ∞ < 𝑥 < +∞ 𝜎𝑥𝑦 (𝑥, 0+ ) = 𝜎𝑥𝑦 0

(5a)

𝐶 (𝑥, + ) 𝑆 (𝑥, − ) 𝜎𝑦𝑦 0 = 𝜎𝑦𝑦 0 = 𝜎(𝑥) = 𝑝1 (𝑥) −𝑎 <𝑥 <𝑎 𝐶 (𝑥, + ) 𝑆 (𝑥, − ) 𝜎𝑥𝑦 0 = 𝜎𝑥𝑦 0 = 𝜏(𝑥) = 𝑝2 (𝑥)

(5b)

AN US

{

CR IP T

𝐶 (𝑥, 𝜎𝑦𝑦 𝑕𝑀 ) = 0 𝐶 (𝑥, 𝜎𝑥𝑦 𝑕𝑀 ) = 0

{

where 𝑝1 (𝑥) and 𝑝2 (𝑥) are known functions.

As the material properties of the FGM coating are affected by many factors, like

M

material design and manufacture process, the variation of the material properties would be complicated and even arbitrary. That means the analytical solution of Eq. (2)

ED

will be very difficult to obtain. In order to solve the problem, the piecewise-exponential model (PE model) is proposed to depict the material properties

PT

of FGM coating. The FGM coating is divided into M layers along the direction of thickness as shown in Fig. 2a. 𝑕𝑛 is the y-value of the interface between n-th layer

CE

and n+1-th layer (𝑛 = 1, 2, … , 𝑀 − 1), while 𝑕0 = 0. The thickness of the n-th layer is 𝑕𝑛 − 𝑕𝑛−1 and the material properties vary exponentially in each layer.

AC

Furthermore, the material properties at the interfaces between adjacent layers are equal to the real ones. By using a series of exponential functions, the material properties of the FGM coating can be approximated accurately. Fig. 2b gives the comparison between the real property and the approximate property.

6

CR IP T

ACCEPTED MANUSCRIPT

AC

CE

PT

ED

M

AN US

(a)

(b)

Fig. 2 The schematics of the layers in the PE Model in the FGM coating. (a) The

layered FGM coating; (b) the variation of the actual property and the approximate property. According to Eq. (1) and the PE model, the shear modulus of the FGM coating

should be described by a series exponential functions as 𝜇𝑛 (𝑦) = 𝜇0𝑛 𝑒 𝜃𝑛 𝑦 7

𝑛 = 1,2, … , 𝑀

(6)

ACCEPTED MANUSCRIPT

where 𝜇0𝑛 and 𝜃𝑛 are the nonhomogeneous parameters of the n-th layer. The real shear modulus of the FGM coating is expressed as 𝜇(𝑦) = 𝑞(𝑦)

(7)

where 𝑞(𝑦) is a continuous function. As the shear moduli at the interfaces in the PE model are equal to the real values, we have 𝜇𝑛 (𝑕𝑛−1 ) = 𝑞(𝑕𝑛−1 ) 𝜇𝑛 (𝑕𝑛 ) = 𝑞(𝑕𝑛 )

𝑛 = 1, 2, … , 𝑀

(8)

CR IP T

{

So the nonhomogeneous parameters 𝜇0𝑛 and 𝜃𝑛 can be solved from Eq. (8) {

𝜃𝑛 =

ln(𝑞(𝑕𝑛 ))−ln(𝑞(𝑕𝑛−1 ))

𝜇0𝑛 =

𝑕𝑛 −𝑕𝑛−1 𝑞(𝑕𝑛 )𝑒 −𝜃𝑛𝑕𝑛

𝑛 = 1, 2, … , 𝑀

(9)

For convenience, the subscript 𝑛 = 0 is used to represent the substrate. As the

AN US

shear modulus is continuous along the thickness and the substrate is homogeneous, the shear modulus of the substrate is

𝜇0 (𝑦) = 𝜇1 (0) = 𝜇01

(10)

In the PE model, the boundary and continuity conditions Eqs. (3-5) can be

M

written respectively as

𝜎𝑀𝑦𝑦 (𝑥, 𝑕𝑀 ) = 0 𝜎𝑀𝑥𝑦 (𝑥, 𝑕𝑀 ) = 0

(11)

𝜎0𝑦𝑦 (𝑥, −𝐻0 ) = 0 𝜎0𝑥𝑦 (𝑥, −𝐻0 ) = 0

(12)

ED

{

CE

PT

{

{

𝜎0𝑦𝑦 (𝑥, 0− ) = 𝜎1𝑦𝑦 (𝑥, 0+ ) − ∞ < 𝑥 < +∞ 𝜎0𝑥𝑦 (𝑥, 0− ) = 𝜎1𝑥𝑦 (𝑥, 0+ )

AC

𝜎0𝑦𝑦 (𝑥, 0− ) = 𝜎1𝑦𝑦 (𝑥, 0+ ) = 𝜎(𝑥) = 𝑝1 (𝑥) { −𝑎 <𝑥 <𝑎 𝜎0𝑥𝑦 (𝑥, 0− ) = 𝜎1𝑥𝑦 (𝑥, 0+ ) = 𝜏(𝑥) = 𝑝2 (𝑥)

(13a)

(13b)

According to Eq. (1) and Eq. (6), the constitutive relations of each coating layer

and the substrate can be expressed as 𝜇𝑛 (𝑦)

𝜎𝑛𝑥𝑥 (𝑥, 𝑦) = 𝑘

𝑛 (𝑥,𝑦)−1

*,𝑘𝑛 (𝑥, 𝑦) + 1-

𝜇𝑛 (𝑦) *,𝑘𝑛 (𝑥, 𝑦) 𝑛 (𝑥,𝑦)−1

𝜎𝑛𝑦𝑦 (𝑥, 𝑦) = 𝑘 {

+ 1-

𝜎𝑛𝑥𝑦 (𝑥, 𝑦) = 𝜇𝑛 (𝑦),

𝜕𝑢𝑛 (𝑥,𝑦) 𝜕𝑥 𝜕𝑣𝑛 (𝑥,𝑦)

𝜕𝑦 𝜕𝑢𝑛 (𝑥,𝑦) 𝜕𝑦

+ ,3 − 𝑘𝑛 (𝑥, 𝑦)+ ,3 − 𝑘𝑛 (𝑥, 𝑦)-

+

𝜕𝑣𝑛 (𝑥,𝑦)

+

𝜕𝑦 𝜕𝑢𝑛 (𝑥,𝑦) 𝜕𝑥

+ 𝑛 = 0, 1, … , 𝑀 (14)

𝜕𝑣(𝑥,𝑦) 𝜕𝑥

-

According to the previous studies (Chen and Erdogan, 1996; Guo and Noda, 8

ACCEPTED MANUSCRIPT

2007), the variation range of the Poisson’s ratio is pretty small in the FGM coating and the variation of the Poisson’s ratio plays little role in the final calculation results. So that it is appropriate to regard the Poisson’s ratio of the FGM coating as a constant. Then we have 𝑘𝑛 (𝑥, 𝑦) = 𝑘 = (3 − 𝜈)/(1 + 𝜈) (for plane stress problem) and 𝑘𝑛 (𝑥, 𝑦) = 𝑘 = 3 − 4𝜈 (for plane strain problem).

𝜕𝜎𝑛𝑥𝑥 𝜕𝑥 {𝜕𝜎 𝑛𝑦𝑦 𝜕𝑦

+ +

𝜕𝜎𝑛𝑥𝑦

CR IP T

The equilibrium equations of each coating layer and the substrate are =0

𝜕𝑦 𝜕𝜎𝑛𝑥𝑦

𝑛 = 0, 1, 2, … , 𝑀

=0

𝜕𝑥

(15)

Substituting the Eq. (14) into the Eq. (15) and applying Fourier transform to x-variable, the ordinary differential equations can be obtained.

𝑎𝑛21 𝑈𝑛 (𝑠, 𝑦) + 𝑎𝑛22 𝑉𝑛 (𝑠, 𝑦) + 𝑎𝑛23

𝜕𝑈𝑛 (𝑠,𝑦) 𝜕𝑦 𝜕𝑈𝑛 (𝑠,𝑦) 𝜕𝑦

𝜕𝑉𝑛 (𝑠,𝑦)

𝜕2 𝑈𝑛 (𝑠,𝑦)

AN US

{

𝑎𝑛11 𝑈𝑛 (𝑠, 𝑦) + 𝑎𝑛12 𝑉𝑛 (𝑠, 𝑦) + 𝑎𝑛13

+ 𝑎𝑛14

+ 𝑎𝑛24

𝜕𝑦 𝜕𝑉𝑛 (𝑠,𝑦) 𝜕𝑦

+ 𝑎𝑛15

+ 𝑎𝑛25

𝜕𝑦 2 𝜕2 𝑉𝑛 (𝑠,𝑦) 𝜕𝑦 2

=0

=0

𝑛 = 0, 1, 2, … , 𝑀(16)

𝑈𝑛 (𝑠, 𝑦) and 𝑉𝑛 (𝑠, 𝑦) are the Fourier transformed from 𝑢𝑛 (𝑥, 𝑦) and 𝑣𝑛 (𝑥, 𝑦) respectively

and

s

is

Fourier

space

variable.

The

coefficients

𝑎𝑛𝑖𝑗

M

(𝑛 = 0, 1, 2, … , 𝑀; 𝑖 = 1, 2; 𝑗 = 1, 2, … , 5) are given in Appendix.

ED

The characteristic roots can be obtained after solving the ordinary differential equations.

(𝑘−3)𝑠2 𝜃𝑛 2

1

PT

𝜆𝑛1,𝑛3 (𝑠) = 2 ,−𝜃𝑛 ∓ √4𝑠 2 + 𝜃𝑛 2 − 4√

1+𝑘 2

𝑛 = 1, 2, … , 𝑀

2

1 (𝑘−3)𝑠 𝜃 𝜆𝑛2,𝑛4 (𝑠) = 2 ,−𝜃𝑛 ∓ √4𝑠 2 + 𝜃𝑛 2 + 4√ 1+𝑘 𝑛 { 𝜆01,02 (𝑠) = −|𝑠| { 𝜆03,04 (𝑠) = |𝑠|

(18)

AC

CE

(17)

Then the expressions of 𝑈𝑛 and 𝑉𝑛 can be obtained. After applying Fourier inverse transform to 𝑈𝑛 (𝑠, 𝑦) and 𝑉𝑛 (𝑠, 𝑦), the displacement and stress fields of the FGM coating can be expressed as 1

+∞

𝑢𝑛 (𝑥, 𝑦) = 2𝜋 ∫−∞ ∑4𝑗=1 𝐸𝑛𝑗 (𝑠)𝐴𝑛𝑗 (𝑠)𝑒 𝜆𝑛𝑗(𝑠)𝑦−i𝑠𝑥 𝑑𝑠 { 𝑛 = 1, 2, … , 𝑀 (19) 1 +∞ 𝑣𝑛 (𝑥, 𝑦) = 2𝜋 ∫−∞ ∑4𝑗=1 𝐹𝑛𝑗 (𝑠)𝐴𝑛𝑗 (𝑠)𝑒 𝜆𝑛𝑗(𝑠)𝑦−i𝑠𝑥 𝑑𝑠

9

ACCEPTED MANUSCRIPT

1

+∞

1

+∞

1

+∞

𝜎𝑛𝑥𝑥 (𝑥, 𝑦) = 2𝜋 𝑒 𝜃𝑛𝑦 ∫−∞ ∑4𝑗=1 𝐵𝑛𝑗 (𝑠)𝐴𝑛𝑗 (𝑠)𝑒 𝜆𝑛𝑗 (𝑠)𝑦−i𝑠𝑥 𝑑𝑠 𝜎𝑛𝑦𝑦 (𝑥, 𝑦) = 2𝜋 𝑒 𝜃𝑛 𝑦 ∫−∞ ∑4𝑗=1 𝐶𝑛𝑗 (𝑠)𝐴𝑛𝑗 (𝑠)𝑒 𝜆𝑛𝑗(𝑠)𝑦−i𝑠𝑥 𝑑𝑠 𝑛 = 1, 2, … , 𝑀(20) 𝜎 (𝑥, 𝑦) = 2𝜋 𝑒 𝜃𝑛𝑦 ∫−∞ ∑4𝑗=1 𝐷𝑛𝑗 (𝑠)𝐴𝑛𝑗 (𝑠)𝑒 𝜆𝑛𝑗 (𝑠)𝑦−i𝑠𝑥 𝑑𝑠 { 𝑛𝑥𝑦 According to Eq. (18), the characteristic roots of the homogeneous substrate are multiple roots. So the displacement and stress fields of the substrate are +∞

1

1

+∞

1

+∞

(21)

CR IP T

𝑢0 (𝑥, 𝑦) = 2𝜋 ∫−∞ ∑4𝑗=1 𝐸0𝑗 (𝑦, 𝑠)𝐴0𝑗 (𝑠)𝑒 𝜆0𝑗(𝑠)𝑦−i𝑠𝑥 𝑑𝑠 { 1 +∞ 𝑣0 (𝑥, 𝑦) = 2𝜋 ∫−∞ ∑4𝑗=1 𝐹0𝑗 (𝑦, 𝑠)𝐴0𝑗 (𝑠)𝑒 𝜆0𝑗(𝑠)𝑦−i𝑠𝑥 𝑑𝑠

𝜎0𝑥𝑥 (𝑥, 𝑦) = 2𝜋 ∫−∞ ∑4𝑗=1 𝐵0𝑗 (𝑦, 𝑠)𝐴0𝑗 (𝑠)𝑒 𝜆0𝑗(𝑠)𝑦−i𝑠𝑥 𝑑𝑠 𝜎0𝑦𝑦 (𝑥, 𝑦) = 2𝜋 ∫−∞ ∑4𝑗=1 𝐶0𝑗 (𝑦, 𝑠)𝐴0𝑗 (𝑠)𝑒 𝜆0𝑗(𝑠)𝑦−i𝑠𝑥 𝑑𝑠 𝜎 (𝑥, 𝑦) = { 0𝑥𝑦

+∞ 4 ∑ 𝐷 (𝑦, 𝑠)𝐴0𝑗 (𝑠)𝑒 𝜆0𝑗(𝑠)𝑦−i𝑠𝑥 ∫ 2𝜋 −∞ 𝑗=1 0𝑗 1

(22)

𝑑𝑠

AN US

In Eqs. (19-22), 𝐸𝑛𝑗 , 𝐹𝑛𝑗 , 𝐵𝑛𝑗 , 𝐶𝑛𝑗 and 𝐷𝑛𝑗 are known functions given in Appendix and 𝐴𝑛𝑗 remain unsolved (𝑛 = 0, 1, 2, … , 𝑀; 𝑗 = 1, … , 4). According to the continuity conditions of stress and displacement fields on the interfaces of adjacent layers in the FGM coating, we have

(23)

ED

M

𝑢𝑛 (𝑥, 𝑕𝑛 ) = 𝑢𝑛+1 (𝑥, 𝑕𝑛 ) 𝑣𝑛 (𝑥, 𝑕𝑛 ) = 𝑣𝑛+1 (𝑥, 𝑕𝑛 ) 𝜎𝑛𝑦𝑦 (𝑥, 𝑕𝑛 ) = 𝜎(𝑛+1)𝑦𝑦 (𝑥, 𝑕𝑛 ) 𝑛 = 1, 2, … , 𝑀 − 1 {𝜎𝑛𝑥𝑦 (𝑥, 𝑕𝑛 ) = 𝜎(𝑛+1)𝑥𝑦 (𝑥, 𝑕𝑛 )

PT

Substituting Eqs. (19-20) into the continuity condition Eq. (23) and then applying Fourier transform to the equations, we have

CE

∑4𝑗=1 𝐸𝑛𝑗 (𝑠)𝐴𝑛𝑗 (𝑠)𝑒 𝜆𝑛𝑗(𝑠)𝑕𝑛 = ∑4𝑗=1 𝐸𝑛+1𝑗 (𝑠)𝐴𝑛+1𝑗 (𝑠)𝑒 𝜆𝑛+1𝑗(𝑠)𝑕𝑛 ∑4𝑗=1 𝐹𝑛𝑗 (𝑠)𝐴𝑛𝑗 (𝑠)𝑒 𝜆𝑛𝑗(𝑠)𝑕𝑛 = ∑4𝑗=1 𝐹𝑛+1𝑗 (𝑠)𝐴𝑛+1𝑗 (𝑠)𝑒 𝜆𝑛+1𝑗(𝑠)𝑕𝑛

∑4𝑗=1 𝐶𝑛𝑗 (𝑠)𝐴𝑛𝑗 (𝑠)𝑒 𝜆𝑛𝑗(𝑠)𝑕𝑛 = ∑4𝑗=1 𝐶𝑛+1𝑗 (𝑠)𝐴𝑛+1𝑗 (𝑠)𝑒 𝜆𝑛+1𝑗(𝑠)𝑕𝑛

𝑛 = 1, 2, … , 𝑀 − 1 (24)

AC

4 4 𝜆 (𝑠)𝑕 𝜆 (𝑠)𝑕 {∑𝑗=1 𝐷𝑛𝑗 (𝑠)𝐴𝑛𝑗 (𝑠)𝑒 𝑛𝑗 𝑛 = ∑𝑗=1 𝐷𝑛+1𝑗 (𝑠)𝐴𝑛+1𝑗 (𝑠)𝑒 𝑛+1𝑗 𝑛

Substituting Eqs. (20, 22) into the boundary conditions Eqs. (11, 12) and then applying Fourier transform to the equations, we have {

∑4𝑗=1 𝐶0𝑗 (−𝐻0 , 𝑠)𝐴0𝑗 (𝑠)𝑒 −𝜆0𝑗(𝑠)𝐻0 = 0 ∑4𝑗=1 𝐷0𝑗 (−𝐻0 , 𝑠)𝐴0𝑗 (𝑠)𝑒 −𝜆0𝑗(𝑠)𝐻0 = 0 {

∑4𝑗=1 𝐶𝑀𝑗 (𝑠)𝐴𝑀𝑗 (𝑠)𝑒 𝜆𝑀𝑗(𝑠)𝑕𝑀 = 0 ∑4𝑗=1 𝐷𝑀𝑗 (𝑠)𝐴𝑀𝑗 (𝑠)𝑒 𝜆𝑀𝑗 (𝑠)𝑕𝑀 = 0 10

(25)

(26)

ACCEPTED MANUSCRIPT

Similarly, the following equations can be obtained from the stress continuity conditions of Eq. (13a). {



+



+

∑4𝑗=1 𝐶0𝑗 (0− , 𝑠)𝐴0𝑗 (𝑠)𝑒 𝜆0𝑗(𝑠)0 = ∑4𝑗=1 𝐶1𝑗 (𝑠)𝐴1𝑗 (𝑠)𝑒 𝜆𝑀𝑗(𝑠)0

∑4𝑗=1 𝐷0𝑗 (0− , 𝑠)𝐴0𝑗 (𝑠)𝑒 𝜆0𝑗(𝑠)0 = ∑4𝑗=1 𝐷1𝑗 (𝑠)𝐴1𝑗 (𝑠)𝑒 𝜆𝑀𝑗(𝑠)0

(27)

For further derivation, it is necessary to introduce a group of auxiliary functions 𝜕

𝑓1 (𝑥) = 𝜕𝑥 ,𝑢1 (𝑥, 0+ ) − 𝑢0 (𝑥, 0− ){ 𝜕 𝑓2 (𝑥) = 𝜕𝑥 ,𝑣1 (𝑥, 0+ ) − 𝑣0 (𝑥, 0− )-

CR IP T

(28)

Considering the following continuity conditions of displacement fields at the coating-substrate interface

𝑢0 (𝑥, 0− ) = 𝑢1 (𝑥, 0+ ) 𝑥 ∉ (−𝑎, 𝑎) 𝑣0 (𝑥, 0− ) = 𝑣1 (𝑥, 0+ )

we have

(29)

AN US

{

𝑓1 (𝑥) = 𝑓2 (𝑥) = 0, 𝑥 ∉ (−𝑎, 𝑎) So that +∞

+𝑎

∫−∞ 𝑓1 (𝑥) 𝑑𝑥 = ∫−𝑎 𝑓1 (𝑥) 𝑑𝑥 = 0 +∞

(31)

+𝑎

∫−∞ 𝑓2 (𝑥) 𝑑𝑥 = ∫−𝑎 𝑓2 (𝑥) 𝑑𝑥 = 0

M

{

(30)

Substituting Eqs. (20, 22) into the auxiliary functions and then applying Fourier

ED

transform to the obtained equations, the following equations can be obtained. 1

𝑎

PT

∑4𝑗=1 𝐸0𝑗 𝐴0𝑗 − ∑4𝑗=1 𝐸1𝑗 𝐴1𝑗 = ∫−𝑎 𝑓1 𝑒 𝑖𝑠𝑥 𝑑𝑥 𝑖𝑠 { 1 𝑎 4 4 ∑𝑗=1 𝐹0𝑗 𝐴0𝑗 − ∑𝑗=1 𝐹1𝑗 𝐴1𝑗 = ∫−𝑎 𝑓2 𝑒 𝑖𝑠𝑥 𝑑𝑥 𝑖𝑠

(32)

CE

Relating Eqs. (24-27, 32), we have a linear system of equations ,𝛹-(4𝑀+4)×(4𝑀+4) ∙ 𝐴(4𝑀+4) = 𝐹(4𝑀+4)

(33)

AC

where

,𝛹-(4𝑀+4)×(4𝑀+4)

,𝜙0 (0)-4×4 0 ⋮ 0 = 0 0 [,𝑋0 (−𝐻0 )-2×4

−,𝜙1 (0)-4×4 ,𝜙1 (𝑕1 )-4×4 ⋮ 0 0 0 0

0 −,𝜙2 (𝑕1 )-4×4 ⋮ 0 0 0 0

⋯ 0 ⋯ 0 ⋱ ⋮ ⋯ ,𝜙𝑀−2 (𝑕𝑀−2 )-4×4 ⋯ 0 ⋯ 0 ⋯ 0

0 0 ⋮ −,𝜙𝑀−1 (𝑕𝑀−2 )-4×4 ,𝜙𝑀−1 (𝑕𝑀−1 )-4×4 0 0

𝐴(4𝑀+4) = ,𝐴01 , … , 𝐴04 , 𝐴11 , … , 𝐴14 , … … , 𝐴𝑀4 -𝑇 1

𝑎

1

𝑎

𝐹(4𝑀+4) = ,i𝑠 ∫−𝑎 𝑓1 𝑒 𝑖𝑠𝑥 𝑑𝑥 , i𝑠 ∫−𝑎 𝑓2 𝑒 𝑖𝑠𝑥 𝑑𝑥 , 0,0, … , … … ,0-𝑇

0 0 ⋮ 0 (34) −,𝜙𝑀 (𝑕𝑀−1 )-4×4 ,𝑋𝑀 (𝑕𝑀 )-2×4 0 ]

(35) (36)

where ,𝜙𝑛 (𝑦)-4×4 (𝑛 = 0, 1, 2, … , 𝑀), ,𝑋0 (−𝐻0 )-2×4 and ,𝑋𝑀 (𝑕𝑀 )-2×4 are given 11

ACCEPTED MANUSCRIPT

in Appendix. After solving Eq. (33), the unknown functions 𝐴0𝑗 can be expressed using auxiliary functions 𝑓1 and 𝑓2 . 𝑎

𝑎

𝐴0𝑗 = 𝑞1𝑗 ∫−𝑎 𝑓1 𝑒 i𝑠𝑥 𝑑𝑥 + 𝑞2𝑗 ∫−𝑎 𝑓2 𝑒 i𝑠𝑥 𝑑𝑥

(37)

where 𝑞1𝑗 = ,𝛹 −1 -𝑗1 /(i𝑠) 𝑞2𝑗 = ,𝛹 −1 -𝑗2 /(i𝑠)

(38)

CR IP T

{

Substituting 𝐴0𝑗 into the loading conditions Eq. (13b), the following group of integral equations can be obtained. 1

+∞

+∞

lim𝑦→0−

+∞ +∞ 1 𝑎 ∫ ,𝑓 (𝑢) ∫−∞ 𝑕21 (𝑦, 𝑠)𝑒 i𝑠(𝑢−𝑥) 𝑑𝑠 +𝑓2 (𝑢) ∫−∞ 𝑕22 (𝑦, 𝑠)𝑒 i𝑠(𝑢−𝑥) 𝑑𝑠-𝑑𝑢 2𝜋 −𝑎 1

AN US

{

𝑎

lim𝑦→0− 2𝜋 ∫−𝑎,𝑓1 (𝑢) ∫−∞ 𝑕11 (𝑦, 𝑠)𝑒 i𝑠(𝑢−𝑥) 𝑑𝑠 +𝑓2 (𝑢) ∫−∞ 𝑕12 (𝑦, 𝑠)𝑒 i𝑠(𝑢−𝑥) 𝑑𝑠-𝑑𝑢 = 𝑝1 (𝑥)

(39)

= 𝑝2 (𝑥)

where 𝑕𝑖𝑗 (𝑦, 𝑠) (𝑖, 𝑗 = 1, 2) are given in Appendix.

In order to obtain the singular integral equations, the asymptotic analyses are conducted. Then we have

2𝜇

01 lim|𝑠|→∞ 𝑕12 (0, 𝑠) = −i 𝑘+1

M

{ 2𝜇01 lim|𝑠|→∞ 𝑕21 (0, 𝑠) = −i 𝑘+1

(40)

ED

A series of mathematic manipulation is conducted after substituting Eq. (40) into Eq. (39). Finally, the singular integral equations can be obtained as 1

𝑎 −4𝜇01 𝑓2 (𝑢)

∫ , 2𝜋 −𝑎

PT

𝑘+1 𝑥−𝑢 𝑎 −4𝜇01 𝑓1 (𝑢) ∫ , 2𝜋 −𝑎 𝑘+1 𝑥−𝑢

{1

+ 𝑘11 (𝑥, 𝑢)𝑓1 (𝑢) + 𝑘12 (𝑥, 𝑢)𝑓2 (𝑢)-𝑑𝑢 = 𝑝1 (𝑥)

+ 𝑘21 (𝑥, 𝑢)𝑓1 (𝑢) + 𝑘22 (𝑥, 𝑢)𝑓2 (𝑢)-𝑑𝑢 = 𝑝2 (𝑥)

(41)

CE

where 𝑘𝑖𝑗 (𝑥, 𝑢), 𝑖, 𝑗 = 1, 2 are given in Appendix. 3. The numerical solutions

AC

In this work, Lobatto-Chebyshev method is used to solve the singular integral

equations. In order to get the standard form of Eq. (41), the variables u and x are normalized 𝑥

𝑢

𝑡 = 𝑎, 𝑤 = 𝑎

(42)

Then Eq. (41) can be transformed into following form 1 1 −4𝜇01 𝑓2 (𝑎𝑤) + 𝑎𝑘11 (𝑎𝑡, 𝑎𝑤)𝑓1 (𝑎𝑤) + 𝑎𝑘12 (𝑎𝑡, 𝑎𝑤)𝑓2 (𝑎𝑤)-𝑑𝑤 ∫ , 2𝜋 −1 𝑘+1 𝑡−𝑤 { 1 1 −4𝜇 𝑓 (𝑎𝑤) + 𝑎𝑘21 (𝑎𝑡, 𝑎𝑤)𝑓1 (𝑎𝑤) + 𝑎𝑘22 (𝑎𝑡, 𝑎𝑤)𝑓2 (𝑎𝑤)-𝑑𝑤 ∫ , 01 1 2𝜋 −1 𝑘+1 𝑡−𝑟

12

= 𝑝1 (𝑎𝑡) = 𝑝2 (𝑎𝑡)

(43)

ACCEPTED MANUSCRIPT

The auxiliary functions 𝑓𝑖 (𝑎𝑤) are written as products of unknown bounded functions 𝑔𝑖 (𝑤) and a weight function

1

.

√(1−𝑤 2 )

𝑓𝑖 (𝑎𝑤) =

𝑔𝑖 (𝑤)

(44)

√(1−𝑤 2 )

As 𝑡, 𝑤 ∈ ,−1,1-, they can be written in discrete forms as 𝑘−1

𝑤𝑘 = cos [𝑛−1 𝜋] , 𝑘 = 1,2, … , 𝑛 2𝑙−1 2𝑛−2

𝜋] , 𝑙 = 1,2, … , 𝑛 − 1

CR IP T

𝑡𝑙 = cos [

(45) (46)

Substituting Eqs. (44-46) into Eq. (43), the integral equations can be transformed into numerical integration form

AN US

1 −4𝜇 𝑔 (𝑤 ) ∑𝑛 * 01 2 𝑘 + 𝑎,𝑔1 (𝑤𝑘 )𝑘11 (𝑡𝑙 , 𝑤𝑘 ) + 𝑔2 (𝑤𝑘 )𝑘12 (𝑡𝑙 , 𝑤𝑘 )-+𝜐𝑘 2(𝑛−1) 𝑘=1 𝑘+1 𝑡𝑙 −𝑤𝑘 { 1 −4𝜇 𝑔 (𝑤 ) ∑𝑛𝑘=1* 01 1 𝑘 + 𝑎,𝑔1 (𝑤𝑘 )𝑘21 (𝑡𝑙 , 𝑤𝑘 ) + 𝑔2 (𝑤𝑘 )𝑘22 (𝑡𝑙 , 𝑤𝑘 )-+𝜐𝑘 2(𝑛−1) 𝑘+1 𝑡𝑙 −𝑤𝑘

= 𝑝1 (𝑡𝑙 )

(47)

= 𝑝2 (𝑡𝑙 )

1

where 𝜐𝑘 = 1, 𝑘 = 2, 3, … , 𝑛 − 1; 𝜐1 = 𝜐𝑛 = 2. Similarly, the same manipulation is applied on Eq. (31), and then a 2𝑛 × 2𝑛 linear system of equations of 𝑔𝑖 (𝑤𝑘 ) (𝑖 = 1, 2; 𝑘 = 1, 2, … , 𝑛) can be obtained after relating the acquired equations with Eq.

M

(47). When using the above method, the amount of the discrete points (𝑛) in Eqs. (45)

ED

and (46) is mainly related to the Young’s modulus ratio, the coating properties variation and the numerical accuracy. Numerical calculations show that usually 30~70

PT

discrete points are enough to obtain accurate results. According to Shbeeb et al. (1999b), the mixed-mode SIFs can be expressed as 2𝜇

𝑎

CE

01 √ 𝐾𝐼 (−𝑎) = lim𝑥→−𝑎 √2(−𝑎 − 𝑥)𝜎0𝑦𝑦 (𝑥, 0− ) = − 𝑘+1 𝑔2 (−𝑎) { 2𝜇01 √𝑎 𝐾𝐼 (𝑎) = lim𝑥→𝑎 √2(𝑥 − 𝑎)𝜎0𝑦𝑦 (𝑥, 0− ) = 𝑘+1 𝑔2 (𝑎)

AC

2𝜇

(48)

√𝑎

01 𝐾𝐼𝐼 (−𝑎) = lim𝑥→−𝑎 √2(−𝑎 − 𝑥)𝜎0𝑥𝑦 (𝑥, 0− ) = − 𝑘+1 𝑔1 (−𝑎) { 2𝜇01 √𝑎 𝐾𝐼𝐼 (𝑎) = lim𝑥→𝑎 √2(𝑥 − 𝑎)𝜎0𝑥𝑦 (𝑥, 0− ) = 𝑘+1 𝑔1 (𝑎)

(49)

Furthermore, the SERRs can be expressed as 𝐺(𝑎) =

𝜋(𝑘+1) 8

,𝐾𝐼2 (𝑎) + 𝐾𝐼𝐼2 (𝑎)-

So the solution of SIFs and SERRs can be obtained after 𝑔𝑖 (𝑤𝑘 ) are solved.

13

(50)

ACCEPTED MANUSCRIPT

4. Results and analyses The following analyses are all set under the plane strain state. There are two kinds of loading conditions in this work: the normal stress loading ( 𝑝1 (𝑥) = 𝜎0 , 𝑝2 (𝑥) = 0) and the shear stress loading (𝑝1 (𝑥) = 0, 𝑝2 (𝑥) = 𝜏0 ). In order to analyze the calculation results conveniently, the obtained SIFs and SERRs are

where

CR IP T

normalized in the form of 𝐾𝑖 (±𝑎)/𝐾0 and 𝐺(±𝑎)/𝐺0 respectively ( 𝑖 = 𝐼, 𝐼𝐼 ),

𝐾0 = 𝑝0 √𝑎, 𝑝0 = 𝜎0 or 𝑝0 = 𝜏0 𝐺0 =

𝜋(𝑘+1) 8𝜇01

𝑝0 2 𝑎

(51) (52)

As the geometry and the properties of the strip are symmetric with respect to

AN US

y-axis and the loadings applied on the crack are uniform, it is found in the calculation results that 𝐾𝐼 (𝑎) = 𝐾𝐼 (−𝑎), 𝐾𝐼𝐼 (𝑎) = −𝐾𝐼𝐼 (−𝑎), 𝐺(𝑎) = 𝐺(−𝑎) when the crack surface is under the normal stress loading and 𝐾𝐼 (𝑎) = −𝐾𝐼 (−𝑎), 𝐾𝐼𝐼 (𝑎) = 𝐾𝐼𝐼 (−𝑎), 𝐺(𝑎) = 𝐺(−𝑎) when the crack surface is under the shear stress loading. For

M

the convenience of analyses, the mode I and mode II SIFs and SERRs shown in the

4.1 Verification

ED

following figures are the 𝐾𝐼 (𝑎), 𝐾𝐼𝐼 (𝑎) and 𝐺(𝑎) respectively.

In order to verify the validity of the PE model and the related algorithm, some

PT

calculating examples are given to for the comparisons with previous studies. It is shown in Fig. 3 that present work agrees well with the results of Chen and

CE

Erdogan (1996). 𝜃𝑎 is the horizontal value. 𝜃 = ln(𝐸1 /𝐸0 )⁄𝑎 , where 𝐸1 and 𝐸0 are the moduli of the top and bottom surfaces of the FGM coating respectively; 𝑎 is

AC

the half length of the crack and 𝑕1 is the thickness of the homogeneous substrate. What is noteworthy is that the FGM coating is divided into three layers in present work to verify the PE model, while one layer is enough for the coating with properties vary exponentially in general analyses.

14

CR IP T

ACCEPTED MANUSCRIPT

Fig.3 Comparisons between the present results and those of Chen and Erdogan (1996)

AN US

In Fig. 4a and Fig. 4b, it is obvious that both mode I and mode II SIFs of the present work are consistent with those of Shbeeb and Binienda (1999), no matter under the normal stress loading or the shear stress loading. In their work, the coating consisting of a FGM layer, whose properties varied exponentially, and a homogeneous

M

layer. 𝜃 and 𝑎 mean the same as Fig.3, and 𝑕1 , 𝑕2 and 𝑕3 are the thicknesses of the homogeneous substrate, the FGM layer and the homogeneous layer, respectively.

ED

The results show that the PE model is able to depict the property variation of the

AC

CE

PT

structure containing homogeneous sublayer.

(a)

15

CR IP T

ACCEPTED MANUSCRIPT

(b)

AN US

Fig 4. Comparisons between the present results and those of Shbeeb and Binienda (1999). (a) mode I SIFs; (b) mode II SIFs.

4.2 Results of Different Coating Property Variation and Structure Geometry Now that the validity of the PE model has been verified, how the coating

M

property variation and the structure geometry affect the fracture behaviors of the interface crack between the FGM coating and the homogeneous substrate can be

ED

investigated. Four kinds of property variation, defined by exponential functions, linear functions, power-law functions with 𝑔 = 2 and 𝑔 = 0.5, respectively, are discussed.

𝑞(𝑥) = 𝑞0 𝑒 𝜃𝑥

CE

written as

PT

When the actual shear modulus is determined by exponential functions, 𝑞(𝑥) is

(53)

When the actual shear modulus is determined by power-law functions, 𝑞(𝑥) is

AC

written as

𝑥

𝑞(𝑥) = 𝑞0 + (𝑞𝑀 − 𝑞0 )(𝑕 ) 𝑔

(54)

𝑀

It is obvious that the shear modulus varies linearly when 𝑔 = 1 in Eq. (54). As the coating-substrate structure is isotropic, the Young’s modulus ratio between the top surface of coating and the interface between coating and substrate 𝐸(𝑕𝑀 )/𝐸(𝑕0 )

is

used

to

depict

the

coating

property

variation,

while

𝐸(𝑕𝑀 )⁄𝐸(𝑕0 ) = 𝑞(𝑕𝑀 )/𝑞(𝑕0 ). Poisson’s ratio is set as 0.3. Considering the actual 16

ACCEPTED MANUSCRIPT

application of FGMs, the cases that modulus ratio 𝐸(𝑕𝑀 )/𝐸(𝑕0 ) varies in the range between 0.1 and 10 are investigated. Firstly, the division mode and the number of coating layers that is enough for calculation convergence should be confirmed. So that the balance between the computing efficiency and accuracy can be reached. For the exponential functions, it is obvious that one layer is enough; for linear functions, the thicknesses of the sub-layers

CR IP T

are all the same with each other. When it comes to the power-law functions, the case is different. As the derivatives of power-law functions are not constants, the zones where the properties vary more drastically should be divided into more sub-layers to ensure the accuracy. So in the present work the ratio between the thicknesses of

AN US

adjacent layers (𝑕𝑛+1 − 𝑕𝑛 )⁄(𝑕𝑛 − 𝑕𝑛−1 ) (𝑛 = 1, 2, … , 𝑀 − 1) is set as 0.9 when 𝑔 = 2, while (𝑕𝑛+1 − 𝑕𝑛 )⁄(𝑕𝑛 − 𝑕𝑛−1 ) = 1.5 when 𝑔 = 0.5. It is shown in Fig. 5a that the results of the division 𝑀 = 6 are almost the same with those of 𝑀 = 8 and 𝑀=10 (𝐻0 is the thickness of the substrate). So it is enough to divide the coating into

M

6 layers for linear functions in present work. From Fig. 5b and Fig. 5c, it can be seen

AC

CE

PT

ED

that 𝑀 = 8 is enough for power-law functions with 𝑔 = 2 and 𝑔 = 0.5 both.

17

M

AN US

CR IP T

ACCEPTED MANUSCRIPT

ED

Fig. 5. Comparisons among the results of different FGM coating divisions under the normal stress loading. (a) Linear functions; (b) power-law functions with 𝑔 = 2; (c)

AC

CE

PT

power-law functions with 𝑔 = 0.5.

18

AN US

CR IP T

ACCEPTED MANUSCRIPT

M

Fig. 6 The variations of mode I SIFs with ln (E(hM )⁄E(h0 )) under the normal stress loading for different coating moduli definitions. (a) exponential functions; (b) linear

AC

CE

PT

ED

functions; (c) power-law functions with 𝑔 = 2;(d) power-law functions with 𝑔 = 0.5.

19

AN US

CR IP T

ACCEPTED MANUSCRIPT

M

Fig. 7 The variations of mode II SIFs with ln (E(hM )⁄E(h0 )) under the normal stress loading for different coating moduli definitions. (a) exponential functions; (b)

𝑔 = 0.5.

AC

CE

PT

ED

linear functions; (c) power-law functions with 𝑔 = 2;(d) power-law functions with

20

AN US

CR IP T

ACCEPTED MANUSCRIPT

M

Fig. 8 The variations of SERRs with ln (E(hM )⁄E(h0 )) under the normal stress loading for different coating moduli definitions. (a) exponential functions; (b) linear

ED

functions; (c) power-law functions with 𝑔 = 2;(d) power-law functions with 𝑔 = 0.5. Figs. 6-8 show the 𝐾𝐼 /𝐾0 , 𝐾𝐼𝐼 /𝐾0 and 𝐺/𝐺0 of the interface crack under the

PT

normal stress loading respectively. In order to investigate the influence of the ratio between crack length and substrate thickness, six different thicknesses of the substrate

CE

are involved. We let 𝑕𝑀 = 𝑎 and 𝐻0 vary from 𝑎 to 100𝑎. The influence of the modulus ratio 𝐸(𝑕𝑀 )⁄𝐸(𝑕0 ) is also taken into consideration. As the modulus ratio

AC

varies from 0.1 to 10, the range of ln (E(hM )⁄E(h0 )) is about (-2.3026, 2.3026). It is obvious that the magnitudes of 𝐾𝐼 /𝐾0 are much larger than those of 𝐾𝐼𝐼 /𝐾0

under the normal stress loading. Considering the expression of SERR, 𝐺/𝐺0 is affected significantly by 𝐾𝐼 /𝐾0 in these cases. So it is reasonable that the variation tendencies of 𝐺/𝐺0 are similar to 𝐾𝐼 /𝐾0 (comparing Fig. 6 and Fig. 8). Both 𝐾𝐼 /𝐾0 and 𝐺/𝐺0 decrease with the increase of ln(E(hM )⁄E(h0 )). Additionally, 𝐾𝐼 /𝐾0 and 𝐺/𝐺0 also decrease with the increase of substrate thickness; when the substrate thickness is more than 5𝑎, 𝐾𝐼 /𝐾0 and 𝐺/𝐺0 almost converge. 21

ACCEPTED MANUSCRIPT

Comparing with 𝐾𝐼 /𝐾0 and 𝐺/𝐺0 , the variation of 𝐾𝐼𝐼 /𝐾0 is more complicated. From Fig. 7, it can be found firstly that 𝐾𝐼𝐼 ⁄𝐾0 = 0 when 𝑕𝑀 = 𝐻0 = 𝑎 and ln(E(hM )⁄E(h0 )) = 0, as both the coating and the substrate are homogeneous and the crack is at the symmetric position of the structure. In some cases the curves of 𝐾𝐼𝐼 /𝐾0 turn from negative to positive, which means the absolute values of 𝐾𝐼𝐼 /𝐾0 decline with the increase of ln(E(hM )⁄E(h0 )) at first and then augment. That also

CR IP T

means the propagation direction of interface crack changes. The crack propagates towards the coating when 𝐾𝐼𝐼 ⁄𝐾0 < 0 and towards the substrate when 𝐾𝐼𝐼 ⁄𝐾0 > 0. It can be found from Figs. 7(a, b, d) that the variation amplitudes of 𝐾𝐼𝐼 /𝐾0 curves shrink visibly with the increase of substrate thickness, while those in Fig. 7c almost

AN US

keep unchanged. Additionally, the variation of 𝐾𝐼𝐼 ⁄𝐾0 is different if the modulus is defined by different functions. The variation rates of 𝐾𝐼𝐼 ⁄𝐾0 reduce with the increase of ln(E(hM )⁄E(h0 )) when the coating modulus is defined by exponential functions (seeing Fig. 7a), while 𝐾𝐼𝐼 ⁄𝐾0 varies more dramatically inversely when the coating

M

modulus is defined by power-law functions with 𝑔 = 2 (seeing Fig. 7c). If the coating modulus is defined by linear functions or power-law functions with 𝑔 = 0.5,

ED

the variation rates of 𝐾𝐼𝐼 ⁄𝐾0 augment with the increase of ln(E(hM )⁄E(h0 )) firstly and then remain almost unchanged or decline.

PT

Figs. 9(a-c) show the comparisons of 𝐾𝐼 /𝐾0 , 𝐾𝐼𝐼 /𝐾0 and 𝐺/𝐺0 among different coating moduli definitions when 𝐻0 /𝑎 = 2 under the normal stress loading.

CE

The variation tendencies keep roughly the same when 𝐻0 changes, so the results of Fig. 9 are representative enough. When the coating modulus is defined by power-law

AC

functions with 𝑔 = 2 and the top surface of the coating is softer than the substrate (ln(E(hM )⁄E(h0 )) < 0), 𝐾𝐼 /𝐾0 and 𝐺/𝐺0 are the lower than those of others; when the top surface of the coating is stiffer than the substrate (ln(E(hM )⁄E(h0 )) > 0), power-law functions with 𝑔 = 0.5 have the lowest 𝐾𝐼 /𝐾0 and 𝐺/𝐺0 . 𝐾𝐼𝐼 /𝐾0 of the four moduli definitions all increase with the increase of ln(E(hM )⁄E(h0 ); and power-law functions with 𝑔 = 2 have the minimum variation range, while power-law functions with 𝑔 = 0.5 have the largest one. 22

AC

CE

PT

ED

M

AN US

CR IP T

ACCEPTED MANUSCRIPT

Fig. 9. Calculation results comparisons under the normal stress loading among different coating moduli definitions when 𝐻0 /𝑎 = 2. (a) 𝐾𝐼 ⁄𝐾0 ; (b) 𝐾𝐼𝐼 ⁄𝐾0 ; (c) 𝐺 ⁄𝐺0 .

23

AN US

CR IP T

ACCEPTED MANUSCRIPT

M

Fig. 10 The variations of mode I SIFs with ln (E(hM )⁄E(h0 )) under the shear stress loading for different coating moduli definitions. (a) exponential functions; (b) linear

AC

CE

PT

ED

functions; (c) power-law functions with 𝑔 = 2;(d) power-law functions with 𝑔 = 0.5.

24

AN US

CR IP T

ACCEPTED MANUSCRIPT

M

Fig. 11 The variations of mode II SIFs with ln (E(hM )⁄E(h0 )) under the shear stress loading for different coating moduli definitions. (a) exponential functions; (b) linear

𝑔 = 0.5.

AC

CE

PT

ED

functions; (c) power-law functions with 𝑔 = 2; (d) power-law functions with

25

AN US

CR IP T

ACCEPTED MANUSCRIPT

M

Fig. 12 The variations of SERRs with ln (E(hM )⁄E(h0 )) under the shear stress loading for different coating moduli definitions. (a) exponential functions; (b) linear

ED

functions; (c) power-law functions with 𝑔 = 2;(d) power-law functions with 𝑔 = 0.5. Figs. 10-12 show the 𝐾𝐼 /𝐾0 , 𝐾𝐼𝐼 /𝐾0 and 𝐺/𝐺0 under the shear stress loading

PT

for the same geometries and material moduli variations as those in Figs. 6-8. Mode II SIFs are dominant under the shear stress loading. So like the relationships between

CE

𝐾𝐼 /𝐾0 and 𝐺/𝐺0 shown in Fig. 6 and Fig. 8 respectively, the variation tendencies of 𝐾𝐼𝐼 /𝐾0 in Fig. 11 and 𝐺/𝐺0 in Fig. 12 are similar to each other. Both 𝐾𝐼𝐼 /𝐾0 and

AC

𝐺/𝐺0 decrease with the increase of ln(E(hM )⁄E(h0 )). From Fig. 12a, it can be seen that 𝐺/𝐺0 almost decreases linearly with the increase of ln(E(hM )⁄E(h0 )) when the coating moduli are defined by exponential functions. For other coating moduli definitions,

𝐾𝐼𝐼 /𝐾0

and

𝐺/𝐺0

decrease

slowly

with

the

increase

of

ln(E(hM )⁄E(h0 )) when ln(E(hM )⁄E(h0 )) < 0, and the decreases accelerate when ln(E(hM )⁄E(h0 )) > 0. Increasing the thickness of substrate can reduce 𝐾𝐼𝐼 /𝐾0 and 𝐺/𝐺0 . But 𝐾𝐼𝐼 /𝐾0 and 𝐺/𝐺0 converge basically when 𝐻0 > 5𝑎. 𝐾𝐼 /𝐾0 in Fig. 10 turn from positive to negative in some cases. That means the 26

ACCEPTED MANUSCRIPT

moving tendency of the upper and lower crack surface turns from opening to pressing against each other with the increase of ln(E(hM )⁄E(h0 )). It can be found that the zero points of 𝐾𝐼 /𝐾0 curves move right with the increase of 𝐻0 . The 𝐾𝐼 /𝐾0 curves keep positive within the range shown in Fig. 10 when 𝐻0 is large enough. Furthermore, the 𝐾𝐼 /𝐾0 curves’ variation amplitudes decrease with the increase of 𝐻0 . That means increasing substrate thickness can reduce the influence of moduli

CR IP T

ratio on mode I SIFs.

The relative locations and variation tendencies of the curves corresponding to particular coating moduli definitions shown in Figs. 13(a-c), corresponding to 𝐾𝐼 /𝐾0 , 𝐾𝐼𝐼 /𝐾0 and 𝐺/𝐺0 respectively, are similar to each other. The relative locations and

AN US

variation tendencies are also remain unchanged roughly when 𝐻0 varies, so only showing the results of 𝐻0 /𝑎 = 2 as examples is reasonable. When the top surface of the coating is softer than the substrate (ln(E(hM )⁄E(h0 )) < 0), power-law functions with 𝑔 = 2 have the lowest 𝐾𝐼 /𝐾0 , 𝐾𝐼𝐼 /𝐾0 and 𝐺/𝐺0 ; when ln(E(hM )⁄E(h0 )) >

AC

CE

PT

ED

M

0, this position is taken by power-law functions with 𝑔 = 0.5.

27

AC

CE

PT

ED

M

AN US

CR IP T

ACCEPTED MANUSCRIPT

Fig. 13. Calculation results comparisons under the shear stress loading among different coating moduli definitions when 𝐻0 /𝑎 = 2. (a) 𝐾𝐼 ⁄𝐾0 ; (b) 𝐾𝐼𝐼 ⁄𝐾0 ; (c) 𝐺 ⁄𝐺0 .

28

AN US

CR IP T

ACCEPTED MANUSCRIPT

M

Fig. 14. The variations of SERRs under the normal stress loading with the increase of coating thickness for different coating moduli definitions. (a) exponential functions;

g=0.5.

AC

CE

PT

ED

(b) linear functions; (c) power-law functions with g=2;(d) power-law functions with

29

AN US

CR IP T

ACCEPTED MANUSCRIPT

M

Fig. 15. The variations of SERRs under the shear stress loading with the increase of

ED

coating thickness for different coating moduli definitions. (a) exponential functions; (b) linear functions; (c) power-law functions with g=2;(d) power-law functions with g=0.5.

PT

Finally, the influence of coating thickness is investigated. We let 𝐻0 /𝑎 = 100

CE

and coating thickness 𝑕𝑀 varies from 0.5𝑎 to 10𝑎 (0.5𝑎, 𝑎, 2𝑎, 5𝑎, 10𝑎). The range of ln(𝑕𝑀 /𝑎) is about (-0.69315, 2.30259). In Fig. 14 and Fig. 15, every curve

AC

corresponds to a modulus ratio 𝐸(𝑕𝑀 )/𝐸(𝑕0 ). When the normal stress loading is applied on the crack surface, it can be seen from Fig. 14 that SERRs decrease with the increase of coating thickness. The decreasing is very dramatic when the coating thickness increases from 0.5𝑎 to 𝑎, and the decreasing amplitudes are over 50%. Then the decreasing amplitudes of SERRs with the increase of 𝑕𝑀 slump when 𝑕𝑀 > 𝑎 (ln(𝑕𝑀 ⁄𝑎) > 0) and the SERRs converge to constants basically when 𝑕𝑀 > 5𝑎 (ln(𝑕𝑀 ⁄𝑎) > 1.60944). The SERRs also decrease with the increase of 𝐸(𝑕𝑀 )⁄𝐸(𝑕0 ). For the shear stress cases (seeing Fig. 15), the SERRs decrease with 30

ACCEPTED MANUSCRIPT

the increase of 𝑕𝑀 in general, but the variations are much less dramatic than the normal stress loading cases, and increasing 𝐸(𝑕𝑀 )⁄𝐸(𝑕0 ) reduces both the SERRs and the variation amplitudes of SERR curves. 5. Conclusions The piecewise-exponential model is developed to solve the mixed-mode interface crack problem of FGM coating-substrate structure with general coating

CR IP T

properties. A series of typical calculations are conducted. The results can be summarized as follows:

(1) Comparing with the results of the previous fracture problems of FGMs with particular properties, the validity and accuracy of the present PE model are

AN US

verified. The divisions of FGM coating for different coating property variations, which take both computational accuracy and efficiency into consideration, are determined through the verification calculations. (2) The influences of several parameters on the fracture behaviors of interface

M

crack are investigated, including the thicknesses of the coating and the substrate, the coating moduli variations and the moduli ratio between the top

ED

surface of the coating and the substrate. These studies are valuable for designing the coating properties variations of FGM coating-substrate

PT

structures.

(3) For the four types of coating moduli definitions, the increase of the moduli

CE

ratio can reduce SERRs under both normal stress loadings and shear stress loadings. If the top surface of the coating is softer than the substrate, the

AC

coating moduli determined by power-law functions with 𝑔 = 2 can help reduce SERRs efficiently; if the top surface of the coating is stiffer than the substrate, the SERRs of the coatings whose moduli determined by power-law functions with 𝑔 = 0.5 are lower than those of others. So if the material compositions of the substrate and the top surface of the coating have been determined, choosing appropriate moduli variations can help reduce SERRs according to the value of the moduli ratio.

(4) When the interface crack length is a constant, increasing the thicknesses of 31

ACCEPTED MANUSCRIPT

the coating and the substrate can both reduce SERRs. But the decreases of SERRs become negligible when the coating or substrate are thick enough. So the rational thicknesses of the coating and the substrate can be obtained, and the unnecessary increase of structure thicknesses for reducing SERRs can also be avoided in the process of designing 6. Acknowledgements

CR IP T

This work is supported by the NSFC (Grant Nos. 11432005 and 11322217). Appendix 𝑎𝑛11 = −

(1+𝑘)𝑠2 𝜇0𝑛 −1+𝑘

𝑎𝑛12 = −i𝑠𝜃𝑛 𝜇0𝑛 𝑎𝑛13 = 𝜃𝑛 𝜇0𝑛

𝑛 = 1, 2, … , 𝑀

2i𝑠𝜇0𝑛

𝑎𝑛21 =

AN US

𝑎𝑛14 = − −1+𝑘 {𝑎𝑛15 = 𝜇0𝑛

(A.1)

i(−3+𝑘)𝑠𝜃𝑛 𝜇0𝑛 −1+𝑘

𝑎𝑛22 = −𝑠 2 𝜇0𝑛 𝑎𝑛23 = −

−1+𝑘 (1+𝑘)𝜃𝑛 𝜇0𝑛

𝑛 = 1, 2, … , 𝑀

(A.2)

M

𝑎𝑛24 =

2i𝑠𝜇0𝑛

−1+𝑘

(1+𝑘)𝑠2 𝜇01

ED

{𝑎𝑛25 =

−1+𝑘 (1+𝑘)𝜇0𝑛

PT

𝑎011 = − 𝑎012 = 0 𝑎013 = 0

CE

(A.3)

2i𝑠𝜇01

𝑎014 = −

AC

−1+𝑘

−1+𝑘

{𝑎015 = 𝜇01 𝑎021 = 0 𝑎022 = −𝑠 2 𝜇01 𝑎023 = − 𝑎024 = 0 {𝑎025 =

𝐸𝑛𝑗 (𝑠) = − { 𝐹𝑛𝑗 (𝑠) = 1

2i𝑠𝜇01

(A.4)

−1+𝑘

(1+𝑘)𝜇01 −1+𝑘

i,(−1+𝑘)𝑠2 −(1+𝑘)𝜆𝑛𝑗 (𝜃𝑛 +𝜆𝑛𝑗 )𝑠,(−3+𝑘)𝜃𝑛 −2𝜆𝑛𝑗 -

32

𝑛 = 1, 2, … , 𝑀

(A.5)

ACCEPTED MANUSCRIPT

𝐵𝑛𝑗 (𝑠) = 𝐶𝑛𝑗 (𝑠) = 𝐷 (𝑠) = { 𝑛𝑗

*−(1+𝑘)𝑠2 +𝜆𝑛𝑗 ,8𝜃𝑛 +(5+𝑘)𝜆𝑛𝑗 -+𝜇0𝑛 (−3+𝑘)𝜃𝑛 −2𝜆𝑛𝑗 ,(−3+𝑘)𝑠2 −(1+𝑘)𝜆𝑛𝑗 2 -𝜇0𝑛

𝑛 = 1, 2, … , 𝑀

(−3+𝑘)𝜃𝑛 −2𝜆𝑛𝑗

(A.6)

i(𝜃𝑛 +𝜆𝑛𝑗 ),−(−3+𝑘)𝑠2 +(1+𝑘)𝜆𝑛𝑗 2 -𝜇0𝑛 𝑠,(−3+𝑘)𝜃𝑛 −2𝜆𝑛𝑗 2𝑠2 𝜇

𝐸0𝑗 (𝑦, 𝑠) = −1+𝑘01 { 𝑗 = 1, 3 2i𝑠√𝑠2 𝜇01 𝐹0𝑗 (𝑦, 𝑠) = ∓ −1+𝑘 2,∓(√𝑠2 +𝑘√𝑠2 )+𝑠2 𝑦-𝜇01

𝐸0𝑗 (𝑦, 𝑠) = −1+𝑘 { 2i𝑠(1∓√𝑠2 𝑦)𝜇01 𝐹0𝑗 (𝑦, 𝑠) = −1+𝑘 𝐵0𝑗 (𝑦, 𝑠) = −

4i𝑠3 𝜇01 2

𝐵0𝑗 (𝑦, 𝑠) =

{𝐷0𝑗 (𝑦, 𝑠) =

2i𝑠,±(5√𝑠2 +𝑘√𝑠2 )−2𝑠2 𝑦-𝜇01 2 −1+𝑘 2i𝑠,∓(√𝑠2 +𝑘√𝑠2 )+2𝑠2 𝑦-𝜇01 2 −1+𝑘

2,(3+𝑘)𝑠2 ∓2(𝑠2 )

ED

𝐸𝑛2 𝑒 𝜆𝑛2𝑕 𝐹𝑛2 𝑒 𝜆𝑛2𝑕 𝐶𝑛2 𝑒 𝜆𝑛2𝑕+𝜃𝑛 𝑕 𝐷𝑛2 𝑒 𝜆𝑛2𝑕+𝜃𝑛 𝑕

𝐸01 (𝑕, 𝑠)𝑒 𝜆01 𝑕 𝐹 (𝑕, 𝑠)𝑒 𝜆01 𝑕 = 01 𝐶01 (𝑕, 𝑠)𝑒 𝜆01 𝑕 [𝐷01 (𝑕, 𝑠)𝑒 𝜆01 𝑕

PT

,𝜙𝑛 (𝑕)-4×4

CE

,𝜙0 (𝑕)-4×4

(A.8a)

−1+𝑘

M

𝐶0𝑗 (𝑦, 𝑠) =

𝑗 = 1, 3

(A.7b)

AN US

−1+𝑘 4(𝑠2 )3⁄2 𝜇01 2

{𝐷0𝑗 (𝑦, 𝑠) = ∓

𝐸𝑛1 𝑒 𝜆𝑛1𝑕 𝐹𝑛1 𝑒 𝜆𝑛1𝑕 = 𝐶𝑛1 𝑒 𝜆𝑛1𝑕+𝜃𝑛𝑕 [𝐷𝑛1 𝑒 𝜆𝑛1 𝑕+𝜃𝑛𝑕

𝑗 = 2, 4

−1+𝑘 4i𝑠3 𝜇01 2

𝐶0𝑗 (𝑦, 𝑠) =

CR IP T

(A.7a)

3⁄2

𝑗 = 2, 4

(A.8b)

𝑦-𝜇01 2

−1+𝑘

𝐸𝑛3 𝑒 𝜆𝑛3 𝑕 𝐹𝑛3 𝑒 𝜆𝑛3𝑕 𝐶𝑛3 𝑒 𝜆𝑛3𝑕+𝜃𝑛 𝑕 𝐷𝑛3 𝑒 𝜆𝑛3𝑕+𝜃𝑛 𝑕

𝐸02 (𝑕, 𝑠)𝑒 𝜆02 𝑕 𝐹02 (𝑕, 𝑠)𝑒 𝜆02 𝑕 𝐶02 (𝑕, 𝑠)𝑒 𝜆02 𝑕 𝐷02 (𝑕, 𝑠)𝑒 𝜆02 𝑕

𝐸𝑛4 𝑒 𝜆𝑛4𝑕 𝐹𝑛4 𝑒 𝜆𝑛4𝑕 𝑛 = 1, 2, … , 𝑀(A.9) 𝐶𝑛4 𝑒 𝜆𝑛4𝑕+𝜃𝑛𝑕 𝐷𝑛4 𝑒 𝜆𝑛4 𝑕+𝜃𝑛𝑕 ]

𝐸03 (𝑕, 𝑠)𝑒 𝜆03 𝑕 𝐹03 (𝑕, 𝑠)𝑒 𝜆03 𝑕 𝐶03 (𝑕, 𝑠)𝑒 𝜆03 𝑕 𝐷03 (𝑕, 𝑠)𝑒 𝜆03 𝑕

𝐸04 (𝑕, 𝑠)𝑒 𝜆04 𝑕 𝐹04 (𝑕, 𝑠)𝑒 𝜆04 𝑕 (A.10) 𝐶04 (𝑕, 𝑠)𝑒 𝜆04 𝑕 𝐷04 (𝑕, 𝑠)𝑒 𝜆04 𝑕 ]

,𝑋0 (−𝐻0 )-2×4 =

𝐶01 (−𝐻0 , 𝑠)𝑒 −𝜆01𝐻0 𝐷01 (−𝐻0 , 𝑠)𝑒 −𝜆01𝐻0

AC

[

𝐶02 (−𝐻0, 𝑠)𝑒 −𝜆02 𝐻0 𝐷02 (−𝐻0 , 𝑠)𝑒 −𝜆02 𝐻0

𝐶 𝑒 𝜆𝑀1 𝑕𝑀 ,𝑋𝑀 (𝑕𝑀 )-2×4 = [ 𝑀1 𝜆 𝑕 𝐷𝑀1 𝑒 𝑀1 𝑀

𝐶03 (−𝐻0 , 𝑠)𝑒 −𝜆03𝐻0 𝐷03 (−𝐻0 , 𝑠)𝑒 −𝜆03𝐻0

𝐶𝑀2 𝑒 𝜆𝑀2 𝑕𝑀 𝐷𝑀2 𝑒 𝜆𝑀2 𝑕𝑀

𝐶𝑀3 𝑒 𝜆𝑀3 𝑕𝑀 𝐷𝑀3 𝑒 𝜆𝑀3 𝑕𝑀

𝐶04 (−𝐻0 , 𝑠)𝑒 −𝜆04 𝐻0 ](A.11) 𝐷04 (−𝐻0 , 𝑠)𝑒 −𝜆04𝐻0

𝐶𝑀4 𝑒 𝜆𝑀4 𝑕𝑀 ] 𝐷𝑀4 𝑒 𝜆𝑀4 𝑕𝑀

(A.12)

𝑕11 (𝑦, 𝑠) = ∑4𝑗=1 𝐶0𝑗 (𝑦, 𝑠)𝑞1𝑗 𝑒 𝜆0𝑗𝑦

(A.13a)

𝑕12 (𝑦, 𝑠) = ∑4𝑗=1 𝐶0𝑗 (𝑦, 𝑠)𝑞2𝑗 𝑒 𝜆0𝑗𝑦

(A.13b)

𝑕21 (𝑦, 𝑠) = ∑4𝑗=1 𝐷0𝑗 (𝑦, 𝑠)𝑞1𝑗 𝑒 𝜆0𝑗𝑦

(A.13c)

33

ACCEPTED MANUSCRIPT

𝑕22 (𝑦, 𝑠) = ∑4𝑗=1 𝐷0𝑗 (𝑦, 𝑠)𝑞2𝑗 𝑒 𝜆0𝑗𝑦

(A.13d)

+∞

𝑘11 (𝑥, 𝑢) = ∫−∞ 𝑕11 (0, 𝑠)𝑒 i𝑠(𝑢−𝑥) 𝑑𝑠 +∞

2𝜇

+∞

2𝜇

(A.14a)

01 𝑘12 (𝑥, 𝑢) = ∫−∞ ,𝑕12 (0, 𝑠) + i 𝑘+1 -𝑒 i𝑠(𝑢−𝑥) 𝑑𝑠 01 𝑘21 (𝑥, 𝑢) = ∫−∞ ,𝑕21 (0, 𝑠) + i 𝑘+1 -𝑒 i𝑠(𝑢−𝑥) 𝑑𝑠

+∞

References

(A.14c) (A.14d)

CR IP T

𝑘22 (𝑥, 𝑢) = ∫−∞ 𝑕22 (0, 𝑠)𝑒 i𝑠(𝑢−𝑥) 𝑑𝑠

(A.14b)

Bai, X.-M., Guo, L.-C., Wang, Z.-H., Zhong, S.-Y., 2013. A dynamic piecewise-exponential model for transient crack problems of functionally graded

AN US

materials with arbitrary mechanical properties. Theor. Appl. Fract. Mech. 66, 41-51.

Binienda, W.K., Arnold, S.M., 1995. Driving force analysis in an infinite anisotropic plate with multiple crack interactions. Int. J. Fract. 71, 213-245.

Chen, Y.F., Erdogan, F., 1996. The interface crack problem for a nonhomogeneous coating bonded to a homogeneous substrate. Journal of the Mechanics and Physics of

M

Solids 44, 771-787.

ED

Choi, H.J., Paulino, G.H., 2010. Interfacial cracking in a graded coating/substrate system loaded by a frictional sliding flat punch. Proc. R. Soc. A Math. Phys. Eng. Sci.

PT

466, 853-880.

Dag, S., Erdogan, F., 2002. A Surface Crack in a Graded Medium Under General

CE

Loading Conditions. Journal of Applied Mechanics 69, 580-588. Delale, F., Erdogan, F., 1983. The Crack Problem for a Nonhomogeneous Plane.

AC

Journal of Applied Mechanics 50, 609-614. Erdogan, F., 1985. CRACK PROBLEM FOR BONDED NONHOMOGENEOUS MATERIALS UNDER ANTIPLANE SHEAR LOADING, American Society of Mechanical Engineers (Paper). Gao, X.W., Zhang, C., Sladek, J., Sladek, V., 2008. Fracture analysis of functionally graded materials by a BEM. Composites Science and Technology 68, 1209-1215. Guo, L.-C., Wu, L., Ma, L., 2004. The interface crack problem under a concentrated load for a functionally graded coating–substrate composite system. Composite 34

ACCEPTED MANUSCRIPT

Structures 63, 397-406. Guo, L., Noda, N., 2014. Investigation methods for thermal shock crack problems of functionally graded materials-Part I: Analytical method. J. Therm. Stresses 37, 292-324. Guo, L., Wang, Z., Noda, N., 2012. A fracture mechanics model for a crack problem

A Math. Phys. Eng. Sci. 468, 2939-2961.

CR IP T

of functionally graded materials with stochastic mechanical properties. Proc. R. Soc.

Guo, L.C., Noda, N., 2007. Modeling method for a crack problem of functionally graded materials with arbitrary properties-piecewise-exponential model. International Journal of Solids and Structures 44, 6768-6790.

AN US

Guo, L.C., Noda, N., Wu, L., 2008. Thermal fracture model for a functionally graded plate with a crack normal to the surfaces and arbitrary thermomechanical properties. Composites Science and Technology 68, 1034-1041.

Huang, G.-Y., Wang, Y.-S., Yu, S.-W., 2004. Fracture analysis of a functionally graded

Structures 41, 731-743.

M

interfacial zone under plane deformation. International Journal of Solids and

ED

Kim, J.H., Paulino, G.H., 2002. Finite element evaluation of mixed mode stress intensity factors in functionally graded materials. Int. J. Numer. Methods Eng. 53,

PT

1903-1935.

Li, M., Tian, Y.L., Wen, P.H., Aliabadi, M.H., 2016. Anti-plane interfacial crack with

CE

functionally graded coating: static and dynamic. Theor. Appl. Fract. Mech. 86, 250-259.

AC

Noda, N., Jin, Z.H., 1993. Steady thermal stresses in an infinite nonhomogeneous elastic solid containing a crack. J. Therm. Stresses 16, 181-196. Noda, N., Jin, Z.H., 1994. A crack in functionally gradient materials under thermal shock. Arch Appl Mech 64, 99-110. Shbeeb, N.I., Binienda, W.K., 1999. Analysis of an interface crack for a functionally graded strip sandwiched between two homogeneous layers of finite thickness. Engineering Fracture Mechanics 64, 693-720. Shbeeb, N.I., Binienda, W.K., Kreider, K.L., 1999a. Analysis of the Driving Forces 35

ACCEPTED MANUSCRIPT

for Multiple Cracks in an Infinite Nonhomogeneous Plate, Part I: Theoretical Analysis. Journal of Applied Mechanics 66, 492-500. Shbeeb, N.I., Binienda, W.K., Kreider, K.L., 1999b. Analysis of the Driving Forces for Multiple Cracks in an Infinite Nonhomogeneous Plate, Part II: Numerical Solutions. Journal of Applied Mechanics 66, 501-506.

of singular integral equations. Q Appl Math 35, 173-183.

CR IP T

Theocaris, P.S., Ioakimidis, N.I., 1977. Numerical integration methods for the solution

Wang, B.L., Mai, Y.W., Noda, N., 2002. Fracture mechanics analysis model for functionally graded materials with arbitrarily distributed properties. Int. J. Fract. 116, 161-177.

AN US

Yan, Z., Jiang, L., 2009a. On a moving dielectric crack in a piezoelectric interface with spatially varying properties. Engineering Fracture Mechanics 76, 560-579. Yan, Z., Jiang, L.Y., 2009b. Study of a propagating finite crack in functionally graded piezoelectric materials considering dielectric medium effect. International Journal of

M

Solids and Structures 46, 1362-1372.

Yildirim, B., Dag, S., Erdogan, F., 2005. Three dimensional fracture analysis of FGM

ED

coatings under thermomechanical loading. Int. J. Fract. 132, 369-395. Zhi-He, J., Naotake, N., 1994. Transient thermal stress intensity factors for a crack in

PT

a semi-infinite plate of a functionally gradient material. International Journal of Solids

AC

CE

and Structures 31, 203-218.

36