Design-oriented dynamic model of deployable fin under time-varying elevated temperature environment

Abstract

Coupling of clearance joint and harsh aerodynamic heating environment is an inevitable nonlinear factor in folding mechanism of the fin of high-speed aircrafts that remarkably modifies natural frequencies and modes of vibration from the initial design state. However, accurately predicting dynamic properties of deployable fin with full consideration of these effects is not common industry practice. A practical semi-analytical model based on Hertz contact theory and ESDU-78035 model is proposed in this study to investigate high-temperature connection stiffness of local hinged–locked mechanisms. Material property degradation and clearance variation caused by thermal expansion are comprehensively considered and quantified in this model. Vibration characteristics of the assembled deployable fin are then solved using finite element method (FEM). The real-time evolutionary process of thermal mode of the fin is discussed. And natural frequencies of fixed-value and time-varying connection stiffness are compared. The simulation results of this study demonstrate that the relative error of structure temperature between the sequential approach and fully coupled simulations is less than 6.98%. The connection stiffness (slope of the load-displacement curve) of the folding mechanism under high temperature conditions decreases by 3.52%, and the variation is mainly caused by the degradation of the elastic modulus of the material, while the clearance change due to the thermal expansion has no significant effect on the slope. The natural frequency of the deployable fin exhibits an inverse correlation with the temperature change trend, and the first three frequencies decrease by 1.67, 7.75, and 16.28 Hz compared to the initial value, respectively.

Nomenclature

a

half-contact width

C _pc

specific heat

d

distance between the spring and the symmetry plane of the inner fin

E

elastic modulus

f

natural frequency

F

force of the outer fin

H _e

boundary layer outer edge enthalpy

H _re

recovery enthalpy

H _w

wall temperature enthalpy

k

slope of the load-displacement curve

K

connection stiffness of the hinged–locked structure

L

axial length of the contact zone

M

moment of the outer fin

M _e

Mach number behind the oblique wave

m _f

fuel consumption rate

m

vehicle mass

P

vertical load

Q _f

heat flux of the fin

Q _b

heat flux of the reference point

R

radius of the curvature

t

time

T

material temperature

V _e

velocity behind the oblique wave

v

flight speed

Greek Symbol

$\upsilon$

Poisson’s ratio

$\delta$

contact deformation

$\chi$

reference factors

$\beta$

oblique wave angle

$\theta_c$

half cone angle

$\rho_c$

density of the material

$\lambda$

thermal conductivity

$\alpha$

angle-of-attack

$\theta$

trajectory inclination angle

1.0 Introduction

Studies have typically focused on the development of next-generation high-speed vehicles with main design specifications of high lift coefficient and ultralong range using enabling technologies, such as large wings or fins. However, modifying the vehicle shape is necessary due to the geometric compatibility of the launch platform. Folding configuration has been developed and extensively utilised in various aeronautic vehicle industries to address contradiction issues between large wing area of vehicle and space constraints of the platform. The division of the wing or fin into separate inboard and outboard segments in the span-wise direction of the out-of-plane deployable configuration presents advantages of low width dimension of the aircraft and thus reduced storage and transportation space [1, 2]. Hinged–locked structure is the core subcomponent of the deployable fin that exhibits typical nonlinear mechanical properties due to the existence of clearance. Local behaviour caused by assembly clearance significantly affects structural dynamic properties of the whole aerofoil. Based on these nonlinear effects, the operation of the high-speed vehicle within the atmosphere over a high Mach number throughout the entire flight further increases the complexity of the problem. The rapid increase in temperature of the body and fin surface due to the presence of extreme aerodynamic heating, especially during downthrust or high-manoeuvering trajectory conditions, further affects the body stiffness and nonlinear mechanical properties of connection structures [3]. Berry et al. [4] provided a review of the research contributions of NASA Langley Research Center (LaRC) to the experimental aerodynamics of slender and winged hypersonic vehicles, and further summarised the aerothermodynamic anomalies and lessons learned in the hypersonic flight programs of the United States. Uyanna et al. [5] reviewed the development and application of thermal protection technology in different periods, and discussed the current challenges, including mass efficient TPS materials, modelling and simulating tools, and sensor measurements systems. Dutch Space Center conducted an in-depth study of the extreme thermal environment during the reentry process of the European Space Agency’s hypersonic reentry test bed (EXPERT) vehicle [6]. The trajectory-based peak heat flux of blunt nose reached 1.69 MW/m², and the maximum temperature of the flaps is up to 2328 K. Furthermore, it is worth noting that Fatemi [6] specifically pointed out that the aerodynamic heating and force loading during the reentry process not only degrade the mechanical properties of the material, but also generate axial and radial gaps between the ceramic heat resistant material nose and metal fuselage under the effect the thermal mismatch. Determining the gap evolution characteristics and implementing effective thermal seal design are crucial.

The traditional design concept holds that the most severe thermal environment is at the aircraft nose cone and leading edge with a local temperature exceeding 1,500 °C. However, the design concept and requirement of thermal protection systems have gradually expanded with the many applications and extensive research on the deployable fin. The implementation of thermal protection systems on the folding mechanism remains challenging due to insufficient high-temperature resistance, cumbersome fabrication process, or limited scalability of existing technologies. The effective prediction of mechanical properties of the connection structure of deployable fin under the ultrahigh-temperature service environment and influence assessment of local connection structure heating on vibration properties of the overall fin in these circumstances have increasingly become potential challenging technical issues faced by engineers [7].

Many experiments have been performed and numerous technologies have been proposed from the extensive investigation of mechanical characteristics of wing and fin under elevated temperature [8]. On the one hand, many theoretical investigations have been performed. Chen et al. [9] developed a physical description of the thermal–structural–vibrational response of low-aspect ratio hypersonic wing, proposed a strategy based on adaptive time-step and hybrid interpolation methods to solve multiphysical coupling variable, and further analysed vibration properties of the wing under sustained aerothermodynamic loads. Shen et al. [10] established a multiphysical field full-time coupling analysis method for an all-movable wing along the trajectory and provided the analysis process of aerothermoelastic stability and transient response. On the other hand, the following experimental investigations have been conducted. NASA established the thermal mechanical coupling experiment test system at an early stage and thoroughly investigated the mechanical characteristics of aerofoil under elevated temperature conditions [11, 12]. Wu [13] and Wang et al. [14] recently proposed a novel experimental method to examine the thermal vibration characteristics of locally fixed wing-like structures in the design of high-speed aircrafts. He et al. [15] explored the nonlinear aerothermoelastic problem of folding fin with preset clearance, carried out a high-temperature thermal modal test on the basis of a thermal environment simulation device, and examined the coupling effect of temperature-dependent elastic modulus and thermal stress on flutter characteristics of folding fin under steady-state temperature field.

These studies demonstrated that spring elements are used in numerous ways to simulate the connection stiffness of hinged–locked structures. The expression of stiffness is mainly based on artificial assumptions, and the stiffness value is tuned through ground vibration test results or model updating technology for a specific system. However, the design-oriented nonlinear mechanical model that can reflect the relationship of geometric dimension, physical characteristics and connection stiffness under ultrahigh temperature requires further investigation. Specifically, aeroelasticity properties or control characteristics of the deployable fin may fail to meet the design requirements under a specific service environment. If the mechanical model describing physical characteristics of the connecting structure is unknown, then the improvement and redesign of the structure can be impossible. In addition to this limitation, the temperature of a real deployable fin is multidimensionally distributed under the condition of manoeuvering flight. Notably, connection stiffness is highly dependent on temperature. However, to the best of the authors’ knowledge, relevant studies typically ignore the influence of high-temperature conditions and clearances on the connection stiffness of the folding mechanism in real high-speed vehicles.

Because of the increased complexity of the problem, the following studies are of great significance to guide the design of folding fin: (1) improving the description and prediction ability of the mechanical theoretical model; (2) calculating the connection stiffness of the hinged–locked structure of the deployable fin, which is exposed to the time-varying aerodynamic heating; and (3) assessing whether frequency characteristics of the deployable fin are still within the effective design range under high-temperature environment given the increased complexity of the problem.

Figure 1. High-speed vehicle and deployable fin with different rotation angle.

Figure 2. (a) Metal-composite integrated deployable fin structure. (b) Partial hinge-locked structure.

Therefore, this study aims to develop a high-efficiency design-oriented mechanical model of out-of-plane deployable fin with hinged–locked connection structure by (1) constructing a semi-analytical model of stiffness characteristics of the local hinged–locked structure under time-varying high temperature to illustrate the influence of temperature rise on the connection stiffness; (2) analysing the temperature-displacement response of the fin surface and clearance variation characteristics of the connecting structure with typical heat flow input; and (3) exploring the vibration characteristics of the deployable fin under extreme thermal environment, followed by obtaining the time-varying frequency of the whole fin system.

2.0 Metal-composite integrated deployable fin

The metal-composite integrated deployable fin for the high-speed vehicle is illustrated in Figs 1 and 2. A new efficient thermal protection method is necessary because the existing metal all-movable fin may be unable to withstand harsh aerodynamic heating conditions. Although heat-resistant composite materials present excellent properties, such as low density, low thermal conductivity, ablation resistance and airflow shear resistance, they fail to meet the mechanical performance requirements of deployable fins due to their limitations, such as low machining accuracy and poor strength and stiffness. The high-performance composite material and metal skeleton are inlaid and connected innovatively in this study to incorporate the satisfactory ablation resistance capability of the composite material and toughness of the metal material comprehensively [16]. The resulting integrated fin can simultaneously achieve enhanced multiple functional properties of thermal protection and mechanical loading, as shown in Fig. 2(a).

The deployable fin specimen mainly consists of the fixed fin segment, movable fin segment, folding mechanism and in-position locking device. The movable fin is a symmetrical trapezoidal section composed of a heat-insulation layer and steel skeleton. The fixed fin is a solid alloy steel structure connected to the aircraft body through a fin shaft. The folding mechanism and in-position locking device between the inner fin and outer fin are illustrated in Fig. 2(b). The inner and outer fin structures are connected by a rotation pair to constrain the translational and rotational degrees of freedom of the outer fin, except for the degree of freedom rotating along the folding axis. The outer fin rotates around the folding shaft under the action of the preload torsion spring in the rotating shaft after the deployment command occurs. The locking pin is inserted into the inner hole under the action of the preload spring to complete the locking process when the outer fin is unfolded to the state that the pin is flush with the inner fin hole.

3.0 Mechanical model

A semi-analytical method is first derived in this section to predict the connection stiffness formed by the relationship of loads and displacements. FEM is then utilised to establish the model of the deployable fin. The hinge-locked connection structure is modeled using a bush element. The mechanical model of the whole deployable fin is finally developed after the connection stiffness extracted from the analytical model is used into the bush element.

It should be noted that this paper focuses on the analysis of the clearance nonlinear connection stiffness (as shown in Fig. 3), including the clearance region and the linear region (or the contact region). The slope of the load-displacement curve in the clearance region is 0 and the slope in the contact region is a finite value of k. The connection stiffness defined in this paper and the corresponding solution method mainly refer to the slope of the linear segment. This is different from the equivalent linear stiffness obtained using the descriptive function method.

Figure 3. Connection stiffness model with clearance.

3.1 Semi-analytical model of hinged–locked structure

The hinged–locked mechanism enables the outer-fin component to rotate downward. The hinged–locked mechanism basically provides free rotation around a folding axis through actuation torsion spring during the morphing process. However, the mechanism can be considered a fixed structure because it may prevent motion by controlling the locking pin. The connection stiffness at the locked state is the main focus of this investigation. The local coordinate system shown in Fig. 4 has been developed as a convenient framework to describe properties of the mechanical structure. The connection stiffness of the hinged–locked structure in six directions is denoted K _x , K _y , K _z , K _θx , K _θy and K _θz . The stiffness under the locked condition is specifically expressed as the ratio of the load increment of the outer fin (ΔF _x , ΔF _y , ΔF _z , ΔM _x , ΔM _y and ΔM _z ) to the corresponding displacement changes (Δx, Δy, Δz, Δθ _x , Δθ _y and Δθ _z ). Notably, the load and structure are considered symmetric or antisymmetric, and coupling effects between different loads are excluded from the following analysis.

Figure 4. Local reference coordinate system.

Figure 5. Analytical model for calculating K _θx . [17].

3.1.1 Rotational stiffness K_θx around the folding axis

The schematic of the hinged–locked structure showed that contributions to the connection stiffness mainly consists of the following parts: (1) deformation caused by interfacial contacts of pins and holes and (2) bending deformation of the locking pin, as demonstrated in Fig. 5 [17]. High deformation in each part corresponds to low rotational stiffness K _θx when the external load ΔM _x is constant. Therefore, each part of the deformation should be calculated separately and then combined to obtain the angle of outer fin and the connection stiffness.

i) Contact deformation

Clearances are reserved during the structural design to ensure the satisfactory operation state of the locking mechanism and avoid the pin from becoming stuck in the nonworking position of the hole. The preset clearance value between the locking pin and the hole wall of the folding mechanism is 1 mm. The vehicle is for one-shot use, thus, without considering the variation of clearance caused by friction and wear. Shape of the locking pin and the hole wall are sufficiently dissimilar, and the interface is nonconforming, thereby implying that the area of contact is remarkably smaller than characteristic dimensions of contacting bodies. The Hertz contact model [18] can be used to calculate the contact deformation as follows:

(1) \begin{align}\frac{1}{{{E_{{\rm{eq}}}}}} = \frac{1}{2}\!\left( {\frac{{1 - v_1^2}}{{{E_1}}} + \frac{{1 - v_2^2}}{{{E_2}}}} \right)\end{align}

(2) \begin{align}\frac{1}{{{R_{{\rm{eq}}}}}} = \frac{1}{{{R_1}}} + \frac{1}{{{R_2}}}\end{align}

(3) \begin{align}a = \sqrt {\frac{{8P{R_{{\rm{eq}}}}}}{{\pi {E_{{\rm{eq}}}}}}} \end{align}

(4) \begin{align}\delta (t) = \frac{{4P\!\left( {1 - v_1^2} \right)}}{{\pi {E_{{\rm{eq}}}}}}\!\left( {\frac{1}{3} + \ln \frac{{4{R_{{\rm{eq}}}}}}{a}} \right)\end{align}

where E _eq and R _eq are the equivalent elastic modulus and equivalent radius, respectively; E ₁ and E ₂ are the elastic modulus of the curvature of the cylinder and plane, respectively; υ ₁ and υ ₂ are the Poisson’s ratio of the curvature of the cylinder and plane, respectively; R ₁ and R ₂ are the radius of the curvature of the cylinder and plane, respectively; a is the half-contact width; δ is the contact deformation; t is the time; and P is the vertical load.

ii) Bending deformation of the locking pin

The contact area is smaller than the pin size and stresses are highly concentrated in contact areas due to the nonconforming contact between the locking pin and the hole wall. The beam model in [19] can be used to simplify the locking pin because local contact distortion exerts an insignificant effect on the macroscopic shape of bending deformation.

Because of multi-points contact between the locking pin and the hole wall, the number of unknown normal and tangential forces is more than the number of equilibrium equations for the system, that is, the mechanical system is non-stationary. For the non-stationary system, the real contact force should minimise the potential energy. Although the equations of Hertz spring model and beam model are analytical, a nonlinear optimisation algorithm is required to solve for the extremes of potential energy, thus we define it as a semi-analytical method. A more detailed discussion can be referred to [17].

3.1.2 Connection stiffness in other directions

The connection stiffness model presented in Fig. 6 is established to solve the translation stiffness values K _x , K _y and K _z as well as the rotational stiffness values K _θy and K _θz . The tensile load carried by the outer fin in the model is transferred through the cooperation of the shaft and hole. The quasi-rigidity assumption is adopted for the body and only the resistance of connections in contact areas is considered while analysing the connection stiffness because bodies provide no or very minimal additional stiffness to the system. The contact stiffness between the shaft and the outer fin hole are denoted k ₁ and k ₂, respectively, and the contact stiffness between shaft and the inner fin hole is denoted k ₃. The overall stiffness of the hinged–locked structure is then represented by a combination of these three parts, k ₁ – k ₃ as the outcome of the analysis.

Figure 6. Contact model of the shaft and hole connection.

The ESDU-78035 model [20] is used to solve the contact stiffness k ₁ – k ₃ because the clearance in the rotation pair is small. The ESDU-78035 model is expressed as follows:

(5) \begin{align}\delta (t) = {F_N}\!\left( {\frac{{{\sigma _i} + {\sigma _j}}}{L}} \right)\left[ {\ln\! \left( {\frac{{4L\left( {{R_i} - {R_j}} \right)}}{{{F_N}\left( {{\sigma _i} + {\sigma _j}} \right)}}} \right) - 1} \right]\end{align}

where σ _i(j) = (1−υ _i(j) ²)/E _i(j); υ _i and υ _j are the Poisson’s ratio of materials of the hole and shaft, respectively; E _i and E _j are the elastic modulus of materials of the hole and shaft, respectively; R _i and R _j are the outer hole and inner shaft radii, respectively; F _N is the total normal load in the contact area; L is the axial length of the contact zone; and t is the time. Dimensions and material properties of different subcomponents along the flight time are used into Equation 5 to obtain the contact stiffness k ₁ – k ₃.

Figure 7. FEM the deployable fin and diagram of the zone partition.

The structural connection stiffness is then calculated under the following assumptions:

1. 1) The translation stiffness K _y and K _z are treated the same, in which springs of k ₁ and k ₂ are parallel and connected in series with k ₃. Thus, the translation connection stiffness K _y (or K _z ) is expressed as follows:

(6) \begin{align}{K_y}\left( {or\,\,{\rm{ }}{K_z}} \right) = \frac{{\left( {{k_1} + {k_2}} \right){k_3}}}{{{k_1} + {k_2} + {k_3}}}\end{align}

1. 2) Assuming that the rotational stiffness K _θy and K _θz are the same and contributions of k ₁ and k ₂ to the torque are considered, then the rotational connection stiffness K _θy (or K _θz ) is expressed as follows:

(7) \begin{align}{K_{\theta y}}\left( {or\,\,{\rm{ }}{K_{\theta z}}} \right) = {k_1}{d_1} + {k_2}{d_2}\end{align}

where d ₁ (or d ₂) is the distance between the spring k ₁ (or k ₂) and the symmetrical plane of the inner fin.

1. 3) Determination of the translation stiffness K _x is indirectly and dominated by the interface between the inner and outer fins. The influence of the body stiffness is significantly greater than that of the contact stiffness due to the large contact area; hence, the contact stiffness K _x is set to 1 × 10¹⁰ N/m as a relatively large value.

3.2 Numerical simulation model

3.2.1 FEM

Finite element analysis is conducted using ABAQUS/Standard to model the complex structure of the inner and outer fins, and quadratic tetrahedral thermal displacement coupling elements (C3D10MT) are assigned (Fig. 7). To verify the mesh irrelevance, a comparison between extremely coarse, coarser, fine and extremely fine is studied based on the variation of temperature of locking pin region versus time. The number of elements in each case is mentioned in Table 1. The variation of temperature versus time for different grid sizes is shown in Fig. 8. The fine mesh with 129,561 elements is chosen based on the study. The metal skeleton and composite heatshield of the outer fin model are established separately and then tied together regardless of the interface layer. Thus, displacement and temperature are considered continuous. Nodes near contact areas of inner and outer fins are rigidly constrained to the control point at the centre of the rotating shaft through MPC. A bush element is established between control points to simulate the connection relationship, and the nonlinear stiffness obtained using previous analytical methods is assigned to the bush element.

Table 1. Number of elements in each case

Table 2. Properties of superalloy steel

Figure 8. Variation of temperature of locking pin region versus time by using different meshes.

The following material properties are assigned to different subcomponents: superalloy steel for the inner fin and skeleton, hardened steel for pins and shaft, and high-performance composite material for heatshield. Degraded material properties are used in the model because of the evident rise of temperature in the entire fin. Tables 2–4 present material properties at elevated temperatures obtained through experimental measurement. The density of the superalloy steel, hardened steel and composite material is set to a constant value of 8,240, 7,780, and 1860 kg/m³, respectively, given that high temperature exerts a negligible influence on the behaviour of density. In addition, mechanical properties of the high-performance composite material in this study are basically unchanged in the analysis temperature range; hence, the constant value of elastic modulus, Poisson’s ratio and thermal expansion coefficient are utilised and emissivity is set to 0.8.

Table 3. Properties of hardened steel

Table 4. Properties of composite thermal insulation materials

3.2.2 Boundary conditions

The root of the all-moveable fin axis of a real aircraft is generally connected to the actuation and transmission system, such as rods, rocker arms and electric motor. The detailed simulation of these complex connections is beyond the scope of this study. Therefore, the simplified fixed displacement boundary condition is utilised through multipoint constraint (MPC) method, which is similar to the resolvent work in Reference [20].

Determining the unsteady heat flux distribution is a fundamental task of heat transfer analysis. Different heat flux and radiation values are applied as input thermal loads or boundaries for various surface zones in FEM to simulate flight heating conditions reasonably. Specifically, the windward and leeward surfaces of the fin are denoted A and B, respectively, and then further divided into wedge (denoted 1) and plane (denoted 2) areas. Therefore, a total of four partitions (A1, A2, B1 and B2) are set, as shown in Fig. 7.

It should be noted that, for the wing or fin of high-speed vehicle, thermal protection structures of the leading-edge region, wedge region and planar region are usually designed separately. Although the distribution of heat flux is not uniform in each region, the same thickness of thermal protection structure is used through each region because of the processing characteristics. Generally, the location with highest heat flow in each region is utilised, ensuring that the designed thickness of heat protection system can meet the requirements along the whole flight profile. Selecting more points and dividing more regions can achieve better performance and obtain detailed thermal evaluation results. However, using existing conservative design methods, that is, only the highest heat flow point in each region is considered during the preliminary design process, is efficient and sufficient in engineering applications. The superiorities and disadvantages of the similar design methodology are specifically discussed in the literature [21].

3.3 Analysis methodology of heat flux

3.3.1 Cold-wall heat flux

Accurately predicting aerodynamic heating flux is the critically important task for the high-speed vehicle design. Generally, based on whether the shape is regular, the surface of the whole vehicle can be divided into large areas and local areas. Large areas include conical and cylindrical parts of aircraft, and the heat flux can be obtained by utilising engineering algorithms conveniently. Furthermore, the estimating formula can be combined with the trajectory to calculate the heat flow at different flight times. The local structural area includes the wing or fin structures. Due to the complex flow fields of local structures, it is often difficult to accurately predict heat flow using empirical formulas. Although wind tunnel test or CFD simulation can be used to obtain the heat flux conditions, obtaining the heat flow through the entire flight process faces the obstruction of computational complexity, time alocation and poor economy.

In this paper, a method (named interference factor method) that combines CFD numerical simulation and engineering algorithms is conducted. A more detailed description of the interference factor analysis method can be found in reference [22]. The analysis process and main calculation formulations are briefly listed as follows:

1. 1) Select a reference point in a large area, and select a point in a local surface, such as different region of fin. Define reference factors χ:

   (8) \begin{align}\chi = {\raise1.5ex\hbox{${{Q_f}}$}\!\! \bigg/\!\! \lower1.5ex\hbox{${{Q_b}}$}}\end{align}

Where: Q _f is the heat flux density of the fin, and Q _b is the heat flux of the reference point.

1. 2) By using STAR-CCM+ software [23], reference factor interpolation tables under different flight conditions, such as flight altitude, Mach number, angle-of-attack and fin yaw can be obtained.

2. 3) The reference enthalpy method [24] is used to calculate the heat flux Q _b of the reference point through the whole flight trajectory. The formula is as follows:

Calculate the oblique wave angle β:

(9) \begin{align}\sin \beta = \sin {\theta _c}{\left( {\frac{{\gamma + 1}}{2} + \frac{1}{{M_\infty ^2{{\sin }^2}{\theta _c}}}} \right)^{\frac{1}{2}}}\end{align}

Where, θ _c is the half cone angle.

Calculate the physical parametres behind the shock wave:

(10) \begin{align}{H_e} = \frac{{\left( {7M_\infty ^2{{\sin }^2}\beta - 1} \right)\left( {M_\infty ^2{{\sin }^2}\beta + 5} \right)}}{{36 \cdot {{\left( {{M_\infty }\sin \beta } \right)}^2}}}{H_\infty }\end{align}

(11) \begin{align}{M_e} = \sqrt {5\left( {\frac{{{T_\infty }\left( {1 + 0.2M_\infty ^2} \right)}}{{{T_e}}} - 1} \right)} \end{align}

(12) \begin{align}{V_e} = {M_e} \cdot 20.0467 \cdot \sqrt {{T_e}} \end{align}

Where, H _e , M _e , and V _e are the static enthalpy, Mach number and velocity behind the oblique wave, respectively.

The Sutherland method is used to obtain viscosity coefficient on the outer edge of boundary layer:

(13) \begin{align}\mu _e^* = 1.458 \times {10^{ - 6}}\frac{{{{\left( {H_e^*} \right)}^{1.5}}}}{{H_e^* + 110.4}}\end{align}

Calculate the reference enthalpy based on the Eckert reference enthalpy method:

(14) \begin{align}H_e^* = 0.28{H_e} + 0.5{H_w} + 0.22{H_{re}}\end{align}

Where, H _w is the wall temperature enthalpy, H _re is the recovery enthalpy, and H _e is the boundary layer outer edge enthalpy.

Calculate the heat flux density at the reference point Q _b

(15) \begin{align}{Q_b} = 0.332\sqrt 3\Pr {^{ - \frac{2}{3}}}{\left( {\frac{{{\rho ^*}{\mu ^*}}}{{{\rho _e}{\mu _e}}}} \right)^{ - 0.5}} \times {\rho _e}{\mu _e}{\left( {0.5{R_{ex}}} \right)^{ - 0.5}} \times \left( {{H_r} - {H_w}} \right)\end{align}

1. 4) Predict the heat flux Q _f of the fin surface using the flowing equation

(16) \begin{align}{Q_f} = \chi \cdot {Q_b}\end{align}

3.3.2 Hot-wall heat flux

Because heat conduction is mainly along the thickness direction, one-dimensional heat transfer model can be used to simplify the analysis. The one-dimensional analysis program is integrated into cold-wall heat flow simulation program, and the outer surface temperature of the structure is iteratively calculated. Thus, the hot-wall heat flux can be obtained.

The one-dimensional heat conduction equation is shown as follows:

(17) \begin{align}{\rho _c}{C_{pc}}\frac{{\partial T}}{{\partial t}} = \frac{\partial }{{\partial y}}\left( {\lambda \frac{{\partial T}}{{\partial y}}} \right)\end{align}

Where, y is the coordinate along the thickness direction, ${\rho _c}$ is the density of the material, ${C_{pc}}$ is the specific heat, and T is the material temperature, and λ is the thermal conductivity.

For the metal-composite integrated structure, the governing equation is:

(18) \begin{align}{\rho _{c,i}}{C_{pc,i}}\frac{{\partial T}}{{\partial t}} = \frac{\partial }{{\partial y}}\left( {{\lambda _i}\frac{{\partial T}}{{\partial y}}} \right)\end{align}

Where, the subscript i represents the material parametres of the i th layer.

For the outer surface, the thermal radiation is considered, and heat transfer balance equation is:

(19) \begin{align}Q - \varepsilon \sigma T_w^4 = - \lambda {\left( {\frac{{\partial T}}{{\partial y}}} \right)_{y = 0}}\end{align}

The adiabatic condition is used for inner surface, and the equation is:

(20) \begin{align}0 = - \lambda {\left( {\frac{{\partial T}}{{\partial y}}} \right)_{y = \delta }}\end{align}

The heat flow equilibrium and temperature conditions of inter-surfaces satisfied:

(21) \begin{align} - {\lambda _1}\left( {\frac{{\partial {T_1}}}{{\partial y}}} \right) = - {\lambda _2}\left( {\frac{{\partial {T_2}}}{{\partial y}}} \right)\end{align}

(22) \begin{align}{T_1} = {T_2}\end{align}

4.0 Verification of sequential approach to generate thermal boundary conditions

A head cone of the vehicle is analysed to discuss the possible loss of accuracy due to adopting the sequential approach to generate thermal boundary conditions instead of fully coupled simulation (CFD/CSD method). The cold-wall and hot-wall heat flux, obtained by using the two methods, are compared. The point is location at the profile x = −1 m from the head of the vehicle, on the windward side.

4.1 Verification of cold-wall heat flux

The following flight condition is used for the simulation: the flight speed is 1000 m/s, the altitude is 10 km, the flow temperature is 223.15 K, the pressure is 26417.3 Pa, the speed of sound is 299.38 m/s, the flight Mach number is 3.34 and the AOA is 0°. In the CFD model of STAR-CCM+, the K-Omega SST turbulence model is selected, and the Coupled Inviscid Flux is set as AUSM. In the cold wall heat flux calculation, the vehicle surface temperature is set to a constant values of 300 K. The head tip of the vehicle is defined as the coordinate origin (x = 0, y = 0, z = 0), and a profile is established at the x = −1 m position of the vehicle head. The heat flux calculation results of the CFD method for this profile are quantified and compared with the results using the engineering heat flux analysis method.

Figures 9(a) and 9(b) present the geometric model used in STAR-CCM+, and the heat flux distribution on the whole surface of the vehicle, respectively. Figure 9(c) shows the heat flux at the profile position by utilising CFD method and reference temperature method. The results show that the heat flux using CFD method at different circumferential angles θ floats at the mean value of 190.14 kW/m², while the results obtained by the engineering method is 197.08 kW/m². The relative error is 3.65%, and it meets the accuracy requirement for the thermal environment design.

Figure 9. (a) Geometry and boundary conditions of CFD model; (b) cold-wall heat flux of the whole vehicle; (c) cold-wall heat flux by using CFD method and reference temperature method.

To further analyse the influence of the angle-of-attack on the calculation accuracy, the models under different AOAs from 0 to 20 degree are established and analysed by CFD method and engineering method, respectively. It should be noted that the results of engineering method represent the windward heat flux (the maximum heat flux value of the profile). The cold-wall heat flux and the relative errors of the two methods are presented in Fig. 10 and detailed data is shown in Table 5. The max relative error is 9.81%, and a good agreement is achieved.

Table 5. Cold-wall heat flux and relative errors between CFD and engineering results at different AOAs

Figure 10. Cold-wall heat flux and relative errors between CFD and engineering methods at different AOAs.

4.2 Verification of hot-wall heat flux

The heat transfer within the structure and radiation of outer surface should be considered in the hot-wall heat flux simulation. The thickness of the thermal insulation materials and steel metal shell is 3 mm and 10 mm, respectively. The flight condition is same with the previous case of cold-wall heat flux simulation. Based on the sequential engineering approach and CFD/CSD coupling method, Figs.11 and 12 present the hot-wall heat flux and temperature response during the flight. The detailed data is shown in Tables 6 and 7, where Q _HE and Q _HCFD represent the hot-wall heat flux by using the engineering method and CFD/CSD method, respectively. T _OE and T _OCFD represent the outer-surface temperature. T _IE and T _ICFD represent the inner-surface temperature.

Figure 11. Hot-wall heat flux results by using CFD/CSD and engineering methods under constant flight conditions.

Table 6. Hot-wall heat flux results of the two methods and relative error

Table 7. Temperature response results of the two methods and relative error

Figure 12. Temperature response of the structure by using CFD/CSD and engineering methods under constant flight conditions.

The results show that:

1. 1) The trends of predicted hot-wall heat flux by utilising the engineering method and CFD/CSD coupling method are in good agreement. The relative errors are 5.23%, 10% and 29.14%, at 50 s, 200 s and the end stage, respectively.

2. 2) The temperature response within the structure calculated by the engineering method and the CFD/CSD method are in good coincident, and the maximum relative error is within 10%, which indicates the correctness and usability of the sequential method.

3. 3) Hot-wall heat flux presents a significant decreasing trend versus time, while the cold-wall heat flux is constant due to the assumption of isothermal wall surface. Therefore, it is necessary to consider the heat conduction within the structure in the thermal boundary calculation.

It should be noted that the hot-wall heat flux at the tip of nose cone is not the highest (as shown in Fig. 11), which is because the thickness of the thermal protection materials is larger than that at other areas in the model.

5.0 Results and discussion

5.1 Thermal response

5.1.1 Flight profile

The heat flux is closely related to flight speed, altitude, aircraft attitude and aircraft shape. The high-speed aircraft is regard as a particle, and a 3-degree of freedom aircraft dynamics model (Equation (23)) is used to solve the flight profile considering gravity, thrust and aerodynamic parametres. And then, the heat flux can be obtained by using the interference factor method as mentioned above.

(23) \begin{align}m\frac{{dv}}{{dt}} &= F\cos \alpha - D - mg\sin \theta \nonumber\\mv\frac{{d\theta }}{{dt}} &= F\sin \alpha + L - mg\cos \theta \nonumber\\\frac{{dx}}{{dt}} &= v\cos \theta \\\frac{{dy}}{{dt}} &= v\sin \theta \nonumber\\\frac{{dm}}{{dt}} &= - {m_f}\nonumber\end{align}

Where, v is the flight speed, α is the angle-of-attack, F is the thrust, m _f is the fuel consumption rate, m is the vehicle mass, θ is the trajectory inclination angle, L is the lift force and D is the drag force.

The angle-of-attack is designed according to the desired goal. The entire flight process of the trajectory can be divided into three stages:

1. 1) After the vehicle is launched, the engine ignites, and the missile undergoes a climb process in which the fuel burns out. Thereafter, the trajectory is a parabolic trajectory similar to a ballistic missile.

2. 2) After re-entering the atmosphere, the projectile does not directly dive to hit the target. The vehicle is pulled up again, achieving a jumping manoeuver. Relying on the lift force without power, the second parabola is completed.

3. 3) Due to the fact that the trajectory undergoes a second takeoff, the speed of the projectile during reentry is somewhat lower than that of the first, and the trajectory tends to flatten under the effect of dense atmospheric lift. Therefore, it is necessary to make the projectile bow and dive after the second reentry.

The variation equation of AOA is:

(24) \begin{align}\alpha \left( t \right) &= - 4{\alpha _{\max ,j}}Z\left( {1 - Z} \right) \nonumber\\[6pt] Z &= {e^{ - {a_j}\left( {t - {t_j}} \right)}},\,j = 1,2,3\end{align}

Where, subscript j represents the number of jumps; α _max,j is the maximum allowable angle-of-attack for the j th jump; and the design value a _i can adjust the speed of the turn; t _i is the jump start time.

The calculated flight profile, AOA, cold-wall heat flux and hot-wall heat flux are shown in the following Fig. 13.

Figure 13. Flight profile and heat flux results during the flight.

Results of heat flux (Figs13(c) and 13(d)) show that: in the initial stage, the heat flux reaches the first peak due to the low altitude and the large AOA. Then the heat flux is relatively small due to the thin air at a higher altitude. When the vehicle performs the second pull-up manoeuver with large AOA, the heat flux increases obviously again. However, when the vehicle performs the final dive, although the altitude is relatively low, and the angle-of-attack attitude is relatively large, the velocity of the end of the flight is relatively small. The heat flux is therefore not as obvious as the first two peaks.

5.1.2 Response of thermal displacement coupling

The heat flux of A1, A2 and B1 is shown in Fig.13(d), and the flux of B2 is 75% of that of B1. The temperature variation of the metal skeleton near the locking pin region versus time is shown in Fig. 14. The temperature rises faster in the initial state due to the peak heat flux and decreases slightly because of the high-temperature outer surface’s thermal radiation at the end of the flight. Two nodes near the locking hole are selected, and the distance change (Fig. 15) in the upper and lower surfaces of the locking hole obtained by subtracting the displacement of the two nodes reveals the clearance variation.

Figure 14. Temperature of metal skeleton near the hole of locking pin versus time.

Figure 15. Clearance change of the locking hole versus time.

5.2 Connection stiffness

5.2.1 Rotational stiffness K_θx

The clearance and temperature obtained by previous analyses are assigned into the semianalytical model of the hinged–locked structure.

First, the reference temperature of the contact zone is determined via the average temperature of the node. Temperature-dependent material properties are used in the Hertz model (Equation 1) to calculate contact deformation.

Second, the temperature rise of the whole locking pin is assumed unclear. Therefore, properties at ambient temperature are used in the beam model to obtain the bending deformation of pin while excluding thermal stress.

Finally, the obtained K _θx values are illustrated in Fig. 16.

The (i) clearance variation due to thermal expansion and (ii) material property degradation, which considers only one factor, are compared to determine the influence of two key factors on the connection stiffness (Fig. 16). The connection stiffness considering both factors K _v, the result of considering the material factor alone K _vm, the result of considering the clearance factor alone K _vc, the reference value K _v, relative error of K _vm (Re _vm) and relative error of K _vc (Re _vc) are listed in Table 8.

Table 8. Values of connection stiffness K_θx at typical moments

Figure 16. Connection stiffness K _θx versus time.

The results in Fig. 16 and Table 8 present the following:

1. 1) The connection stiffness changes nonlinearly with time and shows an inverse correlation with the trend of temperature. The initial state of K _v (1.278 × 10⁵ Nm/rad) is considered the reference, and the value at the end time of heating is 1.233 × 10⁵ Nm/rad, which is 3.52% lower than the initial value.

2. 2) The final K _v is 1.270 × 10⁵ Nm/rad when only the clearance effect is considered. This finding showed a decrease of only 0.59% compared with the initial value. The effect of clearance is insignificant in this analysis because the corresponding expansion of the metal hole at the end is only 0.245 mm, which is smaller than the initial preset clearance of 1 mm.

3. 3) Table 7 also shows that Re _vm and Re _vc are less than 0.73% and 3.41%, respectively, during the heating process. For the connection model with a preset clearance ensures the normal operation of the mechanism, material property degradation due to high temperature is the main cause of the nonlinear variation of the connection stiffness, and the effect of the clearance change caused by thermal expansion is insignificant.

5.2.2 Connection stiffness in other directions

Translational (K _y and K _z ) and rotational (K _θy and K _θz ) stiffnesses are mainly determined by contact properties between the shaft and holes. The temperature of the wall surface and the clearance between the hole and the shaft are obtained and then assigned into Equation (5) to solve the connection stiffness.

Figure 17 shows the variation of K _y , K _z , K _θy , and K _θz and the downward trend in each curve because of the material property degradation. The stiffness of the translational connection is 2.288 × 10⁸ and 2.117 × 10⁸ N/m at the initial and final moments, respectively. If the initial value is used as a reference value, then the final translational stiffness decreases by 7.47%. In addition, the rotational stiffness is approximately two orders of magnitude smaller than the translational stiffness, as determined by d ₁ and d ₂.

Figure 17. Connection stiffness versus time (K _y , K _z , K _θy and K _θz ).

5.3 Vibration

Two analytical approaches exist for the study of the vibration characteristics of a structure with a clearance nonlinear stiffness connection:

1. 1) Consider that the vibration amplitude of folding fin is large. Vibration amplitude of the outer fin includes the gap region and the linear region (Fig. 18(a)). For example, the dynamic response and flutter analysis of folding fin without preload under medium magnitude excitation belongs to this case.

   

   Figure 18. (a) vibration of medium amplitude without preload; (b) vibration of small amplitude with preload.

2. 2) Consider that the folding fin is undergoing a large preset aerodynamic load and a small excitation value. The vibration amplitude of the rudder surface is only in the linear zone (Fig. 18(b)). The folding fin vibrate approximately in the linear stiffness range of small amplitude. For example, while evaluating the natural vibration characteristics of a fin with gap, the outer fin is generally preloaded and then excited by a small amplitude of shaker. The operating range is near the linear segment.

The following thermal modal frequencies discussed in this paper belongs to the second case. The connection stiffness values in the FEM are adopted from the linear segment of the clearance nonlinear stiffness model. The natural frequency is regard as small amplitude vibration in the linear region under preload conditions.

Figure 19. Theoretical prediction of the natural frequency of the deployable fin at typical moments.

5.3.1 Thermal mode

Numerical test cases are developed to investigate the contribution of connection stiffness and material property degradation to the natural frequency of the fin. Two cases are notable, and material property degradation is considered in both cases. The effect of connection stiffness and temperature-dependent material properties are preserved in the first case because the time-varying connection stiffness is used in model. Meanwhile, the constant connection stiffness is used in the second case of interest during the whole heating process (the connection stiffness is maintained equal to the initial value), and the material effect of the fin body is included as the only variable parametre. Figure 19 shows the variation of the first three natural frequencies of the deployable fin. Values at typical moments, including the frequency under variable connection stiffness f _v, frequency under constant connection stiffness f _c and relative error between two stiffness model Re, are listed in Table 9.

Table 9. First three modal frequencies of the deployable fin at typical moments

The results in Fig. 19 and Table 9 indicate the following:

1. 1) Natural frequencies of each order show a downward trend during the heating process. The first three orders of natural frequencies corresponding to the thermal mode are 1.67, 7.75 and 16.28 Hz lower than that of the ambient temperature for the time-varying stiffness model. The first-order mode shape of the deployable fin is an approximate rigid body mode shape with constraints and the second- and third-order mode shapes are self-elastic deformation.

2. 2) The frequency of the time-varying stiffness model is significantly lower than that of the fixed-value spring model because the connection stiffness of the spring decreases at high temperatures.

3. 3) The first-order mode is rigid, the connection stiffness plays a dominant role, the third-order frequency is the elastic deformation of the fin surface and the connection stiffness exerts a marginal effect on the deformation area far from the connection area. Hence, the relative error of the first frequency between the fixed and time-varying stiffness models is 1.46%, while the relative error of the third frequency is only 0.46% at the end of time.

5.3.2 Contribution rate of connection stiffness to frequency

The fixed connection stiffness results only reflect the material property degradation of the fin body, while the time-varying connection stiffness results comprehensively reflect the combined effect of the weakening of the connection stiffness and the variation of material properties. Frequency changes Δf _v and Δf _c relative to their respective initial frequencies are calculated, and the contribution rate index of the connection stiffness change to the frequency is defined as follows:

(25) \begin{align}C = \frac{{\Delta {f_{\rm{v}}} - \Delta {f_{\rm{c}}}}}{{\Delta {f_{\rm{v}}}}} \times 100\% \end{align}

where a large contribution rate index indicates a strong influence of the connection stiffness.

Figure 20 shows the index C at different times. The contribution rate of the connection stiffness to the first-order frequency fluctuates slightly but is maintained at 50%–60% as the heating progresses. The contribution rate to the second- and third-order frequencies is about 30% and 10%, respectively. The results revealed that the change of connection stiffness affects the first-order frequency the most evidently, followed by the second- and third-order frequencies.

Figure 20. Contribution rate of the connection stiffness on frequency at different times.

6.0 Conclusions

1. 1) The connection stiffness K _θx based on Hertz theory and the beam model expresses the connection properties of rotation around the folding axis. The connection stiffness in other directions (containing the translation K _y and K _z and rotational K _θy and K _θz stiffness values) is also presented using the ESDU-78035 model and the combination relationship of each contact area. The theoretical model proposed in this study can innovatively reveal the influence of temperature on connection stiffness because of its advantages of clear mechanical mechanism, strong physical interpretability and high computational efficiency.

2. 2) The final connection stiffness of hinged–locked mechanisms decreased by about 3.52% of the initial value under elevated temperatures. The comparative analysis of connection stiffness showed that the variation of clearance caused by thermal expansion exerts a marginal effect because the connection stiffness (slope of the load-displacement curve) of hinged–locked mechanisms with preset clearance under high temperature is mainly determined by material property degradation. Therefore, the connection stiffness at elevated temperatures in the design can be estimated by multiplying the connection stiffness under ambient temperature and the correction coefficient of the material elastic modulus.

3. 3) The first three orders of thermal natural frequencies of the fin decrease by 2.25%, 3.45% and 3.73% compared with that at ambient temperature. The comparison between the time-varying spring and fixed-value stiffness models demonstrated that the first-order frequency is dominantly affected by the connection stiffness while the third-order frequency is mainly affected by material property degradation of the fin body.

The method and results in this study can provide a theoretical basis for the analysis of thermal vibration characteristics of the deployable fin for high-speed vehicles. The framework of this investigation provides an effective design and evaluation criterion for heat insulation systems of deployable mechanisms and offers a reference for explaining anomalies of thermal modal experiments. The proposed method can be utilised in the analysis of similar deployable fins or folding wings.