Abstract
Hypersonic aircraft have been rapidly developed in recent years both theoretically and experimentally. Aerothermoelastic simulation is very challenging due to its inherent complexity, but physical tests are a workable approach. Flutter tests with variable speed are a popular alternative to hypersonic tests which provide nonstationary structural response data. This paper proposes a nonstationarity assessment method based on energy distribution in the time-frequency domain. The proposed method reveals the nonstationarity level corresponding to the appropriate modal identification algorithm or flutter boundary prediction (FBP) method. Several classic flutter criteria are utilized to build a hypersonic aircraft FBP framework. Numerical simulation and experimental applications demonstrate the effectiveness and feasibility of the proposed method, which facilitates accurate flutter predictions for the subcritical turbulence response during hypersonic flutter flight.
1. Introduction
Flutter quality is a significant issue in the development of hypersonic aircraft. Problematic aspects of hypersonic flight include extreme aerodynamic heating and complex hypersonic flow field effects [1–3]. Numerical simulation is an economic (and popular) research method, but the environmental conditions of physical flight experiments are rarely simulated. The causes of nonlinearity and complexity in hypersonic flows, such as real gas effects, boundary layer-shock wave interactions, changes in material properties due to aerodynamic heating, and nonlinear aeroelastic interactions between the flow field and structures [4, 5], as well as the monolithic coupling of aerodynamic thermal characteristics, aerodynamic force, elastic force, and inertia, are highly complex. The partitioned simulation method is currently the most often applied approach to hypersonic aircraft research, but it unfortunately requires that weak coupling characteristics be ignored per the limited capability of existing software programs.
Many hypersonic aerothermoelastic experimental studies have been conducted over the past several decades [6–9]. The flutter behavior of various wing and rudder models has been determined (e.g., airfoil shape, airfoil thickness, attack angle, leading-edge bluntness, and wing-body interference) [10–13]. Existing test and data processing methods are a workable foundation for further theoretical development hypersonic aircraft simulation techniques. “Flutter” is a phenomenon characterized by system instability due to structural modal coupling [14–16] as evidenced by the structural response of undamped oscillation. In a hypersonic flutter test, structural responses directly reflect the process of modal change and energy alternation; the responses provide useful indicators of various flutter trends. Structural response signal analysis can also help to further improve both simulation techniques as well as innovative structure designs. There is, to this effect, urgent demand for effective hypersonic flutter flight test data processing methods.
Flutter flight tests are applied to conventional aircraft by launching several flight stages with different airspeeds to approximate the aerodynamic force as a stationary excitation. Several algorithms for modal parameter identification or flutter boundary prediction (FBP) have been established based on the stationary stochastic process theory; they include fast Fourier transform (FFT) [17, 18], random decrement technique (RDT) [19], natural excitation technique combined with the Eigen system realization algorithm (NExT-ERA) [20], time series analysis based on the auto-regressive (AR) model [21], stochastic subspace identification (SSI) methods [22], and Hilbert-Huang Transform method [23].
Hypersonic aircraft have extremely short acceleration periods and very high Mach numbers. Flutter tests, however, require an efficient, effective progression of variable airspeeds. The structural response signal from such tests are nonstationary, and the solutions must be distinguished per the stationary situations. When there is low nonstationarity, the algorithms discussed above may be still feasible; otherwise, new methods for nonstationary signal processing are necessary for accurate estimation. There is not yet a standardized definition of nonstationarity—or any standard method for assessing it—due to inherent deficiencies in the mathematical theory.
This paper proposes a novel assessment theory for nonstationarity accompanied by classic flutter boundary prediction methods. A data processing framework for structural response was designed accordingly, and numerical simulations and experiments validate the proposed method and framework.
2. Nonstationarity Level Assessment
The nonstationarity level of the target signal directly affects the signal processing method and accuracy. There are several important aspects of nonstationarity, as discussed below.
(a) Degree of Cyclostationarity [24]. Cyclostationarity is defined by the distance between the autocorrelation function of periodic stationary signals and the generalized stationary autocorrelation function :That is, it is the integral of the least square error between normalized and according to the delay time .
(b) Hilbert Time-Frequency Spectrum [25]. The Hilbert spectrum of the signal is first analyzed on the time-frequency plane; then a marginal spectrum is obtained:
The nonstationarity level is defined as follows:where nonstationarity is defined as a function of frequency, because certain frequency components in the signal may be nonstationary while others remain stationary.
(c) Correlation Integral [26]. Consider a time series with the following correlation integral:where N is the total number of samples, is the norm space, R is the reference distance, and H is the unit step function.
For the stationary signal , its distance sequence is the normal distribution. The mathematical expectation of is and the standard deviation is . The nonstationarity of the assessed signal with the reference distance is defined as follows:
The first of the above definitions relates to periodic nonstationary signals without generality; the second is computationally complex and does not have clear boundaries for stationary signals; the third applies only to specific engineering applications, which is hardly generalizable to the analysis of normal nonstationary signals.
For any signal , the time-frequency distribution based on the Hermite window function iswhere is the order Hermite function, defined asPer the definition of “stationary”, the first moment and secondary moment are independent of time:where is the time delay and is the mathematical expectation.
The local time-frequency distribution must be equal to the mean of the global distribution at any time point. Their relationship iswhere is the time corresponding to the local time spectrum, is the frequency of the local time spectrum, and is the expectation of distribution.
For a nonstationary case, the energy distribution in the time-frequency domain is at and the entire energy distribution of signal is constituted as a set . If is a nonstationary signal, its time-frequency power spectral density will vary over time. A quantized value is proposed here to represent energy fluctuations in the time-frequency domain. The variance of indicates the level of nonstationarity.where is the mathematical expectation; is the calculation method for the distance between the energy distributions, which is defined as follows:Equation (12) can be used to describe the deviation of the energy distribution of the measured signal with that of the whole time series using the distance concept in the time-frequency domain, which is consistent with the level of nonstationarity. This definition is widely applicable because it is based on the basic theory of time-frequency analysis.
A flow chart of the entire nonstationarity level assessment process is provided in Figure 1.
3. Hypersonic Flutter Analysis Framework
The feasibility of common flutter analysis methods for hypersonic flutter tests is discussed in this section. There are two main structural response analysis methods as applied to flutter tests: structural modal damping and system stability. The former extrapolates the flutter critical velocity from the subcritical response by fitting the dangerous modal damping; the latter obtains the flutter boundary by establishing a dynamic model and calculating the stability parameters.
(a) Structural Modal Damping. The concept of modal damping in flutter test analysis is not limited to static structural damping but applies to mixtures of structural damping, aerodynamic damping, material damping, and all the physical quantities which suppress structural vibration. Hypersonic aircraft structural flutter is a kind of self-excited vibration with aerodynamic force as the primary excitation source. This mechanism allows the structural intrinsic modals to change with time over the course of the test. Technically, modal damping is difficult to generalize or estimate reasonably based on the observed vibration response of the structure—especially in the case of multimodal coupling, complex structures, and structural or material nonlinearity. Flutter modal damping is a relatively small value that is highly sensitive to observation noise, calculation error, computer word length, and other effectually uncontrollable factors.
There are a few existing methods for estimating flutter modal damping from structural responses. The half-power bandwidth method, for example, is derived from a single-degree-of-freedom mass spring system that is usually only applicable to single-mode or multimode decoupled systems. Under dense mode or low signal-to-noise ratio (SNR) conditions, the frequency spectrum is irregular, peaks are overlapped, or an abundance of burrs significantly affects the effectiveness of this method. The random decrement technology is also derived from the single-degree-of-freedom mass spring system wherein any response signal can be split into three parts: random response, initial velocity response, and initial displacement response. When the random response and initial velocity response are partially removed, only the initial displacement response is retained as a free decay signal or unit impulse response of the system. Corresponding damping estimates can be obtained via empirical formula or envelope fitting. This method is only fully applicable to systems with single degrees of freedom and requires a set quantity of observed samples. The Miramand method is based on the transfer function theory. When the excitation signal is determinable statistically, all five modal parameters can be obtained through rational function fitting and the calculation of residues. Generally, this method is suitable for both single-mode and multimode situations. Prony methods are based on complex exponential function fitting—there are many deformation algorithms which apply to these methods, for example, the matrix pencil method, which can be used for multimodal signal analysis to directly obtain the frequency and damping parameters.
Traditional damping-based methods mainly work within the frequency domain under the integral transformation concept. The demands for stationarity and length in the data may not be suited to hypersonic flutter analysis. The progression of variable airspeeds in a hypersonic flutter test does not satisfy the assumption of a stationary stochastic process, which is the theoretical basis of these damping methodologies. Some researchers [27, 28] have successfully solved modal identification for nonstationary measured signals by using time-domain modeling methods such as ARMA or TVAR, which are mathematically equivalent to the system stability criterion discussed here.
(b) System Stability Criterion. Considering the structure response as the response of a dynamic system, the stability criterion transforms the flutter problem into a system instability problem; i.e., a discrete system which can be used investigates whether the poles of the system transfer function all lie within the unit circle. As opposed to modal damping methods, system stability methods are defined by the coefficients of the eigenpolynomial of the system. The procedure can be decomposed into two steps. The system is first modeled per the structure responses from a flutter test; then the relationships between the flutter criterion and system coefficients are established.
Generally, the ARMA model can be formed by mixing an AR model representing the response series and a moving average (MA) model as the white noise excitation. The difference equation of the ARMA model iswhere is the structural response signal for the aircraft, is the noise (described as aerodynamic excitation), and is the order of the system, which can be calculated from the AIC criterion or based on a priori information such as the actual number of system modals. Based on the measured data, the coefficients with in (14) are the key factors in evaluating system stability.
There are many system stability criteria with different forms but mathematical equivalence. The most commonly used are Jury’s criterion and Lyapunov’s criterion.
Jury’s criterion defines the square matrix , as follows:where the stability of the model system should meet the following conditions:
Lyapunov’s criterion is defined by the matrix :When The condition of Lyapunov’s criterion for system stability is that the matrix is positive definite, which meanswhere is the principal minor of the matrix .
The system stability criteria above are based on stationary stochastic process theory. To evaluate nonstationary data from hypersonic flutter tests with progressive variable airspeed, and to characterize the time-varying modal parameters for coupling in the test, the coefficients in (14) must also vary with time. Consider the following time-varying parameter model:where the coefficients and are also the functions of sampling point . Generally, the atmospheric turbulence excitation can be approximated as Gaussian white noise, so (21) can be simplified as a time-varying autoregressive (TVAR) model:When the dynamic data has been modeled based on the measured signals, the system stability criterion can be calculated by substituting the time-varying coefficients into (17) and (20) and then extrapolated to obtain a flutter boundary prediction.
The system stability criteria and time-varying model do not necessitate the assumption of a stationary stochastic process. They are feasible for nonstationary situations caused by hypersonic flutter tests with variable airspeeds. Unlike damping-based methods, however, they lack a clear physical meaning to guide structure design and their accuracy is directly affected by dynamic modeling quality. It is important to select an appropriate modeling algorithm.
(c) Hypersonic Flutter Analysis Framework. According to the characteristics of the hypersonic flutter test and the flutter analysis methods discussed above, we drew a block diagram of the flutter prediction system (Figure 2) which includes nonstationarity level assessment, measured data preprocessing, modal damping estimation, dynamic data modeling, system stability calculation, curve fitting and extrapolation, basic I/O interfacing, and other necessary modules. Preprocessing includes statistical analysis and digital filtering. The modal damping is estimated within a stationary situation to reflect actual engineering requirements. The system stability calculations are designed for the rest of the nonstationary case. The dynamic data modeling module contains several time-varying modeling algorithms based on the adaptive filter theory and Bayes method.
The parameters necessary upon system initiation depend on the characteristics and mechanism of the test aircraft. The nonstationarity of the measured signal is assessed online; then appropriate algorithms are selected accordingly. In the case of nonstationarity, dynamic data modeling is performed after the signal preprocessing; the flutter criteria of the system stability are estimated simultaneously, the confidences of which are according to the level of nonstationarity as determined. Prediction results are then secured by curve fitting and extrapolation with the criteria against the Mach number.
4. Data Simulation
(a) Nonstationary Simulation. To validate our nonstationarity level assessment based on the energy distribution in the time-frequency domain, we generated two nonstationary signals with variable amplitude and frequency. The amplitude of signals does not affect the nonstationarity level, so the rate of variation speed of the amplitude simulated was considered to be the independent variable. We investigated the start frequencies, termination frequencies, and rates of speed variation in a variable frequency case. The results are shown in Figures 3 and 4.
In Figure 3, the abscissa represents the change rate of the signal amplitude from the beginning to the end, while the ordinate represents the nonstationarity as calculated by the proposed method. The middle point of the curve corresponds to the stationary signal, which has the minimum nonstationarity level. The left side of the minimum point corresponds to the signal with amplitude that decreases over time; the right side shows the signal with amplitude that increases over time. The near-linear symmetry on both sides of the curve indicates that the nonstationarity brought about by changes in amplitude and other factors is consistent and has a direct proportion to the change rate.
In Figure 4, the axis X represents the initial frequency of the signal, the axis Y is the terminal frequency, and the axis Z is the nonstationarity as calculated. The minimum value is a line at the bottom of the tridimensional plane of face, where the starting frequency and ending frequency are equal. The nonstationarity gradually rises on both sides of the plane and exhibits strong linearity and symmetry characteristics. To this effect, the proposed assessment method has good indicative significance.
(b) Flutter Test Simulation with Variable Speed. According to the general mechanism of structural flutter, that is, the coupling of two-order structural modes, the structural modes change over time throughout the flutter process. We selected a numerical simulation to simulate the dynamic process of said structural response and the signals of the flutter test with variable speed. represent the frequencies and damping factors of the two-degree-of-freedom structural modes. If the modal frequencies of the dynamic system gradually converge while the modal damping decreases at a certain slope, the dynamic system changes from stable to unstable. For example, consider two frequencies start from , , converge, and end at , , where n is the length of the samples; the modal damping descends from , to , . The modal frequencies of the dynamic system converge as shown in Figure 5, while the modal damping decreases at a respective slope. The sample index serves as an indicator of airspeed acceleration due to the linear relation between the sample index and the modal coupling. The corresponding design flutter point is 4,451. Flutter tests are generally conducted under subcritical conditions, so we selected 4,096 samples in the front as structural response signals under subcritical conditions (Figure 6).
We used the proposed hypersonic flutter monitoring and analysis system to process this signal. Figure 7 shows Jury’s criterion with curve fitting and extrapolations as an example. The stability of the dynamic system during flutter is clearly reflected in the trend in the figure.
More prediction results from the different criteria are listed in Table 1 to further demonstrate the accuracy of the proposed framework and the system, which effectively meet the requirements of hypersonic flutter tests (especially the performance of F- (3)).
5. Application to Flutter Flight Tests with Variable Speed
We used data from a wind tunnel test and transonic flight flutter test (as hypersonic flutter test data is relatively scarce) to verify the effectiveness and feasibility of the proposed algorithm.
(a) Wind Tunnel Test with Progressive Variable Speed. Compared to a flutter test with actual aircraft, the FBP result in a wind tunnel test is more closely consistent with the design value. The data can also be verified at lower risk. In our test, the sampling rate was 128 Hz and the acquisition time was 162 seconds. The airspeed was gradually and smoothly increased until flutter occurred. The torsion signals acquired from the strain gauge are shown in Figure 8 by comparison against the variable airspeed shown in Figure 9. The flutter consistently occurred at about 90 seconds when the airspeed was around 28 m/s.
We first assessed nonstationarity as shown in Figure 10. Fluctuations in nonstationarity were fairly smooth prior to the flutter but dramatically increased after the vibration diverged.
Next, signals of the front 60 seconds were selected as the subcritical condition for flutter prediction with a maximum speed is 22m/s to further test the proposed framework. The results of system stability criteria F- (3) at each time point and the FBP from curve fitting and extrapolation are shown in Figure 11. The FBP result is 27.17, which is not only consistent with the physical experiment, but also indicative of a monotonous descending tendency at the early stage.
(b) Transonic Flight Flutter Test with Progressive Variable Speed. The selected test data came from a specific type of aircraft transonic flutter flight test which includes both subsonic and transonic periods. Each measured signal was collected from the sensor of an accelerator located on the middle of the left airfoil, as shown in Figure 12, with variable airspeed and flight altitude as shown in Figures 13 and 14, respectively. The sampling rate was 256 Hz and the acquisition time was 132 seconds.
The nonstationarity levels calculated via the proposed method are shown in Figure 15. The flight speed is far from the flutter velocity, so the signal has lower stationarity than the wind tunnel test signal but is affected by unsteady aerodynamic force or other flight environment characteristics. As shown in Figures 13 and 14, during the 40th second to 60th second time period, the aircraft altimeter and speedometer spiked simultaneously—when the aircraft broke through the sound barrier—accompanied by an abrupt change in nonstationarity level. Above the sound barrier, this effect no longer exists; the nonstationarity level drops back down. This indicates that the proposed assessment method is sensitive to the conditions of measured signals.
We selected orthogonal basis functions method for TVAR modeling in subsequent data processing per the nonstationarity levels as calculated. The system stability criteria at each time point and the FBP results from curve fitting and extrapolation are shown in Figure 16. The maximum airspeed of the flight test must be limited far below the flight envelope, especially for variable airspeed situations, to ensure safety. As opposed to wind tunnel tests, there is no true value for comparison. However, the consistent and monotonous descending tendency validates the effectiveness of the proposed framework.
6. Conclusion
A flutter analysis methodology based on structural responses of hypersonic flutter tests is proposed in this paper as a response to the inherent drawbacks to existing aerothermoelastic numerical simulations. The fundamental problems with nonstationarity assessment theory and flutter prediction criteria are discussed; then a novel signal processing framework is introduced based on hypersonic engineering. Simulations and experimental results indicate that the proposed assessment method accurately indicates nonstationarity; the proposed flutter analysis framework for hypersonic flutter test also successfully predicts the flutter boundary forming the subcritical turbulence response in a variable airspeed scenario.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This work is supported by the Aeronautical Science Foundation of China (20151353017), National Natural Science Foundation of China (Grant no. 11302175), and the Fundamental Research Funds for the Central Universities (3102017HQZZ013).