Robust Stability Clearance of Flight Control Law Based on Global Sensitivity Analysis

To validate the robust stability of the flight control system of hypersonic flight vehicle, which suffers from a large number of parametrical uncertainties, a new clearance framework based on structural singular value (μ) theory and global uncertainty sensitivity analysis (SA) is proposed. In this framework, SA serves as the preprocess of uncertain model to be analysed to help engineers to determine which uncertainties affect the stability of the closed loop system more slightly. By ignoring these unimportant uncertainties, the calculation of μ can be simplified. Instead of analysing the effect of uncertainties on μwhich involves solving optimal problems repeatedly, a simpler stability analysis function which represents the effect of uncertainties on closed loop poles is proposed. Based on this stability analysis function, Sobol’s method, the most widely used global SA method, is extended and applied to the new clearance framework due to its suitability for system with strong nonlinearity and input factors varying in large interval, as well as input factors subjecting to random distributions. In this method, the sensitive indices can be estimated via Monte Carlo simulation conveniently. An example is given to illustrate the efficiency of the proposed method.


Introduction
In the past decades, clearance of flight control system has been paid great attention by the air force of many countries.Rational clearance improves not only the reliability and safety of flight control system, but also feedback valuable information to designers.Clearance methods split neatly into two types based on the different clearance principle.One is the so-called analytical-model-based method (AMBM) [1] and another is simulation-data-based method (SDBM) [2,3].AMBM analyses stable margin, robust stability and some other criteria based on advanced control theories which include (1) structural singular value () analysis [4,5]; (2) polynomial-based clearance [6]; (3) v-gap analysis; (4) bifurcation and continuation method [7]; and (5) optimizationbased clearance [8].AMBM has been improved and extended to more criteria and more types of air vehicles [9][10][11][12][13][14]. SDBM demands high fidelity simulation platform and combines experiment design techniques with decision sciences.Generally speaking, AMBM is suitable for early stages of the whole controller designing process, whereas SDBM is suitable for late stages.These two types of methods are complementary.
For hypersonic flight vehicle (HFV) whose maximum flight speed may reach to about 20 Mach, the principal criterion to validate is robust stability of its attitude control system because of the existence of uncertain coefficients.The uncertainties are mainly caused by the following three factors.First, the aerodynamic coefficients of HFV, obtained via wind tunnel test and computational fluid mechanics tools, often suffer from much more serious deviations than those of supersonic and subsonic vehicles.Second, it is still unknown how the external environment affects the HFV due to lack of flight experience.Finally, erosion and corrosion may damage the aerodynamic configuration of the vehicle and hence change the aerodynamic coefficients.Many advanced methods [15][16][17] have the potential to improve the ability and performance of the attitude control system, which may also lead to more complicated controllers and require more rational clearance.
Among these model-based clearance methods mentioned above,  analysis is the favourable method for robust stability criterion due to its nonconservatism.Calculation of  with pure real uncertainties [18] had always been a challenge until the optimization-based approach [19] was proposed.Hybrid

Robust stability analysis
Clearance result display Uncertainty sensitivity analysis optimization [8], which combines sequential quadratic programming (SQP) with intelligent optimization, can produce an exact estimation of .However, the computational burden is still heavy due to iteratively sweeping frequencies.In particular, when the number of uncertainties is big, it will take a long time to get  values at a series of frequency points.

Clearance result display
In order to reduce the computational burden of , engineers usually ignore some unimportant uncertain factors according to their experience.How to find those minor uncertain factors of the HFV when we have little experience?The answer is sensitivity analysis (SA).SA methods can be broken down into local SA [20,21] and global SA [22][23][24][25].Global SA has been applied widely in industrial production, ecological engineering, hydrology systems, and so on [26] due to its two main advantages.First, it is suitable for input factors varying in broad range and system model with strong nonlinearity.Second, it can consider input factors obeying certain random distributions.
In this work, we proposed a new clearance framework for robust stability clearance of the control system of HFV based on  theory and global uncertainty SA (GUSA), which is shown in Figure 1.The GUSA serves as the preprocess of the  analysis to help engineers to determine those minor uncertainties.Because measuring the direct effect of uncertainties on  involves computing  repeatedly, a novel indirect GA index based on our newly defined stability analysis function is adopted.This function represents the effect of uncertainties on dominant closed loop poles.With this new GA index, we can extend the widely used Sobol's method [23] and apply it to analyse the effect of uncertainties on stability of the flight control system.It is worth noting that the framework is suitable for not only  analysis but also other robust analysis approaches.Furthermore, GUSA can provide valuable information about the effect of aerodynamic coefficients on stability for subsequent design and validation.
The paper is organized as follows.In Section 2, the uncertain model of the attitude control system of HFV is established and  framework is presented.The Sobol's method is extended and applied to GUSA in Section 3. Section 4 gives an example to illustrate the efficiency of the proposed method.In Section 5, some conclusions are stated.

Problem Description
The motion equations of the hypersonic flight vehicle are given as follows [27] whose nominal values are given in the form of look-up tables.  , ,   ,   , and   are uncertain coefficients varying in some certain intervals around their nominal values.
C P Linearizing ( 1) and ( 2) along the nominal flight trajectory according to the methods introduced in [4,18], the uncertain model P can be written as where is the control input, and   ∈ R 6×6 and   ∈ R 6×3 are uncertain matrices in the affine form as follows: where  is the number of normalized uncertainties after combination and linearization and   ∈ [−1, 1] are uncertain parameters that represent the uncertainties of these coefficients.The uncertain model ( 4) is not unique via different modelling methods.Because of different designing and manufacturing environment for different vehicles, uncertainties may obey certain random distribution.Consider for instance, the normal distribution   ∼ (,  2 ) or the uniform distribution   ∼ (, ).Suppose that a state feedback controller C, which can stabilize the nominal system of (3), is where   ∈ R ×1 is the state vector of controller,  = [      ]  is the command attitude angles to be tracked, Then, the closed loop system, as shown in Figure 2, can be written as where It is clear that the stability of the system (7), suffering from affine interval uncertainties   , is determined by the system matrix  0 +Δ.In  framework, the linear fraction representation (LFR) of ( 7), denoted as F(, Δ), should be modelled first according to the approach presented in [4,18].Consider where  represents known part of the system with input  and output , Δ = diag{     } represents uncertainties of diagonal structure, rank(⋅) is the rank of a matrix, and   ∈ R (+)×  and   ∈ R   ×(+) are obtained by full rank decomposition.Consider At a specific frequency point , structural singular value of F(, Δ) (9) is defined as where ‖ ⋅ ‖ ∞ is the infinite norm of a matrix.Define  * = max   Δ (()).Since ‖Δ‖ ∞ ≤ 1, the sufficient condition of robust stability is max{‖Δ‖ ∞ } < 1/ * ; that is,  * < 1.Furthermore, the system stays stable even if the uncertain intervals extend to For simplicity, we denote  Δ (()) as  hereinafter.
According to definition (11), if the problem min has an optimal solution  * , then the estimation of  is 1/ * .By combining global intelligent algorithm (such as particle swarm optimization (PSO) and differential evolution (DE) [28]) with local SQP, this problem can be iteratively solved [8].In order to get the estimation of  * = max  , we need to solve the optimal problem (12) at a series of frequency points.When the number of uncertainties is large, the dimension of the searching space is also large.In order to cover the whole space uniformly, the number of particles of PSO (or analogy of other intelligent algorithms) should exponentially increase with respect to the number of uncertainties.Then, the computational burden will increase exponentially.Therefore, our aim is to find those   s which affect the stability dramatically and ignore remaining unimportant   s when calculating  * .If we measure the effect of uncertainties on  * directly, we need to calculate  * repeatedly.So we introduce an indirect stability analysis function.
As we know, a system is stable if and only if all of its closed loop poles lie in the left half of the complex plan.Therefore, we define the nonlinear stability analysis function as follows: where   (⋅) is the th eigenvalue of a matrix and Re[⋅] is the real component of a complex.The effect of   on the stability of the closed loop system can be represented by (13) and SA can be invoked.It is clear that the stability analysis function ( 13) is also suitable for other robust stability clearance method.
All local SA methods need to calculate the differential or derivative; that is, they are based on the analytical model.However, the nonlinear function (13) does not satisfy this requirement.So we will use the global Sobol's method to deal with this problem in the upcoming sections.

Sobol's Method. Consider the nonlinear model in the form
where x = [ 1 ,  2 , . . .,   ] are  independent input factors and  is the scalar output.If x is defined over Ω = {x | 0 ≤   ≤ 1,  = 1, 2, . . ., }, then the function can be decomposed as the following form [23]: where   =   (  ),   =   (  ,   ), and so on.The uniqueness of ( 15) is guaranteed by where 1 ≤ These functions in (15) are obtained by and similarly for higher orders, where And the partial variance is defined as Then the decomposition (15) leads to the following decomposition of the variance of The normalized sensitivity index is defined as where Supposing X ∈ R × is the matrix of  groups of  input factors where and P and Q are two samples of X generated by MC.Then we define that is the combined matrix with column  from Q and all other  − 1 columns from P, where P  and Q  are samples of X  .Similarly, we define The estimation below can be obtained [24,25] as follows: where p  , q  , (P () Q )  , and (Q () P )  are, respectively, the th row of P, Q, P ()  Q , and Q () P .

Uncertainties Sensitivity Analysis Based on Sobol's Method.
It is worth noting that the input factors of the nonlinear model considered in Sobol's method are limited to 0 ≤   ≤ 1.
However, the normalized uncertainties of the flight control system is   ∈ [−1, 1].One manner is to renormalize the uncertainties to the interval [0, 1].Here we will introduce an easier method.
Considering the nonlinear model ( 14), if its domain of definition is we can get the below theorem.
Remark.In fact, Sobol's method can be extend to any nonlinear function defined over where  −  = 1.
The uncertain matrix Δ of the system matrix in (7) can be rewritten as where   =   /2 ∈ [−0.5, 0.5] are the normalized uncertainties.It is clear that the uncertain space of   is Π defined in (30).Substituting ( 7) and (36) into the nonlinear function ( 13), we can get where  = ( 1 ,  2 , . . .,   ),  ∈ Π.The Sobol's method can be applied to analyse the effect of uncertainty terms   on the real part of the closed loop poles.The sample matrices of  can be generated randomly according to the distribution of .
Generally, the whole clearance process can be summarized into the following.
Step 1: establish the uncertain model of the attitude control system and normalize the uncertainties to [−0.5, 0.5].
Step 2: construct stability analysis function as (37) according to the specific uncertain model.
Step 3: generate the sample matrices  and  of  groups of uncertainties .
Step 4: calculate output samples of the stability analysis functions (37) according to MC simulations with input samples  and .
Step 5: calculate the SA index according to (29) and ignore the minor uncertainties with relatively small SA index.
Step 6: calculate the singular structure value according to (12).

Examples
In this section, a numerical example is given to reveal the effectiveness of the proposed method.Consider the nominal flight trajectory of a HFV depicted as Figure 3, from which 30 operation points for robust stability clearance have been selected.At first, we will illustrate the application of the new clearance process with one single point.Then, we will analyse all these operation points together.
Because the pitch motion has no coupling with yaw and roll motion, the controller can be designed separately.For convenience, a low order controller that can stabilize the nominal system has been designed as follows: where  1 = −99.9845(40) Then the system matrix of the closed loop system is where Calculate the global sensitivity indices of the function (37) with respect to   =   /2 ∼ (−0.5, 0.5) according to the steps presented in Section 3. In order to illustrate the consistency and repeatability of the proposed method, 10 MC simulations have been conducted.The obtained global sensitivity indices   are depicted in Figures 4 and 5.It is evident that (1) the results of 10 simulations are nearly consistent; (2)  1 ,  3 , and  4 are extremely close to 0 and obviously smaller than other indices; and (3) when the sample number  increases from 20000 to 50000, the results converge accordingly.Varying   in the interval [−1, 1] singly and fixing   = 0,  ̸ = , the distribution of the dominant closed loop poles is shown as Figure 6.This can directly validate the SA results.
Since   ≈ 0,  = 1, 3, 4, we can ignore   ,  = 1, 3, 4 in the analysis of robust criterion. is calculated via hybrid optimization method at 136 frequency points  = [0.01,0.02, . . ., 0.09, 0.1, 0.2, . . ., 3, 4, 5, . . ., 100] in the MAT-LAB environment.Contrastive results are shown in Figures 7 and 8.When only considering the remaining 5 major uncertainties (  ,  = 2, 5, 6, 7, 8), the  curve is very close to the situation of considering all 8 uncertainties.The computational time, however, decreases from 3.3387 s to 2.2903 s when the number of uncertainties decreases from 8 to 5. When taking into account all the operation points together, the proposed clearance frame can save a large amount of time.Because  * ≈ 1.1 > 1, the controller cannot stabilize the system in the whole uncertain space.However, we can get the valuable information that the approximate allowable uncertain space is |  | < 1/1.1 ≈ 0.9 and feed it back to designers.
The global sensitivity indices   of all operation points are plotted in Figure 9.Note that   ≈ 0,  = 1, 3, 4, holds for   is relative index.It means that we can order the effect of all uncertainties based on the GUSA index, but we cannot find a threshold value to ignore unimportant uncertainties.There is a tradeoff between precision and computational time when ignoring unimportant uncertainties.For instance, if we need to reduce more computational time, we can ignore one more uncertainty  6 .

Conclusion
In this paper, a new  framework based on GUSA is proposed.GUSA serves as the preprocessing step of  analysis and obviously release the computational burden.Based on our  proposed stability analysis function, Sobol's method is extended and applied to determine the major uncertainties because of its maturity and flexibility.An example is given to illustrate the efficiency of the proposed method.In the potential applications, this work provides an easy-to-follow procedure for engineers to digest and expand.

Figure 1 :
Figure 1: The new clearance process.

Figure 2 :
Figure 2: Closed loop system of the attitude control of HFV.

Figure 7 :
Figure 7: Structural singular value of the uncertain systems.

Figure 8 :
Figure 8: Structural singular value in low frequency interval of the uncertain systems.

Figure 9 :
Figure 9: GSA of all operation points.

Figure 10 :
Figure 10: Computational time of  for all operation points.
are the so-called first-order sensitivity indices and   and   1  2 ,...,  are the second-order and the th-order indices, respectively.  1  2 ,...,  indicates sensitivity of  to the interaction among   1 ,   2 , . . .,    .So the total sensitivity indices   for the input   is defined as the sum of all indices relating to it.