# Modeling of Thermal and Mechanical Behavior of ZrB2-SiC Ceramics after High Temperature Oxidation

Authors: Jun Wei; Lokeswarappa R. Dharani; K. Chandrashekhara; Gregory E. Hilmas; William G. Fahrenholtz
Journal: Journal of Ceramics
 (2014)

---

## Abstract

The effects of oxidation on heat transfer and mechanical behavior of ZrB2-SiC ceramics at high temperature are modeled using a micromechanics based finite element model. The model recognizes that when exposed to high temperature in air ZrB2-SiC oxidizes into ZrO2, SiO2, and SiC-depleted ZrB2 layer. A steady-state heat transfer analysis was conducted at first and that is followed by a thermal stress analysis. A “global-local modeling” technique is used combining finite element with infinite element for thermal stress analysis. A theoretical formulation is developed for calculating the thermal conductivity of liquid phase SiO2. All other temperature dependent thermal and mechanical properties were obtained from published literature. Thermal stress concentrations occur near the pore due to the geometric discontinuity and material properties mismatch between the ceramic matrix and the new products. The predicted results indicate the development of thermal stresses in the SiO2 and ZrO2 layers and high residual stresses in the SiC-depleted ZrB2 layer.


---

## Body

## 1. Introduction

Ultrahigh temperature ceramics (UHTCs) such as zirconium diboride and hafnium diboride (ZrB2 and HfB2) have been proposed for thermal protection of hypersonic aerospace vehicles, which may be exposed to temperatures above 1500°C in oxidizing environments. These materials are chemically and physically stable above 1600°C and have melting points above 3000°C [1]. In particular, ZrB2 because of its lower theoretical density is attractive for aerospace applications [2]. Exposure of solid zirconium boride (ZrB2 (s)) to air at elevated temperatures results in its oxidation to solid zirconia (ZrO2 (s)) and liquid boria (B2O3 (l)). The oxidation resistance of ZrB2 (s) can be improved by adding SiC (s) to promote the formation of a silica-rich scale. At high temperature, above 1100°C, SiC (s) oxidizes by reaction to form SiO2 (l) which has a lower volatility and a higher melting point and viscosity compared with B2O3 (l) [3–5]. Based on the experimental observations, Fahrenholtz [3] proposed a reaction sequence for the formation of the SiC-depleted layer during the oxidation of ZrB2-SiC at 1500°C in air. The oxide scales that form on ZrB2-SiC consist of an outer layer of SiO2, a middle layer of porous ZrO2, sometimes filled with SiO2, and a layer of SiC-depleted ZrB2 adjoining the unoxidized ZrB2-SiC at around 1500°C [2–6]. Parthasarathy et al. [7] developed a chemical reaction model for the oxidation of ZrB2-SiC ceramics to predict the thicknesses of the above three new productions. For temperatures below ~1600°C, an external glassy SiO2 layer forms and completely fills in pores of the porous ZrO2 scale whereas at higher temperatures, the glassy scale recedes due to evaporation of SiO2 (l) so that it only partially fills the pores in the ZrO2 layer.The region of particular interest, from a mechanical perspective, is the interface between the pores and the corner of the pores in the ZrO2 scale. The pore itself may or may not be filled with liquid SiO2 (l). The interface therefore consists of three materials (ZrO2 scale, solid/liquid SiO2, and SiC-depleted ZrB2 layer) of significantly different thermal and mechanical properties. This thermomechanical mismatch and geometric discontinuity would lead to residual stresses, and additional stress concentrations during the cool-down process from the processing temperature, thereby leading to potential cracking. Some researchers [8–12], having performed furnace oxidation and high velocity thermal shock tests on ZrB2-SiC, have indeed shown cracking in the ZrO2 scale.The purpose of this study is to develop a thermal and mechanical simulation model for ZrB2-SiC ceramics after oxidation. A steady-state heat transfer analysis was conducted using finite element analysis (FEA) modeling. An adpative remeshing technique is employed in both heat transfer and thermal stress analysis. A “global-local modeling” technique is used to combine finite element with infinite element for the thermal stress and the stress concentration analysis near a pore. Temperature, thermal, and residual stress distributions will be presented.


## 2. FEA Model and Simulation Procedure for ZrB2-SiC after Oxidation

To simplify the problem, the ZrO2 scale was assumed to be of uniform thickness with regularly distributed pores. The pores were assumed to be straight, columnar in structure without tortuosity. A cylindrical representative volume unit (CRVU) was constructed and further treated as a two-dimensional (2D, pseudo-3D) axisymmetric problem subjected to local heating as shown in Figure 1. It is assumed that the body is stress free prior to heating.Figure 1
A schematic of a model with cylindrical representative volume unit (CRVU) subjected to local heating.The oxide scale that forms on ZrB2-SiC consists of an outer layer of SiO2, a middle layer of porous ZrO2, and a layer of SiC-depleted ZrB2 next to the unoxidized ZrB2-SiC [3]. Based on the chemical oxidation models [7] and the experimental observation [3], the FEA models for ZrB2-SiC ceramic after oxidation at high temperature were created as shown in Figure 2. An adaptive remeshing zone was created to cover the SiO2, ZrO2, ZrB2 (SiC-depleted) and part of the ZrB2-SiC base near the pore. The temperature dependent dimensions of the ZrO2 scale (crystalline oxide), glassy SiO2, and SiC-depleted ZrB2 layer in ZrB2-20 vol% SiC were obtained from the chemical oxidation model [7]. The unoxidized ZrB2-SiC ceramic is treated as a macroscale continuous solid with properties of a predetermined ratio of 4 : 1 of ZrB2 to SiC (ZrB2-20 vol% SiC).FEA models for ZrB2-SiC ceramic after oxidation. (a) Lower temperature (external layer of SiO2), and (b) higher temperature (partial evaporation of the glassy SiO2).
(a)
(b)The heat conduction equation for an axisymmetric problem can be expressed as(1)ρc∂T∂t=kr∂∂rr∂T∂r+k∂2T∂z2,
where t is the time, T is the temperature, r and z are polar axis and longitudinal axis, ρ is the mass density, c is the specific heat, and k is the thermal conductivity. The thermoelastic model is given by
(2)σ=Dε-αΔT,
where [D] is the elasticity matrix, α the coefficient of thermal expansion, ΔT the temperature increment, σ the stress vector, and ε the strain vector.The heat flux condition is given by(3)q=-βT-T0,
where, q is the heat flux, β is the surface film coefficient, and T0 is the sink temperature. All simulations were conducted using ABAQUS finite element code.The temperature dependent thermal and mechanical properties of the solid phases, needed for the heat transfer and mechanical analyses, can be found in the literature or databases [13, 14]. However, the thermal conductivity of liquid phase of SiO2 (l) and the elastic constants cannot be found. As such, the temperature dependence of the thermal conductivity of liquid SiO2 (l) and the elastic constants have to be predicted based on thermodynamics and some available test data. The predictive methods used for calculating the above properties are outlined in the next section. The cylindrical representative volume unit with equivalent pore diameter was treated as a 2D axisymmetric model (pseudo-3D). The modeling involves a steady-state heat transfer analysis representing local heat-up to calculate the temperature distribution and then a transient heat transfer analysis for 30 minutes representing a cool-down event to calculate the residual temperature distribution. The resulting temperature distributions were then applied to a thermomechanical finite element model to calculate the thermal stress distribution in the material. Adaptive remeshing technique was employed for the heat transfer analysis to improve accuracy. A “global-local modeling”, along with the adaptive remeshing technique, is used to combine finite element with infinite element for thermal stress analysis. The procedure is summarized in Figure 3.Figure 3
Flow diagram showing the procedure and steps for thermal and mechanical analyses.


## 3. Thermal and Mechanical Properties

As mentioned earlier, the thermal conductivity and elastic constants of liquid phases of SiO2 are not readily available. In an earlier work [15], the authors developed a method for calculating the thermal conductivity of liquid SiO2 at a given temperature. The following thermal conductivity equation for a liquid by Hirschfelder et al. [16] is used in the method [15]:
(4)k=2.8kBAρTM2/3cpcV1/2cl.In the above equation,kB is Boltzmann’s constant, n is molecules per unit volume for the liquid, A is Avogadro’s number, M is molar mass, and ρ(T) is the temperature dependent bulk density of the liquid, cp and cV are the specific heats at constant pressure and at constant volume, respectively, and cl is speed of sound in the liquid. The temperature dependent specific heat at constant pressure, cp, for liquid SiO2 was reported in [14]. Then, the specific heat at constant volume, cV, is calculated using the following relationship [17]:
(5)cV=cp2cp+α2Tcl2,
where α is the coefficient of thermal expansion. The temperature dependence of density and coefficient of thermal expansion of liquid SiO2 were given in [18]. The speed of sound in liquids of SiO2 was found in [19]. With all the needed parameters, the thermal conductivity of liquid SiO2 was calculated using (4).Using the temperature dependent values for density and speed of sound in liquid of SiO2, the bulk and shear moduli (K,G) of liquid SiO2 were calculated using the Newton-Laplace equation [20–24]:
(6)K=ρCL2;G=ρCT2,
where CL and CT are the sound velocity of longitudinal and transverse wave, respectively. To simplify the problem, the temperature dependent elastic properties were used for the liquid phase of SiO2 instead of the viscous properties because the stress state in liquid phase of SiO2 was not of interest in the present study.


## 4. Results and Discussions

A 2D (pseudo-3D) 4-node linear axisymmetric heat transfer quadrilateral element was used in the thermal analysis. Heat flux was used as an error indicator variable to control the adaptive remeshing rule [25]. The dimensions of the new products after oxidation were taken from the chemical oxidation model [7].Two steps were used in the heat transfer analyses. The first step in the heat transfer analysis was a steady-state analysis representing local heating at the top surface to calculate the temperature distribution. The second step was a transient heat transfer analysis for 30 minutes representing a cooling event to room temperature to predict residual temperature distribution. The surface heating temperature was set asTheat (1780 K or 2240 K) during heating and 293 K during cooling. Outside the local heating area, the sink temperature was set at 293 K. The initial temperature of the material was 293 K. The heating, cooling and sink temperature conditions are summarized in Figure 4. The surface film coefficient, β, was set as 2500 W/(m2·K) during the heating representing a high speed fluid flow and 100 W/(m2·K) during cooling assuming a cooler fluid flow next to a solid boundary in air. The surface film coefficient was set as 100 W/(m2·K) at all other boundaries during both heating and cooling.Figure 4
The heating, cooling, and sink temperature conditions used in the analysis.


### 4.1. Results of Heat Transfer Analysis

The calculated temperature distributions in the body after surface heating temperatures of 1780 K and 2240 K are shown in Figures5 and 6, respectively. In the following results, the temperatures shown in parenthesis correspond to the case of 2240 K. The maximum temperature at the top surface of the outer SiO2 layer is 1168 K (1492 K) which is less than the applied heating temperature of 1780 K (2240 K). This is due to the effect of the surface film coefficient on heat transfer between a fluid and a solid and the thermal conduction at the boundaries. The temperatures at the interface between the outer SiO2 layer and the ZrO2, and at the interface between the oxide scale and ZrB2, are 1160 K (1432 K) and 1148 K (1404 K), respectively. The temperature at the bottom surface is 1124 K (1370 K) which is much less than the heating temperature applied at the top surface. The temperatures at locations shown in Figures 5 and 6 were also calculated for different heating temperatures. The predicted temperatures in various materials and at the interfaces are linearly dependent on heating temperature, except in ZrO2 layer. This deviation could be due to the increase in ZrO2 layer thickness accompanied by a decrease in SiO2 layer thickness. Figure 7 shows the predicted heat flux distribution of ZrB2-SiC after steady-state analysis for heating to 2240 K. It is seen that a heat flux concentration occurs at the pore corner due to the geometric discontinuity and thermal conductivity mismatch.Figure 5
Predicted temperature distribution of ZrB2-SiC after steady-state analysis at heating to 1780 K.Figure 6
Predicted temperature distribution of ZrB2-SiC after steady-state analysis for heating to 2240 K.Figure 7
Predicted heat flux distribution of ZrB2-SiC after steady-state analysis for heating to 2240 K.


### 4.2. Results of Thermal Stress Analysis

In the thermal stress analysis, the layout of infinite elements and finite elements, as well as the displacement constraints for the stress analysis shown in Figure8, are used. The distribution of maximum principal stresses for the steady-state heating at 1780 K is shown in Figure 9(a). The maximum principal stress distribution after cooling from 1780 K is shown in Figure 9(b). The maximum value of the maximum principal stresses of 946 MPa occurs at the top surface of SiO2 layer. The temperature at this location is about 1166 K (Figure 5), which is below the glass melting point. The brittle glassy SiO2 is sensitive to tensile stress with an average tensile strength of 364±57 MPa [26]. Therefore, a tensile stress of 946 MPa may induce cracking in the SiO2 layer. The maximum value of the maximum principal stresses of 568 MPa occurs at the upper corner of the pore in the ZrO2 layer and is less than the flexural strength (900 MPa) of ZrO2 [27]. The maximum stress in the ZrB2 is 451 MPa and occurs near the lower corner of the pore. This is higher than the measured bend strength of ZrB2 [28]. The largest principal stress in the ZrB2-SiC is 191 MPa, located near the lower corner of the pore. The flexural strength of ZrB2-SiC is 1000 MPa [29–31]. For the cool-down case, the residual stresses in Figure 9(b) vary between 153 and 360 MPa.Figure 8
The infinite element, the finite element, and the boundary conditions for the thermal stress analysis.Maximum principal stresses distribution in the enlarged area near pore after steady-state thermal analysis at 1780 K.
(a)
Heating-up
(b)
Cool-downThe distribution in the maximum principal stresses near the pore for steady-state heating to 2240 K and cool-down from 2240 K to 293 K are shown in Figures10(a) and 10(b), respectively. The stresses for heating vary between 182 and 2702 MPa while the corresponding stresses for cool-down vary between 281 and 529 MPa. Once again, the maximum values occur in the SiO2 (2702 MPa) and ZrO2 (2224 MPa) near the pore. This may initiate tensile cracking. These results are consistent with the experimental observations by Levine et al. [8] who performed furnace oxidation and high velocity thermal shock of ZrB2 + 20 vol.% SiC ceramic tests. Their results show that both pores and cracks appeared in the ZrO2 when oxidized in air at 1927°C for ten 10-min cycles. The highest maximum principal stress near the lower corner of the pore in the ZrB2 layer is 657 MPa shown in Figure 10(a) which is higher than the bending strength of the ZrB2 as indicated above. For the cool-down case, the maximum principal stresses shown in Figure 10(b) are less than the respective material strengths except in the ZrB2 layer.Maximum principal stresses distribution in the enlarged area near pore after steady-state thermal analysis at 2240 K.
(a)
Heating-up
(b)
Cool-downThe results were also obtained for additional temperatures. The variation of the maximum principal stress at indicated locations in Figures9 and 10 with temperature T is shown in Figure 11(a) for heating up to T and in Figure 11(b) for cool-down from T to 293 K. The thermal stress at heating temperature T in the ZrB2-SiC matrix is relatively small and does not vary much with heating temperature T (Figure 11(a)). The residual stress at 293 K (Figure 11(b)) in the ZrB2-SiC decreases with the increasing heating temperature T. Watts et al. [32] measured thermal residual stresses in ZrB2-30 vol% SiC composites using neutron diffraction. Their results indicated that stresses begin to accumulate at about 1673 K during cool-down from the processing temperature of 2172 K. The stress increased to an average compressive stress of 880 MPa in the SiC phase and to an average tensile stress of 450 MPa in the ZrB2 phase. By using the rule of mixtures for 34 vol% SiC, the stress in the SiC (880 MPa) converts to an equivalent stress of 453 MPa which is very close to the measured stress of 450 MPa [32].Thermal and residual maximum principal stress at indicated locations for different heating temperatures.
(a)
Thermal stress
(b)
Residual stress


## 4.1. Results of Heat Transfer Analysis

The calculated temperature distributions in the body after surface heating temperatures of 1780 K and 2240 K are shown in Figures5 and 6, respectively. In the following results, the temperatures shown in parenthesis correspond to the case of 2240 K. The maximum temperature at the top surface of the outer SiO2 layer is 1168 K (1492 K) which is less than the applied heating temperature of 1780 K (2240 K). This is due to the effect of the surface film coefficient on heat transfer between a fluid and a solid and the thermal conduction at the boundaries. The temperatures at the interface between the outer SiO2 layer and the ZrO2, and at the interface between the oxide scale and ZrB2, are 1160 K (1432 K) and 1148 K (1404 K), respectively. The temperature at the bottom surface is 1124 K (1370 K) which is much less than the heating temperature applied at the top surface. The temperatures at locations shown in Figures 5 and 6 were also calculated for different heating temperatures. The predicted temperatures in various materials and at the interfaces are linearly dependent on heating temperature, except in ZrO2 layer. This deviation could be due to the increase in ZrO2 layer thickness accompanied by a decrease in SiO2 layer thickness. Figure 7 shows the predicted heat flux distribution of ZrB2-SiC after steady-state analysis for heating to 2240 K. It is seen that a heat flux concentration occurs at the pore corner due to the geometric discontinuity and thermal conductivity mismatch.Figure 5
Predicted temperature distribution of ZrB2-SiC after steady-state analysis at heating to 1780 K.Figure 6
Predicted temperature distribution of ZrB2-SiC after steady-state analysis for heating to 2240 K.Figure 7
Predicted heat flux distribution of ZrB2-SiC after steady-state analysis for heating to 2240 K.


## 4.2. Results of Thermal Stress Analysis

In the thermal stress analysis, the layout of infinite elements and finite elements, as well as the displacement constraints for the stress analysis shown in Figure8, are used. The distribution of maximum principal stresses for the steady-state heating at 1780 K is shown in Figure 9(a). The maximum principal stress distribution after cooling from 1780 K is shown in Figure 9(b). The maximum value of the maximum principal stresses of 946 MPa occurs at the top surface of SiO2 layer. The temperature at this location is about 1166 K (Figure 5), which is below the glass melting point. The brittle glassy SiO2 is sensitive to tensile stress with an average tensile strength of 364±57 MPa [26]. Therefore, a tensile stress of 946 MPa may induce cracking in the SiO2 layer. The maximum value of the maximum principal stresses of 568 MPa occurs at the upper corner of the pore in the ZrO2 layer and is less than the flexural strength (900 MPa) of ZrO2 [27]. The maximum stress in the ZrB2 is 451 MPa and occurs near the lower corner of the pore. This is higher than the measured bend strength of ZrB2 [28]. The largest principal stress in the ZrB2-SiC is 191 MPa, located near the lower corner of the pore. The flexural strength of ZrB2-SiC is 1000 MPa [29–31]. For the cool-down case, the residual stresses in Figure 9(b) vary between 153 and 360 MPa.Figure 8
The infinite element, the finite element, and the boundary conditions for the thermal stress analysis.Maximum principal stresses distribution in the enlarged area near pore after steady-state thermal analysis at 1780 K.
(a)
Heating-up
(b)
Cool-downThe distribution in the maximum principal stresses near the pore for steady-state heating to 2240 K and cool-down from 2240 K to 293 K are shown in Figures10(a) and 10(b), respectively. The stresses for heating vary between 182 and 2702 MPa while the corresponding stresses for cool-down vary between 281 and 529 MPa. Once again, the maximum values occur in the SiO2 (2702 MPa) and ZrO2 (2224 MPa) near the pore. This may initiate tensile cracking. These results are consistent with the experimental observations by Levine et al. [8] who performed furnace oxidation and high velocity thermal shock of ZrB2 + 20 vol.% SiC ceramic tests. Their results show that both pores and cracks appeared in the ZrO2 when oxidized in air at 1927°C for ten 10-min cycles. The highest maximum principal stress near the lower corner of the pore in the ZrB2 layer is 657 MPa shown in Figure 10(a) which is higher than the bending strength of the ZrB2 as indicated above. For the cool-down case, the maximum principal stresses shown in Figure 10(b) are less than the respective material strengths except in the ZrB2 layer.Maximum principal stresses distribution in the enlarged area near pore after steady-state thermal analysis at 2240 K.
(a)
Heating-up
(b)
Cool-downThe results were also obtained for additional temperatures. The variation of the maximum principal stress at indicated locations in Figures9 and 10 with temperature T is shown in Figure 11(a) for heating up to T and in Figure 11(b) for cool-down from T to 293 K. The thermal stress at heating temperature T in the ZrB2-SiC matrix is relatively small and does not vary much with heating temperature T (Figure 11(a)). The residual stress at 293 K (Figure 11(b)) in the ZrB2-SiC decreases with the increasing heating temperature T. Watts et al. [32] measured thermal residual stresses in ZrB2-30 vol% SiC composites using neutron diffraction. Their results indicated that stresses begin to accumulate at about 1673 K during cool-down from the processing temperature of 2172 K. The stress increased to an average compressive stress of 880 MPa in the SiC phase and to an average tensile stress of 450 MPa in the ZrB2 phase. By using the rule of mixtures for 34 vol% SiC, the stress in the SiC (880 MPa) converts to an equivalent stress of 453 MPa which is very close to the measured stress of 450 MPa [32].Thermal and residual maximum principal stress at indicated locations for different heating temperatures.
(a)
Thermal stress
(b)
Residual stress


## 5. Conclusion

A “global-local modeling” technique is used combining finite element with infinite element for thermal stress analysis for the oxidation effects on heat transfer and mechanical behavior of ZrB2-SiC ceramics at high temperature. Thermal conductivity was calculated for the liquid phase of SiO2 based on a theoretical formulation. The predicted temperature at the top surface of the outer SiO2 layer is less than the applied heating temperature due to the surface film coefficient effect on the heat transfer between a fluid and a solid and the thermal conduction at the boundaries. An increase in ZrO2 layer thickness, accompanied by a decrease in SiO2 layer thickness, during oxidation will affect heat transfer in the body. Heat flux concentration occurs at the pore corner due to the geometric discontinuity and the material property mismatch. Thermal and residual stress concentrations occur near the pore due to geometric discontinuity and the material properties mismatch between the ceramic matrix and the new products. Thermal stresses in the surface oxide layers consisting of SiO2 and ZrO2, are higher than their respective materials strengths. Thermal and residual stresses in the layer of a new oxidation product of SiC-depleted ZrB2 layer for both heating and cooling cases are higher than the material strength. Therefore, it is expected that damage may initiate in the layers of new oxidation products.



---

*Source: 169748-2014-11-11.xml*