# Heat Transfer and Ablation Prediction of Carbon/Carbon Composites in a Hypersonic Environment Using Fluid-Thermal-Ablation Multiphysical Coupling

Authors: Xuewen Sun; Haibo Yang; Tao Mi
Journal: International Journal of Aerospace Engineering
 (2020)

---

## Abstract

Carbon/carbon composites are usually used as a thermal protection material in the nose cap and leading edge of hypersonic vehicles. In order to predict the thermal and ablation response of a carbon/carbon model in a hypersonic aerothermal environment, a multiphysical coupling model is established taking into account thermochemical nonequilibrium of a flow field, heat transfer, and ablation of a material. A mesh movement algorithm is implemented to track the ablation recession. The flow field distribution and ablation recession are studied. The results show that the fluid-thermal-ablation coupling model can effectively predict the thermal and ablation response of the material. The temperature and heat flux in the stationary region of the carbon/carbon model change significantly with time. As time goes on, the wall temperature increases and the heat flux decreases. The ablation in the stagnation area is more serious than in the lateral area. The shape of the material changes, and the radius of the leading edge increases after ablation. The fluid-thermal-ablation coupling model can be used to provide reference for the design of a thermal protection system.


---

## Body

## 1. Introduction

With the development of hypersonic technology, the requirements of thermal protection materials are becoming higher and higher [1, 2]. Carbon/carbon composites are widely used in the thermal protection system of aircraft because of their high latent heat and good strength at high temperature [3]. The flow field around the aircraft affects the heat transfer and ablation of carbon/carbon composites, changing the thickness and shape of the thermal protection layer. The ablation of carbon/carbon composites changes the temperature and component distribution of the flow field, and then it will affect the ablation of carbon/carbon composites conversely [4, 5]. There is a strong coupling between a hypersonic flow field and carbon/carbon composite ablation. Therefore, it is important to establish a fluid-thermal-ablation coupling model and predict the temperature distribution and ablation response of carbon/carbon composites.Prediction of the ablation response of a thermal protection material is one of the key problems in the development of hypersonic vehicles. In recent years, relevant scholars have carried out some researches on it. Martin and Boyd established a coupling model between flow and a thermal protection material, mainly analysed the temperature response of a thermal protection material in a hypersonic environment, and validated the numerical method with IRV-2 [6]. They cross-studied the response of a thermal protection material in the rocket nozzle environment, established the coupling model of a flow field and a carbon/phenolic material, and studied the flow field environment and material response [7]. Mortensen and Zhong considered the ablation and real gas effect on the surface of a material, established a thermochemical nonequilibrium model of the flow field, and analysed the effect of ablation on the hypersonic boundary layer [8]. Chen et al. developed a fully implicit ablation response program for predicting the ablation response of a material. The program was used to simulate the flow of pyrolysis gas, ablation, and shape change of a thermal protection material [9, 10]. Kumar established a numerical coupling model between ablation and a flow field and analysed the interaction between pyrolysis gas and an external flow field [11]. Li et al. established a nonlinear pyrolysis model for carbonized composites and predicted ablation response [12]. Candler et al. studied the oxidation and ablation process of carbon matrix composites. The finite rate reaction model and chemical equilibrium gas-solid reaction model were compared to predict the material response, and the spherical cone model was used to verify the prediction [13]. Qin et al. established a multiscale thermochemical ablation model of carbon/carbon composites. The ablation response of a material was predicted considering the reaction of matrix and fiber [14]. Yin et al. [15], Meng et al. [16], and Chen [17] calculated the ablation response of carbon/carbon composites and predicted the shape change and temperature distribution of a carbon/carbon model. The existing studies for the response of a thermal protection material usually predict the material response by one-way coupling, without considering the effects of material temperature and ablation shape change on the external flow field.In this paper, a fluid-thermal-ablation coupling model of carbon/carbon composites in a hypersonic environment is established considering the thermochemical nonequilibrium of a flow field, heat transfer, and ablation of material. The temperature distribution, ablation recession, and shape change of the carbon/carbon model at different times are analysed.


## 2. Governing Equations

The numerical simulation of carbon/carbon composite ablation needs to consider the physical and chemical reactions, including aerodynamic heating, thermal response, and ablation boundary movement of a material. In this paper, Fluent is used to calculate the external flow field, and the aerodynamic thermal load can be obtained. Abaqus is used to calculate the transient thermal response of the material. Considering the finite rate ablation model, the ablation recession of carbon/carbon composites is predicted. The surface of the material model is tracked by a mesh moving strategy. Finally, Mpcci is used to realize the data exchange between a flow field and the material model.


### 2.1. Governing Equations of Aerodynamic Flow

For the thermochemical nonequilibrium flow, the Navier-Stokes equations with a chemical source term are expressed as follows:(1)∂Q∂t+∂E∂x+∂F∂y−∂Ev∂x+∂Fv∂y=S,Q=ρ1⋮ρnsρuρvρEρeve,E=ρ1u⋮ρnsuρu2+pρuvρHuρeveu,F=ρ1v⋮ρnsvρuvρv2+pρHvρevev,S=ω̇1⋮ω̇ns000ω̇ve,Ev=qx1⋮qxnsτxxτxyuτxx+vτxy+qx+qvex+ρ∑i=1nsDimhi∂Yi∂xqvex+ρ∑i=molDimevei∂Yi∂x,Fv=qy1⋮qynsτyxτyyuτyx+vτyy+qy+qvey+ρ∑i=1nsDimhi∂Yi∂yqvex+ρ∑i=molDimevei∂Yi∂y,where Q is the conservative variable. E and F are inviscid fluxes in the x and y directions, respectively. Ev and Fv are viscous fluxes in the x and y directions, respectively. S is the source term reflecting the effect of flow field chemical nonequilibrium. u and v are the velocities in the x and y directions, respectively. p, ρ, and Yi are the pressure, density, and mass fraction of components of the gas, respectively. qxi and qyi are the diffusivity of components in the x and y directions, respectively. Dim and hi are the diffusion coefficient and the absolute enthalpy per unit mass of components, respectively. qx, qy, qvex, and qvey are translational-rotational heat fluxes and vibrational heat fluxes in the x and y directions, respectively. ω̇i and ω̇ve are mass generation rate and vibration energy, respectively.The Gupta chemical kinetic model is used to calculate the component change caused by chemical reactions. The thermodynamic nonequilibrium in this paper is characterized by the Park two-temperature model. The specific chemical reaction model is referred to in Reference [18].


### 2.2. Fluid-Solid Interface and Ablation Rate Model

The energy conservation between the flow field and the material model needs to be satisfied. Figure1 shows the energy transfer at the gas-solid interface, which satisfies the equation as follows:
(2)q=qconv+qrad‐in+ṁchcs−ṁhcs−qrad‐out,where q is the heat flux into the material model. qconv is the aerodynamic heat flux. qrad‐in is the radiation heat flux by an external flow field. ṁchcs is the heat flux generated by the chemical reaction on the material surface. ṁhcs is the heat flux carried away by the ablation loss. qrad‐out is the radiation heat flux dissipation of the material to a flow field.Figure 1
Energy transfer at the gas-solid interface.The finite rate ablation model is used to calculate the surface ablation of carbon/carbon composites. The mass loss caused by mechanical erosion and surface melting is not considered. Surface ablation of carbon/carbon composites mainly includes oxidation and nitridation of carbon. The surface chemical reaction mechanism and mass loss rate are as follows [16]:
(1)
The oxidation of carbon:O+Cs⟶COThe mass loss rate is(3)ṁ1=ρCOkTω2πmOβOMCMO.
(2)
The oxidation of carbon:O2+2Cs⟶2COThe mass loss rate is(4)ṁ2=2ρCO2kTω2πmO2βO2MCMO2.
(3)
The nitridation of carbon:N+Cs⟶CNThe mass loss rate is(5)ṁ3=ρCNkTω2πmNβNMCMN,where ρ is the density of gas. Ci is the mass fraction of the inflow components. k is the Boltzmann constant. Tw is the wall temperature. mi is the mass of components. βi is the surface chemical reaction efficiency. Mi is the molar mass of components.
(6)βO=0.63exp−1160/Tω,βO2=0.5,βN=0.3.The total mass loss rate of the material is(7)ṁ=ṁ1+ṁ2+ṁ3.


### 2.3. Material Thermal Response and Ablation Boundary Tracking

The heat conduction will affect the temperature distribution, oxidation properties, and ablation of a material. Therefore, the heat conduction of a material must be considered in order to predict the ablation response of the material. The heat conduction in the material follows the Fourier heat equation and energy conservation. The governing equation of heat conduction in the material can be written as follows:(8)ρscp∂T∂t−∂∂xiλ∂T∂xi=q,where t is the time. T is the temperature. ρs is the material density. cp is the specific heat capacity of the material. λ is the heat conductivity of the material. q is the heat flux.The ablation of carbon/carbon composites is a dynamic process, and the model shape changes constantly. In order to predict the ablation response, it is necessary to track the location of the ablation surface. The linear ablation rateṠ of the material can be calculated by
(9)ṁ1+ṁ2+ṁ3=ρsṠ.The movement of nodes on a material surface can be calculated by the ablation rate and time stepΔt. The moving direction of nodes is perpendicular to the surface, and the displacement of nodes is as follows.
(10)δxδy=ṠΔtnxṠΔtny.The Abaqus mesh motion algorithm is used to simulate ablation surface movement, and the user-defined subroutine umeshmotion is used to update the location of ablation surface nodes. The ALE (Arbitrary Lagrangian-Eulerian) mesh adaptive technology is used to remesh the internal grids of the model, avoiding the excessive deformity of the grids, and finally realizing the simulation of the ablation surface recession process.


### 2.4. Fluid-Thermal-Ablation Coupling Strategy

The coupling model of hypersonic flow and carbon/carbon composites is realized by a zoning method. According to the physical space, it is divided into fluid and solid domains, as shown in Figure2. The coupling between the fluid domain and the solid domain is both a physical and a chemical process in which hypersonic aerodynamic heating, heat conduction, and ablation interact through the coupling interface. Ωf and Ωs are the fluid domain and the solid domain, respectively. Γ is a coupled interface, where a solid provides wall temperature and displacement boundaries for a fluid, and the fluid provides aerothermal loads for the solid. The data of a flow field and the material model can be exchanged repeatedly on the coupling interface, and the data exchange between two domains is realized by Mpcci. Figure 3 shows the data exchange between mismatched grids.Figure 2
Fluid-solid interaction.Figure 3
Data exchange between mismatched grids.In the numerical calculation of material ablation in a hypersonic environment, different solvers are used for the solid and fluid regions. From the physical aspect, the time for hypersonic flow to reach steadiness is much smaller than the time for transient structural heat transfer and ablation. The time scale discontinuity reflected in the numerical calculation is that the fluid solver can allow an extremely much smaller time step size than the solid solver, which presents a computational challenge in practical engineering applications. Loosely coupling analysis strategy can avoid the detailed transient aerodynamic flow analysis that usually requires significant computational effort. The accuracy and efficiency of the coupling procedure mainly depend on the coupling time step size. In this paper, we refer to the research results of the literature [19], and the time step size is 0.01 s.The fluid-thermal-ablation coupling process is shown in Figure4. The specific process is as follows: (1) At the initial time t0, the wall temperature and position of the carbon/carbon model are transferred to the flow field through the coupling interface as the boundary conditions for flow field calculation. (2) Based on the boundary conditions, the flow field is calculated to obtain the wall heat flux. (3) The heat transfers from the flow field to the solid domain through the coupling interface. (4) Based on the boundary conditions, the heat conduction in the material model from time t0 to time t1 is calculated, and the temperature distribution and ablation of the material at time t1 are obtained. (5) The wall temperature and position at time t1 is transferred to the fluid domain as the boundary condition of the flow field. (6) At this point, the calculation of a coupling time step is completed, and then the next coupling time step is continued, and so on until the whole coupling time is over.Figure 4
Fluid-thermal-ablation loosely coupling schematic.


### 2.5. Verification of Numerical Model for Thermochemical Nonequilibrium Flow

A cylinder model is considered for validating a numerical model for thermochemical nonequilibrium flow. Its experimental study was conducted in the pure air of DLR HEG [20]. The free stream velocity was 4776 m/s, the pressure was 687 Pa, the temperature was 694 K, and the wall temperature was 300 K. The initial mass fraction of N2, O2, N, O, and NO were 0.23, 0.77, 0, 0, and 0, respectively. The diameter of the cylinder model was 90 mm, and the length was 380 mm. The computational grids and boundary conditions of the model are shown in Figure 5.Figure 5
The computational grids.The calorimetric perfect gas model and the five-component (N2, O2, NO, O, and N) thermochemical nonequilibrium gas model are established, respectively. Figure 6 is the distribution of Mach number along a stationary line with different gas models. The detachment distance of a shock wave using the calorimetric perfect gas model is 15.1 mm, and using the thermochemical nonequilibrium model is 11.8 mm, which is more consistent with the experimental value of 11.9 mm. Figure 7 is the comparison contour between the numerical and experimental results. It can be seen that the numerical result using the thermochemical nonequilibrium gas model is in good agreement with the experimental result, which also verifies the reliability of the thermochemical nonequilibrium model.Figure 6
Distribution of Mach number along a stationary line.Figure 7
Comparison between the numerical and experimental results.
(a)
(b)


## 2.1. Governing Equations of Aerodynamic Flow

For the thermochemical nonequilibrium flow, the Navier-Stokes equations with a chemical source term are expressed as follows:(1)∂Q∂t+∂E∂x+∂F∂y−∂Ev∂x+∂Fv∂y=S,Q=ρ1⋮ρnsρuρvρEρeve,E=ρ1u⋮ρnsuρu2+pρuvρHuρeveu,F=ρ1v⋮ρnsvρuvρv2+pρHvρevev,S=ω̇1⋮ω̇ns000ω̇ve,Ev=qx1⋮qxnsτxxτxyuτxx+vτxy+qx+qvex+ρ∑i=1nsDimhi∂Yi∂xqvex+ρ∑i=molDimevei∂Yi∂x,Fv=qy1⋮qynsτyxτyyuτyx+vτyy+qy+qvey+ρ∑i=1nsDimhi∂Yi∂yqvex+ρ∑i=molDimevei∂Yi∂y,where Q is the conservative variable. E and F are inviscid fluxes in the x and y directions, respectively. Ev and Fv are viscous fluxes in the x and y directions, respectively. S is the source term reflecting the effect of flow field chemical nonequilibrium. u and v are the velocities in the x and y directions, respectively. p, ρ, and Yi are the pressure, density, and mass fraction of components of the gas, respectively. qxi and qyi are the diffusivity of components in the x and y directions, respectively. Dim and hi are the diffusion coefficient and the absolute enthalpy per unit mass of components, respectively. qx, qy, qvex, and qvey are translational-rotational heat fluxes and vibrational heat fluxes in the x and y directions, respectively. ω̇i and ω̇ve are mass generation rate and vibration energy, respectively.The Gupta chemical kinetic model is used to calculate the component change caused by chemical reactions. The thermodynamic nonequilibrium in this paper is characterized by the Park two-temperature model. The specific chemical reaction model is referred to in Reference [18].


## 2.2. Fluid-Solid Interface and Ablation Rate Model

The energy conservation between the flow field and the material model needs to be satisfied. Figure1 shows the energy transfer at the gas-solid interface, which satisfies the equation as follows:
(2)q=qconv+qrad‐in+ṁchcs−ṁhcs−qrad‐out,where q is the heat flux into the material model. qconv is the aerodynamic heat flux. qrad‐in is the radiation heat flux by an external flow field. ṁchcs is the heat flux generated by the chemical reaction on the material surface. ṁhcs is the heat flux carried away by the ablation loss. qrad‐out is the radiation heat flux dissipation of the material to a flow field.Figure 1
Energy transfer at the gas-solid interface.The finite rate ablation model is used to calculate the surface ablation of carbon/carbon composites. The mass loss caused by mechanical erosion and surface melting is not considered. Surface ablation of carbon/carbon composites mainly includes oxidation and nitridation of carbon. The surface chemical reaction mechanism and mass loss rate are as follows [16]:
(1)
The oxidation of carbon:O+Cs⟶COThe mass loss rate is(3)ṁ1=ρCOkTω2πmOβOMCMO.
(2)
The oxidation of carbon:O2+2Cs⟶2COThe mass loss rate is(4)ṁ2=2ρCO2kTω2πmO2βO2MCMO2.
(3)
The nitridation of carbon:N+Cs⟶CNThe mass loss rate is(5)ṁ3=ρCNkTω2πmNβNMCMN,where ρ is the density of gas. Ci is the mass fraction of the inflow components. k is the Boltzmann constant. Tw is the wall temperature. mi is the mass of components. βi is the surface chemical reaction efficiency. Mi is the molar mass of components.
(6)βO=0.63exp−1160/Tω,βO2=0.5,βN=0.3.The total mass loss rate of the material is(7)ṁ=ṁ1+ṁ2+ṁ3.


## 2.3. Material Thermal Response and Ablation Boundary Tracking

The heat conduction will affect the temperature distribution, oxidation properties, and ablation of a material. Therefore, the heat conduction of a material must be considered in order to predict the ablation response of the material. The heat conduction in the material follows the Fourier heat equation and energy conservation. The governing equation of heat conduction in the material can be written as follows:(8)ρscp∂T∂t−∂∂xiλ∂T∂xi=q,where t is the time. T is the temperature. ρs is the material density. cp is the specific heat capacity of the material. λ is the heat conductivity of the material. q is the heat flux.The ablation of carbon/carbon composites is a dynamic process, and the model shape changes constantly. In order to predict the ablation response, it is necessary to track the location of the ablation surface. The linear ablation rateṠ of the material can be calculated by
(9)ṁ1+ṁ2+ṁ3=ρsṠ.The movement of nodes on a material surface can be calculated by the ablation rate and time stepΔt. The moving direction of nodes is perpendicular to the surface, and the displacement of nodes is as follows.
(10)δxδy=ṠΔtnxṠΔtny.The Abaqus mesh motion algorithm is used to simulate ablation surface movement, and the user-defined subroutine umeshmotion is used to update the location of ablation surface nodes. The ALE (Arbitrary Lagrangian-Eulerian) mesh adaptive technology is used to remesh the internal grids of the model, avoiding the excessive deformity of the grids, and finally realizing the simulation of the ablation surface recession process.


## 2.4. Fluid-Thermal-Ablation Coupling Strategy

The coupling model of hypersonic flow and carbon/carbon composites is realized by a zoning method. According to the physical space, it is divided into fluid and solid domains, as shown in Figure2. The coupling between the fluid domain and the solid domain is both a physical and a chemical process in which hypersonic aerodynamic heating, heat conduction, and ablation interact through the coupling interface. Ωf and Ωs are the fluid domain and the solid domain, respectively. Γ is a coupled interface, where a solid provides wall temperature and displacement boundaries for a fluid, and the fluid provides aerothermal loads for the solid. The data of a flow field and the material model can be exchanged repeatedly on the coupling interface, and the data exchange between two domains is realized by Mpcci. Figure 3 shows the data exchange between mismatched grids.Figure 2
Fluid-solid interaction.Figure 3
Data exchange between mismatched grids.In the numerical calculation of material ablation in a hypersonic environment, different solvers are used for the solid and fluid regions. From the physical aspect, the time for hypersonic flow to reach steadiness is much smaller than the time for transient structural heat transfer and ablation. The time scale discontinuity reflected in the numerical calculation is that the fluid solver can allow an extremely much smaller time step size than the solid solver, which presents a computational challenge in practical engineering applications. Loosely coupling analysis strategy can avoid the detailed transient aerodynamic flow analysis that usually requires significant computational effort. The accuracy and efficiency of the coupling procedure mainly depend on the coupling time step size. In this paper, we refer to the research results of the literature [19], and the time step size is 0.01 s.The fluid-thermal-ablation coupling process is shown in Figure4. The specific process is as follows: (1) At the initial time t0, the wall temperature and position of the carbon/carbon model are transferred to the flow field through the coupling interface as the boundary conditions for flow field calculation. (2) Based on the boundary conditions, the flow field is calculated to obtain the wall heat flux. (3) The heat transfers from the flow field to the solid domain through the coupling interface. (4) Based on the boundary conditions, the heat conduction in the material model from time t0 to time t1 is calculated, and the temperature distribution and ablation of the material at time t1 are obtained. (5) The wall temperature and position at time t1 is transferred to the fluid domain as the boundary condition of the flow field. (6) At this point, the calculation of a coupling time step is completed, and then the next coupling time step is continued, and so on until the whole coupling time is over.Figure 4
Fluid-thermal-ablation loosely coupling schematic.


## 2.5. Verification of Numerical Model for Thermochemical Nonequilibrium Flow

A cylinder model is considered for validating a numerical model for thermochemical nonequilibrium flow. Its experimental study was conducted in the pure air of DLR HEG [20]. The free stream velocity was 4776 m/s, the pressure was 687 Pa, the temperature was 694 K, and the wall temperature was 300 K. The initial mass fraction of N2, O2, N, O, and NO were 0.23, 0.77, 0, 0, and 0, respectively. The diameter of the cylinder model was 90 mm, and the length was 380 mm. The computational grids and boundary conditions of the model are shown in Figure 5.Figure 5
The computational grids.The calorimetric perfect gas model and the five-component (N2, O2, NO, O, and N) thermochemical nonequilibrium gas model are established, respectively. Figure 6 is the distribution of Mach number along a stationary line with different gas models. The detachment distance of a shock wave using the calorimetric perfect gas model is 15.1 mm, and using the thermochemical nonequilibrium model is 11.8 mm, which is more consistent with the experimental value of 11.9 mm. Figure 7 is the comparison contour between the numerical and experimental results. It can be seen that the numerical result using the thermochemical nonequilibrium gas model is in good agreement with the experimental result, which also verifies the reliability of the thermochemical nonequilibrium model.Figure 6
Distribution of Mach number along a stationary line.Figure 7
Comparison between the numerical and experimental results.
(a)
(b)


## 3. Physical and Numerical Models

### 3.1. Carbon/Carbon Leading Edge

The ablation response of a carbon/carbon wedge-shaped leading edge is studied. The half-angle of the model is 10°, and the radius of the leading edge is 0.025 m. The geometry of the leading edge is shown in Figure8. The leading edge suffers a hypersonic flow of Mach number 8 at an altitude of 20 km. The static temperature is 217 K, and the static pressure is 5475 Pa. The mass fractions of N2 and O2 are 79% and 21%, respectively. The properties of carbon/carbon composites are shown in Table 1.Figure 8
Geometry of the leading edge.Table 1
Properties of carbon/carbon composites.
Density (kg/m3)
Specific heat (J/(kg·K))
Thermal conductivity (W/(m·K))
Radiation coefficient
Elastic modulus (Gpa)
Poisson’s ratio
1800
840
15
0.8
69
0.3


### 3.2. Fluid-Thermal-Ablation Coupling Model

Fluent is used to establish the numerical simulation model of a flow field. Half of the model is established as the leading edge is symmetrical. The grids near the wall in the flow field are refined to meet the requirement of heat flux simulation. Abaqus is used to establish the carbon/carbon composite model, and the number of grids is 2355. The simulation models of the flow field and the carbon/carbon leading edge are shown in Figures9 and 10, respectively. At the initial time, the temperature of the leading edge is 300 K, the time increment Δt is 0.01 s, and the total time of coupling is 30 s. The temperature distribution, surface ablation rate, and ablation shape of the material model are analysed.Figure 9
Grids of the flow field.Figure 10
Grids of the carbon/carbon leading edge.


### 3.3. Grid Independence Analysis

The quality of the flow field grids is the key factor affecting the prediction of aerodynamic heating. The grid independence study is made with three different groups. The groups of the flow grids are shown in Table2.Table 2
Different groups of grids.
Item
Grid spacing of the first layer near the wall (m)
Maximum grid spacing (m)
Total cell numbers
Coarse (unrefined)
0.001
0.005
1856
Medium
0.0005
0.0034
5251
Fine
0.0002
0.0019
10502Figure11 shows the change of heat fluxes at the stagnation point with different grids; it can be seen that the accuracy of the medium grids is sufficient.Figure 11
Stagnation heat fluxes with different grids.


## 3.1. Carbon/Carbon Leading Edge

The ablation response of a carbon/carbon wedge-shaped leading edge is studied. The half-angle of the model is 10°, and the radius of the leading edge is 0.025 m. The geometry of the leading edge is shown in Figure8. The leading edge suffers a hypersonic flow of Mach number 8 at an altitude of 20 km. The static temperature is 217 K, and the static pressure is 5475 Pa. The mass fractions of N2 and O2 are 79% and 21%, respectively. The properties of carbon/carbon composites are shown in Table 1.Figure 8
Geometry of the leading edge.Table 1
Properties of carbon/carbon composites.
Density (kg/m3)
Specific heat (J/(kg·K))
Thermal conductivity (W/(m·K))
Radiation coefficient
Elastic modulus (Gpa)
Poisson’s ratio
1800
840
15
0.8
69
0.3


## 3.2. Fluid-Thermal-Ablation Coupling Model

Fluent is used to establish the numerical simulation model of a flow field. Half of the model is established as the leading edge is symmetrical. The grids near the wall in the flow field are refined to meet the requirement of heat flux simulation. Abaqus is used to establish the carbon/carbon composite model, and the number of grids is 2355. The simulation models of the flow field and the carbon/carbon leading edge are shown in Figures9 and 10, respectively. At the initial time, the temperature of the leading edge is 300 K, the time increment Δt is 0.01 s, and the total time of coupling is 30 s. The temperature distribution, surface ablation rate, and ablation shape of the material model are analysed.Figure 9
Grids of the flow field.Figure 10
Grids of the carbon/carbon leading edge.


## 3.3. Grid Independence Analysis

The quality of the flow field grids is the key factor affecting the prediction of aerodynamic heating. The grid independence study is made with three different groups. The groups of the flow grids are shown in Table2.Table 2
Different groups of grids.
Item
Grid spacing of the first layer near the wall (m)
Maximum grid spacing (m)
Total cell numbers
Coarse (unrefined)
0.001
0.005
1856
Medium
0.0005
0.0034
5251
Fine
0.0002
0.0019
10502Figure11 shows the change of heat fluxes at the stagnation point with different grids; it can be seen that the accuracy of the medium grids is sufficient.Figure 11
Stagnation heat fluxes with different grids.


## 4. Results and Discussion

The distribution of the flow field and the carbon/carbon composite responses at different times are obtained. Figure12 is the distribution of the Mach number and temperature at 20 s. It can be seen that the shock wave is formed at the front of the leading edge, the Mach number near the stagnation point is close to 0, and the maximum temperature of the flow field at the stagnation point can reach 2860 K.Figure 12
Distribution of flow field at 20 s.
(a)
Mach number
(b)
TemperatureFigure13 is the ablation morphology of the carbon/carbon leading edge at 20 s compared with the initial state. It can be seen that the material model has an obvious recession, and the ablation in the stagnation area is more serious than the lateral area. This phenomenon can also be seen from the distribution of the recession depth along the surface in Figure 14. The main reason is that the temperature around the stagnation region is higher and the chemical reaction is more active, which leads to a faster ablation rate.Figure 13
Ablation morphology of the leading edge.Figure 14
Recession depth distribution along the surface.Figure15 is the heat flux distribution along the surface at different times. It can be seen that the maximum heat flux occurs in the stagnation region. As time goes on, the heat flux in the stagnation area decreases gradually. This is because the heat flux is related to the temperature gradient near the wall. The larger the temperature gradient is, the larger the heat flux is. At the initial time, the temperature of the material surface is low, and the temperature gradient near the surface is large, which leads to a large heat flux. With the heat transferring to the material model, the temperature difference between the material model and the flow field decreases gradually which leads to the decrease of heat flux.Figure 15
Heat flux distribution along the surface.Figure16 shows the temperature distribution along the surface. It can be seen that the temperature in the stagnation area rises rapidly as time goes on, while the temperature far from the stagnation area rises slowly. At the same time, the change of the wall heat flux in Figure 15 is verified.Figure 16
Temperature distribution along the surface.Figure17 is the temperature distribution of the carbon/carbon leading edge at different times. It can be seen that under the aerodynamic thermal load, the high temperature region of the material always concentrates near the stationary region, and the temperature gradient is large. As time goes on, heat gradually transfers to the interior of the model.Figure 17
Temperature distribution of the carbon/carbon leading edge at different times.
(a)
1 s
(b)
10 s
(c)
20 s
(d)
30 sFigure18 shows the change of the stagnation point temperature with time. As the aerodynamic heating of the hypersonic flow, the temperature at the stagnation point rises rapidly in the initial stage, then slows down gradually. Finally, the temperature at the stagnation point tends to be stable. This is due to the low temperature of the carbon/carbon composites and the high heat flux in the initial stage, which leads to the rapid increase of the temperature in the stagnation zone. As time goes on, the temperature of carbon/carbon composites increases gradually, resulting in the decrease of the heat flux in the stagnation area.Figure 18
The temperature at the stagnation point.Figure19 is the ablation morphology of the carbon/carbon leading edge at different times. As time goes on, the recession depth of the material model increases gradually, and the radius of the leading edge increases gradually. The front of the leading edge becomes blunt due to the ablation, and this leads to a reduction in the aerodynamic heating.Figure 19
Ablation morphology of the leading edge.Figure20 is the ablation rate at the stagnation point. In the initial stage, the ablation rate is small. As time goes on, the temperature of the material rises and the chemical reaction rate accelerates, which results in the rapid increase of the ablation rate of the material. The heat flux on the material surface is gradually stable, and the surface temperature tends to be stable, which leads to the stable ablation rate of the material.Figure 20
Ablation rate at the stagnation point.


## 5. Conclusions

In this paper, a new model for heat transfer and ablation response of carbon/carbon composites has been proposed. The conclusions are as follows:(1)
Considering flow field thermochemical nonequilibrium, heat transfer, and ablation of material, a fluid-thermal-ablation coupling model of a hypersonic flow field and carbon/carbon leading edge is established. The heat transfer and ablation response of the carbon/carbon model are predicted
(2)
In the initial stage, the heat flux in the stagnation area is the largest. As time goes on, the temperature of the carbon/carbon model increases gradually, the temperature gradient in the stagnation area decreases, and the heat flux decreases
(3)
The carbon/carbon ablation is serious in the stagnation area than in other regions. The shape of the carbon/carbon model changes after ablation, and the radius of the leading edge increasesThe fluid-thermal-ablation coupling model can be used to provide some reference for the design of thermal protection systems.



---

*Source: 9232684-2020-02-03.xml*