Random input parameters.

## Abstract

Thermal protecti.on system (TPS) is one of the most important subsystems of hypersonic vehicles which are subjected to severe aerodynamic heating. The reliability and structural integrity of TPS are crucial to the structural safety and integrity of hypersonic aircrafts. During the design and service stages of the TPS, there are numerous inevitable uncertainties in aerothermal environment, material properties, manufacture and assembly errors, analysis modeling errors, etc., which have great impact on the reliability assessment of the TPS. In this chapter, the probabilistic heat transfer for TPS is presented to achieve light-weight TPS with high reliability. The uncertainties and their modeling and sampling methods are introduced at first. Then the finite element model and the precise integration method for transient probabilistic heat transfer of TPS are given. Based on the theoretical and numerical methods, the probabilistic design and reliability assessment procedure for TPS are finally formed. A typical multilayer ceramic composite material TPS is used as example to show the effectiveness of the proposed method.

### Keywords

- thermal protection system
- uncertainty
- probability method
- reliability assessment
- precise time integration

## 1. Introduction

When the vehicle flies in the atmosphere with hypersonic speed, the air is strongly compressed and rubbed, which converts kinetic energy into heat, resulting in a sharp increase of the ambient air temperature. The air temperature in stagnation point of the vehicle can reach 400°C when Ma ≥ 3.As the flight Mach number further increases, a strong bow shock will be formed in front of the head. The gas temperature in the shock layer is very high, causing the dissociation and ionization of air molecules, generating thermal radiation effects. Elevated temperatures can easily lead to damage to crew members and instruments in the cabin and adversely affect the structural integrity and reliability of the projectile. Meanwhile, the heat radiation of the shock wave layer and the optical window is also the infrared interference source of the optical system. Therefore, thermal protection is always a major issue for the design of hypersonic vehicles.

There are different types of TPSs according to the aerothermal conditions. NASA has classified the TPS into three categories, namely the passive, semi-passive, and active thermal protection system. The structure and technology of the active TPS are more complicated, which is the key and a difficult point for the research of the reusable hypersonic vehicle in the future. Currently, ablative and non-ablative passive (semi-passive) TPSs are being used extensively in hypersonic vehicles, both of which are typical multilayer TPSs. The non-ablative passive or reusable TPSs mainly include the rigid ceramic tile TPS, the metallic TPS and the ceramic matrix composite (CMC) shingle TPS.

The ablative thermal protection structure uses carbon, silicon, phenolic, and other heat-resisting ablative composites to cover the surface of vehicles. The material undergoes complex physical and chemical changes such as decomposition, melting, evaporation, sublimation, and erosion. The heat is taken away by the consumption of the surface material [1], and the working mechanism is shown in Figure 1. Under the action of heat flow, the resin matrix and other inside materials are pyrolyzed after absorbing heat. The generated gas flocks to the boundary layer to prevent heat from transmitting into the inside material, causing a heat blocking effect. At the same time, the resin and other carbonized products caused by heating are deposited on the surface to form a “carbonization layer,” which acts as an outward radiation.

Rigid ceramic tile TPS consists of ceramic tiles (high temperature and emissivity coating), nomex flexible strain isolation pad (SIP), room temperature vulcanizing (RTV) silicone, and gap filler, and the multilayer structure is shown in Figure 2. The latest generation of this kind of TPS is the toughened uni-piece fibrous reinforced oxidation-resistant composite (TUFROC), which is used in the X-37B spacecraft. This new ceramic composite structure not only withstands the heat generated during reentry but also solves the bottleneck problems of ceramic tiles in the high-temperature environment, such as thermal cracking and oxidation. TUFROC comprises a treated carbonaceous cap composed of refractory oxidation-resistant ceramic carbon insulation (ROCCI), which provides dimensional stability to the outer mold line, while the fibrous base material provides maximum thermal insulation for the vehicle structure. This structure can be reused, and its structure is shown in Figure 3.

The research on metallic TPS started in the late 1950s to early 1960s. Metallic TPS has the advantages of light quality, durability, strong maneuverability, and low cost. At present, the main structure of metallic TPS is superalloy honeycomb TPS panel, which can be used for windward or leeward surfaces of hypersonic vehicles with 600–1200°C. The structure of an improved metallic TPS design, adaptable robust metallic operable reusable (ARMOR) insulation, is shown in Figure 4. In this structure, the honeycomb sandwich panels are simultaneously subjected to aerodynamic heating, pneumatic pressure, and sound pressure. In order to solve the connection difficulties, the outer honeycomb sandwich panel is connected to support structure at each corner of the panel box by the metal bracket. The pressure load is transmitted through the four flexible supports and the beam structure of the panel box, not only avoiding thermal interference but also to allow free thermal expansion of the outer surface [4, 5].

The concept of CMC shingle TPS was originally preferred by HERMES. Compared with the traditional rigid ceramic tile TPS, it has the advantages of strong impact resistance, easy loading and unloading, easy to monitor, and so on. Compared with metallic TPS, it has the advantages of high carrying capacity, high temperature resistance, producing larger size components, and lower maintenance costs. As shown in Figure 5, this kind of TPS consists of CMC shingle, sealing materials, internal insulation materials, and connectors. The structural elements of CMC shingle TPS can be divided into two categories, namely, the element to withstand mechanical load (mechanical enclosure, fasteners, etc.) and element to withstand heat load (internal insulation, insulation felt, etc.). In the design stage, the mechanical and thermal protection functions can be treated separately so that the best material can be selected in each field.

For a variety of TPSs, design and maintenance is a comprehensive system issue that includes the acquisition of aerodynamic thermal, determination of materials and sizing, numerical and experimental verification, flight tests, health monitoring and management, maintenance, and other aspects. Various uncertainties exist in all stages, such as ballistic deviation, deviation of material property, model errors, data measurement errors, system assembly deviations, external loading conditions, and other unknown uncertainties. How to fully consider the influence of these uncertainties in the process of design and use and to set system margin and safety factor reasonably is the key issue for the reliability and structural integrity of TPS. Therefore, in the following sections, the treatment of various uncertainties and the probabilistic design method for TPS will be systematically introduced.

## 2. Uncertainties in TPS and their analysis methods

### 2.1. Classification and acquisition of uncertainties

There are two kinds of uncertainties in the design and service stage of TPS. One is based on the nature of the uncertainty and the other is based on the source of the uncertainty. According to the nature of uncertainties, they can be divided into the following three categories [7]:

Stochastic variability: It means that uncertainty can be clearly recognized but cannot be reduced by means of experiment or analysis, such as natural fluctuations in the atmospheric environment and so on. If the quantitative change law of this kind of uncertainty can be obtained or it can be parameterized, it can be analyzed by Monte Carlo simulation.

Structural uncertainty: It mainly refers to the errors and uncertainties caused by mathematical models and numerical discretization, such as the simplification of physical processes and boundary conditions, temporal and spatial discretization errors, and so on. These kinds of uncertainties require reliable validation and verification of the model to control their effect on the design and application of TPS.

Parametric uncertainty: It mainly refers to the uncertainties of the various input parameters required in the analysis and design, such as the catalytic and emissivity of the wall, the rate of chemical reaction, transition criteria, and so on.

According to the sources of the uncertainties, they can be roughly divided into the following categories:

Trajectory uncertainties: It is mostly due to navigation uncertainty and orbital maneuver uncertainty, which can be modeled by the overshoot and undershoot of trajectory design.

Aerodynamic heating uncertainties: In hypersonic flight, complex ionization and chemical reaction occur in the boundary layer. Besides, uncertainties in flow conditions and boundary conditions are also factors causing aerodynamic uncertainties. The uncertainty of aerodynamic heating is one of the key factors in the size design of TPS.

Material uncertainties: For non-ablative TPS, the deviation of parameters, such as thermal conductivity, heat capacity, surface emission coefficients, and other thermal physical parameters, is mainly considered. For ablative TPS, the abovementioned uncertainties must be considered both before and after ablation. In general, material uncertainty is obtained by experiment. Material uncertainty is also one of the major factors in the size design of TPS.

Uncertainties of processing and assembly error: The error of processing and assembly refers to the geometric uncertainty, such as thickness error, geometric error, gap in assembly process and so on.

Analysis model uncertainties: This kind of uncertainty refers to the error introduced in the description, modeling, and numerical calculation process of TPS. The uncertainty of the numerical calculation part can be controlled and estimated by the validation and verification of the model.

Uncertainties of ground test: This type of uncertainty mainly refers to the ground test equipment that cannot accurately reproduce the true working process of aerodynamic thermal environment. In general, the impact of this type of uncertainty can be eliminated by reasonably setting the thickness margin of the TPS.

Other unknown uncertainties: Other uncertainties not included above.

As mentioned above, there are many uncertainties in the design, analysis, and use of TPS. Due to its complexity, it is almost impossible to accurately obtain the degree and distribution of uncertainty through a large amount of test data. Therefore, in probability analysis, it is often necessary to combine various methods to estimate the uncertainty of the input parameters. In general, there are six ways to get the uncertainty of the system parameters [8]. These ways are:

A combination of a variety of data sources to estimate the mean and standard deviation of the parameters. This is the only approach that is statistically rigorous and should be used wherever possible.

Performing an independent assessment of the data source. This method can be used when only a single data source is available.

Using recommended data in relevant reviews or reports.

Using the uncertain data given or cited in the original data.

Adopting the similarity rules.

Relying on the expert judgments.

### 2.2. Sampling methods and surrogate models

#### 2.2.1. Sampling methods

In the probability and reliability analysis of TPS, the sampling procedure is always inevitable. Therefore, how to obtain the statistical characteristics of the system accurately with a few samples is the key to improve the overall analysis efficiency. We introduce three commonly used sampling methods in probability analysis, namely the Monte Carlo sampling method, the Latin hypercube sampling method, and the adaptive importance sampling method.

(1) Monte Carlo sampling method.

Monte Carlo sampling can be directly used to simulate the real process of various projects. When the input uncertainty can be quantified, this method can simulate any process of the TPS in practice. A significant disadvantage of this method is that there is no memory feature for the sampling process, which may lead to a concentration or agglomeration of sampling points around a local region of the sample space. If the sample data are distributed in the entire input parameter space unevenly, those concentrated data points in the simulation cycle are equivalent to double counting, which does not provide effective reference value any more.

(2) Latin hypercube sampling method.

Latin hypercube sampling method has the feature of sampling memory, which can effectively overcome the problem of low efficiency caused by data concentration in Monte Carlo sampling method. Besides, it forces that the sampling points must be discretely distributed throughout the sampling space. Therefore, in general, for the same precision results of same problem, Latin hypercube sampling method can reduce the simulation cycle number by 20–40% than direct Monte Carlo sampling method. It needs to be pointed out is that the number of simulation cycles is determined by the problem itself to achieve a desired accuracy or confidence level.

In Latin hypercube sampling method, the range of each random variable is divided into several intervals. A representative value is randomly selected from each interval of a random variable. These representative values form a sample of this random variable, which ensures that all possible values are extracted during the sampling process, called the memory feature in the sampling process.

(3) Adaptive importance sampling method.

In the analysis of system probability properties, the distribution of output parameters in the whole sample space is concerned. Therefore, we should avoid sample concentration in a certain area when sampling. But in the analysis of system failure probability or reliability, more attention should be paid to the samples near the failure critical state. However, the failure probability of the TPS is very low. If the direct Monte Carlo or the Latin hypercube method is used, the samples near the failure state are sparse and thus the final failure probability estimation is not accurate enough. In this case, important sampling methods can be used to change the sampling center of gravity in order to increase the number of samples in the failure state.

#### 2.2.2. Surrogate models

In the analysis of probability and reliability of TPS, the direct calculation using the sampling method introduced in the previous article consumes a large amount of computational resources. However, it is an effective way to improve computational efficiency by replacing the complex relationship between random input and stochastic output using the existing data to construct the surrogate model. At present, the surrogate models commonly used in the engineering field include response surface method (RSM), Kriging model, and so on.

(1) Response surface method.

The basic idea of RSM is to use a proper function to approximate a function that cannot be explicitly expressed. The main purpose is to reduce computation as much as possible with acceptable accuracy. To achieve this goal, the probabilistic features of the original state surface must be described as accurately as possible, especially near the design point. RSM assumes that the influence from random input variables can be expressed using mathematical functions. In general, the function is a polynomial of one or two orders, such as:

where
_{,} and

For TPS, although it is unrealistic to accurately replace the real function by the response surface function, sensitivity analysis based on the response surface can be easily performed to provide guidance for subsequent parameter designs. In practice, the RSM is often used to solve the probability and sensitivity near the checking point where the performance function is easier to fit. Meanwhile, the simplified form of the response surface also helps to reduce the computational efforts of system analysis when analyzing the sensitivity of a TPS with multiple random variables.

Although RSM requires fewer samples than the Monte Carlo method, it is not as versatile as Monte Carlo simulation techniques. In order to use RSM, the function relationship among random variables must be smooth and continuous. If this case is not satisfied, this method cannot be used. Besides, the number of required cycles depends on the number of random input variables. If there are many random variables, fitting the response surface requires a large number of samples, which lead to the loss of efficiency. As an unbiased estimate, the response surface function value calculated at each sample point has an error, even at the design points.

(2) Kriging model.

The Kriging model method is a statistical prediction method based on the stochastic process. It can find the optimal, linear, and unbiased estimates of regional variables with the smoothing effect and the least statistical variance of estimation. The Kriging model method has the characteristics of local estimation, which makes it easier to obtain the ideal fitting effect in solving the problem of higher nonlinearity.

The Kriging model is composed of the global model and the local deviation, which can be expressed as

where

where

where
* x* is:

where
* x* and sampling point

The variance of the global model is estimated as:

The relevant parameters

When

## 3. Structural response analysis considering uncertainties

When designing and analyzing a TPS, the impact of uncertainties on the structural responses should be taken into account due to divergence of material property parameters and geometry uncertainties from manufacture errors. These uncertainties must be properly introduced into the probabilistic finite element model.

### 3.1. Characterization of uncertainties and probabilistic finite element model

Parametric uncertainties of TPS structures, once acquired, should be used for building up a finite element model using commercial software packages like ANSYS and ABAQUS. Taking a CMC shingle TPS, for example, it is a typical multi-plate structure, where the top part is a C/C composite cover plate with a radiation layer, followed by an insulation felt, and a high-temperature and a low-temperature insulating layer. And there is a strain isolation pad between the low-temperature insulating layer and the vehicle surface. The finite element model is shown in Figure 6. Here, the uncertainties in material properties and geometries are considered.

Firstly, we consider the uncertainties of thermal conductivity and specific thermal capacity. These parameters are nonlinear dependent on the working temperature. Here, fitted polynomial functions are applied to relate the parameter and the temperature, where only the constant component varies to represent the uncertainties. For example, a 4-order polynomial is used to fit 10 experimental data of thermal conductivity shown in Figure 7 (a) and (b) demonstrates its upper and lower bounds. Such characterization can be generalized for other nonlinear parameters that vary with temperature.

### 3.2. Precise time integration for transient thermal conduction

The temperature distribution can be obtained through transient analysis once the finite element model is built. The analysis efficiency is one of the most crucial factors since massive samples and calculations have to be carried out during the probabilistic analysis of a system. Usually, the direct integration is used for the transient heat conduction analysis, which shows good generality but high dependency of precision and stability on time step length. To improve the precision, efficiency, and stability of transient thermal analysis, the precise time integration (PTI), which is based on the exponential integrator, is applied here.

The discretization equation of thermal conduction is

where

where

The analytical transient solution can be written as

where

Generally, we only need to carry out those calculations at the exact moments and do not need to start from

where

Noting that

and Eq. (11) can be rewritten as

where

If the inhomogeneous terms

and

Now we get an equivalent homogeneous format as

And the solutions can be calculated from

The equality above expands the dimensions, but avoids matrices inversions, which makes the algorithm more efficient. The key to this algorithm is how we acquire the matrix exponential

The precise time integration for transient thermal conduction has some impressive advantages, including the insensitivity of time step length, unconditional stability, and high efficiency. Figure 8 shows a comparison between PTI algorithm and the conventional time difference method for a transient heat transfer problem with analytical solutions. It can be easily observed that the precise integration method is far more accurate than the conventional time difference method used in ANSYS .

## 4. Probabilistic design and reliability assessment method for thermal protection system

The main failure mode of TPS for thermal reliability assessment is the overtemperature of bondline. Based on the theories and methods introduced in the previous sections, the procedures for probabilistic design and reliability evaluation are shown in Figure 9. The three main modules are introduced in detail as the following.

(1) Deterministic design.

Deterministic design, also known as nominal design, is often used to determine the initial geometry of TPS. First, thermal load at different locations needs be predicted numerically or experimentally according to the aerothermal environment using the initial configuration. Secondly, the proper material types and structural forms of TPS according to estimated thermal load need to be determined. We can finally get an initial thickness

where

(2) Probabilistic design.

Firstly, the source, range, and distribution of the system uncertainties are determined, and the parameters are reasonably set using the methods in Section 2. Then, the probability structural thermal response is analyzed by the methods described in Section 2 and 3. Finally, the probability characteristics of the structure temperature field are obtained. The highest value of bondline overtemperature is usually taken as a random output variable, which is:

where

(3) Reliability assessment.

The thermal reliability of a TPS can be defined as:

where

In the probabilistic design, it is necessary to set a reasonable range of reliability. If the system meets the reliability requirements, the design process ends, or else the system parameters need to be adjusted to carry out the reliability evaluation again.

## 5. Application of probabilistic design in thermal protection system

The probabilistic design of a typical CMC shingle TPS is used as an example to show the effectiveness of the design procedure in the previous section. Here, the uncertainty of the material properties and thickness is considered.

### 5.1. Model uncertainty

The typical composite TPS investigated here is shown in Figure 6 of Section 3.1, assuming that 18 random input variables follow the truncated normal distribution. The standard deviation is 3% of the mean value, and the range of the value is mean ± 5%. The random input parameters are listed in Table 1. Taking H1 (thickness of CMC cover) as an example, the probability density and cumulative probability distribution curve are shown in Figure 10.

Variable | Physical meaning | Mean | Std. dev. | Lower bound | Upper bound |
---|---|---|---|---|---|

H1 | Thickness of CMC cover | 4.0E-3 | 1.33E-4 | 3.8E-3 | 4.2E-3 |

H2 | Thickness of fiber cloth | 3.0E-3 | 1.0E-4 | 2.85E-3 | 3.15E-3 |

H3 | Thickness of high-temp insulation layer | 2.0E-2 | 6.67E-4 | 1.90E-2 | 2.10E-2 |

H4 | Thickness of low-temp insulation layer | 4.5E-2 | 1.5E-3 | 4.275E-2 | 4.725E-2 |

H5 | Thickness of strain isolation pad | 3.0E-3 | 1.0E-4 | 2.85E-3 | 3.15E-3 |

H6 | Thickness of skin | 4.0E-3 | 1.33E-4 | 3.8E-3 | 4.20E-3 |

A4 | Sp.ht. of CMC cover | 815.4 | 27.18 | 774.63 | 856.17 |

B4 | Sp.ht. of fiber cloth | 1000.0 | 33.33 | 950.0 | 1050.0 |

D5 | Sp.ht. of high-temp insulation layer | 900.30 | 30.10 | 855.28 | 945.32 |

E2 | Sp.ht. of low-temp insulation layer | 990.0 | 33.00 | 940.50 | 1039.5 |

F4 | Sp.ht. of strain isolation pad | 1000.0 | 33.33 | 950.0 | 1050.0 |

G2 | Sp.ht. of skin | 840.0 | 28.00 | 798.0 | 882.0 |

L5 | Heat conductivity of CMC cover | 47.4 | 1.59 | 45.32 | 50.09 |

M3 | Heat conductivity of fiber cloth | 0.23 | 7.67E-3 | 0.2185 | 0.2415 |

O7 | Heat conductivity of high-temp insulation layer | 6.125E-2 | 2.04E-3 | 5.82E-2 | 6..43E-2 |

P6 | Heat conductivity of low-temp insulation layer | 5.195E-2 | 1.73E-3 | 4.94E-2 | 5.45E-2 |

Q4 | Heat conductivity of strain isolation pad | 6.10E-2 | 2.03E-3 | 5.80E-2 | 6.41E-2 |

S2 | Heat conductivity of skin | 0.65 | 2.16E-2 | 0.6175 | 0.6825 |

### 5.2. Probability characteristics and sensitivity analysis

For a given heat flux, a nonlinear transient thermal analysis is conducted considering uncertainties listed in Table 1. The highest bondline temperature during the computational interval is selected as the random output parameter. In order to reduce the computational load, Latin hypercube sampling method is used. Then, the response surface is fitted with the calculated samples. Next, the probability characteristic of random output is obtained by sampling again on the response surface. In order to ensure the accuracy of the fitting, the number of samplings for fitting the response surface needs to be at least twice the number of samples, that is, ((18 + 1) × (18 + 2)/2) × 2 = 380 times.

Sampling 10,000 times on response surface, the histogram and the cumulative probability distribution are shown in Figure 11. The results show that the maximum bondline temperature is basically in accordance with the normal distribution with a mean of 77.83°C, a standard deviation of 6.765, a maximum of 104.23°C, and a minimum of 59.31°C. Based on the results of the probabilistic analysis, the reliability can be evaluated using the method described in Section 4.

The sensitivity of the random output to the random input parameters can be obtained by probability analysis, as shown in Figure 12. The top five input variables that have great impact on the output variables of the TPS are H4, E2, H3, P6, and H6, respectively, namely the thickness of low-temperature insulation layer, the specific heat capacity of the low-temperature insulation layer, thickness of high-temperature insulation layer, heat conductivity of low-temperature layer, and the thickness of the vehicle structure. Among them, the thickness of the low-temperature insulation layer has the greatest impact on the maximum bondline temperature with negative correlation. In the first five parameters, only P6 is positively correlated.

According to the sensitivity analysis result, in order to improve the reliability of the investigated TPS, manufacturing errors and the deviation of specific heat of low-temperature insulation layer should be strictly controlled.

## 6. Conclusion

With the rapid improvement of experimental and computational technologies, the uncertainties of thermal protection system are more easily acquired and modeled, which makes the probabilistic design of highly reliable TPS possible. In this chapter, the types and sources of uncertainties during the design and service stages of the thermal protection system are summarized. The mathematical treatment and finite element modeling methods of various uncertainties are introduced. The dimensional expanding precise time integration method is used to improve the numerical accuracy and efficiency of transient heat conduction analysis. Based on the theoretical and numerical methods, the probabilistic design and analysis framework is finally formed and is validated via an example of a typical CMC shingle TPS.

The following is recommended for the future research: (1) theoretical and experimental methods for uncertainty acquisition of TPS; (2) probabilistic modeling and fast analysis method of structural thermal response considering uncertainty; (3) the influence of structural cracks and defects on the reliability of TPS; and (4) reliability and risk assessment method of TPS based on probability analysis.

## Acknowledgments

This work is supported by the National Natural Science Foundation of China (Grant No. 11502037).

## References

- 1.
Meseguer J, Pérez-Grande I, Sanz-Andrés A. Thermal protection systems. In: Meseguer J, Pérez-Grande I, Sanz-Andrés A, editors. Spacecraft Thermal Control. 1st ed. Holand: Academic Press, Elsevier; 2012. pp. 305-325. DOI: 10.1533/9780857096081.305 - 2.
Rodriguez AC, Snapp CG. Orbiter thermal protection system lessons learned. In: AIAA SPACE 2011 Conference & Exposition; 27–29 September 2011; Long Beach, California. United States: AIAA; 2011 - 3.
Stewart DA, Leiser DB. Lightweight TUFROC TPS for hypersonic vehicles. In: 14th AIAA/AHI Space Planes and Hypersonic Systems and Technologies Conference; 6–9 November 2006; Canberra, Australia. America: AIAA; 2006 - 4.
Blosser ML. Fundamental modeling and thermal performance issues for metallic thermal protection system concept. Journal of Spacecraft and Rockets. 2004; 41 (2):195-206. DOI: 10.2514/1.9182 - 5.
Poteet CC, Blosser ML. Improving metallic thermal protection system hypervelocity impact resistance through numerical simulations. Journal of spacecraft and rockets. 2004; 41 (2):221-231. DOI: 10.2514/1.9193 - 6.
Pichon T, Barreteau R, Soyris P, et al. CMC thermal protection system for future reusable launch vehicles: Generic shingle technological maturation and tests. Acta Astronautica. 2009; 65 (1):165-176. DOI: 10.1016/j.actaastro.2009.01.035 - 7.
Glass DE. Ceramic matrix composite (CMC) thermal protection systems (TPS) and hot structures for hypersonic vehicles. In: 15th AIAA Space Planes and Hypersonic Systems and Technologies Conference; 28 April-1 May 2008, Dayton, Ohio. United States: AIAA; 2005 - 8.
Wright MJ, Bose D, Chen YK. Probabilistic modeling of aerothermal and thermal protection material response uncertainties. AIAA Journal. 2007; 45 (2):399-410. DOI: 10.2514/1.26018 - 9.
Zhong WX, Williams FW. A precise time step integration method. Proceedings of the Institution of Mechanical Engineers, Part C: Journal of Mechanical Engineering Science. 1994; 208 (6):427-430. DOI: 10.1243/PIME_PROC_1994_208_148_02 - 10.
Lin J, Shen W, Williams FW. A high precision direct integration scheme for structures subjected to transient dynamic loading. Computers & structures. 1995; 56 (1):113-120. DOI: 10.1016/0045-7949(94)00537-D