The thermal shock resistance of solids is analysed for a plate subjected to a sudden temperature change under the framework of hyperbolic, non-Fourier heat conduction. The closed form solution for the temperature field and the associated thermal stress are obtained for the plate without cracking. The transient thermal stress intensity factors are obtained through a weight function method. The maximum thermal shock temperature that the plate can sustain without catastrophic failure is obtained according to the two distinct criteria: (i) maximum local tensile stress criterion and (ii) maximum stress intensity factor criterion. The difference between the non-Fourier solutions and the classical Fourier solution is discussed. The traditional Fourier heat conduction considerably overestimates the thermal shock resistance of the solid. This confirms the fact that introduction of the non-Fourier heat conduction model is essential in the evaluation of thermal shock resistance of solids.
Thermal shock resistance is one of the most important performance parameters in solids for high temperature environments which cause thermal stresses and risks for thermal shock damage. Examples are as varied as energy conversion systems, electronic devices and cutting tools . A common measure of thermal shock resistance is the maximum jump in surface temperature which a material can sustain without cracking. Early studies on the thermal shock resistance of materials use the criterion that the maximum thermal stress that appears on the materials under thermal shock equals or exceeds the strength limit of the materials. This is the well-known stress-based criterion. In a series of investigations, Hasselman [2–5] introduced thermal shock resistance parameters by means of comparing the thermal shock behaviour of ceramic materials in terms of their physical and mechanical properties. After the initiation of the fracture mechanics concept, the thermal shock resistance of materials has been determined by the fracture mechanics-based criterion as well as the stress-based criterion. The motivation is that cracks or defects are commonly observed during thermal shock [6,7]. Li et al.  modelled the thermal shock of ultra-high temperature ceramics under active cooling. Rangaraj & Kokini  investigated the fracture behaviour of single-layer zirconia-bond coat alloy composite coatings under thermal shock. Zhou & Kokini  analysed the effect of surface pre-crack on the fracture of thermal barrier coatings under thermal shock. Zhou et al.  studied the thermal shock resistance of ceramics through biomimetically inspired nanofins. It is well accepted that the thermal shock resistance of solids ceramic is strongly affected by factors such as the heat conductivity, the geometric shape and the size of the sample, which govern the temperature gradient, the crack density and the duration of thermal stresses [6,12,13]. The influence of surface cracking, crack depth and shape on the structural integrity of the reactor pressure vessel during pressurized thermal shock was studied . In addition, the material non-homogeneity, temperature dependence of the material properties and multiple cracking also have considerable influence on the thermal shock resistance behaviour of solids [15–17].
Traditionally, analysis and prediction of thermal shock resistance of solids have been based on classical Fourier heat conduction, which has very different physical bases from those in non-Fourier approaches. According to the classical Fourier heat conduction law, it is assumed that the speed of heat propagation in a body is infinite. The body will be affected by the boundary condition or initial condition at the instant. However, in reality, the speed of heat propagation in a body is always finite, making a thermal response behave like a wave . Laser heating and nanotechnology have created the problems of heat conduction within very small time and length scales. For example, a carefully controlled incident beam can be used to heat up a very small area at a rate of up to 180 K s−1 for a few nanoseconds . In such situations, researchers have reported that the predictions by Fourier heat conduction do not agree well with experimental observations. Maurer & Thompson  observed that the surface temperature of a slab taken immediately after a sudden thermal shock is 300 K higher than that predicted by Fourier law. The disagreement between Fourier prediction and the experimental observation is rooted in the unrealistic propagation speed of the thermal signal adopted by Fourier law. To better describe the wave-like behaviour of heat conduction, instead of using Fourier law, the hyperbolic equation, which takes finite heat travelling speed into account, was presented as a constitutive equation that was coupled with the local energy balance [21,22]. Since then, considerable effort has been devoted to the solutions for the non-Fourier heat conduction equation for the cases of one-dimensional heat media [23–28], semi-infinite media  and layered composite media [30,31]. Further, temperature dependence of the thermal conductivity, specific heat and thermal diffusivity was also studied . It is well established that hyperbolic, non-Fourier theories of heat conduction are required for certain problems of practical interest in contemporary engineering. Wang et al.  studied the non-Fourier heat conductions in nanomaterials based on the thermomass theory and suggested that the inertial effect of high-rate heat and the interactions between heat and surface in confined nanostructures dominate the non-Fourier heat conduction in nanomaterials. In fact, Fourier's law has met great challenges in heat conductions in many situations, such as ultra-small scales (both temporal and spatial scales) , heat conduction in bio-materials , large heat flux and quick heating or cooling [36,37]. In particular, under high heating rate, the non-Fourier effect is expected to have a pronounced influence on the stress state around the crack border [38–40]. In this situation, without considering the non-Fourier effect, error and misinterpretation of the predicted temperature and thermal stress may result. A better understanding of thermal shock resistance behaviour of materials associated with non-Fourier heat conduction is critical. However, there are still open questions that need to be addressed regarding the thermal shock of solids associated with non-Fourier heat conduction.
This paper studies the transient temperature and thermal stress for a crack in a plate under the framework of non-Fourier heat conduction. A thermal shock condition is applied on the surfaces of the plate. Time-varying temperature field and stress field for the un-cracked media are established in a closed form. Transient crack tip stress intensity factors are obtained exactly. The thermal shock resistance of the materials is given. The results show that if the non-Fourier effect is ignored, then thermal shock resistance of the material will be dramatically overestimated. This is true especially for small plate thicknesses. Therefore, application of non-Fourier heat conduction is essential for correct evaluation of thermal shock resistance of material under small length scale.
2. Hyperbolic, non-Fourier heat conduction
Suppose in a coordinate system x (whose components are xi, where i=1, 2, 3) there is a body occupying space Ω, which is surrounded by surface S. The temperature inside the body may vary from point to point and with time. Let T(x,t) denote this temperature which is assumed to be a continuous functions of the coordinates xi and time t. If there is no internal heat generation, the heat flow is controlled by the following conduction equation: 2.1where ρ(x) is the mass density, c(x) the specific heat. Hereafter, the summations over the indices i and j will be assumed when appearing twice in an equation.
One of the non-Fourier approaches is the wave theory of heat conduction, which uses the relaxation behaviour to describe pulsed heat transport at short time intervals (see [23,24] for a review of the literature). For this situation, the governing heat-transfer equation becomes 2.2in which λ is the thermal conductivity, qi are the components of the heat flux vector q and ∂T/∂xj the temperature gradients, τq is the thermal relaxation time, which is related to the collision frequency of the molecules within the energy carrier. The thermal relaxation time (τq) can be calculated using the thermal wave speed (C) and the thermal diffusivity (κ), where κ=λ/ρc, as τq∝κ/C2. If τq is much smaller than the time interval for a particular transient process, equation (2.2) can be expanded using a first-order Taylor series expansion to yield 2.3Substitution of equation (2.3), together with its derivative with respect to xi, into equation (2.1) gives the governing equation in terms of temperature: 2.4The hyperbolic heat conduction equation for the non-Fourier heat conduction law, equation (2.4), is referred to as the thermal wave equation [23,24].
3. Thermal stress associated with hyperbolic heat conduction
Consider the problem of a one-dimensional temperature rise in a plate along its thickness direction as shown in figure 1. Denote b as the thickness of the plate (therefore, x∈[0,b]). In this case, equation (2.4) can be re-written as 3.1where , and is a length parameter. If there is a sudden constant temperature change T0 applied at x=0 and x=b surfaces of the plate, the temperature field can be solved from the method of variable separating, to give 3.2where is the non-dimensional thickness of the plate, and 3.3where λ1m and λ2m are the roots of the following characteristic equation: 3.4and 3.5The temperature is symmetric with respect to . Through extensive numerical exercises, it is observed that equation (3.2) allows the following results: 3.6However, a rigorous mathematical proof of equation (3.6) from equation (3.2) is difficult and is beyond the scope of this paper.
Denote the strain of the plate perpendicular to the through-thickness direction as εyy. The compatibility condition that needs to be satisfied is d2εyy/dx2, giving a linear distribution of εyy along the x direction εyy=Ax+B, where A and B are constants. Because of symmetric heating, the plate stretches uniformly but does not bend. Therefore, the constant A must be zero so that εyy is independent of all spatial dimensions including x, and depends only on time t. Suppose that the plate is very thin in the thickness direction as well as in the z direction so that the plane stress condition holds. Therefore, σxx=0 and σzz=0, and the constitutive equation becomes σyy=E(εyy−αT). In addition, the plate expands freely so that there is vanishing axial force: . It follows that the temperature-induced stress associated with the temperature distribution is 3.7As a comparison, the temperature field based on the classical Fourier heat conduction will also be given here. In this case, the governing equation is ∂T/∂t=κ∇2T. In current one-dimensional case, after normalization, it becomes . The solution is 3.8The thermal stress associated with classical Fourier temperature distribution is 3.9
The solution of the problem is presented in non-dimensional form. This avoids having to identify the range of parameters where the effect substantially becomes important. Specially, a characteristic length parameter l0 of the medium is identified to normalize its scale. We expect that the length scale of the problem is small enough so that the classical Fourier heat conduction model fails but large enough so that the hyperbolic heat conduction-based continuum approach applies (otherwise, the definitions of thermal flow and crack have no meaning). The length scale for this could be at the micrometre scale. For example, the thermodynamic properties at atmospheric conditions (25°C) are : the thermal diffusivity λp/[ρc]p=71.4 m2 s−1, [ρc]p=2854.2 kJ K m−3 and [τq]p=6.35×10−12 s for platinum, and thermal diffusivity , [ρc]g=1620.6 kJ K m−3 and [τq]g=1×10−10 s for quartz glass. The subscripts p and g stand, respectively, the platinum and quartz glass. Accordingly, the characteristic length parameter (l0) is l0=21.3×10−6 m for platinum and l0=11.7×10−6 m quartz glass. These are at the micrometre scale, which is long enough for continuum mechanics to be applicable.
To describe the results in a more convenient way, a reference stress σ0 is introduced as σ0=EαT0. Variation in the dimensionless thermal stress with the position in the plate is plotted in figure 2 for selected value of time for . Owing to symmetry, the stress is only displayed for . For comparison, classic solutions are also shown. Obviously, the non-Fourier solution reduces to the classical Fourier solution for sufficiently large values of time. It can be seen that if T0>0, the region near the surface of the plate experiences a compressive stress while a tensile zone is developed at the centre of the plate. The compressive stress is largest at the surface of the plate, and the tensile stress is largest at the centre of the plate. The largest stresses as a function of time are plotted in figure 3 for selected values of plate thickness. The values of the stress based on classical Fourier heat conduction are also shown. It can be seen that the magnitude of the maximum compressive stress is σ0=EαT0 which is attained at the surface of the plate at the beginning of the heat shock (i.e. at ). The maximum tensile stress σ0 based on the classic heat conduction has no relationship with the thickness of the plate. By observing the numerical results, it is found that σ0 is σ0≈0.309EαT0, which appears at . These expressions are in non-dimensional forms and are valid for any materials. For non-classic heat conduction, the maximum tensile stress depends strongly on the thickness of the plate. It appears approximately at if is not too large. However, if the thickness of the plate is large enough, the tensile stress may reach its maximum value at a higher heating time. For sufficiently large values of plate thickness and heating time, there is no difference between the non-Fourier solution and the classical Fourier solution.
Another fact observed from figure 3 is that the stress based on the non-Fourier model displays discontinuity at a certain time for a given location , or at a certain location for a specified time . This reflects the wave behaviour of the thermal stress produced by the non-Fourier temperature variation, which also displays the wave behaviour as shown in figure 2. On the contrary, thermal stress field calculated by the classical Fourier heat conduction is a smooth, continuous function of time and space .
In the case of a cold shock condition (T0<0), at any location of the plate, the stress will change its sign from the heating condition. Accordingly, the surface of the plate will be tensile and the centre will become compressive. The magnitudes of the thermal stress subjected to a cooling will be the same as those in a heating condition.
4. Cold shock cracking of the plate
The quenching-strength test is a common laboratory approach to determine the thermal shock resistance of ceramics [6,41]. In this experiment, the material specimen is heated to a specific temperate T0, and then put in a cold liquid of room temperature (reference temperature, denoted as zero). Under a cold shock, the maximum tensile stress is attained at the surface of the plate while the compressive stress is largest at the centre of the plate. For most materials, their tensile strength is less than their compressive strength. Thus, for a material plate under a cold shock, its thermal shock resistance will be controlled by the tensile stress region at the surface of the plate. As the surface of the plate undergoes a transient tensile stress under a cold shock, it is then expected that a crack may initiate at the edge of the plate under the cooling condition. The crack grows unstably until its tip reaches the compressive region.
Consider again the infinite plate of figure 1 subjected to a cold shock. The plate may contain a number of large cracks on the scale of the plate thickness, then it is expected that cracking will commence from the ‘worst flaw’. A rational definition of ‘worst flaw’ is the one which has the largest transient mode I stress intensity factor . The problem can be synthesized from the general solution for an edge crack of length c in the plate as shown in figure 4 ( is the non-dimensional crack length). As the crack plane is normal to the upper and lower surfaces of the plate, it does not perturb the transient temperature distribution in this arrangement. Therefore, determination of the temperature distribution and the resulting thermal stress for the un-cracked plate would be quite straightforward. The stress given by equation (3.7), with opposite sign, becomes the crack surface tractions in the quasi-static, mixed boundary-value problem. Define the mode I stress intensity factor K1 according to . By using weight function, the thermal stress intensity factor K1 for an edge crack in a plate can be obtained by numerical integration : 4.1where , , and F1 is a non-dimensional weight function, given in Appendix A. Here, is the stress given by equation (3.7).
In the following analysis, the thermal stress intensity factors will be normalized by the reference stress intensity factor, which is or , where m stands for metre. In figure 5, the normalized thermal stress intensity factor K1 is plotted as a function of thermal shock time for different values of crack length . The normalized thickness of the plate is considered to be . For any given , the stress intensity factor increases with time from the initial zero value, displays a peak and then decreases as the thermal shock time continues to increase. When the thermal shock time reaches infinity (the steady state), the stress intensity factor vanishes. This means that the maximum thermal stress intensity factor occurs only at a finite thermal shock time (the transient state). For a larger crack length, the thermal shock time for the stress intensity factor to reach its peak value is higher. In addition, the thermal stress intensity factor predicted by the non-Fourier heat conduction is a discontinuous function of time. The peak value of K1 appears approximately at for the crack length region considered (note that both and are non-dimensional). Unlike non-Fourier heat conduction, classical Fourier heat conduction gives smooth and continuous curves for thermal stress intensity factors.
Figure 6 plots the stress intensity factor K1 as a function of crack length for selected values of thermal shock time . It can be seen that K1 increases with crack length, to a peak value, and then decreases as the crack length further increases. Moreover, the results clearly show that the peak value of K1 predicted by the non-Fourier heat conduction is considerably higher than the corresponding value based on the classical Fourier heat conduction. In other words, the thermal stress intensity factors are less if non-Fourier heat conduction is adopted.
Obviously, from figure 6, it can be seen that for a crack of a given initial length, the temporal behaviour is alike for all values of thermal shock time. The envelope of this family of curves is the locus of the highest stress intensity factors reached for every crack depth at any thermal shock time. Envelops for the values of plate thickness , 10 and 20 are displayed in figure 7, where the reference thermal stress intensity factor is so that the normalized peak value of the stress intensity factor does not depend on the plate thickness for the classical Fourier heat conduction. For each curve, the peak value of the normalized stress intensity factor Kpeak increases with crack length, displays a maximum value and then decreases with increasing crack length. The value of Kpeak for non-Fourier heat conduction has a strong dependency on the plate thickness and is always higher than that calculated by classical Fourier heat conduction. The difference between the non-Fourier Kpeak and classical Fourier Kpeak is significant when the thickness of the plate is small. In other words, Fourier solution holds only when the thickness of the plate is sufficiently large.
5. Hot shock cracking of the plate
The above analysis shows that under the heating condition, the plate has the largest tensile stress at its centre. It is then expected that a crack may initiate at the plate centre under a hot shock. The problem can be synthesized from the solution for a centre crack of length 2a in a plate. The surfaces of the plate are subjected to a sudden temperature increase T0. Similar to the case of the edge crack problem, by using the weight function, the thermal stress intensity factor K1 can be obtained by numerical integration (; figure 8): 5.1where , , , F2 is a non-dimensional weight function given in appendix B, and is the stress given by equation (3.7).
In figures 9 and 10, the thermal stress intensity factor K1 is plotted against time and crack length for . K1 is normalized such that it does not depend on the crack length and plate thickness. It can be seen that the thermal stress intensity factor can increase or decrease with time and crack length. For a larger crack length, the thermal shock time needed for K1 to reach its peak value is shorter. This is owing to the fact that the crack tips for a larger crack are closer to the surfaces of the plate where the thermal shock load is applied. Similar to the cold shock, the envelope of the curves is the locus of the peak stress intensity factors reached for every crack length at any thermal shock time. Envelops for the values of plate thickness 10 and 20 are displayed in figure 11, where the reference thermal stress intensity factor is . Once again, the normalized peak value of the stress intensity factor does not depend on the plate thickness for the classical Fourier heat conduction. However, the thermal stress intensity factor from the non-Fourier heat transport has a very different behaviour. Similar to the case of edge crack in a plate under a cold shock, the plate thickness has a profound influence on the normalized thermal stress intensity factor if non-Fourier heat conduction is adopted. The thinner the plate, the higher the normalized thermal stress intensity factor is for non-Fourier heat conduction.
6. Thermal shock resistance of the plate
Thermal shock resistance is a major issue in the design of materials for high or lower temperature applications. Generally, there are two distinct criteria for the determination of thermal shock resistance: the stress-based criterion and the fracture mechanics-based criterion.
(a) Fracture mechanics-based criterion
For each curve of figure 7 (cold shock) and 11 (hot shock), the peak value of the stress intensity factor increases with crack length, displays a maximum value and then decreases with increasing crack length. The maximum values of the peak stress intensity factor, , are dependent on the plate thickness and are listed in tables 1 and 2, respectively, for cold shock and hot shock. For comparison, for classical Fourier heat conduction is also given. Note that does not depend on the plate thickness for classical Fourier heat conduction but has a strong dependence on the plate thickness for non-Fourier heat conduction. In addition, classical Fourier heat conduction considerably underestimates the thermal stress intensity factor, especially for small plate thicknesses. Based on the fracture mechanics criterion, the material will fail if the condition is satisfied, where KIc is the fracture toughness of the material. Thus, the critical thermal shock temperature TfractureC (i.e. the thermal shock resistance) of the material can be sustained without catastrophic failure can be found from 6.1As is a complicated function of , equation (6.1) can only be solved numerically to obtain TfractureC. The results of TfractureC for selected values of plate thickness are listed in tables 1 and 2. The normalized thermal shock resistance based on the fracture mechanics criterion, , is an increasing function of plate thickness. Classical Fourier heat conduction overestimates the thermal shock resistance, especially for small plate thicknesses. When the thickness of the plate is sufficiently large, prediction of the thermal shock resistance by non-Fourier heat conduction approaches that by classical Fourier heat conduction.
(b) Stress-based criterion
The stress-based criterion is such that the maximum thermal stress that appears inside the plate is equal to the strength limit of the material (σc for compression and σb for tension). As σb is smaller than σc for most materials, only tensile stress needs to be considered. According to figure 3, the maximum tensile stress for the cold shock and hot shock can be obtained. The temperature sustainable by TstressC predicted by the stress-based criterion is given by equalling and σb. The results of TstressC are given in tables 1 and 2 for selected values of plate thickness. For cold shock, the surfaces of the plate experience the maximum thermal stress which is the same for any value of plate thickness, whether non-Fourier heat or classical Fourier heat conduction is adopted. Accordingly, the cold shock resistance TstressC of the material is the same. Under a hot shock, the maximum tensile stress at the centre of the plate is different for different values of plate thickness if non-Fourier heat conduction is adopted. As a result, based on non-Fourier heat conduction, the hot shock resistance of the material TstressC depends on the thickness of the plate. From the results, it can be seen that TstressC/(σb/[Eα]) is an increasing function of plate thickness. The classical Fourier solution holds when the thickness of the plate is large enough.
Finally, the true thermal shock resistance of the material Tc should be such that, under a thermal shock, the material will not fail either by ‘less of fracture toughness’ or by ‘less of strength limit’. Therefore, 6.2Information given in tables 1 and 2 may be useful to researchers for the design of thermally resistive materials under strong thermal shock environments.
Analytical models for non-Fourier heat conduction and the associated thermal fracture mechanics for a material plate subjected to a thermal shock on its surfaces were established. The temperature field and thermal stresses for the un-cracked plate were obtained in closed form. The thermal shock stress intensity factor histories were obtained numerically by the weight function method. The maximum temperature that the materials can sustain without catastrophic failure was investigated according to the maximum local tensile criterion and the maximum stress intensity factor criterion. Dependency of the thermal shock resistance of the material on the thickness of the plate was identified.
Based on the non-Fourier heat conduction analysis, it is found that the thermal shock resistance of the materials is very sensitive to the size of the materials. Thermal stress intensity factors are considerably higher for small plate thicknesses than for large plate thicknesses. The thermal shock resistance of the materials is weakened significantly for large material sizes. Classical Fourier heat conduction overestimates the thermal shock resistance of the materials.
This research was supported by the National Science Foundation of China (project no. 11172081) and Shenzhen Research Innovation Foundation, China (project no. JCYJ20120613150312764).
A1where , A2 A3 A4 and A5
- Received December 29, 2012.
- Accepted February 13, 2013.
- © 2013 The Author(s) Published by the Royal Society. All rights reserved.