Acessibilidade / Reportar erro

One-Dimensional Thermal Analysis Model for Charring Ablative Materials

ABSTRACT:

This paper presents a one-dimensional model for the analysis of the charring ablative materials used in spacecraft thermal protection systems. The numerical method is based on an implicit finite difference formulation of the governing equations written for a system of mobile coordinates that accounts for the possible presence of surface recession. The maximum allowable operating temperature for the adhesive layer of the junction between the heat shield and the substructure is used as a design parameter for determining the minimum heat shield thickness. A case study on the re-entry of the Stardust capsule is presented. The model proposed as a useful dimensioning tool for the preliminary design phase of the heat shields of spacecraft entering the atmosphere. The model was validated through a survey of the literature related to the dimensioning of thermal shields, but based on numeric programs of highly representative industrial standards.

KEYWORDS:
Thermal protection system; Ablative materials; Thermal analysis

INTRODUCTION

The term ablation refers to the process of removing a material surface through vaporization, chemical reactions, and/or erosion. In the aerospace field, ablative materials are mainly used in the manufacture of heat shields. This type of thermal protection mechanism dissipates the high entering heat fluxes and the corresponding thermal loads via a phase change in the material, which results in the loss of ablative material mass. The use of ablative materials has been a classical approach in the design of Thermal Protection Systems (TPSs) for over 60 years in a wide range of applications. To date, all NASA planetary probes have been equipped with ablative heat shields. From a phenomenological perspective, ablative materials can be divided into the following two classes: “non-charring”, in which no chemical reaction occurs inside the material (e.g., carbon-carbon or silica), and “charring”, in which the material pyrolyzes (e.g., phenolic compounds, plastic resins, and ablative structural ceramics). The charring materials are the most commonly used materials for atmospheric re-entry. Most ablative materials consist of a fibrous composite with organic resins as binders. Although each type of material shows specific behaviors, a common phenomenology can be highlighted. Overall, the ablation process involves a wide range of physical and chemical phenomena, most of which mutually interact. During the atmospheric hypersonic flight of a spacecraft, a bow shock forms, increasing the temperature near the surface of the vehicle and creating an interaction between the bow shock and the boundary layer. The viscous flow in the boundary layer, in turn, increases the wall temperature, and the heat is transferred to the heat shield by energized particles through radiation and convection. The heat is thus transmitted by conduction from the outer surface of the shield to the entire underlying coating layer. Then, when the virgin ablative material has undergone adequate heating, a change of state begins. In a charring ablative material, the heated resin undergoes a decomposition known as pyrolysis, which generates gaseous products, primarily hydrocarbons. This resin pyrolysis also produces a carbonized porous residue - the char - that settles on the reinforcing fiber of the composite. This process is usually endothermic. The developed pyrolysis gases press towards the underlying impermeable ablative virgin material, generating a pressure increase that then pushes these gases through the porous structure of the char until their exit from the free front with injection into the boundary layer. During this transit, the gas mixture can undergo further chemical reactions, such as cracking processes, resulting in the formation of smaller molecules and reactions that generate carbon residues. These residues settle in the preexisting char, thereby reducing the porosity and modifying the thermal conductivity of the charred layer. All these reactions must be strongly endothermic to achieve effective thermal protection and heat removal. The gases that permeate through the char zone remove additional energy by convection, attenuating the conduction of heat to the underlying reaction zone. The pyrolysis gases, once they have reached the surface of the material, are injected into the boundary layer. This outgoing heat flux results in a significant beneficial blocking effect on convective transfer from the free current to the body. In addition, the pyrolysis gases can chemically react with the gases of the boundary layer, resulting in a net heating effect on the surface. In addition, to a lesser degree, the composition of the ablation products can change the amount of radiative heating and its spectral distribution. Moreover, the chemical reactions between the material surface and the chemical species present in the boundary layer may deteriorate the surface material, leading to surface recession, and further chemical reactions can affect the outer surface of the char. For example, the carbon residue can undergo oxidation processes (exothermic) via the surrounding fluid in atmospheres with suitable compositions. Clearly, the interactions of the ablative materials and their products with the surrounding environmental gases are much more complex than described here, as are all the mechanisms involved in determining the final entering heat flux. Ultimately, ablation functions as a thermal barrier by dissipating part of the entering thermal energy through “sacrificing” the coating material; in this manner, the ablated material absorbs a considerable amount of heat while creating a barrier effect through the flow of outgoing gases generated during the chemical and physical transformations of the material. Finally, mechanical erosion phenomena that are capable of causing substantial removal of the thermal protection material occur consistently. In addition, in the case of a melting material, the melting of the surface layer influences the whole heat transfer process. Consequently, all these phenomena reduce the initial thickness of the ablative layer. This surface recession involves a dual effect: the influence of heat transfer to the substrate and the modification of the body shape, causing variations in the aerothermal field. Figure 1 presents a schematic of a charring ablative TPS with the exchanges and phenomena involved.

Figure 1
Schematic of charring ablative TPS.

Understanding the set of physical and chemical phenomena previously described is essential for developing an appropriate model simulation of the ablation process, the complexity of which is evident from the foregoing discussion. Furthermore, many factors involved in this process - for example, determining the material properties and the actual operating conditions - are sources of variability and are affected by estimation or measurement uncertainties. Alternatively, ablative materials can be classified based on the behavior of their surface layer into melting and non-melting types (Nathan and Bindu 2005Nathan RP, Bindu S (2005) Low temperature ablative heat shield for re-entry vehicles, AIAA Paper 2005-5063.). In the first type, which consists of thermoplastics, the surface liquid layer resulting from the fusion of the material is removed immediately after forming, resulting in the surface exposed to the heat flux always being new. This behavior generally exhibits poor efficiency. The second type - non-melting materials - is also divided into two subtypes: high-temperature ablators (HTAs) and low-temperature ablators (LTAs). Examples of HTAs include carbon-carbon and carbon-silicon carbide (ceramic matrix composites). Three-dimensional carbon-carbon composites are mainly used in the nose tip and any leading edges of spacecraft wings. Because of their strong mechanical characteristics, which are preserved at high temperatures, carbon-carbon composites ensure protection over long periods. These composites oxidize at temperatures above 1,100 K. In contrast, in LTAs, chemical ablation is preceded by mechanical ablation; in general, considerable deterioration of the mechanical properties occurs with increasing temperature. In this case, thermal ablation (sublimation) becomes appreciable above 3,000 K. LTA materials are usually thermosetting plastics that undergo a carbonization process. Typically, LTA materials are used for re-entry missions, spacecraft, or ballistic missiles, which experience high heat fluxes for a short duration (< 5,700 W/cm2). The plastics are reinforced with fibers (e.g., carbon epoxy or carbon phenolic) oriented along preferential directions. Phenolic materials, polyamides, and polybenzimidazole are considered suitable because of their ability to form char and their high pyrolysis heat values.

CHRONOLOGICAL SURVEY AND STATE OF THE ART OF NUMERICAL MODELS FOR ABLATIVE TPSs

Since the beginning of space research, the treatment of the thermal problem related to ablative materials has proved to be crucial. Many models of increasing complexity and representativity have been introduced, starting with the Landau (1950)Landau HG (1950) Heat conduction in a melting solid, quarterly of applied mathematics 8(1):81-94. one-dimensional model and followed by the models developed by Roberts (1958)Roberts L (1958) A theoretical study of stagnation-point ablation, NASA TN 4392., Moyer and Rindal (1958)Moyer CB, Rindal RA (1968) An analysis of the coupled chemically reacting boundary layer and charring ablator - Part II - Finite difference solution for the in-depth response of charring materials considering surface chemical and energy balances, NASA CR-1061., and Swann et al. (1965)Swann RT, Pittman CM, Smith JC (1965) One-dimensional numerical analysis of the transient response of thermal protection systems, NASA TN D-2976.. Multidimensional representations were introduced since 1970; examples include those of Curry (1974)Curry DM (1974) Two-dimensional analysis of heat and mass transfer in porous media using the strongly implicit procedure, NASA TN D-7608. and Pittman and Howser (1972)Pittman CM, Howser LM (1972) Numerical analysis of the transient response of an axisymmetric ablative char layer considering internal flow effects, NASA TN D-6895.. One important aspect is modeling the phenomena that occur on the free surface of the material via the interactions between the surface ablation and the surrounding fluid. Since 1992, significant progress has been made regarding the coupling problem of the ablative material response using solvers for the viscous shock layer under more or less wide hypotheses, such as the following: presence or absence of surface recession; ablation products modeled as single or multiple species; gaseous products that react or do not react with the external flow; presence or absence of induced turbulence in the boundary layer resulting from the injection of pyrolysis gases; and thermochemical non-equilibrium in the boundary layer. Despite the high representativity of the available models, the first results acquired using actual data, such as those obtained by comparing the information provided by the Galileo spacecraft on its entry into the Jovian atmosphere, showed substantial differences in the thickness of the actual char layer relative to the theoretical and experimental forecasts made before the mission. To date, the main conclusions that can be drawn from the experience of several decades of using ablative heat shields are as follows:

  • Despite the enormous amount of research in this area, the design methodology for ablative heat shields and, above all, the models on which these designs are based, still require considerable improvement.

  • A very limited number of materials have been used; thus, the development and testing of new materials are required.

  • The mass fraction of the TPS can reach very high values, so investigations aiming to optimize the heat shield mass are necessary.

  • The uncertainties inherent to this problem, such as uncertainty regarding the material properties under real flight conditions, may affect the accuracy of forecasts.

  • Following this introduction, the next section describes the model and governing equations, along with notes on discretization, nodal schemes and nodal equations. A case study and a summary are presented in the subsequent sections.

THE ABLATIVE THERMOPHYSICAL MODEL

The model used here is one-dimensional and specific for charring materials. It is suitable for the conceptual development phase of a vehicle and its mission. At this phase, the analysis can be limited to the stagnation point, namely, the point with the maximum entering heat flux, which can be considered the dimensioning element. The TPS is modeled as a layer of ablative material glued on the underlying structure by a thin adhesive layer. The maximum operating temperature allowable for this layer, which ensures its sealing, is used as a design parameter for determining the minimum heat shield thickness. The ablative layer is subdivided into three distinct areas, depending on the conditions that occur during the running time: the mature char, the reaction zone and the virgin material. In the adopted model, the thickness of the bond-line is considered negligible compared to the rest of the material and thus does not substantially contribute to the thermal insulation.

MODEL ASSUMPTIONS

In the governing equations of the problem presented below, the following significant modeling assumptions are made:

  • The ablative material decomposes from the virgin state to porous char in a well-defined area called the reaction zone.

  • The thermal hysteresis phenomena involving the properties of the char are absent, particularly regarding the thermal conductivity.

  • The reaction zone is defined by two significant temperatures. The beginning of the pyrolysis reaction is identified by the temperature Tabl, and full carbonization occurs at the temperature Tchar. These temperature limits are defined for each type of material by evaluating the thermogravimetric experimental diagrams. The generation of the pyrolysis gases occurs exclusively in this area.

  • The gases generated in the reaction zone pass through the upper layers without loss of pressure and inject themselves into the overlying boundary layer. This step is considered instantaneous; therefore, the residence time of the gases in the char is null.

  • A local thermal equilibrium exists between the char and the gases that reside in its pores.

  • After their formation, the gases are not subject to further chemical reactions with the other materials or the environment.

  • The interface with the internal environment of the spacecraft can be adiabatic or radiative and/or convective, and the invariability of its temperature can be enforced.

  • The possible influence of the thermal stress on the material characteristics is neglected.

GOVERNING EQUATIONS

This section presents the set of assumed equations (Matting 1970Matting FW (1970) Analysis of charring ablation with description of associated computing program, NASA TN D-6085.; Clark 1973Clark RK (1973) An analysis of a charring ablator with thermal nonequilibrium, chemical kinetics, and mass transfer. NASA TN D-7180.; Chen and Milos 1999Chen Y-K, Milos FS (1999) Ablation and thermal analysis response program for spacecraft heatshield analysis. Journal of Spacecraft and Rockets. 36(3):475-483., 2000Chen Y-K, Milos FS (2000) Two-dimensional implicit thermal response and ablation program for charring materials on hypersonic space vehicles. AIAA Paper 2000-0206., 2005Chen Y-K, Milos FS (2005) Three-dimensional ablation and thermal response simulation system. AIAA Paper 2005-5064.). This set of equations is solved numerically, as described later. The general scheme adopted, for which the equations are written, is shown in Fig. 1.

The internal energy balance equation

The internal energy balance is expressed according to the following equation:

(1) ρ c p T t x = x k T x q ˙ R t + H d h ¯ ρ t y + S ˙ ρ c p T x t + m ˙ g h d x t

where ρ = density (kg/m3), cp = specific heat at constant pressure (J/(kg·K)), T = temperature (K), k = thermal conductivity (W(m·K)), R = internal radiative heat flux (W/m2), Hd = pyrolysis enthalpy (J/kg), = partial heat of charring (J/kg) - defined in Eq. 6, Ṡ = char recession rate (m/s), g = pyrolysis gas mass flow rate (kg/(m2·s)), y = fixed coordinate system (m), x = mobile coordinate system, y-S (m), t = time (s).

In Eq. 1, the mobile coordinate system x moves with the surface during recession, whereas the y coordinate system is fixed; at the initial instant, the two systems are coincident. From first to last, the individual terms of Eq. 1 represent the accumulation rate of thermal energy, the net rate of the thermal energy that is transferred by conduction and internal radiation, the energy consumption during pyrolysis, the convective rate caused by the movement of the mobile coordinate system, and the convective rate caused by the pyrolysis gases, respectively. The term corresponding to surface recession is an additional term arising from the adopted modeling approach. The transmission of heat through porous materials involves numerous energy transport mechanisms, including both conduction through solid and gaseous phases and radiation through the porous structure (internal radiative heat flux). Appropriately modeling the simultaneous conduction and internal radiation is numerically complex and requires the separate determination of many individual optical and thermal properties (Marschall et al. 2001Marschall J, Maddren J, Parks J (2001) Internal radiation transport and effective thermal conductivity of fibrous ceramic insulations, AIAA Paper 2001-2822.). Consequently, engineering investigations, as in the present case, are generally limited to pure conduction using an “effective” thermal conductivity that depends on temperature and pressure, assuming an internally opaque body. The local specific heat is a function of temperature for both the virgin ablative material and the material in the char state. In the reaction zone, for partially charred material (ρc < ρ < ρv ), the specific heat is expressed by the rule of mixtures, where the subscript “v” refers to the properties of the virgin state and the subscript “c” refers to the char. The expression is as follows:

(2) c p = τ c p v + 1 τ c p c

where the weight ratio τ defined by Eq. 3 is the mass fraction of virgin material with respect to the total of virgin material and char, which is used to obtain the correct local density:

(3) τ = 1 ρ c / ρ / 1 ρ c / p v

The thermal conductivity k, which is a function of temperature and pressure, is weighted in a similar manner:

(4) k = τ k v + 1 τ k c

According to the hypothesis that defines the reaction zone based on specific temperatures, the following linear relationship, for continuity, is assumed for the density:

(5) ρ = ρ v ρ c T T abl T abl T char

The enthalpy of pyrolysis gas hd is also a function of the temperature and pressure, while , the partial heat of charring, is defined as follows:

(6) h ¯ = p v h v ρ c h c ρ v ρ c

Internal decomposition equation

Regarding pyrolysis decomposition, a global model based on the Arrhenius relationship is adopted, and the coefficients are related to the complete plastic ablative composite, not to the individual components. The relation is expressed here through a formulation that is less common but more appropriate for introducing the data available in the literature (Williams and Curry 1992Williams SD, Curry DM (1992) Thermal protection materials - Thermophysical property data, NASA RP 1289.):

(7) ρ t = K ρ ρ c ρ v ρ c n e B / T

where K = collision frequency factor (kg/(m3·s)), n = decomposition reaction order and B = activation temperature (K).

Internal mass balance equation

The internal decomposition transforms part of the solid into pyrolysis gases, and because of the hypothesis of one-dimensional quasi-static flow and the impermeability of the interface with the virgin material zone, the mass flow of the pyrolysis gases is linked to this decomposition by the following relation:

(8) m ˙ g y = ρ t

Surface energy balance equation

The conditions at the ablating free surface are determined by the convective and radiative heat transfer and by the thermochemical interactions that result in surface exchange with the hot gases of the boundary layer. The surface energy balance equation is expressed by the following relationship:

(9) Q ˙ in = q ˙ c , blow + q ˙ rad + q ˙ comb F σ ε T w 4 T 4

where in = net total heat flux at the surface (W/m2), qc,blow = net hot wall convective heat flux (W/m2), rad = entering radiative heat flux (W/m2), comb = combustion heat flux (W/m2), F = exterior view factor, σ = Stefan-Boltzmann constant (W/(m2·K4)), ε = surface emissivity, TW = wall temperature, (K), T = freestream temperature (K).

In Eq. 9, the aerodynamic heating is split, as usual, into two parts intended to be handled in different manners. Here, this differentiation is even more important because the fraction of the convective origin can be significantly reduced by injecting pyrolysis gases into the boundary layer. Instead, this phenomenon generally has no appreciable effect on the contribution of the radiative origin. The consequences for the aerodynamic heating caused by this mass transfer have been extensively studied in the literature. To quantify this “blocking effect”, this work uses a second-order model, which is a function of the mass flow of the outgoing gases, known as “transpiration theory” (Swann et al. 1965Swann RT, Pittman CM, Smith JC (1965) One-dimensional numerical analysis of the transient response of thermal protection systems, NASA TN D-2976.).

For a spacecraft that moves with velocity V in an atmosphere with specific heat cp atm, the wall enthalpy is:

(10) h w = c p atm T w

and the total enthalpy is:

(11) h 0 = c p atm T + V 2 2

The following coefficient is defined:

(12) A = h 0 q ˙ c , w α c m ˙ c + α g m ˙ g

where c,w = cold wall convective heat flux (W/m2), c = char removal rate (kg/(m2·s)); and where the coefficients αc and αg are used to differentiate the molecular weight of gases in the boundary layer from that of the injected pyrolysis gases. The coefficient αc also takes into account the part of the char that is mechanically removed rather than sublimated.

The hot wall convective heat flux (i.e., the cold wall convective heat flux multiplied by the hot wall reduction factor) is defined as follows:

(13) q ˙ con = q ˙ c , w 1 h w h 0

The final expression assumed for the net hot wall convective heat flux, that is, inclusive of the effect of escaping gas block, is given by the following:

(14) A < 2 . 2 5 q ˙ c , b l o w = 1 0 . 7 2 4 A + 0 . 1 3 A 2 q ˙ con A 2 . 2 5 q ˙ c , b l o w = 0 . 0 4 q ˙ con

Equation 14 shows that a minimum value for the blocking effect occurs (in this case, the ratio between the two fluxes is set equal to 0.04) when the parameter assumes values above 2.25, i.e., when the second-order approximation begins to differentiate itself appreciably from more detailed theoretical models in practice.

Regarding the combustion of the ablation products in the boundary layer, complete combustion is assumed according to the following relationship:

(15) q ˙ comb = m ˙ c Δ H c

where ΔHc = heat of combustion per unit mass (J/kg).

Surface recession

Multiple investigations in the literature have evaluated the removal of char via chemical, thermal, or mechanical phenomena or a combination thereof; these works have led to the definition of a number of correlation relationships that depend, in general, on the specific material involved (see, e.g., Swann et al. 1965Swann RT, Pittman CM, Smith JC (1965) One-dimensional numerical analysis of the transient response of thermal protection systems, NASA TN D-2976.). In this paper, a set of experimental data concerning the rate of the surface recession S. of various materials is used; this data set can be described as either a function of the surface temperature or a function of the total entering heat flux. As a result of char removal, the surface moves relative to the fixed coordinate system. The distance between the starting position of the surface and its current position, i.e., the lost thickness or total surface recession, is given by:

(16) S = 0 t S ˙ dt

SCHEMATIC OF THE ABLATIVE THERMAL MODEL

The implementation phase of the numerical model for the treatment of TPS ablative materials has been previously addressed by the author (Mazzaracchio and Marchetti 2010aMazzaracchio A, Marchetti M (2010a) A probabilistic sizing tool and Monte Carlo analysis for entry vehicle ablative thermal protection systems. Acta Astronautica 66(5-6):821-835. doi: 10.1016/j.actaastro.2009.08.033
https://doi.org/10.1016/j.actaastro.2009...
) and has been applied in various works (Mazzaracchio and Marchetti 2010bMazzaracchio A, Marchetti M (2010b) Coupled aeroassisted orbital plane change manoeuvre and thermal protection system optimisation, IAC-10-C2.7.4, 61st International Astronautical Congress, September 27-October 1; Prague, Czech Republic., 2011Mazzaracchio A, Marchetti M (2011) Spacecraft aerodynamics and heat shield characteristics impact on optimal aeroassisted coplanar orbital transfer, IAC-11-C2.7.9, 62nd International Astronautical Congress, October 3-7; Cape Town, South Africa.; Mazzaracchio 2013Mazzaracchio A (2013) Heat shield mass minimization for an aerocapture mission to Neptune. International Journal of Aerospace Innovations 5(3-4):83-93. doi: 10.1260/1757-2258.5.3-4.83
https://doi.org/10.1260/1757-2258.5.3-4....
, 2015Mazzaracchio A (2015) Flight-path angle guidance for aerogravity-assist maneuvers on hyperbolic trajectories. Journal of Guidance, Control, and Dynamics 38(2):2015:238-248. doi: 10.2514/1.G000670
https://doi.org/10.2514/1.G000670...
, 2016Mazzaracchio A (2016) Coupled versus uncoupled optimal solution for thermal and dynamic problems in spacecraft atmospheric flight. World Journal of Engineering 13(1):53-60. doi: 10.1108/WJE-02-2016-006
https://doi.org/10.1108/WJE-02-2016-006...
). All the considerations, assumptions and schematizations reported here are reprocessed from these studies. The industrial standard numerical programs commonly considered as a reference, when available, are generally very laborious, mainly because of the large amount of input data required. Moreover, in general, these thermal analysis programs require a preliminary execution of other surface chemistry software to generate the part of the input data that represent the phenomena localized on the surface and to quantify the enthalpy of the pyrolysis gases. In turn, these chemistry simulation codes require large amounts of input data and knowledge of the composition of both materials and pyrolysis gases. In addition, codes with high representativity require the development of a thermal model that must be compared and calibrated with experimental data, resulting, for example, from test campaigns with arcjets on the ablative materials used. These data sets must be related to the surface chemistry under a wide range of appropriate conditions. Generally, the absence or poor definition of the items described above prevents obtaining high-fidelity codes and results in problems with the accuracy and convergence of the solution. Examples of such codes are: the “Aerotherm's Charring Material Thermal Response and Ablation Program” (CMA) of the Acurex Corporation and the NASA Ames Research Center's “Fully Implied Ablation and Thermal” (FIAT). Less complex codes and faster methods, like the “virtual ablation method” (Ko et al. 2007Ko L, Gong WL, Quinn RD (2007) Reentry thermal analysis of a generic crew exploration vehicle structure, NASA, TM-214607.), which are characterized, however, by significant representativity, are useful when a high-fidelity model is not available or not usable for some reason. Specifically, these methods can be utilized during the preliminary design phase, when a large number of studies should be quickly finalized. The above issue is just one of the reasons that led to the development and implementation of the numerical model presented here. The numerical method used here is based on an implicit finite difference formulation (backward time, centered space) of the governing equations described earlier. The assumption of such a scheme solves a system of equations at each time step. Therefore, from the perspective of code development and the computational effort required, this scheme is more expensive than a similar achievable explicit method (e.g., using the forward time differences). However, the enormous advantage is that in the thermal field, using parabolic partial differential equations, the adopted implicit scheme is unconditionally stable in time and distance (therefore, the magnitude of the time step is not limited by a convergence criterion), unlike an explicit approach. Figure 2 shows a section of the schematic adopted, from which the convention adopted to distinguish the different layers and various states and conditions of the material may be deduced, inter alia. The free surface may regress with the formation of a reaction zone and a carbonized layer. The ablative part of the TPS is originally composed of a layer of virgin ablative material superimposed on one or more layers of materials, eventually different from each other, which compose the substructure and can be spaced with gaps to allow circulation of the cooling fluid. In accordance with the possible presence of surface recession, the equations are written for a system of mobile coordinates for which the upper free face of the material constitutes the moving surface, as previously described. Following this assumption, the ablative material layer is divided into a predefined and fixed number of nodes spaced from each other by ΔX, the value of which thus depends on the current position of the upper face. Using this measure, the recession surface is modeled continuously, thereby eliminating the need for the cancellation or “condensation” of nodes. However, this approach leads to terms arising from the time derivative of the spatial step ΔX. The procedure for dimensioning the TPS starts with an estimated initial value for the thickness provided by the user; subsequently, the minimum thickness determination is performed via a “targeted” percentage reduction of the thickness, thus using an iterative scheme based on the secant method with tolerance on the resulting thickness. The choice of the secant method allows for smooth convergence through 10-15 iterations and, in this case, is preferable to a Newton method, which, although characterized by greater convergence speed, has difficulty estimating the derivative.

Figure 2
Schematic of the implemented thermal model.

The equations for each node i obtained by moving all terms with unknown temperatures T' (i.e., those at the time t + dt) to the left-hand side and those with known temperatures to the right-hand side are as follows:

(17) A i T i 1 + B i T i + C i T i + 1 = D i

Ultimately, at each iteration, the following tridiagonal system must be solved on N nodes of the discretization:

(18) B 1 T 1 + C 1 T 2 = D 1 A 2 T 1 + B 2 T 2 + C 2 T 3 = D 2 A 3 T 2 + B 3 T 3 + C 3 T 4 = D 3 A N T N 1 + B N T N = D N

Some notes on the realized model are as follows:

  • The reaction zone remains, even when a cooling phase begins, usually towards the end of the mission. Thus, the temperature drops below the temperature that indicates the start of ablation Tabl. In this context, the reaction zone can be described as an area where the carbonization is not complete but that no longer participates in the generation of pyrolysis gases. In the absence of experimental data for this hybrid situation, the properties of the virgin material were assigned. Conversely, the thermal capacity of the pyrolysis gases was evaluated over the entire thickness, which is more or less completely carbonized.

  • In the mature char zone, the formation of gases does not occur; thus, the following is valid: gi = gi + 1.

  • The continuity of the thermal properties of the ablative material as the state varies (mature char/reaction zone/ablative virgin material) allows calculation of their nodal values only once, unlike the interface nodes between different materials, for which the node must be considered as belonging to the first material and then to the next material.

  • In the case of missions in non-oxidizing atmospheres, such as Neptune, the entering heat flux resulting from oxidative phenomena involving leaking gases is annulled.

  • There is the possibility of imposing temperature invariability on the last node corresponding to the inner face cab side (boundary condition). In this case, the coefficients DN-1 and DN are recalculated; this implies that the order of the system listed in Eq. 18 should be lowered by a degree. Therefore, at each integration step, the temperature of the penultimate node is calculated directly, while the temperature of the last node is assigned to the constant value again.

  • The surface recession is subject to the condition that the surface temperature is greater than Tchar. This situation implies that recession occurs for only the mature char state. The values obtained in the simulation show that the beginning and end of the recession occur at temperatures much greater than Tchar; therefore, this condition is always satisfied.

  • The resulting various ΔX values, which were calculated starting from the thickness and number of nodes present in the different layers, must all be of the same order of magnitude. Therefore, it is recommended to choose “thickness/number nodes” ratios that are comparable. Otherwise, instability phenomena may occur because of discontinuities in the value of ΔX.

The hypotheses taken on the discretization and the schematizations assumed for each node and the corresponding nodal equations that comprise the system are illustrated in detail in the next subsection.

KEY ASSUMPTION ON THERMAL MODEL

This section describes some of the key assumptions concerning the discretization by finite differences applied in the thermal model.

Linearized radiation

The radiation phenomena, the terms of which depend on the 4th power of T'i, may be directly taken into account in an implicit method only after their linearization. For this purpose, the classic expression was adopted:

(19) T i 4 = T i + Δ T 4 = T i 4 1 + Δ T T i 4 T i 4 1 + 4 Δ T T i = 4 T i 3 T i 3 T i 4

Moving surface (surface recession)

Surface recession involves introducing additional terms to the energy balance equations. The effect of the time derivative on the thermal capacitive term is as follows:

(20) d dt Δ X ρ c p T = Δ X ρ c p dT dt + ρ c p T d Δ X dt

where

(21) d Δ X dt = d dt S init NP 1 = S ˙ NP 1

Radiation between two parallel flat surfaces

In the presence of a gap with radiation, the classic relations are adopted:

Form factor:

(22) F 12 = 1

Equivalent emissivity:

(23) ε eq = 1 ε 1 + 1 ε 2 1 1

Equivalent thermal conductivity (electrical analogy)

The adoption of the electrical analogy leads to the following expression for thermal conductivity (Fig. 3):

Figure 3
Equivalent thermal conductivity.

(24) k i 1 , i = 1 a Δ X k i 1 + b Δ X k i

if

(25) a = b = Δ X 2 k i 1 , i = 1 1 2 k i 1 + 1 2 k i

thus:

(26) X k T X = 1 1 2 k i 1 + 1 2 k i T i 1 T i Δ X 1 1 2 k i + 1 2 k i + 1 T i T i + 1 Δ X Δ X

Nodal schemes and nodal equations

Figures 4-16 show the schemes used (Curry 1965Curry DM (1965) An analysis of a charring ablation thermal protection system, NASA TN D-3150.) for each node, depending on the state and the conditions to which it may be subjected. These figures also show the contributions of each phenomenon considered in the composition of the corresponding coefficients of the various equations of Eq. 17. The contributions marked with a red arrow provide energy to the node, whereas contributions marked with a blue arrow remove energy.

Figure 4
Nodal scheme - first node: front surface.

Figure 5
Nodal scheme - state 1: mature char zone.

Figure 6
Nodal scheme - state 2: mature char/reaction zone interface.

Figure 7
Nodal scheme - state 3: reaction zone.

Figure 8
Nodal scheme - state 4: reaction zone/ablative virgin material interface.

Figure 9
Nodal scheme - state 5: ablative virgin material.

Figure 10
Nodal scheme - state 6: ablative virgin material/1st substructure material interface.

Figure 11
Nodal scheme - state 7: Jth substructure material.

Figure 12
Nodal scheme - state 8.1: Jth/Jth+1 substructure material without gap.

Figure 13
Nodal scheme - state 8.2: last node interface Jth substructure material with gap (radiation + convection).

Figure 14
Nodal scheme - state 8.3: first node interface Jth substructure material with gap (radiation + convection).

Figure 15
Nodal scheme - state 9.1: last node Jth substructure material with adiabatic surface.

Figure 16
Nodal scheme - state 9.2: last node Jth substructure material with radiation and convection/cabin interior.

Each equation of Eq. 17 that constitutes the system of equations in Eq. 18 belongs to one of 13 types, as described in Eqs. 27-39. These equations refer to the corresponding number of schematizations assumed and described in Figs. 4-16.

  • Front surface node equation:

    (27) m ˙ g 1 c p abl 1 + S ˙ ρ 1 c p 1 + ρ 1 c p 1 Δ X 1 2 Δ t + 1 Δ X 1 2 k 1 + Δ X 1 2 k 2 1 2 ρ 1 c p 1 S ˙ NP 1 T 1 + + m ˙ g abl 2 c p 2 + 1 Δ X 1 2 k 1 + Δ X 1 2 k 2 + ρ 2 c p 2 S ˙ NP 1 1 2 NP 1 T 2 = Q in + ρ ˙ g 2 Δ X 1 2 H d 2 ρ 1 c p 1 Δ X t 2 Δ t T 1

    where the entering heat flux is given by Eq. 9.

  • State 1 equation:

    (28) 1 Δ X 1 2 k i = 1 + Δ X 1 2 k i T i = 1 m ˙ g i c p abl i + ρ i c p i S ˙ NP i + 1 2 NP 1 + 1 Δ X 1 2 k i 1 + Δ X 1 2 k i + 1 Δ X 1 2 k i + Δ X 1 2 k i + 1 + ρ i c p i Δ X 1 Δ t ρ i c p i S ˙ NP 1 T i + m ˙ g i + 1 c p abl i + 1 + 1 Δ X 1 2 k i + Δ X 1 2 k i + 1 + ρ i + 1 c p i + 1 S ˙ NP i 1 2 NP 1 T i + 1 = ρ i c p i Δ X 1 Δ t T i

  • State 2 equation:

    (29) 1 Δ X 1 2 k t 1 + Δ X 1 2 k i T i 1 m ˙ g i c p abl i + ρ i c p i S ˙ NP i + 1 2 NP 1 + 1 Δ X 1 2 k t 1 + Δ X 1 2 k i + 1 Δ X 1 2 k i + Δ X 1 2 k i + 1 + ρ i c p i S ˙ NP 1 T i + m ˙ g i + 1 c p abl i + 1 + 1 Δ X 1 2 k i + Δ X 1 2 k i + 1 + ρ i + 1 c p i + 1 S ˙ NP i 1 2 NP 1 T i + 1 = ρ i c p i Δ X 1 Δ t T i + ρ ˙ g i + 1 Δ X 1 H d i

  • State 3 equation:

    (30) 1 Δ X 1 2 k i = 1 + Δ X 1 2 k i T i 1 m ˙ g i c p abl i + ρ i c p i S ˙ NP i + 1 2 NP 1 + 1 Δ X 1 2 k i = 1 + Δ X 1 2 k i + 1 Δ X 1 2 k i + Δ X 1 2 k i + 1 + p i c p i S ˙ NP 1 T i + m ˙ g i + 1 c p abl i + 1 + 1 Δ X 1 2 k i + Δ X 1 2 k i + 1 + ρ i + 1 c p i + 1 S ˙ NP i 1 2 NP 1 T i + 1 = ρ i c p i Δ X i Δ t T i ρ ˙ g i Δ X i H d i ρ ˙ g i + 1 Δ X 1 H d i + 1

  • State 4 equation:

    (31) 1 Δ X 1 2 k i 1 + Δ X 1 2 k i T i 1 m ˙ g i c p abl i + ρ i c p i S ˙ NP i + 1 2 NP 1 + 1 Δ X 1 2 k i 1 + Δ X 1 2 k i + 1 Δ X 1 2 k i + Δ X 1 2 k i + 1 + ρ i c p i Δ X 1 Δ t ρ i c p i S ˙ NP 1 T i + 1 Δ X 1 2 k i + Δ X 1 2 k i + 1 + ρ i + 1 c p i + 1 S ˙ NP i 1 2 NP 1 T i + 1 = ρ i c p i Δ X 1 Δ t T i m ˙ g i Δ X 1 H d i

  • State 5 equation:

    (32) 1 Δ X 1 2 k i = 1 + Δ X 1 2 k i T i 1 ρ i c p i S ˙ NP i + 1 2 NP 1 + 1 Δ X i 2 k i 1 + Δ X 1 2 k i + 1 Δ X 1 2 k i + Δ X 1 2 k t + 1 + ρ i c p i Δ X 1 Δ t ρ i c p i S ˙ NP 1 T i + 1 Δ X 1 2 k i + Δ X 1 2 k i + 1 + ρ i + 1 c p i + 1 S ˙ NP i 1 2 NP 1 T i + 1 = ρ i c pi Δ X 1 Δ t T i

  • State 6 equation:

    (33) 1 Δ X 1 2 k i = 1 , 1 + Δ X 1 2 k i , 1 T i 1 1 Δ X 1 2 k i 1 , 1 + Δ X 1 2 k i , 1 + 1 Δ X 2 2 k i , 2 + Δ X 2 2 k i + 1 , 2 + ρ i , 1 c p i , 1 Δ X 1 2 Δ t + ρ i , 2 c p i , 2 Δ X 2 2 Δ t T i + 1 Δ X 2 2 k i , 2 + Δ X 2 2 k i + 1 , 2 T i + 1 = ρ i , 1 c p i , 1 Δ X 1 2 Δ t + ρ i , 2 c p i , 2 Δ X 2 2 Δ t T i

  • State 7 equation:

    (34) 1 Δ X j 2 k i 1 , j + Δ X j 2 k i , j T i 1 1 Δ X j 2 k i 1 , j + Δ X j 2 k i , j + 1 Δ X j 2 k i , j + Δ X j 2 k i + 1 , j + ρ i , j c p i , j Δ X j Δ t T i + 1 Δ X j 2 k i , j + Δ X j 2 k i + 1 , j T i + 1 = ρ i , j c p i , j Δ X j Δ t T i

  • State 8.1 equation:

    (35) 1 Δ X j 2 k i 1 , j + Δ X j 2 k i , j T i = 1 1 Δ X j 2 k i 1 , j + Δ X j 2 k i , j + 1 Δ X j + 1 2 k i , j + 1 + Δ X j + 1 2 k i + 1 , j + 1 + ρ i , j c p i , j Δ X j 2 Δ t + ρ i , j + 1 c p i , j + 1 Δ X j + 1 2 Δ t T i + 1 Δ X j + 1 2 k i , j + 1 + Δ X j + 1 2 k i + 1 , j + 1 T i + 1 = ρ i , j c p i , j Δ X j 2 Δ t + ρ i , j + 1 c p i , j + 1 Δ X j + 1 2 Δ t T i

  • State 8.2 equation:

    (36) 1 Δ X j 2 k i 1 , j + Δ X j 2 k i , j T i 1 1 Δ X j 2 k i 1 , j + Δ X j 2 k i , j + h j + 4 σ T i 3 1 ε j + 1 ε j + 1 1 + ρ i , j c p i , j Δ X j 2 Δ t T t + h j + 4 σ T i + 1 3 1 ε j + 1 ε j + 1 1 T i + 1 = ρ i , j c p i , j Δ X j 2 Δ t T i + 3 σ 1 ε j + 1 ε j + 1 1 T i + 1 4 T i 4

  • State 8.3 equation:

    (37) h j 1 + 4 σ T i 1 3 1 ε j 1 + 1 ε j 1 T i 1 1 Δ X j 2 k i , j + Δ X j 2 k i + 1 , j + h j 1 + 4 σ T i 3 1 ε j 1 + 1 ε j 1 + ρ i , j c p i , j Δ X j 2 Δ t T i + 1 Δ X j 2 k i , j + Δ X j 2 k i + 1 , j T i + 1 = ρ i , j c p i , j Δ X j 2 Δ t T i 3 σ 1 ε j 1 + 1 ε j 1 T i 4 T i 1 4

  • State 9.1 equation:

    (38) 1 Δ X j 2 k i 1 , j + Δ X j 2 k i , j T i 1 1 Δ X j 2 k i 1 , j + Δ X j 2 k i , j + ρ i , j c p i , j Δ X j 2 Δ t T i = ρ i , j c p i , j Δ X j 2 Δ t T i

  • State 9.2 equation:

    (39) 1 Δ X j 2 k i 1 , j + Δ X j 2 k i , j T i 1 1 Δ X j 2 k 1 1 , j + Δ X j 2 k i , j + h env + F env σ 4 T i 3 + ρ i , j c p i , j Δ X j 2 Δ t T i = p i , j c p i , j Δ X j 2 Δ t T i h env T env F env σ 3 T i 4 T env 4

CASE STUDY

The validation of the developed model was performed in Mazzaracchio and Marchetti (2010a)Mazzaracchio A, Marchetti M (2010a) A probabilistic sizing tool and Monte Carlo analysis for entry vehicle ablative thermal protection systems. Acta Astronautica 66(5-6):821-835. doi: 10.1016/j.actaastro.2009.08.033
https://doi.org/10.1016/j.actaastro.2009...
via two case studies by comparing the findings with the literature results related to the dimensioning of thermal shields based on highly representative industrial standard numeric programs. In particular, the verification generated extremely positive results and was performed by comparing the values found by Chen et al. (2006)Chen Y-K, Squire T, Laub B, Wright M (2006) Monte Carlo analysis for spacecraft thermal protection system design. AIAA Paper 2006-2951. using the FIAT code for the TPSs of the Stardust and Mars Exploration Rovers probes. The Stardust capsule's re-entry under nominal conditions is presented in Figs. 17-20. The temporal trends of the main kinematic characteristics of the flight (altitude vs. time and velocity vs. time, respectively Figs. 17 and 18) and the convective and radiative components of the heat fluxes (Fig. 19) together with the temperature trend within the shield are presented, highlighting the zones where the ablative material is subdivided (Fig. 20). The pyrolysis reaction starts almost immediately and is followed by the beginning of the charring of the ablative layer (t = 6.5 s). At 20 s, the initial slow surface regression begins, and its intensity gradually increases. At t = 53.5 s, the maximum peak of the entering heat flux, which is the sum of the two convective and radiative components, is reached. At this time, the surface temperature reaches its maximum value of 3,740 K. Starting from a time of approximately 75.5 s, the maximum temperature of the shield no longer occurs on the surface: instead, from this moment on, the temperature trend has a maximum at points inside the ablative layer. At t = 83 s, the surface regression of the char layer stops, with a total thickness of consumed ablation material of 7.1 mm. The run then proceeds until 100 s, at which point the TPS mission is considered to be complete and the entering heat flow ceases.

Figure 17
Stardust capsule re-entry: altitude vs. flight time.

Figure 18
Stardust capsule re-entry: velocity vs. flight time.

Figure 19
Stardust capsule re-entry: entering heat fluxes vs. flight time.

Figure 20
Stardust capsule re-entry: TPS temperature trend.

As for the superficial regression, its value measured at the capsule return was about 15 mm. However, this value takes into account the total erosion during the flight. Conversely, Covington et al. (2004)Covington MA, Heinemann JM, Goldstein HE, Chen Y-K, Terrazas-Salinas I, Balboni JA, Olejniczak J, Martinez ER (2004) Performance of a low density ablative heat shield material. 37th AIAA Thermophysics Conference, Paper No. AIAA-2004-2273, June 28-July 1, 2004; Portland, Oregon., for values of maximum entering heat flux of 1,200 W/cm2, i.e., for nominal conditions equal to those hypothesized here, the estimate of the thickness of lost ablative material is about 10 mm: value extensively within of the interval at 1 s with respect to the result obtained by Mazzaracchio et al. (2010a)Mazzaracchio A, Marchetti M (2010a) A probabilistic sizing tool and Monte Carlo analysis for entry vehicle ablative thermal protection systems. Acta Astronautica 66(5-6):821-835. doi: 10.1016/j.actaastro.2009.08.033
https://doi.org/10.1016/j.actaastro.2009...
using this ablative model with a Monte Carlo simulation.

CONCLUSIONS

A numerical program was presented to solve the one-dimensional thermal problem for charring ablative materials. The results of the Stardust return capsule were presented as a case study. The results obtained show excellent agreement with similar literature data obtained using industrial standard programs with high representativity like the NASA Ames Research Center's “Fully Implied Ablation and Thermal” (FIAT) code for the TPSs.

REFERENCES

  • Chen Y-K, Milos FS (1999) Ablation and thermal analysis response program for spacecraft heatshield analysis. Journal of Spacecraft and Rockets. 36(3):475-483.
  • Chen Y-K, Milos FS (2000) Two-dimensional implicit thermal response and ablation program for charring materials on hypersonic space vehicles. AIAA Paper 2000-0206.
  • Chen Y-K, Milos FS (2005) Three-dimensional ablation and thermal response simulation system. AIAA Paper 2005-5064.
  • Chen Y-K, Squire T, Laub B, Wright M (2006) Monte Carlo analysis for spacecraft thermal protection system design. AIAA Paper 2006-2951.
  • Clark RK (1973) An analysis of a charring ablator with thermal nonequilibrium, chemical kinetics, and mass transfer. NASA TN D-7180.
  • Covington MA, Heinemann JM, Goldstein HE, Chen Y-K, Terrazas-Salinas I, Balboni JA, Olejniczak J, Martinez ER (2004) Performance of a low density ablative heat shield material. 37th AIAA Thermophysics Conference, Paper No. AIAA-2004-2273, June 28-July 1, 2004; Portland, Oregon.
  • Curry DM (1965) An analysis of a charring ablation thermal protection system, NASA TN D-3150.
  • Curry DM (1974) Two-dimensional analysis of heat and mass transfer in porous media using the strongly implicit procedure, NASA TN D-7608.
  • Ko L, Gong WL, Quinn RD (2007) Reentry thermal analysis of a generic crew exploration vehicle structure, NASA, TM-214607.
  • Landau HG (1950) Heat conduction in a melting solid, quarterly of applied mathematics 8(1):81-94.
  • Matting FW (1970) Analysis of charring ablation with description of associated computing program, NASA TN D-6085.
  • Marschall J, Maddren J, Parks J (2001) Internal radiation transport and effective thermal conductivity of fibrous ceramic insulations, AIAA Paper 2001-2822.
  • Mazzaracchio A, Marchetti M (2010a) A probabilistic sizing tool and Monte Carlo analysis for entry vehicle ablative thermal protection systems. Acta Astronautica 66(5-6):821-835. doi: 10.1016/j.actaastro.2009.08.033
    » https://doi.org/10.1016/j.actaastro.2009.08.033
  • Mazzaracchio A, Marchetti M (2010b) Coupled aeroassisted orbital plane change manoeuvre and thermal protection system optimisation, IAC-10-C2.7.4, 61st International Astronautical Congress, September 27-October 1; Prague, Czech Republic.
  • Mazzaracchio A, Marchetti M (2011) Spacecraft aerodynamics and heat shield characteristics impact on optimal aeroassisted coplanar orbital transfer, IAC-11-C2.7.9, 62nd International Astronautical Congress, October 3-7; Cape Town, South Africa.
  • Mazzaracchio A (2013) Heat shield mass minimization for an aerocapture mission to Neptune. International Journal of Aerospace Innovations 5(3-4):83-93. doi: 10.1260/1757-2258.5.3-4.83
    » https://doi.org/10.1260/1757-2258.5.3-4.83
  • Mazzaracchio A (2015) Flight-path angle guidance for aerogravity-assist maneuvers on hyperbolic trajectories. Journal of Guidance, Control, and Dynamics 38(2):2015:238-248. doi: 10.2514/1.G000670
    » https://doi.org/10.2514/1.G000670
  • Mazzaracchio A (2016) Coupled versus uncoupled optimal solution for thermal and dynamic problems in spacecraft atmospheric flight. World Journal of Engineering 13(1):53-60. doi: 10.1108/WJE-02-2016-006
    » https://doi.org/10.1108/WJE-02-2016-006
  • Moyer CB, Rindal RA (1968) An analysis of the coupled chemically reacting boundary layer and charring ablator - Part II - Finite difference solution for the in-depth response of charring materials considering surface chemical and energy balances, NASA CR-1061.
  • Nathan RP, Bindu S (2005) Low temperature ablative heat shield for re-entry vehicles, AIAA Paper 2005-5063.
  • Pittman CM, Howser LM (1972) Numerical analysis of the transient response of an axisymmetric ablative char layer considering internal flow effects, NASA TN D-6895.
  • Roberts L (1958) A theoretical study of stagnation-point ablation, NASA TN 4392.
  • Swann RT, Pittman CM, Smith JC (1965) One-dimensional numerical analysis of the transient response of thermal protection systems, NASA TN D-2976.
  • Williams SD, Curry DM (1992) Thermal protection materials - Thermophysical property data, NASA RP 1289.

Edited by

Section Editor: Marcia Mantelli

Publication Dates

  • Publication in this collection
    2018

History

  • Received
    03 Aug 2017
  • Accepted
    12 Jan 2018
Departamento de Ciência e Tecnologia Aeroespacial Instituto de Aeronáutica e Espaço. Praça Marechal do Ar Eduardo Gomes, 50. Vila das Acácias, CEP: 12 228-901, tel (55) 12 99162 5609 - São José dos Campos - SP - Brazil
E-mail: submission.jatm@gmail.com