Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction

Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction

Engineering Fracture Mechanics xxx (2016) xxx–xxx Contents lists available at ScienceDirect Engineering Fracture Mechanics journal homepage: www.els...

886KB Sizes 2 Downloads 106 Views

Engineering Fracture Mechanics xxx (2016) xxx–xxx

Contents lists available at ScienceDirect

Engineering Fracture Mechanics journal homepage: www.elsevier.com/locate/engfracmech

Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction Xue-Yang Zhang, Xian-Fang Li ⇑ School of Civil Engineering, Central South University, Changsha 410075, PR China

a r t i c l e

i n f o

Article history: Received 25 May 2016 Received in revised form 26 October 2016 Accepted 27 November 2016 Available online xxxx Keywords: Time–fractional heat conduction Cracked plate Heat shock Thermal stress intensity factor Super-diffusion

a b s t r a c t A time-fractional heat conduction equation is applied to analyze the thermal shock fracture problem of a cracked plate. For a thermoelastic plate subjected to heat shock at its surfaces, analytical temperature field and thermal stresses are obtained by using Laplace transform and finite sine transform under the assumption that crack and deformation do not alter the temperature field. With this solution, thermal stress intensity factors at the crack tips are numerically calculated through the weight function method for both an edge and a center crack, respectively. The influences of fractional order describing super-diffusion, normal diffusion, and sub-diffusion on the thermal stress intensity factors are discussed. Thermoelastic fields and the thermal stress intensity factors exhibit pronounced wave–like propagation characteristics for super-diffusion or strong diffusion, and have a similar trend to normal diffusion for sub-diffusion or weak diffusion. The most dangerous crack length and position are discussed for cold and hot shock. The classical thermal stress intensity factors can be recovered from the present results only setting the fractional order to unity. Ó 2016 Published by Elsevier Ltd.

1. Introduction The thermal shock fracture problems of a cracked plate are frequently observed in many engineering applications. Great progress has been made within the framework of the classical Fourier heat conduction theory and a large number of papers have been published. For example, a thermal shock problem of an elastic strip with an edge crack was studied by Nied [1], who considered the influence of surface convection during cooling on singular stress field near the crack tip and showed a great dependence of the thermal stress intensity factors on the Biot number. Lu and Fleck [2] further investigated two cases of thermal shock resistance, a cold shock resistance for an edge crack and a hot shock resistance for a central crack, respectively. It should be noted that according to the classical Fourier heat conduction theory, the speed of heat propagation in a body is infinite. It means that a thermal disturbance at any position in a body can be felt at any other point instantly, which is physically impossible. To overcome this disadvantage, a hyperbolic heat conduction theory was developed by taking finite speed of heat propagation into account. Using hyperbolic heat conduction equation, Chang and Wang [3] obtained the analytical solutions of temperature fields and thermal stresses and analyzed the crack propagation behavior in a semi-infinite brittle medium subjected to heat shock. Wang and Li [4] evaluated the thermal shock resistance of a plate subjected to a sudden temperature change under the framework of hyperbolic, non-Fourier heat conduction. Wang et al. [5] studied the thermal shock fracture of a semi-infinite medium based on the dual-phase-lag heat conduction model.

⇑ Corresponding author. E-mail address: [email protected] (X.-F. Li). http://dx.doi.org/10.1016/j.engfracmech.2016.11.033 0013-7944/Ó 2016 Published by Elsevier Ltd.

Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

2

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

Nomenclature a c c E Ea;b ðzÞ F 1;2 F ; F 1 h Hðt Þ K0 KI L t T0 T 1;2

a e

C

j k h

r

thermal diffusivity coefficient crack length non-dimensional crack length Young’s modulus the Mittag–Leffler function weight functions Fourier sine transform and its inversion layer thickness heaviside unit step function stress intensity factor reference value stress intensity factor Laplace transform time initial temperature applied temperature at the top and bottom surfaces fractional order strain Gamma function non-dimensional time thermal expansion coefficient temperature change stress

In recent years, as a natural extension of the classical differential and integral calculus, fractional or non-integer order calculus has been applied to various fields. Mainardi [6] provided a fractional calculus approach to simulate the processes of some basics phenomena in physics, such as relaxation, diffusion, oscillations and wave propagation. Many researchers applied fractional order differential equations as a useful approach to describe anomalous diffusion in complex systems [7]. According to the value of fractional order, anomalous diffusion can be characterized by the time-fractional diffusionwave equation and be divided into sub-diffusion (weak diffusion), normal diffusion, super-diffusion (strong diffusion) and ballistic diffusion [8]. Sub-diffusive transport occurs in widely different systems such as fractals, glasses, porous and random media. Super-diffusion has been observed in porous glasses, biological systems, inhomogeneous rocks, etc. Povstenko [9] described a time-fractional diffusion-wave equation accounting for memory effects and then applied time-fractional heat conduction equation to deal with thermoelastic problems of a medium. In [10], Povstenko considered an infinite medium with a spherical cavity in the framework of a time-fractional heat conduction theory and showed that the values of non-dimensional speed of maximum propagation v a decrease and approach 1 with the fractional order a increasing from the value 1 to 2. In [11] the radial diffusion in a cylinder is described using time-fractional diffusion equation and the distinction between sub-diffusion and super-diffusion is discussed by numerical analysis. The fractional heat conduction theory has been extended by Sherief et al. [12], Youssef [13]. In particular, Ezzat et al. [14,15] applied the theory of fractional order to treat electro- and magneto-thermoelasticity problems. Obraczka and Kowalski used fractional order heat transfer in ceramic materials [16]. In fired-clay ceramics, Wilson et al. found that fractional order diffusion obeys (time)1/4 law [17,18]. Li and Wang [19] showed that the thermal conductivity depends on the system size L and fractional order power. More applications can be seen in [20]. The mathematical properties of the time-fractional diffusion-wave equation which interpolates the parabolic, classical Fourier heat conduction equation and the wave equation are discussed in [21,22]. The aim of this paper is to apply time-fractional heat conduction equation to deal with fracture of a cracked plate under thermal shocks at its material surfaces. The thermal shock problems for both edge cracks and center cracks are analyzed under the assumption that crack and deformation do not alter the temperature field. Analytic temperature change and thermal stress are obtained by using the finite sine transform and Laplace transform. The thermal stress intensity factors are obtained through the weight function method. The variation of the thermal stress intensity factor as a function of thermal shock time and crack length are depicted graphically. In addition, the effect of different fractional orders reflecting superdiffusion, normal diffusion, and sub-diffusion on thermal stresses and the thermal stress intensity factors are compared and analyzed.

2. Formulation of the problem Consider a cracked elastic plate of thickness h occupying the region 0 6 x 6 h; 1 < y < 1, and subjected to thermal shock at the top and bottom surfaces. Here a Cartesian coordinate system is chosen such that the thickness-wise direction is denoted as the x-axis. An edge or internal crack lies at the segment c 6 x  1 or ðh  cÞ=2 6 x 6 ðh þ cÞ=2, as shown in Fig. 1 Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

3

(a) and (b), respectively. For convenience, plane stress problems are assumed, or the plate width is sufficiently small as compared to its thickness and length. Within the framework of uncoupled thermoelasticity in which the strain and crack do not alter the temperature field, the problem can be solved by decomposing the problem into three parts. We first determine the temperature field in the plate without crack, then evaluate thermal stresses resulting from change in temperature, and finally determine the stress intensity factors near the crack tips. 3. Temperature field Assume that an initial temperature of the elastic plate is denoted as T 0 , and temperature changes, T 1  T 0 and T 2  T 0 , are suddenly applied at the top and bottom surfaces of the plate, respectively. Boundary conditions are stated as

hðh; y; tÞ ¼ ðT 1  T 0 ÞHðtÞ

ð1Þ

hð0; y; t Þ ¼ ðT 2  T 0 ÞHðt Þ

ð2Þ

where Hðt Þ is the Heaviside unit step function, hðx; y; tÞ ¼ T ðx; y; tÞ  T 0 represents change in temperature in the plate and we denote h1 ¼ T 1  T 0 , h2 ¼ T 2  T 0 hereafter. In the present analysis, we only consider that the temperature change alters stress distribution and conversely, elastic deformation as well as crack does not give rise to the temperature change. Due to the above assumption, one infers that the temperature change h is independent of y. Thus, one suffices to determine the transient temperature distribution along its thickness direction. To this end, a generalized heat conduction equation is invoked. That is, the temperature change h should satisfy the following time-fractional heat conduction equation [9]

@ a hðx; t Þ @ 2 hðx; t Þ ¼a ; @t a @x2

0 < a < 2;

ð3Þ

where a is the thermal diffusivity coefficient with the dimension ½a = m2sa [7], and the Caputo fractional derivative

@ a hðx; t Þ 1 ¼ @t a Cðn  aÞ

Z 0

t

ðt  sÞna1

@ n hðx; sÞ ds; n  1 < a < n; @ sn

ð4Þ

is used. Here CðÞ is the Gamma function, n is the integer part of a þ 1. The fractional order parameter a and n have the following relations [8]:

0 < a < 1; ðn ¼ 1Þ for sub-diffusion;

a ¼ 1; for normal diffusion; 1 < a < 2; ðn ¼ 2Þ for super-diffusion: Note that Eq. (3) for a ¼ 1 reduces to the classical Fourier heat conduction equation. Specially, for two limiting cases, Eq. (3) corresponds to steady temperature field if taking a ¼ 0, and ballistic diffusion (or propagating heat waves without attenuation) if taking a ¼ 2 [23]. In addition, initial conditions are furnished below

(a)

(b)

Fig. 1. A plate of thickness h with a crack subjected to thermal shocks T 1 and T 2 on its top and bottom surfaces, respectively; (a) an edge crack of length c, (b) a center crack of length c.

Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

4

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

hðx; 0Þ ¼ 0; 0 < x < h; for 0 < a < 2;

ð5Þ

@hðx; 0Þ ¼ 0; 0 < x < h for 1 < a < 2: @t

ð6Þ

Eq. (6) is not required for sub-diffusion and normal diffusion. In other words, if 0 < a < 1, only one initial condition (5) is needed, whereas for 0 < a < 2, both (5) and (6) are required. Boundary conditions on fractional order or anomalous diffusion are interpreted from thermal equilibrium [24]. Similar to [25], the solution of time-fractional heat conduction Eq. (3) can be obtained by using the finite sine transform with respect to the spatial coordinate x and the Laplace transform with respect to time t. To achieve this, we appeal to the finite sine transform and inverse sine transform in the domain 0 6 x 6 h, which are defined by

F fhðx; t Þg  ~hðnk ; t Þ ¼ F 1

n

Z

h

hðx; t Þ sin ðxnk Þdx;

ð7Þ

0

1 o X ~hðnk ; t Þ  hðx; t Þ ¼ 2 ~hðn ; t Þ sin ðxnk Þ; h k¼1 k

ð8Þ

respectively, where nk ¼ kp=h; k ¼ 1; 2; . . .. With the help of the following property

( ) h i @ 2 hðx; t Þ ¼ n2k ~hðnk ; t Þ þ nk hð0; t Þ  ð1Þk hðh; t Þ ; F 2 @x

ð9Þ

performing the finite sine transform of both sides of Eq. (3) and initial conditions (5) and (6), bearing the boundary conditions (1) and (2) in mind we obtain a h i d h~ðnk ; t Þ þ an2k ~hðnk ; tÞ ¼ nk Hðt Þ h2  ð1Þk h1 ; a dt

ð10Þ

subject to

~hðnk ; 0Þ ¼ 0; for 0 < a < 2;  d~hðnk ; t Þ  ¼ 0; for 1 < a < 2: dt 

ð11Þ ð12Þ

t¼0

Since the Laplace transform of the Caputo derivative of ~ hðnk ; tÞ has the following property

 ( )  n1 m~ X @ a ~hðnk ; t Þ a ~ a1m @ hðnk ; t Þ ¼ s h ðnk ; sÞ  L s  a m  @t @t m¼0

; n  1 < a < n:

ð13Þ

t¼0

~ stands for the Laplace transform of ~ where s is the Laplace transform parameter and h h, i.e. ~ h ¼ Lf~ hg, we perform the Laplace transform to both sides of Eq. (10) subject to the initial conditions (11) and (12) and have

~h ðnk ; sÞ ¼

h i s1 nk h2  ð1Þk h1 sa þ an2k

ð14Þ

Now, taking the inverse Laplace transform of Eq. (14) yields

(

~hðnk ; t Þ ¼ L1

sa

s1 þ

an2k

)

h i i  h nk h2  ð1Þk h1 ¼ t a nk Ea;aþ1 an2k ta h2  ð1Þk h1 ;

ð15Þ

where Ea;b ðzÞ is the Mittag–Leffler function in two parameters a and b, namely

Ea;b ðzÞ ¼

1 X j¼0

zj ; a > 0; b > 0; z 2 C: Cðaj þ bÞ

ð16Þ

Then we perform the inverse finite sine transform to both sides of Eq. (15) and finally obtain the temperature change as follows

hðx; t Þ ¼

1 i  h 2 X t a n Ea;aþ1 an2k t a h2  ð1Þk h1 sin ðxnk Þ: h k¼1 k

ð17Þ

Furthermore, inserting Eq. (16) into Eq. (17), after some collection and algebra we get

hðx; t Þ ¼ h1 

1 1  2ðh1  h2 Þ X   4h1 X sin ðxnk Þ  sin ðxnk Þ  Ea an2k t a  1  Ea an2k t a nk h n h k¼1;3... k k¼1

ð18Þ

Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

5

where Ea ðzÞ is the Mittag–Leffler function of Ea;b ðzÞ with b ¼ 1, and the identity 1 2 X ½1  cos ðnk hÞ sin ðxnk Þ ¼1 h k¼1 nk

ð19Þ

has been used. If introducing the following dimensionless parameters

pffiffiffi

gk ¼ kp; j ¼ ata=2 =h; x ¼ x=h;

ð20Þ

one can express h in (18) as

hðx; t Þ ¼ h1  4h1

1 X sin ðxgk Þ

gk

k¼1;3...

1 X     sin ðxgk Þ  Ea j2 g2k  2ðh1  h2 Þ 1  Ea j2 g2k k¼1

gk

ð21Þ

In what follows, we consider two special cases. Firstly, if the temperature loading on both the top and bottom surfaces is identical, viz., T 2 ¼ T 1 (or h2 ¼ h1 ), from Eq. (21) we have

hðx; t Þ ¼ h1  4h1

1 X sin ðxgk Þ

gk

k¼1;3...

  Ea j2 g2k :

ð22Þ

In the particular case of fractional order a ¼ 1, the Mittag–Leffler function reduces to E1 ðxÞ ¼ ex , and the solution of timefractional temperature field collapses to the classical temperature field given in [4]:

hðx; t Þ ¼ h1  4h1

1 X sin ðxgk Þ

gk

k¼1;3...

ej

2 g2 k

:

ð23Þ

Secondly, if letting the temperature on the bottom surface remain unchanged, viz. T 2 ¼ T 0 ðh2 ¼ 0Þ, from Eq. (21) we have

hðx; t Þ ¼ h1  4h1

1 X sin ðxgk Þ

gk

k¼1;3...

1 X     sin ðxgk Þ  Ea j2 g2k  2h1 1  Ea j2 g2k : k¼1

gk

ð24Þ

4. Thermal stresses Heat flow within a plate induces temperature change and then transient stress response takes place. If we assume that the plate under consideration is very thin in the z direction and infinite in the y direction, then the problem is in plane stress state [26]. Suppose the plate is free of traction at both surfaces x ¼ 0 and x ¼ h. Therefore, rxx ¼ 0, rzz ¼ 0; rxy ¼ 0, and all non-vanishing field quantities are independent of y and z. Thus, the compatibility condition simplifies to 2

d

eyy

dx

2

¼ 0;

ð25Þ

which gives

eyy ¼ Ax þ B

ð26Þ

where A and B are the unknown to be determined from the following boundary conditions

Z

h

ryy ðx; tÞdx ¼ 0;

ð27Þ

ryy ðx; tÞxdx ¼ 0:

ð28Þ

0

Z

h

0

Therefore, with the following constitutive equation for thermoelastic media





ryy ðx; tÞ ¼ E eyy  khðx; tÞ :

ð29Þ

where E and k are Young’s modulus and the thermal expansion coefficient of the plate, respectively, we get the desired transient thermal stress in terms of dimensionless variables below

rðx; tÞ ¼ 4Ekh1

1 X

   2 sin ðxgk Þ Ea j2 g2k  2

k¼1;3;...

gk

gk

  1 X 1  Ea j2 g2k 2 1 þ cos ðgk Þ ðh2  h1 Þ þ ð2  6xÞ  sin ðxgk Þ þ 2Ek k¼1

gk

gk

gk

ð30Þ

For the case of equal temperature loading at two surfaces T 2 ¼ T 1 ðh2 ¼ h1 Þ, the transient thermal stress becomes

Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

6

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

rðx; tÞ ¼ 4Ekh1

   2 sin ðxgk Þ : Ea j2 g2k  2

1 X

ð31Þ

gk

gk

k¼1;3;...

As a check, we take a ¼ 1, meaning normal diffusion, and the above-obtained transient thermal stress is identical to the previous result [4]. Similarly, when two surfaces are imposed different temperature changes with an assumption of unchanged temperature at the bottom surface without loss of generality, from Eq. (24) one gets a transient thermal stress as

rðx; tÞ ¼ 4Ekh1

   2 sin ðxgk Þ Ea j2 g2k  2

1 X

gk

k¼1;3;...

 2Ekh1

1 X

  1  Ea j2 g2k 2

gk

k¼1

gk

gk

þ ð2  6xÞ

1 þ cos ðgk Þ

gk

 sin ðxgk Þ :

ð32Þ

5. Stress intensity factors Usually, the position of crack embedded in plate is arbitrary. To capture its essential feature of the influence of crack, it is sufficient to analyze two typical cases, an edge crack and a central crack since they are idealized to highly simplify a plate containing an isolated mode I crack of length c and reflect the two limiting cases of a flaw at arbitrary position. 5.1. A plate with an edge crack Consider an infinite plate, as shown in Fig. 1(a), under heat shock. Since the crack plane is normal to the surface of the plate, it does not perturb transient temperature distribution. Therefore, the temperature and the resulting thermal stress fields without crack can be straightly used for the cracked plate. Due to the assumption of traction-free crack surfaces, we must superpose a negative stress field to the crack surface. Thus the singular elastic field for a plate with an edge crack is only induced by the negative of the stress r given by (30). Therefore, using the appropriate weight function, the stress intensity factor K I near the edge crack tip can be calculated by the following integral which can be easily found in Tada’s prestigious and widely used handbook [27]:

pffiffiffi Z 1 2 h 1 rðx; tÞF 1 ðx; cÞ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi dx K I ¼ pffiffiffiffiffiffi 3=2 pc ð1  cÞ 1c 1  ½ð1  xÞ=c2

ð33Þ

Here, c ¼ c=h;  x ¼ x=h; F 1 ð x; cÞ is a non-dimensional weight function defined by Eq. (A1) in Appendix A. It is noted that the above integral computation is effective for a positive tensile stress since a negative compressive stress does not give rise to crack opening but closing. Of course, from another point of view, a negative stress intensity factor may be understood as a shield effect to prevent crack from advancing. 5.2. A plate with a central crack For a plate with a central crack (see Fig. 1(b)), with the obtained results, similar to the case of an edge crack, with the aid of the weight functions the stress intensity factors K I at crack tips A and B are given by the following integrals [27]:

KI ¼

pffiffiffipffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Z h tan ðpc=2Þ

1=2þc=2

1=2c=2

 rðx; tÞF 2 ðx; cÞ sin ðpðx  1=2ÞÞ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 þ dx sin ðpc=2Þ 1  ½cos ðpc=2Þ= cos ðpðx  1=2ÞÞ2

ð34Þ

 rðx; tÞF 2 ðx; cÞ sin ðpðx  1=2ÞÞ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1  dx sin ðpc=2Þ 1  ½cos ðpc=2Þ= cos ðpðx  1=2ÞÞ2

ð35Þ

for tip A, and

KI ¼

pffiffiffipffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Z h tan ðpc=2Þ

1=2þc=2

1=2c=2

for tip B, where c ¼ c=h;  x ¼ x=h, and F 2 ð x; cÞis a non-dimensional weight function given by Eq. (A2) in Appendix A. In particular, for the case of equal temperature loads T 2 ¼ T 1 ðh2 ¼ h1 Þ, the stress intensity factors for two tips A and B become

pffiffiffipffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Z K I ¼ 2 h tan ðpc=2Þ

1=2þc=2

1=2

rðx; tÞF 2 ðx; cÞ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi dx: 1  ½cos ðpc=2Þ= cos ðpðx  1=2ÞÞ2

ð36Þ

Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

7

6. Numerical results and discussion In this section, an illustrative example is given and numerical computations are carried out. Here as a representative, we only present numerical results for a cracked plate subject to the same heat shock at both surfaces, i.e. T 2 ¼ T 1 or h2 ¼ h1 .

Fig. 2. Evolution of dimensionless temperature change h=h1 as a function of the thickness  x for various fractional orders: (a) a ¼ 1:75; (b) a ¼ 1; (c) a ¼ 0:25.

Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

8

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

Fig. 3. Dimensionless thermal stress distribution

r=r0 along the thickness x direction for various fractional orders: (a) a ¼ 1:75; (b) a ¼ 1; (c) a ¼ 0:25.

Firstly, transient temperature distributions are shown for different values of fractional order a ¼ 1:75; 1; 0:25 in Fig. 2 (a)–(c), respectively. According to Fujita [21], the points taking its maximum propagate with a finite speed. With the fractional order parameter a increasing from 1 to 2, the speed drops from infinity to the (finite) wave speed of wave equation. Here it should be noted that the parameter j ¼ 0:5 implies that the wave front just reaches the center point, while it does not Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

9

yet reach the center point for j < 0:5 and passes through the center point for j > 0:5, respectively. From Fig. 2(a), one sees that the transient temperature change curves for super-diffusion exhibit wave-like characteristics for three different dimensionless times j ¼ 0:25, j ¼ 0:5 and j ¼ 0:75. The surface temperature remains unchanged, as expected. The temperature at the center position varies dramatically, which remains the initial temperature for j ¼ 0:25 since the wave front does not

Fig. 4. Evolution of the normalized stress intensity factors K I =K 0 for an edge crack against the parameter a ¼ 1; (c) a ¼ 0:25.

j for various fractional orders: (a) a ¼ 1:75; (b)

Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

10

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

arrive, then rises above the surface temperature for j ¼ 0:5 due to arrival of wave front and finally the temperature at the whole plate exceeds the surface temperature for j ¼ 0:75. By inspecting Fig. 2(b) and (c), we see that the temperature change at the center position exhibits a similar trend for a ¼ 1; 0:25. In particular, wave-like characteristics disappear. The temperature change for sub-diffusion is in proximity to that for normal diffusion (a ¼ 1). With the above-obtained temperature change in plate, the response evolution of the non-dimensional thermal stress rðx; tÞ=r0 along the thickness direction, where r0 ¼ Ekh1 , is depicted in Fig. 3(a)–(c) for different values of fractional order. It is evident from Fig. 3(a)–(c) that the transient thermal stress in plate has remarkably wave behavior for super-diffusion with 1 < a < 2, whereas it has a gradual tendency to its steady state for sub-diffusion with a  1. In particular, when subjected to hot shock (h1 > 0), we find that for sub-diffusion and normal diffusion (0 < a  1), the region close to the surfaces experiences a transient compressive stress while the region near the center is in a state of tensile stress. However, due to the wave behavior of super-diffusion, the thermal stresses for 1 < a < 2 become tensile and compressive, in an alternative manner, at the centre and surface zone. Similarly, under the action of cold shock (h1 < 0), the stress at any location of the plate will change its sign from the heating condition. Accordingly the tensile zone will appear at the region close to the surface and the centre will become compressive for sub-diffusion. The magnitude of thermal stress will not change. It implies that in the case of a  1, a central crack embedded in a plate obeying sub-diffusion heat conduction is the most dangerous for hot shock and an edge crack is the most dangerous for cold shock, respectively. This conclusion coincides with that in [2] for normal diffusion heat conduction. However, in the case of 1 < a < 2, although a central crack and an edge crack are still dangerous, they do not relate directly hot or cold heat shock. In analyzing crack problems, the stress intensity factor as an important fracture parameter describes the intensity of elastic stress fields near the crack tip. From Eqs. (33)–(36), it is easily seen that in addition to time-dependent feature, the obtained thermal stress intensity factors depend not only on the values of normalized crack length but also on the distribupffiffiffiffiffiffi tion of thermal stresses. Here, stress intensity factors are normalized by K 0 , i.e. K I =K 0 , where K 0 ¼  phEkh1 > 0 for cold pffiffiffiffiffiffi shock and K 0 ¼ phEkh1 > 0 for hot shock, respectively, c being the length of an edge crack or a central crack. For the case of an edge crack in a plate subjected to cold shock, Fig. 4(a)–(c) shows the normalized thermal stress intensity factors K I =K 0 against j for different values of a. In Fig. 4(a) with 1 < a < 2, due to the wave-like behavior for super-diffusion, the stress intensity factors change their signs with the dimensionless time parameter j increasing. Note that here a negative stress intensity factor has no physical meaning and only implies that it impedes crack growth in theory [28,29]. In practice, when a negative stress intensity factor occurs, opposing faces of a pre-existing crack near the crack tip come into contact. For the case of a  1, although each curve of K I =K 0 displays a peak value after a finite shock time and decreases to zero as j ! 1 for fractional order a  1, the wave-like response behavior observed in Fig. 4(a) disappears. The maximum stress intensity factor occurs only at a finite thermal shock time ðj < 0:3Þ, which depends on the crack length and diffusion coefficient. In practice, of much interest is the maximum thermal stress intensity factor. Fig. 5 shows the evolution of maximum nor=K 0 against the normalized crack length c for various values of a. The curves of malized thermal stress intensity factors K max I the normalized stress intensity factors first grow to a maximum value, and then gradually have a reduction with crack length c increasing. This indicates that for cold shock, the dangerous crack length is shorter cracks, rather than longer cracks. This conclusion agrees with that given in [1]. The reason is due to the very sharp decrease in the transient tensile thermal stress away from the cooled surface. By inspecting Fig. 3(a)–(c), in the case of h1 < 0 one has r1 < 0. Thus tensile stress occurs near the cold surface, whereas compressive stress occurs in the center zone. Consequently, for a short edge crack loaded by tensile stress gives rise to a positive thermal stress intensity factor. With crack length rising, the crack surface enters a region of compressive stress, which cancels the contribution of partial tensile stress. Thus, the thermal stress intensity factors grad-

Fig. 5. Dependence of the maximum normalized stress intensity factors K max =K 0 on the dimensionless edge crack length c. I

Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

11

ually become smaller. From Fig. 5, the position of peak value of thermal stress intensity factors is also dependent on the fractional order a. For the case of a central crack embedded in a plate subjected to hot shock, Fig. 6(a)–(c) shows the normalized stress intensity factors K I =K 0 against the parameter j for different values of a. Similar to the trend of an edge crack, the curves first climb

Fig. 6. Evolution of the normalized stress intensity factors K I =K 0 for a center crack against the parameter a ¼ 1; (c) a ¼ 0:25.

j for various fractional orders: (a) a ¼ 1:75; (b)

Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

12

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

Fig. 7. Dependence of the maximum normalized stress intensity factors K max =K 0 on the dimensionless center crack length c. I

to the maximum values and then drop down with increasing parameter j for a  1 (Fig. 6(b) and (c)), while the wave-like behavior is again seen in Fig. 6(a) for super-diffusion a > 1. =K 0 are also plotted against the normalized crack length c for a In Fig. 7, the maximum thermal stress intensity factors K max I center crack. Compared to the maximum stress intensity factor of an edge crack, that of a center crack of equal length is slightly smaller in magnitude. This is what we expect since the effect of free boundary increases the stress intensity factors of an edge crack. It is worth noting that for a center crack, the whole crack length is c, not 2c. In other words, for hot shock, center cracks are prone to advance, and moreover, the most dangerous crack length for central cracks is that its length is half of the layer thickness or more. 7. Conclusions Using time-fractional heat conduction, thermal shock fracture for a thermoelastic plate was studied. When subjected to thermal shock on the top and bottom surfaces, analytical expressions for the temperature and thermal stresses in an uncracked plate were obtained for two typical hot/cold shock cases of the same or different thermal shock at the top and bottom surfaces. The thermal stress intensity factors were obtained numerically by the weight function method. For sub-diffusion, the temperature change, thermal stresses, and thermal stress intensity factors have a similar trend to those for the classical Fourier heat conduction. For super-diffusion, they exhibit a wave-like propagation property. This is due to fact that the super-diffusion equation (1 < a < 2) interpolates the classical heat conduction equation (a ¼ 1) and the classical wave equation (a ¼ 2). A central crack embedded in a plate is the most dangerous for hot shock and an edge crack is the most dangerous for cold shock, respectively. Acknowledgements This work was supported by the National Natural Science Foundation of China (Nos. U1534207 and 11672336). Appendix A F 1 ð x; cÞ used in Eq. (33) is shown as follow [27]:



2

3 1  x 1  x 1  x F 1 ðx; cÞ ¼ f 1 ðcÞ þ f 2 ðcÞ þ f 4 ðcÞ ; þ f 3 ðcÞ c c c

ðA1Þ

where

f 1 ðcÞ ¼ 0:46 þ 3:06c þ 0:84ð1  cÞ þ 0:66c2 ð1  cÞ ; 5

2

f 2 ðcÞ ¼ 3:52c2 ; f 3 ðcÞ ¼ 6:17  28:22c þ 34:54c2  14:39c3  ð1  cÞ1:5  5:88ð1  cÞ5  2:64c2 ð1  cÞ2 ;

Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033

X.-Y. Zhang, X.-F. Li / Engineering Fracture Mechanics xxx (2016) xxx–xxx

13

f 4 ðcÞ ¼ 6:63 þ 25:16c  31:04c2 þ 14:41c3 þ 2ð1  cÞ1:5 þ 5:04ð1  cÞ5 þ 1:98c2 ð1  cÞ2 ; x; cÞ in Eqs. (34)–(36) reads [27] and F 2 ð

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi



2

p p 2x  1 1 : F 2 ðx; cÞ ¼ 1 þ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  1 1  cos c c 2 p2  4

ðA2Þ

References [1] Nied HF. Thermal shock fracture in an edge-cracked plate. J Therm Stress 1983;6(2–4):217–29. [2] Lu T, Fleck N. The thermal shock resistance of solids. Acta Mater 1998;46(13):4755–68. [3] Chang DM, Wang BL. Transient thermal fracture and crack growth behavior in brittle media based on non-Fourier heat conduction. Eng Fract Mech 2012;94:29–36. [4] Wang BL, Li JE. Thermal shock resistance of solids associated with hyperbolic heat conduction theory. Proc R Soc A Math Phys Eng Sci 2013;469 (2153):20120754. [5] Wang B, Li JE, Yang C. Thermal shock fracture mechanics analysis of a semi-infinite medium based on the dual-phase-lag heat conduction model. Proc R Soc A Math Phys Eng Sci 2015;471(2174):20140595. [6] Mainardi F. Fractional relaxation-oscillation and fractional diffusion-wave phenomena. Chaos Solit Fract 1996;7(9):1461–77. [7] Metzler R, Klafter J. The random walk’s guide to anomalous diffusion: a fractional dynamics approach. Phys Rep 2000;339(1):1–77. [8] Kimmich R. Strange kinetics, porous media, and NMR. Chem Phys 2002;284(1):253–85. [9] Povstenko YZ. Fractional heat conduction equation and associated thermal stress. J Therm Stress 2004;28(1):83–102. [10] Povstenko YZ. Fractional heat conduction equation and associated thermal stresses in an infinite solid with spherical cavity. Q J Mech Appl Math 2008;61(4):523–47. [11] Povstenko YZ. Fractional radial diffusion in a cylinder. J Mol Liq 2008;137(1–3):46–50. [12] Sherief HH, El-Sayed A, El-Latief AA. Fractional order theory of thermoelasticity. Int J Solids Struct 2010;47(2):269–75. [13] Youssef HM. Theory of fractional order generalized thermoelasticity. J Heat Transfer 2010;132(6):061301. [14] Ezzat MA, El Karamany AS. Theory of fractional order in electro-thermoelasticity. Eur J Mech A Solids 2011;30(4):491–500. [15] Ezzat MA. Magneto-thermoelasticity with thermoelectric properties and fractional derivative heat transfer. Physica B 2011;406(1):30–5. [16] Obraczka A, Kowalski J. Heat transfer modeling in ceramic materials using fractional order equations. In: Advances in the theory and applications of non-integer order systems. Springer; 2013. p. 221–9. [17] Wilson MA, Hoff WD, Hall C, McKay B, Hiley A. Kinetics of moisture expansion in fired clay ceramics: a (time)1/4 law. Phys Rev Lett 2003;90 (12):125503. [18] Wilson MA, Carter MA, Hall C, Hoff WD, Ince C, Savage SD, et al. Dating fired-clay ceramics using long-term power law rehydroxylation kinetics. Proc R Soc A Math Phys Eng Sci 2009;465(2108):2407–15. [19] Li B, Wang J. Anomalous heat conduction and anomalous diffusion in one-dimensional systems. Phys Rev Lett 2003;91(4):044301. [20] Povstenko Y. Fractional thermoelasticity, vol. 219. Springer; 2015. [21] Fujita Y. Integrodifferential equation which interpolates the heat equation and the wave equation. Osaka J Math 1990;27:309–21. [22] Luchko Y, Mainardi F, Povstenko Y. Propagation speed of the maximum of the fundamental solution to the fractional diffusion-wave equation. Comput Math Appl 2013;66(5):774–84. [23] Joseph DD, Preziosi L. Heat waves. Rev Mod Phys 1989;61(1):41–73. [24] Korabel N, Barkai E. Boundary conditions of normal and anomalous diffusion from thermal equilibrium. Phys Rev E 2011;83(5):051113. [25] Agrawal OP. Solution for a fractional diffusion-wave equation defined in a bounded domain. Nonlinear Dyn 2002;29(1–4):145–55. [26] Erdogan F, Wu BH. Crack problems in FGM layers under thermal stresses. J Therm Stress 1996;19(3):237–65. [27] Tada H, Paris PC, Irwin GR. The stress analysis of cracks handbook. 3rd ed. ASME; 2000. [28] Li X-F, Liu G-L, Lee KY. Effects of T-stresses on fracture initiation for a closed crack in compression with frictional crack faces. Int J Fract 2009;160 (1):19–30. [29] Li X-F, Lee KY, Tang G-J. Kink angle and fracture load for an angled crack subjected to far-field compressive loading. Eng Fract Mech 2012;82:172–84.

Please cite this article in press as: Zhang X-Y, Li X-F. Thermal shock fracture of a cracked thermoelastic plate based on time–fractional heat conduction. Engng Fract Mech (2016), http://dx.doi.org/10.1016/j.engfracmech.2016.11.033