Abstract
The metallic thermal protection system (MPTS) is a key technology for reducing the cost of reusable launch vehicles, offering the combination of increased durability and competitive weights when compared with other systems. A two-stage optimization strategy is proposed in this paper to improve thermal and mechanical performance of MTPS while minimizing its weight. Sensitivity analysis is conducted to determine the influence of various design variables, and these variables are classified into two groups according to the results. Three main parameters are optimized in the first-stage optimization, and others are optimized in the second stage. Combining the response surface method with genetic algorithm ensures the computational efficiency of this strategy. Optimum results show that the thermal insulation capability of MTPS is improved by 20% with little mass cost, validating the feasibility of this optimization strategy.
1. Introduction
Hypersonic vehicles have to withstand extremely high aerodynamic heating and pressure loads during flight, which makes thermal protection system a key factor for flight security. Considering the maintenance costs and turnaround time between flights, the metallic thermal protection system is now widely used on middle temperature district of hypersonic vehicles. A large number of research works analyzing the performance of a particular design configuration have been done in recent years. However, there is little information available to indicate how to improve the performance of the MTPS.
Blosser [1, 2] derived an approximate equation to predict the thickness and mass of insulation required to limit the maximum temperature reached by an insulated structure. The result was compared with a one-dimension finite element heat transfer model. Xie [3] built a simplified heat transfer model for a multilayer insulation system based on which the globally convergent method of moving asymptotes method was utilized for searching optimal solutions. Results show that the usage of multilayer insulation materials can save more than 17% weight compared with single-layer insulation. Ji [4] built a two-dimension numerical model to investigate the combined radiation and conduction heat transfer in multilayer insulations. It turned out that an optimal number of insulation layers for best thermal performance exists. However, these works focus merely on the insulation, ignoring the influence of other components on MTPS performance. Poteet [5] developed a coupled iterative analysis process to perform sensitivity studies and to size the MTPS system preliminarily. A special one-dimensional finite element thermal analysis model was established and then combined with static defection and failure mode analysis to obtain minimum weight configurations. However, it is impossible to simulate detailed temperature distribution according to this model. Blosser [6] built a two-dimensional model of MTPS and analyzed the key factors that govern the performance of MTPS. The heat capacity of the underlying structure was found to have the greatest effects on thermal performance. Ng [7, 8] analyzed the thermomechanical behaviour of the damaged thermal protection system and examined the effect of damage location using a three-dimensional finite element model. Zhang [9] studied the effects of different material combinations on performance of pyramidal core sandwich panels. They also investigated some specific design variables which are significant for thermomechanical properties of the structure and eventually to get a good balance between insulation performance, load-bearing capacity, and weight. Wang [10] studied the thermal protection system which incorporates graded insulation materials and multilayer ceramic matrix composite cellular sandwich panels to avoid drawbacks of traditional structure. With specific graded insulation materials, the system can achieve better insulation performance. The buckling strength can be significantly improved, and the thermal short effects can be suppressed. These studies provide the relative importance of modelling assumptions and design parameters but not with the way to search for the optimum results. Chen [11] proposed an innovative concept of structural efficiency to access thermal-mechanical performance of the lattice-core thermal protection system. It can help develop the thermal protection system with better thermal and mechanical performance at the same time. Jiang [12] optimized the layout of insulations of five different integrated thermal protection systems (ITPS) on the basis of three-dimensional models. Zhao [13] established two different configurations of ITPS and optimized them with simulated annealing algorithm. Ji [14] worked on analyzing and optimizing the lightweight corrugated-core sandwich integrated thermal protection system. They digged some lightening holes and applied different materials to the structure to reducing structural weight. Also, by changing the hole size and web dip angle, a modified structure was obtained. Meng [15] optimized the TPS for the solar array of stratospheric airships. Taking thermal protection effects, structure weight, and output power as objective functions, nondominated sorting genetic algorithm was developed to solve this multiobjective optimization problem. Yang [16] established a topological optimization method for ITPS to reduce equivalent efficiency and to maintain structural stiffness. With the modification of layout of materials, the new design can provide a good balance between thermal insulation and load-bearing capacity. Compared to ITPS, MTPS is much more complex. It will bring about considerable computational burden if optimization is conducted directly on the three-dimensional model. Effects of different parameters on structural vary far from each other; in this case, optimizing all parameters simultaneously is unacceptable. In order to avoid rounding error, variables are classified into two groups and optimized in an iterative two-stage optimization process.
This paper aims to develop a feasible strategy for MTPS optimization, improving the thermal and mechanical performance effectively. Equivalent models are built and verified by detailed three-dimensional model to keep the balance between computational cost and accuracy.
2. Methodology
2.1. Finite Differential Method
The MTPS studied in this paper consists of lightweight fibrous insulation encapsulated between two honeycomb sandwich panels. The analysis of large-scale sandwich panel structures can be done with detailed FE models, but computational resources required are limiting. By using the finite differential method, the effective thermal properties can be estimated so that the size and complexity of the FE model decrease dramatically and thus the computational cost can be reduced, while acceptable level of accuracy maintained.
The sandwich panel considered in this study is composed of two metallic face sheets perfectly bonded to a metallic honeycomb core. Analysis can be conducted on a hexagonal cell due to the periodic and symmetric configuration of honeycomb panels. The following section outlines the methods for the estimation of equivalent thermal conductance coefficient of honeycomb core. Geometric parameters are graphically presented in Figure 1.
To utilize the finite differential method, the metallic part of the hexagonal cell is divided into seven elements along the longitudinal axes (the direction of z axis), while the air inside remains integral. Two assumptions are made here: (1) thermal convection inside the cavity is ignored. (2) in-plane heat transfer is ignored and temperature distribution in each element is uniform. The heat flows through the honeycomb core can be defined as follows:
Heat which flows though adjacent elements is as follows:where is thermal conductivity of material and is the height of element.
Heat flow through the air iswhere is thermal conductivity of gas inside honeycomb core which in this case is air. Although this part of heat is relatively small, it still cannot be ignored.
Cavity radiation is essential in heat transfer through honeycomb panels. Radiation between element m and element n can be defined as
is radiation view factor between the above two elements, which could be expressed aswhere is the angle between radiation and normal direction of element surface and l is the center-to-center spacing between two elements. In this study, the honeycomb cell is divided into seven elements along the thickness direction; hence, three kinds of view factors should be taken into account, as shown in Figure 2. A simplified method proposed by Buschman and Pittman [17] is adopted to calculate these view factors.
The first type is radiative exchange between the face sheets, which is defined aswhere R is the ratio between circumradius of the honeycomb cell and honeycomb core height.
The second type includes radiative exchange between the face sheets and the honeycomb core elements. The view factor between face sheet and ith element of the honeycomb core is expressed aswhere and .
The third type, view factors between honeycomb core elements, consists of two different conditions: view factor between element adjacent to face sheets and element in the middle, view factor between middle elements, and view factor between elements and themselves. These view factors are shown as follows:where , , and , .
Considering boundary conditions shown in Figure 1, thermal equilibrium equation for all seven elements can be achieved:
Temperature for all seven elements can be obtained by solving equations above, and then effective thermal conductivity can be calculated by solving Equation (1). Furthermore, taking thermal load that MTPS bears into account, material nonlinearity cannot be ignored. An iterative method is proposed in order to solve this problem, and the flow diagram for the overall procedure based on the finite differential method can be seen in Figure 3.
2.2. Response Surface Method
Even though the equivalent model is utilized to replace the detailed honeycomb panel model, embedding the overall MTPS finite element model into optimization codes directly will bring about enormous computational effort due to its huge number of cells. In this case, the response surface method is adopted to establish global approximation function of several thermal and mechanical responses. The accuracy of optimization depends on the accuracy of response surface, which strongly relies on the choice of sample points while the fitting method is determined.
Specific to the MTPS analyzed in this paper, three geometric dimensions are selected as design variables according to sensitivity analysis. These variables includes the thickness of fibrous insulation , thickness of outer honeycomb panel , and thickness of encapsulated panel . The quadratic polynomial is used to fit the response surface, generating ten undetermined coefficients. Chosen responses arewhere is the maximum temperature of the substrate structure, is the maximum out-plane deformation of the outer surface of MTPS, and is the maximum out-plane deformation difference of the outer surface of MTPS. , , and (i = 0, 1, …, 9) are coefficient matrixes fit by the least square method.
The accuracy of response surface increases with the increase in the number of identification points; however, the computational effort increase as well. In order to reduce computational costs while maintaining accuracy, the progressive strategy proposed by Carrere [18] is introduced. Begin with 7 × N + 1 sample points, increase the number of sample points by 80% every time according to residual error and confirm the final number while the accuracy is acceptable. N is the number of design variables and the residue error is defined as:
where is the approximation result and is the finite element simulation result. All samples must be chosen highly randomly, so that the Latin Hypercube Sampling (LHS) described by McKay [19] is adopted. Compared with the Monte-Carlo method, LHS is designed to build the input distribution by less sampling.
The number of sample points is increased by two times in this study. When the number of points comes to 52, the average residual error for , , , is 0.21%, 13.78%, and 9.86%, respectively, and the error decrease is less than 10% with further point number increase.
2.3. Two-Stage Optimization Strategy
A two-stage optimization strategy combining RSM and genetic algorithm (GA) is proposed in this study since it could be more effective when it comes to a complex system like MTPS. As a foundation, a convincible finite element model is established to conduct thermo-mechanical coupling analysis. Sensitivity analysis based on the parametric model is carried out after that, indicating that structure responds differently to different variables. Quantifying the effects of different variables on structural responses provides a basis for choosing design variables in the following optimizations and proving the validity of the optimum result.
Variables having greater influence on structural responses are selected to be optimized in the first-stage optimization. In this stage, significant performance improvement is expected to be obtained, that is reducing MTPS mass as much as possible while maintaining its thermal and mechanical performance. The main function of the MTPS is to protect the vehicle structural from overheating and also to provide an acceptable aerodynamic surface to prevent premature transition to turbulent flow during the atmospheric reentry. , and are chosen as the constraints in this case. Considering computational time, all these responses are approximated by RMS described in section 2.2. Optimum result of the first-stage optimization serves as the initial design model for the second-stage optimization, in which the outer honeycomb panel is optimized in detail using genetic algorithm. Multiple objective functions exist in this stage, thus Pareto Genetic Algorithm proposed by Horn [20] is introduced to solve this problem. Pareto optimality is a state of allocation of resources in which it is impossible to make any one individual better off without making at least one individual worse off. Tradeoffs can be within Pareto frontier and consideration of every parameter is not necessary.
The overall performance of MTPS is not taken into account in the second-stage optimization, so that the optimum result of this stage should be compared with the first-stage optimum result. Error between two stages is utilized to determine whether iterative optimization is necessary. When the error is no more than 10%, the final optimum configuration can be confirmed. Figure 4 shows the framework of this optimization strategy in detail.
2.4. Finite Element Modelling
The superalloy honeycomb MTPS panel shown in Figure 5 [21] is studied in this paper. It is designed to be mechanically attached directly to a smooth, continuous substructure. The outer surface is comprised of a foil-gage Inconel 617 honeycomb sandwich, and the inner surface is a titanium honeycomb sandwich with part of face sheet and core removed to save weight. Beaded, foil-gage, Inconel 617 sheets form the sides to encapsulate saffil insulation. The parameter of the panel coated Nomex felt pad that prevents hot gas flow beneath the panel provides preload to the mechanical fasteners and helps damp out panel vibrations. Considering the fact that the insulation is not assumed to bear any aerodynamic pressure, thermal model and structural model should be established separately. Results calculated by these 3-D finite element models are utilized to verify the accuracy of approximate functions obtained by the response surface method.
In this section, an equivalent model for honeycomb sandwich panel is first built. The outer and inner surface of MTPS are honeycomb sandwich panels. If detailed finite element models of these sandwich panels are built in thermal-structural coupling analysis, it will increase computational burden enormously. The honeycomb sandwich panels will be replaced by these equivalent models as parts of the overall MTPS when conducting heat transfer analysis and thermal stress analysis.
2.5. Equivalent Model for Honeycomb Sandwich Panel
In order to alleviate the computational burden in the future work, an equivalent model for the detailed sandwich panel is necessary. Detailed honeycomb sandwich panel is shown in Figure 6. The effective thermal conductivity and the equivalent stiffness matrix are calculated. Based on these equivalent properties, an anisotropic panel can be built to substitute the original model in the finite element analysis. This part of work is conducted to validate the accuracy of this equivalent model.
The honeycomb sandwich panel used in the experiment is composed of the top face sheet, the honeycomb core, and the bottom face sheet. Materials used for the honeycomb sandwich panel are shown in Table 1. A detailed finite element model for one honeycomb cell is built to simulate the experiment. The boundary conditions are shown in Figure 1.
Transient heat transfer analysis is conducted on this honeycomb cell. The total height of this cell is 7.82 mm, with the core having a height of 7.5 mm and the sheets 0.16 mm for each. The thickness of the honeycomb core is 0.076 mm. The honeycomb core circumcircle radius is 6 mm. As the applied thermal shock rate is 30°C/s, we calculated the bottom face sheet temperature for 30 seconds. The comparison between experimental results and finite element analysis results is shown in Table 2.
According to the comparison between experimental data and numerical analysis results, the detailed FE model is accurate enough. Thus, the detailed FE model can be used as the standard to verify the accuracy of equivalent model. Effective thermal conductivity can be obtained by finite differential method proposed in last section, as well as the detailed finite element analysis. Results are shown in Table 3.
As the temperature changes, the thermal properties of materials change as well. As shown in Table 3, the effective thermal conductivity varies far from each other with or without the material nonlinearity considered. In this case, when this nonlinearity is taken into account, a more accurate result can be obtained. By comparing the temperature distribution and the calculated effective thermal properties calculated by finite differential method and finite element method, errors between these two methods are acceptable. The finite differential method is accurate in calculating the effective thermal conductivity of the honeycomb panel, but it helps reduce the computational cost rapidly compared to the finite element method.
As for the mechanical properties of the honeycomb panel, we analyzed a honeycomb core cell, as shown in Figure 7.
(a)
(b)
When is small enough, we have
The subscript ‘s’ represents the properties for material. Other subscripts 1, 2, and 3 represent honeycomb core properties in x, y, and z directions, respectively. The equivalent stiffness matrix can be written as
The overall honeycomb sandwich panel dimension is . For the honeycomb core, , r = 3 mm, and H = 5 mm. Thickness of the top and bottom face sheets are both 0.2 mm. Two short edges of the panel are clamped. Material properties are shown in Table 4.
Both the detailed finite element model and the equivalent model are built for modal analysis. Results are shown in Table 5. The first-order modal shape is shown in Figure 8.
(a)
(b)
Errors between two models are around 5%, which means the equivalent model is accurate in global stiffness. The equivalent model can substitute the detailed model in the following analysis.
3. Heat Transfer Analysis Model
Heat transfer analysis is conducted to obtain the time-depend temperature distribution in the whole system and to find the maximum structural temperature by searching the recorded temperature history for every structural node. The temperature distribution is then applied to the stress analysis model to realize the coupled thermal and structural analysis, and the maximum structural temperature will be served as one of three constraints in the first-stage optimization. Materials are shown in Figure 5. More details about this model are shown in Figure 9. and are heights of the honeycomb core. The cores are encapsulated by two metallic sheets. These face sheets share the same thickness, denoted as .
Initial geometric parameters are shown in Table 6.
Figure 10 shows the thermal finite element model, and contacts between different components are supposed to be perfect in this case. This finite element model is established in ABAQUS to carry out heat transfer analysis, and the linear heat transfer solid element DC3D8 is used in this model. This analysis is applied to a MTPS panel located on the windward surface of the space shuttle orbiter during atmospheric entry. Absorbed heat flux is balanced by the heat radiated from the surface at any instant in time on the top surface of MTPS panel, and these radiation equilibrium temperatures are depicted in Figure 11. All other surfaces are assumed to be adiabatic, which means that all input heat is reserved inside the MTPS structure. While natural convection exists on the bottom surface of substrate structure and heat transfers between adjacent MTPS modules, this kind of simplification may lead to a conservative optimum design. Despite all this, this simplified model can help improve the computational efficiency and proved to be acceptable for the preliminary design phase.
Transient heat transfer analysis results at nodes A to D are plotted in Figure 12. Comparing the temperature between node D and analysis result in [23], it is found that they have the same trend. In our work, we assume all faces, beside the top face, are adiabatic. The time when the maximum temperature occurs is earlier than reference. The maximum temperatures, which is one of the key indicators in the future work, are almost the same.
The substructure temperature distribution at time 800 s is shown in Figure 13. It is obvious that temperature at node F is highest due to the connecting bolts. Maximum temperature on substructure is 809.7 K, which is beyond the temperature substructure can bear. It is necessary to optimize MTPS to protect the substructure.
3.1. Thermal Stress Analysis Model
The outer honeycomb panel of MTPS is fully exposed to atmosphere and thus bears relatively high aerodynamic pressure and severe temperature gradient. Carrying out thermal stress analysis is rather important for MTPS to maintain the structural integrity during flight. A structural finite element model is used to calculate the deflection of the outer honeycomb sandwich panel. Saffil insulation, Nomex felt, and underlying structure are deleted compared with thermal finite model shown. Uniform pressure loading read from Figure 11 is applied on the sandwich outer surface, and nodal temperatures extracted from heat transfer analysis are applied over the entire model at specific times of interest.
Boundary conditions of MTPS depend on its attachment methods to substructures. In this case, the model is constrained at the bottom surface of four bolts with all degrees of freedoms fixed. All three rotational freedoms of edges on the outer surface are fixed to simulate the actual interaction between modules in the MTPS panel array. The thermal stress analysis model is shown in Figure 14. This finite element model is established in ABAQUS to carry out thermal stress analysis. Solid element C3D82 is used here.
4. Results and Discussion
4.1. Sensitivity Analysis
Sensitivity analysis is performed to examine the effects of design variables on MTPS performance. These variables are as follows: the thickness of insulation ; the thickness of outer honeycomb panel ; the thickness of sides ; the thickness of outer honeycomb core ; the radius of outer honeycomb core r. The initial values and boundaries of these variables are listed in Table 7, and other structural parameters are hold constant throughout this study. The equivalent mechanical model of honeycomb sandwich panel is applicable mealy when ratio between and is no more than .
Three responses are selected to describe the structural performance: the maximum temperature of the substrate structure , the maximum out-plane deformation of the outer surface of MTPS , and the maximum out-plane deformation difference of the outer surface of MTPS . Reducing weight is always the main purpose of any aircraft design and thus effect on MTPS mass is ought to be taken into account. Parts of the sensitivity analysis results are plotted in Figure 15.
(a)
(b)
(c)
(d)
Quantified effects of individual design variables are listed in Table 8, in which all variables are assumed to be increased fivefold compared with their original values. Some conclusions could be drawn from the sensitivity analysis:(i)The insulation thickness has the greatest effect on maximum temperature, which means is the most important factor for improving MTPS thermal performance. However, since the insulation is a non-load-carrying structure, it has a definite disadvantage for the structural mechanical performance. Moreover, the thicker the insulation, the higher the encapsulated sides. The growing height of encapsulated sides will lead to a rapid increase in total mass.(ii)The outer honeycomb sandwich panel is exposed to atmosphere and is the key factor for maintaining the aerodynamic surface of MTPS. Results show that MTPS mechanical performance is extremely sensitive to . Increasing the thickness of honeycomb panel leads to the enhancement of thermal expansion effect. As a result, deformation of the outer surface increases, yet within an acceptable range. On the other hand, increasing height brings about larger bending stiffness and thus decreases deformation difference efficiently.(iii)The encapsulated side plays an important role in transferring force from the outer surface to substructure. However, this side can easily lead to heat shorts because of metal conduction. All responses are very sensitive to , especially the deformation difference. The outer honeycomb panel is too weak compared with thick side, and thicker side can enhance the effect of thermal mismatch. Conversely, if the side is too thin, it is highly possible that MTPS cannot stay integral during flight. In conclusion, must be carefully selected in the design process.(iv)Compared with other design variables, and r have less effect on each response. Their influence on structural deformation is directly related to the stiffness of the outer honeycomb panel. Since the order of magnitude varied far from all design variables, multistage optimization must be taken to improve optimization quality.
According to the sensitivity analysis, design domains for variables , , , , and r are set to [15 mm, 35 mm], [5 mm, 13 mm], [0.5 mm, 2.0 mm], [0.05 mm, 0.3 mm], and [1 mm, 5 mm], respectively.
4.2. First-Stage Optimization
In the first-stage optimization, response surface method is used to alleviate the computational burden. 52 sample points are utilized, and the residual for all samples based on the response surface method is plotted in Figure 16. Fitting coefficients for response surface functions are listed in Table 9.
The first-stage optimization problem is defined aswhere , and are design variables restricted to [15 mm, 35 mm], [5 mm, 13 mm], and [0.5 mm, 2.0 mm], respectively. is the total height of MTPS. W is the function of MTPS mass, described as
Optimum results are shown in Table 10.
The optimum result obtained by response surface fits the detailed finite model well with all errors less than 10%, verifying the accuracy and feasibility for applying global approximate response surface model to optimization procedure. Substructure overheating is the major problem for initial design, so the main purpose of optimization is to improving thermal insulation capability while maintaining the mechanical performance. According to the sensitivity analysis, this target could be reached by increasing insulation thickness and decreasing side thickness. Optimum results show that the insulation thickness is increased from 20.00 mm to 34.3 mm while side thickness decreased from 1.0 mm to the lower boundary value 0.5 mm, which is consistent with the prediction. Final maximum temperature is 647.03 K, indicating that optimization works well in this stage. However, MTPS mass rises from 0.337 kg to 0.36 kg, slightly beyond the set limit by 2.8%. No further improvement space remains in this stage, and adjustment will be made in the following section.
4.3. Second-Stage Optimization
The first-stage optimized model is utilized as the initial design model for the second-stage optimization. This stage is conducted to optimize the outer honeycomb panel in detail so that all responses exist in the first-stage optimization could not be used here. This multiobjective optimization is described in Table 11.
The detailed optimization model is defined aswhere is the reciprocal of equivalent stiffness matrix, is the equivalent thermal conductivity, and is the surface density.
Points in the Pareto frontier are listed in Table 12.
Mechanical performance of the honeycomb panel is relatively more sensitive to and r compared with its thermal performance, thus more attention should be paid to and . Four solutions with higher out-plane stiffness and lower surface density compared with initial design are chosen to compensate for the exceeding weight in the first-stage optimal result. These variables are applied in the complete MTPS analysis model, and results are shown in Table 13.
A configuration with higher equivalent stiffness and lower density can be obtained by decreasing the thickness and the radius of the honeycomb core. Solution No. 25 is selected as the final result of the second-stage optimization.
Thermo-mechanical analysis for the optimal model is conducted in ABAQUS as well. Results of the two-stage optimization are listed in Table 14. Maximum structural temperature drops from 809.70 K to 648.45 K, preventing the substrate structure from overheating. Maximum outer surface deformation and deformation difference are increased slightly. However, the deformations are 2.9% and 0.49% of the total MTPS height, respectively, still maintaining the mechanical performance of the structure. This optimum configuration turns out to combine thermal and mechanical performance perfectly, achieving the goal of optimization.
5. Conclusions
A two-stage optimization strategy suitable for preliminary design of MTPS is described, which allows estimation of structural design variables and MTPS weight at selected vehicle locations. Although this study is performed for a specific MTPS concept and a particular aerodynamic load profile, this strategy remains practicable for different MTPS configurations and load conditions. Main conclusions of this paper are as follows:(i)MTPS thermal performance is most sensitive to insulation thickness while its mechanical performance is more sensitive to outer honeycomb panel thickness and side thickness. and r show relatively less effects on structural responses and they are not supposed be optimized simultaneously with other variables.(ii)Response surfaces fit well in the first-stage optimization, and MTPS thermal performance is highly improved in this stage. Compared with thermal insulation performance, mechanical property requirements are easier to be satisfied for MTPS. Variables having stronger effects on responses are determined in this stage, and further improvement can be obtained only by reducing side thickness, which has met the lower boundary already.(iii)The second-stage optimization is applied to the outer honeycomb panel. Pareto frontier is constructed in this stage, and the optimum result is chosen from these Pareto-optimum solutions. By adjusting parameters of honeycomb core, the maximum temperature of substructure and mass of MTPS are further reduced.
Data Availability
The data used to support the findings of this study are included within the article.
Disclosure
This paper was presented in APISAT Technical Program.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
The authors would like to express their gratitude to all those who helped them during the writing of this thesis. A special acknowledgement should be shown to Prof. Li Shu, the supervisor of Qiuyi Xu, for his professional guidance and patience. The authors also wish to extend their thanks to Qiuyi Xu’s friends in the lab, Wang Yan and Meng Yang, for their encouragement and great support. This work was supported by the Industry-University Cooperation Project of Aviation Industry Corporation of China.