Competing effects of surface catalysis and ablation in hypersonic reentry aerothermodynamic environment

Zhiling CUI , Jin ZHAO ,,c, Guice YAO ,*, Jun ZHANG , Zhihui LI ,,Zhigong TANG , Dongsheng WEN ,c

a School of Aeronautical Science and Engineering, Beihang University, Beijing, 100083, China

b China Aerodynamics Research and Development Center, Mianyang, China

c Ningbo Institute of Technology - Beihang University, Ningbo, Zhejiang 315100, China

KEYWORDS Ablation fluid dynamics;Gas-solid reaction;Graphene;Hypersonic reentry;Molecular dynamics simulation;Recombination effects

Abstract Under hypersonic flow conditions, the complicated gas-graphene interactions including surface catalysis and surface ablation would occur concurrently and intervene together with the thermodynamic response induced by spacecraft reentry.In this work, the competing effects of surface heterogeneous catalytic recombination and ablation characteristics at elevated temperatures are investigated using the Reactive Molecular Dynamics (RMD) simulation method. A Gas-Surface Interaction (GSI) model is established to simulate the collisions of hyper-enthalpy atomic oxygen on graphene films in the temperature range of 500-2500 K.A critical temperature Tc around 900 K is identified to distinguish the graphene responses into two parts: at T ＜Tc, the heterogeneous surface catalysis dominates, while the surface ablation plays a leading role at T ＞Tc. Contradicting to the traditional Arrhenius expression that the recombination coefficient increases with the increase of surface temperature,the value is found to be relatively uniform at T ＜Tc but declines sharply as the surface temperature increases further due to the competing ablation effect. The occurrence of surface ablation decreases the amounts of active sites on the graphene surface for oxygen adsorption, leading to reduced recombination coefficient from both Langmuir-Hinshelwood(L-H) and Eley-Rideal (E-R) mechanisms. It suggests that the traditional Computational Fluid Dynamics (CFD) simulation method, which relies on the Arrhenius-type catalysis model, would result in large discrepancies in predicting aerodynamic heat for carbon-based materials during reentry into strong aerodynamic thermal environment.

1. Introduction

Large-scale spacecraft in orbit faces the problems of deorbiting at the end of life and disintegrates due to tremendous aerodynamic/thermal environment during reentering the earth.1-3The high-temperature thermo-chemical non-equilibrium gas flows produced by spacecraft reentering dense atmospheric surrounding will create cumulative effect of ablative pyrolysis,retreating and deformation, and even thermal ablation damage.4,5Therefore, Thermal Protection System (TPS) has been developed as one of the key techniques for the design of hypersonic vehicles due to the severe heating environment encountering during the hypersonic flight.6-10Ablative TPS materials are used as sacrificial materials to mitigate the incoming heat through phase changes and mass loss, whereas for those non-ablative TPS materials on reusable hypersonic vehicles, catalytic surface heat transfer that occurs at the gassurface interface has been considered as a non-negligible phenomenon for the large amounts of heat flux on the TPS surface.11,12The recombination of dissociated atomic oxygen(O) or atomic nitrogen (N) due to the high temperature effect on their gaseous compounds is a typical surface catalytic exothermic reaction, where the TPS surface acts as the catalyst.13-17Generally, a weak catalytic surface could cause a large decrement in surface heat flux as it prevents the recombination of the dissociated atoms on the surface; hence reduced exothermic heat releases. In most of the cases, the surface catalytic recombination and ablation not only occur simultaneously, but also affect each other. Understanding the dual effects of catalytic recombination and ablation is essential for accurate heat transfer prediction for hypersonic reentry flow.

Surface catalytic and ablative characteristics of various TPS materials have been investigated by researchers with great interests for improving both performance and reliability of future hypersonic aircraft.18,19It has been found that the erosion process modeling and simulation of dissociated gases impinging on TPS materials are significant, which can cause severe failures for reentry vehicles in serious aerothermodynamic environments.20-23A catalytic recombination coefficient γ,which is defined as the ratio of the recombined atomic flux to the impinging atomic flux,24,25is usually used to quantify the effect of surface catalysis. Traditional numerical strategies,e.g. Computational Fluid Dynamics (CFD) models, in hypersonic flow always employ a constant coefficient assumption to predict the aerodynamic heating, termed as non-catalytic (i.e.,γ=0),fully-catalytic(i.e.,γ=1),or a finite catalytic assumption,which results in large difference in the estimation of heat flux. It has been shown that the predicted aerodynamic heat could differ in 3-4 times in the stagnation point under hypersonic reentry flow conditions by using either a non-catalytic or fully-catalytic assumption.26Several experimental studies have shown that the recombination coefficient is affected strongly by the wall temperature.However,the results vary among different research groups. The uncertainties and discrepancies of various experimental facilities, the nature of applied testing techniques and several other unavoidable factors, such as different composition and microstructure of the TPS surface,surface roughness, gas composition and so forth,14,15,27-29are all potentially accountable, as well as the neglection of possible surface oxidation and ablation phenomena.

Recently, an atomistic-scale numerical technique using Reactive Force Field (ReaxFF) potential based on classical Molecular Dynamics (MD) simulation method, also known as Reactive Molecular Dynamics(RMD)method,has brought valuable nanoscale visualization and thermo-chemistry insight during the dissociated flow eroding process,30-33The heterogeneous catalysis recombination coefficients for silica surfaces have been successfully simulated by RMD by a few researchers.34-37Yang et al.36used MD methods to study the role of argon in the catalytic recombination of oxygen atoms on the silica surface, applying the potential energy surface and the atomic surface recombination model analysis. It was found that oxygen was more likely to form chemical bonds on top of silicon atoms, while argon had little effect on the oxygen atom recombination on the silica surface. MacKay et al.37combined MD and Monte Carlo methods to calculate the recombination rate of hydrogen atoms on the silica surface.By calculating the reaction rate through the Langmuir kinetic model, it was estimated that the recombination coefficient γ could be in a wide range of 10-4-100when the gas-phase radical density varied within 1012-1016cm-3and the temperature within 10-2000 K. Sˇimonka et al.38performed large-scale RMD simulations on different crystal faces of Silicon Carbide(SiC) to study the orientation dependence of the early oxidation process of SiC. Their results show that the oxidation rate of the C-face is the highest,followed by the m-,a-and Si-face.In addition,the difference of the growth rates among different surfaces decreased over time. Such differences in the adsorption capability associated with different surface characteristics clearly show their influences on the subsequent surface reactions.

Carbon-based nanomaterials, such as graphene, have been recently considered as one of the most promising materials in a wide range of applications because of its outstanding electronic, optical, thermal and mechanical properties.39-41Srinivasan et al.42simulated the etching process of high-enthalpy atomic oxygen hitting the diamond surface using RMD. The results showed that various functional groups such as ethers,peroxides, oxygen radicals and dioxanes could be formed on the surface, which was consistent with previous experiments and calculations based on the first principle. Poovathingal et al.43conducted a large-scale RMD study on high-enthalpy atomic oxygen oxidation of Highly Oriented Pyrolytic Graphite (HOPG) based on the ReaxFF force field. Their results predicted that oxygen coverage would occur first, followed by the carbon removal reaction, and finally wide and deep voids formed, similar to the experimental observation. The simulation also predicted that the most abundant product species was O2, followed by CO2, and CO was the least product.Cui et al.44used the RMD method to investigate the surface catalysis and ablation behaviors under the impingement of high-enthalpy atomic oxygen on graphene surfaces of seven different configurations. It revealed the importance of the surface configuration,and the graphene surface with multiple layers in a zigzag configuration showed the best resistance to the surface ablation. In addition, it also showed the adsorption of oxygen atoms and the generation of oxygen molecules first during the ablation process, indicating that the surface catalysis and ablation processes were intervened together. But in terms of heterogeneous catalytic recombination and microablation effect for graphene materials, unfortunately, this is still not well understood. The mechanism insight analyses within a large temperature range to evaluate various chemical behaviors and reaction rates quantitatively are still limited.

Different from silica-based TPS, both the effects of surface catalysis and ablation are present at elevated temperatures for carbon-based materials.As the temperature increases,it would be expected that impinging atoms on the surface in a highenthalpy flow could participate in a few reactions such as surface adsorption, surface recombination, and surface oxidation or surface ablation,depending on the temperature and surface characteristics. The same material may exhibit different ablation mechanisms in different thermal environments, especially for short time and high heat flux environment such as hypersonic reentry, where all the reactions take place on the deformation and retreat of the surface, including pyrolysis,ablation and formation of liquid layer.21,45The revelation of the competing effects could shed light to those complicated phenomena including numerical forecast of disintegration of the end-of-life spacecraft during the process of falling reentry.46

As we have demonstrated from our previous work,44the adsorption of oxygen atoms and the formation of oxygen molecules affect significantly the subsequent ablation process.In most of the experimental or numerical studies as reviewed above,however,the investigations have been focused on either surface ablation under high temperatures where the catalysis reactions have been neglected or contained inside the ablation effect, or surface catalysis without the influence of surface ablation. In order to reveal the intervened mechanisms of surface catalysis and ablation, and their competing effects under hypersonic reentry flow conditions, this work will examine detailed surface reactions under a range of temperatures.RMD simulation approach is used via a Gas-Solid Interface(GSI)model to investigate the heterogeneous chemical kinetics of dissociated oxygen impinging on graphene surfaces. The fundamental aero-thermochemical characteristics of graphene materials are presented, with a focus on the recombination coefficients of heterogeneous oxygen atom in the presence of different surface reaction kinetics. A microscopic view of the competing effects of surface catalysis and ablation is revealed.

2. Reactive molecular dynamics simulation details

Molecular dynamics simulations based on the ReaxFF reactive force field,47known as RMD, were carried out to investigate the successive hyperthermal oxygen atom collisions on a single-layer grapheme.48The system energy Esysin ReaxFF is calculated according to the sum of multiple energy terms:

The valence angle energy Evaland torsion angle energy Etorsare bond order functions, and their energies tend to be zero when the bond breaks. The non-bond interactions of the two atoms, including the intermolecular van der Waals potential EvdWaalsand electrostatic interaction potential ECoulomb, are obtained by calculating the forces between atomic pairs. The energy compensation term of the lone pair electron, such as Elp, over-coordination energy Eover, under-coordination energy Eunder,compensation function Epenand four-body conjugate energy Econj, can be used to describe the breaking and formation of chemical bonds in different chemical environments. The force field parameter ReaxFFC/H/Odeveloped by Chenoweth et al.49is used in this work. The parameter is parameterized according to the training set47including atomic charge, bond length, valence and twist energy, heat of formation and reaction energy of various hydrocarbons.

The GSI model was built up, as shown in Fig. 1, which involved the impinging atomic oxygen as the gas phase and a single-layer graphene as the solid slab. For the solid phase,the single-layer graphene structure, containing 864 carbon atoms in total, had a periodical cell dimension of 51.16 A°× 44.30 A° in the x-y plane, and the direction normal to the single-layer graphene (z direction) had a cell dimension of 240 A° . For the gas phase, the impingement of oxygen atoms was initiated 5 A° above the surface while in-plane coordinates of the oxygen atom were chosen randomly. The translational energy of each oxygen atoms was given as 5.0 eV with the direction normal to the surface.48Each oxygen atom trajectory was initiated every 1.0 ps,and the total simulation lasted for a maximum duration of 500 ps.

In all the simulations, the periodical boundary conditions were applied in both × and y directions, while reflective boundary conditions were imposed in the z direction. The established GSI systems were initially energy minimized, followed by the temperature equilibrium with Number-Volumn-Temperature (NVT) canonical ensemble using a Berendsen thermostat. Aimed at investigating the wall temperature effect on the graphene response behaviors,the graphene temperature in the proposed GSI model was equilibrated to eight objective surface temperatures T at 500,700,900,1100,1300,1500,2000 and 2500 K, respectively. The simulation setup has been carefully examined and validated,where more details can be found in our previous work.34

3. Results and discussion

3.1. Validation

Beyond the careful examination of the applied force field and simulation setup as in our previous work,34the temperature effect on the physical structure of graphene sheet after equilibrium is analyzed in this section as an additional validation study.

The Radial Distribution Function (RDF) gC-C(r) profiles for graphene sheet are characterized to investigate the spatial correlations of carbon atoms for graphene films under various equilibrated temperatures, as depicted in Fig. 1. It can be observed that both the peak positions and their corresponding RDF values are comparable with published experimental and simulation result.50,51No significant variation of the gC-C(r)curves for graphene is found with the increase of the surface temperature from 500 K to 2000 K. Moreover, the first three typical peaks that are located in the short range at around 1.44, 2.45 and 2.80 A° , respectively, can be well explained by the length of the carbon-carbon(C-C)bond and the diagonal distance in the same carbon hexagon. Then, the following three peaks at around 3.75,4.24 and 4.99 A° indicate the strong spatial correlation between two adjacent carbon hexagons.The first valley can be found near the radial distance of 1.85 A° ,which indicates the length of the C-C bond broken.52Based on the equilibrium structure of graphene film at various temperatures, the temperature-dependent chemical kinetics process of graphene surface by the collisions of atomic oxygen is investigated in the following sections.

3.2. Interface evolution of graphene film by hyper-enthalpy collisions

After ensuring that the equilibrated surface temperatures of graphene film are achieved as shown in Fig. 2(a), the continuous bombardment processes of atomic oxygen are analyzed.The dynamic evolution of both solid phase and gas-phase composition at various surface temperatures are analyzed in this section.

The oxidation of graphene film during the hyper-enthalpy oxygen collisions can be observed for all the cases, as illustrated by the oxygen density contours at the gas-solid interface at T = 1100 K in Fig. 2(b). It can be found that partial continuous impinging oxygen atoms favor the multi-coordinate adsorption sites with the evolution of time,saturating the maximum number of active adsorption sites on the graphene surface after around 300 ps. Random oxygenated functional groups on its basal planes indicate the oxidation of graphene surface to form Graphene Oxide (GO). The adsorbed oxygen atoms are believed to participate in the further surface catalysis and surface ablation reactions of graphene. Though it has similar morphology to graphene oxidation occurrence, it is found that the chemical kinetics and transport phenomenon vary significantly with the surface temperature of graphene film through analyzing the gasified composition during the oxygen collisions.

Taking four cases with the graphene surface temperatures of 500, 700, 2000 and 2500 K as illustrations, the quantitative calculation of each gas-phase composition evolution can be found in Fig. 3 as analysis here. (A) For cases with surface temperature T = 500 and 700 K, the number of O2molecules dominates as the majority among all the gasified components after hyper-enthalpy collisions, and increases continuously with time, while the CO2and CO molecules are also found to be generated with a relatively low amount during the collision process,as presented in Figs.3(a)and(b).(B)In contrast,as the graphene surface temperature increases to 2000 and 2500 K, the increasing rate of CO2and CO yields accelerates apparently, while the yield rate of O2component decreases,as shown in Figs.3(c)and(d),revealing that the ablative reaction on the graphene surface gradually intensifies, while the surface catalytic effect gradually weakens. The above difference of gasified composition with the variation of graphene surface temperatures indicates that,beyond the oxidation reaction on the graphene surface, both surface catalysis and surface ablation events concurrently occur at the interface during the hyper-enthalpy collisions. In essence, the concurrent catalytic and ablative reaction rates not only influence one another,but also are temperature-dependent significantly.To deepen our understanding of the chemical kinetics mechanisms, the detailed surface ablation and heterogeneous catalysis behaviors depending on surface temperature are presented in Section 3.3 and 3.4, respectively.

3.3. Temperature-dependent heterogeneous catalysis recombination analysis

The number of molecular oxygens O2recombining from the dissociated atomic oxygen with the variation of graphene surface temperature are calculated as depicted in Fig.4,where the graphene surface acts as a catalyst to this exothermic heterogeneous surface catalytic recombination reaction:

It can be found that the surface catalytic effect of the graphene surface is dependent on its surface temperature obviously. As the increase of graphene temperature, the recombined amount of oxygen molecules is almost monotonously decreasing.

As illustrated and visualized with the snapshots through RMD simulations in Fig.4,the combination of oxygen atoms into oxygen molecules can be divided into the following two steps. First, oxygen atoms adsorb on the surface of graphene to form oxygen-containing functional groups;second,the oxygen atoms combine to form oxygen molecules through the Eley-Rideal(E-R)or the Langmuir-Hinshelwood(L-H)mechanisms under the catalysis of the graphene surface. Some O2molecules are formed by one adsorbed oxygen atom on the graphene surface recombining with another atomic oxygen from the incoming particle flux due to the collisions, known as the E-R mechanism. There are also some O2molecules formed by the recombination of two adsorbed atomic oxygen binding onto the graphene surface due to the diffusion,known as the L-H mechanism.

With the obtaining of recombined oxygen molecule number, the evolution of catalytic recombination coefficient as a function of time can be correspondingly calculated as shown in Fig. 5. The catalytic recombination coefficient γO2is calculated by the following formula:

It can be seen that the heterogeneous catalytic recombination coefficients for all the cases are becoming equilibrium after around 300 ps,which is in accordance with the saturation point of oxygen adsorption on the active sites of graphene surfaces. However, differing from the traditional Arrhenius expression, which shows an increasing trend of γO2as a function of surface temperature,15the obtained catalytic recombination coefficient γO2for graphene shows a different tendency. The equilibrated γO2decreases with the increase of graphene surface temperature. This is due to the concurrent competing surface ablation reaction occurring when graphene surface temperature increases from 500 K to 2500 K, as illustrated in Fig. 5 with the side-view snapshots of graphene surface morphology. The participation of more impinging atoms at elevated temperature leads to reduced recombination coefficient.

To address the effect of concurrent competing surface ablation on the heterogeneous catalysis with the variation of surface temperature, the calculated values of lnγ, i.e., by averaging the catalytic recombination coefficient between 450 and 500 ps from Fig. 5, as a function of the temperature inverse for graphene film are shown in Fig. 6. Together with Fig. 5, it can be observed clearly that when the temperature increases from 500 K to 900 K, the average catalytic recombination coefficient varies in a limited extent around 0.28, while the heterogeneous catalytic recombination coefficient declines with the increasing temperature beyond 900 K, as a critical point, till 2500 K. This is because the surface ablation rate of graphene is significantly accelerated as the temperature increases continuously.Obvious voids can be intuitively visualized through the top-view snapshots of graphene surface with the evolution of time in Fig. 6. The increased surface ablation would result in two issues:(A)reduce the maximum number of active sites on graphene surface for oxygen adsorption before heterogeneous catalysis; (B) consume some of impinging atomic oxygen due to the combination with carbon atoms on the graphene surface, leading to the production of CO and CO2molecules. Therefore, the number of recombined oxygen molecules declines,resulting in a decreasing trend of γO2at elevated surface temperature.

This phenomenon implies that the traditional Arrhenius fitting53to predict the finite-rate catalytic recombination coefficient is conditional and highly material dependent. It may be suitable for silica-based TPS,but for graphene-based materials where there is a competing surface ablation phenomenon occurring simultaneously, the traditional Arrhenius linear fitting becomes invalid. The neglection of such competing reactions shall result in big mistake in predicting the exothermic heat associated with the recombination. The detailed concurrent surface ablation effect is presented in the next section to reveal its interrelations with the surface catalysis.

3.4. Temperature-dependent surface ablation characteristic of graphene

It can be calculated that for surface temperature increasing from 500 to 900 K,the ablation percentage of graphene surface varies a little around 0.015, while when the surface temperature continuously ascends from 900 K, the ablation reaction rate accelerates significantly. For the case with graphene surface temperature at 2500 K, the surface ablation percentage of graphene is around 0.17 at the end of oxygen bombardment calculations. This intensified surface ablation phenomenon is directly related to the decline of heterogeneous catalytic recombination coefficient profiles as the graphene surface temperature increases. The surface ablation rates kTof graphene can be further quantified from the previous carbon atom (C)lost ratio analysis with the following expression:

As the description of Arrhenius fitting law,this surface rate constant kTcan be written as a function of temperature given by

The calculated values of lnkTas a function of temperature inverse are presented in Fig. 8 and compared with lnγ profiles to observe the relationship between the surface ablation rate kTand heterogeneous catalytic recombination coefficient γTin the complex chemical kinetics reaction process. The linear fitting formulas of ln kTand lnγTas a function of temperature inverse are as follows:

The presence of competing reactions between surface ablation and heterogeneous catalysis can be identified clearly. The top views of graphene surface at different temperatures at the end of oxygen bombardment are also qualitatively compared in Fig. 8. A critical temperature around Tc= 900 °C can be regarded as a demarcation of the dominance of surface catalysis and surface ablation.

Due to the surface ablation reaction process, the carbon atoms on the graphene surface are combined with oxygen atoms to form CO and CO2molecules, leaving voids on the graphene surface.The variation of CO and CO2gasified molecules with the evolution of time at different surface temperatures is shown in Fig. 9, associated with the snapshot illustrations of detailed CO and CO2formation process:

The exact kinetic rates of above two chemical reactions can be further examined by analyzing the generation rates of CO2and CO molecules during the surface ablation process, as shown in Fig. 10. The linear fitting formulas of lnkCO2and ln kCOas a function of temperature inverse are as follows

It can be observed that when surface temperature is below 900 K,the majority of gasified composition during the surface ablation reaction is CO2molecules, while CO becomes the majority when T is greater than 900 K.This is also related with the amount of O2caused by the heterogeneous catalysis reaction rates.Higher heterogeneous catalysis reaction rate when T is below 900 K can contribute to more recombined oxygen molecules, promoting the possible concurrent oxidation reaction: C + O2→CO2.

When the messengers came to fetch her she was terribly frightened, for she knew that it was her wicked stepmother who in this way was aiming at her life

4. Conclusions

To advance the understanding of graphene-based thermal environment materials,this work conducted a detailed reactive molecular dynamics simulation to investigate the complicated surface reactions during atomic oxygen impingement onto a graphene surface at elevated temperatures, with the focus on the competing effects of surface heterogeneous catalytic recombination and ablation involved in hypersonic reentry aerothermodynamic environment. A gas-surface interaction model was established to simulate the collisions of hyperenthalpy atomic oxygen on the graphene film in the temperature range of 500-2500 K. The results can be summarized as follows:

(1) A strong nonlinear temperature-dependent surface chemical response of graphene film is observed due to the atomic oxygen collisions, and the surface temperature influences significantly the reaction types and kinetics.

(2) A critical surface temperature Tcaround 900 K is identified to distinguish the dominant phenomenon of surface reactions: at T ＜Tc, the heterogeneous surface catalysis dominates, while the surface ablation plays a leading role at T ＞Tc.

(3) Contradicting to the traditional Arrhenius expression that the recombination coefficient increases with the increase of surface temperature, the value is found to be relatively uniform at T ＜Tcbut declines sharply as the surface temperature increases further.Traditional CFD simulation method based on Arrhenius expression to predict the aerodynamic heat for carbon-based materials may result in large discrepancies.

(4) In contrast,the kinetic rate of surface ablation shows an opposite effect to that of the surface recombination coefficient, revealing a strong competition between surface ablation and heterogeneous catalysis. Specially, as the ablation time goes on, irregular network deformation appears on the surface of the material, and the ablation surface retreats so that finally the structure materials are destructed.

(5) The occurrence of surface ablation decreases the amounts of active sites on the graphene surface for oxygen adsorption, leading to reduced recombination coefficient from both L-H and E-R mechanisms.

(6) A constitutive chemical non-equilibrium, finite-rate reactive dynamics model to consider the competing oxidation, heterogeneous surface catalysis and ablation behaviors concurrently is essential for high-accuracy CFD solutions of aerodynamics heating in hypersonic reentry flows, which provides an effective way to simulate microscopic defects evolution, molecular bond breaking and failure during the reentry process.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgements

This research was supported by the project of the Manned Space Engineering Technology (No. ZS2020103001), the National Natural Science Foundation of China (No.52006004),National Numerical Wind Tunnel Project of China(No.NNW2018-ZT3A05)and the Open Fund of Key Laboratory of Icing and Anti/De-icing (No. IADL20190102).