Robust aerodynamic parameter estimation method based on variational Bayesian inference

Abstract

The invention belongs to the technical field of aerospace, and particularly relates to a robust aerodynamic parameter estimation method based on variational Bayesian inference. A traditional scheme, such as a smoothing method of a smoother, for estimating the optimal performance of aerodynamic parameters is often based on the assumption that measurement noise meets Gaussian distribution. However, in practical applications, noise distribution often deviates from Gaussian hypothesis, for example, heavy tail density may occur in applications involving flicker noise, air turbulence and the like. The factors can cause the Gaussian hypothesis of noise measurement to be invalid, so that the performance can be greatly reduced when aerodynamic parameters are estimated based on a filter or smoother method. According to the method, the measurement noise is modeled by using the zero-mean Pearson type VII distribution, the novel, simple and easy-to-implement aircraft aerodynamic parameter estimation method is provided by means of a variational Bayesian inference method, and the requirements of the aerospace field on high-precision estimated values of aircraft aerodynamic parameters are met.

Description

Technical Field

[0001] The invention belongs to the field of aerospace technology, and in particular relates to a robust aerodynamic parameter estimation method based on variational Bayesian inference. Background Art

[0002] In the field of aerospace, accurate aircraft dynamics models are critical for designing control laws, verifying performance, and realizing flight simulation. Such models are generally expressed by some coefficients or parameters, and the determination of these parameters requires reference to different flight conditions. Among the many parameters, the aerodynamic coefficients, which are closely related to aerodynamic forces and moments, are particularly critical. The specific values ​​of these coefficients are mostly derived from wind tunnel test data or computational fluid dynamics (CFD) technology. However, wind tunnel tests and CFD simulations often cannot fully reproduce the actual flight environment. Therefore, in the aviation industry, flight experiments are usually conducted and aerodynamic coefficients are estimated based on experimental data. Common aerodynamic coefficient estimation methods include filtering methods and smoothing methods. However, these two methods can only obtain the best performance when estimating aerodynamic parameters under the assumption that the measurement noise conforms to the Gaussian distribution. However, in actual applications, due to various factors, the noise distribution often deviates from the Gaussian assumption. For example, heavy-tail density may appear in applications involving flicker noise, air turbulence, etc. These factors will invalidate the Gaussian assumption of the measurement noise, which may lead to a significant decline in the performance of estimating aerodynamic parameters based on filter or smoother methods. Summary of the invention

[0003] The purpose of the present invention is to propose a robust aerodynamic parameter estimation method based on variational Bayesian inference, realize adaptive adjustment of measurement data noise variance, and obtain accurate estimation values ​​of aerodynamic parameters. Specifically, a novel, simple and easy-to-implement aircraft aerodynamic parameter estimation method is proposed by using a zero-mean Pearson type VII distribution to model measurement noise, and with the help of a variational Bayesian inference method, to meet the demand for high-precision estimation values ​​of aircraft aerodynamic parameters in the aerospace field.

[0004] The technical solution of the present invention is:

[0005] Consider a general parameterized discrete system described by the following state-space model (SSM)

[0006] x t =f(x t-1 ,θ)+w t

[0007] y t =h(x t ,θ)+v t

[0008] Where t represents the sampling time, is the state quantity of the aircraft dynamic system, Represents about x t The observation results, is the unknown parameter vector that specifies the process map f(·) and the measurement map h(·), w t The noise is composed of all sensor noises, which has a mean of 0 and a variance of the diagonal matrix Q t Gaussian distribution, v t is the measurement noise, w t and v t Are independent of each other.

[0009] Measurement noise v t It is generally assumed to obey Gaussian distribution. However, in most applications, due to unexpected interference and failure of sensors, the measurement noise is no longer Gaussian distributed, and its distribution exhibits heavy-tailed characteristics. In order to better illustrate this property, a heavy-tailed distribution is used instead of the Gaussian distribution used for frequency to simulate the measurement noise. Specifically, in this work, a zero-mean Pearson type VII distribution is adopted, that is,

[0010] v t ~P VII (0,R,α t )

[0011] in, is the scale matrix, α t is the degrees of freedom. The density of the Pearson type VII distribution is

[0012]

[0013] Here, Γ(·) represents the gamma function.

[0014] The purpose of the present invention is to Estimate θ and the statistics of the two noises, i.e.,

[0015] In order to avoid the complex nonlinear optimization problem for the unknown parameters θ, a state augmentation method is adopted and the unknown parameters θ are regarded as random variables. Specifically, θ is modeled as a Gaussian random walk, i.e.

[0016] θ t =θ t-1 +p t

[0017] in, is a Gaussian white noise. Define an augmented state as The augmented process noise is The augmented SSM is

[0018] z t =g(zt-1 )+q t

[0019] y t =h(z t )+v t

[0020] Among them, g(z t-1 ) is defined as

[0021]

[0022] The estimation of θ comes from the smoothed parameter state. However, since the model involves the Pearson type VII distribution, the common smoothing algorithm cannot be used to directly solve the estimation problem of unknown parameters. To meet this challenge, another expression of the Pearson type VII distribution is used. Specifically, it can be expressed as a mixture of weighted Gaussian distributions:

[0023]

[0024] in is a gamma distribution, where a t >0 and b t > 0 is its parameter. In addition, a non-informative Jeffrey prior is imposed on R and Q, i.e.,

[0025]

[0026] The posterior distribution of all variables to be estimated can be decomposed into:

[0027]

[0028] Where Z is defined as

[0029] Due to the complexity of the distribution involved, it is difficult to directly solve the augmented state by solving the maximum likelihood probability problem, and then estimate the unknown parameters. To address this challenge, an approximate Bayesian inference method, namely the variational Bayes (VB) method, is used. Specifically, in VB, a variational distribution q(Z,{ω t},Q,R), so q(Z,{ω t},Q,R) and p(Z,{ω t},Q,R|Y) is the smallest, that is,

[0030]

[0031] where, for simplicity, · denotes the corresponding variable, is the set of all possible solutions of q(·). In addition, the mean field approximation is applied, where the variational distribution can be decomposed into

[0032]

[0033] Based on the update rule of the VB method, there are the following update steps for each variational distribution:

[0034]

[0035] In the formula, is the complete distribution of the system identification problem under consideration, given by

[0036]

[0037] Obviously, the variational distributions are coupled to each other. Therefore, these variational distributions can be updated in an iterative manner, that is, only one parameter is updated while keeping the other parameters unchanged. In the following, the update expression of each variational distribution is derived.

[0038] 1. Update of q(Z)

[0039] First, let’s discuss the update of q(Z), which can be simplified to

[0040]

[0041] Note that logp(z t |z t-1 ,Q) and logp(y t |z t ,R,ω t ) can be expressed as

[0042]

[0043] in, C1 and C2 are constants.

[0044] Substituting the above expression, we can get

[0045]

[0046] Essentially, q(Z) can be viewed as the joint density function of the SSM, which can be efficiently solved in the Gaussian approximate smoothing framework. Typical solutions include volumetric Kalman smoothing, unscented Kalman smoothing, and so on.

[0047] 2. Update of q(Q)

[0048] q(Q) follows an inverse Wishart distribution, that is,

[0049]

[0050] Among them, Ψ is given by

[0051] Ψ=∫(z t -g(z t-1 )) T (z t -g(z t-1 ))q(z t ,z t-1 )dz t dz t-1

[0052] Furthermore, the expected value of Q is calculated as

[0053]

[0054] 3. Update of q(R)

[0055] q(R) also follows an inverse Wishart distribution, that is,

[0056]

[0057] Among them, Φ is given by

[0058]

[0059] Furthermore, the expected value of R is calculated as

[0060]

[0061] 4.q(ω t )

[0062] q(ω t ) is a gamma distribution, that is, where a t ,b t is updated to,

[0063]

[0064] Among them, B is given as

[0065]

[0066] Implementing a smoothing-type algorithm on the augmented SSM will yield an estimated augmented state, namely Then we can get the estimated values ​​of the randomized unknown parameters at different times. use Then we can finally get the estimated value of the unknown parameter by taking the average value, that is,

[0067]

[0068] The beneficial effect of the present invention is that the aerodynamic parameter estimation method proposed in the present invention can provide a robust estimation of aerodynamic parameters when the measurement noise has a heavy-tailed distribution due to wild value disturbances, etc., thereby providing a guarantee for the establishment of a high-precision model of the aircraft system. BRIEF DESCRIPTION OF THE DRAWINGS

[0069] Figure 1 The statistical results of the estimated values ​​of the method proposed by the present invention;

[0070] Figure 2 The NMSE characteristics of the estimated value of the method proposed by the present invention;

[0071] Figure 3 The convergence characteristics of the method proposed in the present invention. DETAILED DESCRIPTION

[0072] The present invention is described in detail below in conjunction with the accompanying drawings and simulation examples. Taking the aerodynamic parameter identification of the longitudinal short-period motion of the HFB-320 aircraft as an example, the practicality of the present invention is demonstrated by comparing two existing solutions within the EM framework. The first comparison algorithm combines UKS with EM (EM-UKS). The second comparison algorithm uses a combination of a robust nonlinear smoother and EM (EM-RS). In EM-RS, a Huber cost function with a parameter of 1.35 is introduced to improve the state estimation accuracy, thereby helping to improve the system identification performance. The number of iterations of these methods is set to 50.

[0073] The longitudinal motion dynamics model of HFB-320 is:

[0074]

[0075] The drag, lift, and pitching moment coefficients are modeled as

[0076]

[0077] Among them, δ e is the elevator deflection degree, and V0 is the reference speed.

[0078] The measurement at time t is It includes true airspeed, angle of attack, pitch attitude, pitch rate, pitch acceleration, acceleration along the x-axis, and acceleration along the y-axis. t is represented as

[0079]

[0080] Among them, ω t To measure noise, longitudinal and vertical force coefficients C X , CZ for

[0081] C X =C L sinα t -C D cosα t

[0082] C Z =-C L cosα t -C D sinα t

[0083] The parameter vector θ to be estimated is given by In the simulation, the initial values ​​of all unknown parameters are set to 0, and the initial values ​​of Q and R are both unit matrices. The continuous dynamic model of HFB-320 is discretized using the Runge-Kutta fourth-order method. The values ​​of all geometric quantities (e.g., mass, moment of inertia, etc.) as well as the measured data are taken from published software packages. These measured values ​​are considered perfect, which means that there are no outliers. In order to simulate outliers, the measured values ​​used in the simulation are generated as

[0084]

[0085] in, represents the perfect data recorded, and is the additional outlier perturbation, which is given by

[0086]

[0087] Where λ is the pollution ratio, is with Correlated noise covariance.

[0088] The present invention simulates 100 rounds of Monte Carlo random experiments to collect results. The estimated error bars of the unknown parameters are as follows: Figure 1 As shown in Figure 2. It can be observed that both the robust algorithms (i.e., the proposed method and EM-RS) can provide more accurate estimates than EM-UKS. In addition, the proposed method provides lower uncertainty estimates compared to EM-RS. As λ increases, the uncertainty of the estimation results obtained by EM-RS and EM-UKS increases significantly with the increase of λ, while the uncertainty of the proposed method remains stable.

[0089] The present invention also evaluates the normalized mean square error (NMSE) of the estimated parameters, defined as

[0090]

[0091] where θ iis the true value of the i-th component of the unknown parameter θ, is the kth round of Monte Carlo analysis of θ i The estimated results.

[0092] Figure 2 The NMSE of the aerodynamic coefficients estimated by different methods is plotted as a function of λ. It can be seen that the proposed method maintains a low NMSE at different λ values, verifying the accuracy and robustness of the proposed method. Although the NMSE of the three methods will increase with the increase of λ due to the increase in the number of outliers, the proposed method is less sensitive to λ.

[0093] Figure 3 The convergence of the proposed method is shown. Note that only the convergence of four aerodynamic coefficients is shown in the figure, and the other parameters are omitted due to the similarity of the convergence. The simulation results show that the method can effectively converge to the true value, and the convergence is achieved in the first 10 to 15 iterations. This shows that the estimation algorithm is robust and can accurately determine the aerodynamic coefficients after a relatively small number of iterations.

[0094] In summary, a robust system identification method is proposed. Specifically, a Pearson type VII distribution is used to capture the heavy-tailed nature of noise associated with measurements affected by outliers. The longitudinal aerodynamic parameters of the HFB-320 are estimated in simulations to verify the effectiveness of the proposed method. Numerical results show that the proposed method can reliably estimate unknown aerodynamic parameters and outperform existing state-of-the-art solutions.

Claims

1. A robust aerodynamic parameter estimation method based on variational Bayesian inference, characterized in that: The following steps are involved: S1. Establish the state space model as: x t =f(x t-1 ,θ)+w t y t =h(x t ,θ)+v t Where t represents the sampling time, is the state quantity of the aircraft dynamic system, Represents about x t The observation results are: n and m are the dimensions of the system state and observation respectively. is the unknown parameter vector that specifies the process map f(·) and the measurement map h(·), w t is the noise composed of all sensor noises, which follows a Gaussian distribution with a mean of 0 and a variance of the diagonal matrix Q. t is the measurement noise, w t and v t are independent of each other; A zero-mean Pearson type VII distribution is used to simulate the measurement noise, namely: v t ~P VII (0,R,a t ) in, is the scale matrix, α t is the degree of freedom; definition Set the parameter estimates to estimate θ from Y and the statistics of the two noises S2. Model θ as a Gaussian random walk: i t =θ t-1 +p t in, is a Gaussian white noise; define an enhancement state as The enhanced process noise is The enhanced state-space model is reinterpreted as: z t =g(z t-1 )+q t y t =h(z t )+v t Among them, g(z t-1 ) is defined as S3. The expression defining the Pearson type VII distribution is: Among them, ω t is an auxiliary variable used to write the Pearson type VII distribution as with the auxiliary variable ω t The finite Gaussian distribution of is a gamma distribution, a t >0 and b t > 0, impose a non-informative Jeffrey prior on R and Q, namely: S4. Decompose the posterior distribution of all variables to be estimated into: Where Z is defined as Using variational Bayesian inference, we find a variational distribution q(Z,{ω t },Q,R), so that q(Z,{ω t },Q,R) and p(Z,{ω t The kl divergence between},Q,R|Y) is the smallest, that is: where · represents the corresponding variable, is the set of all possible solutions of q(·); applying the mean field approximation, the variational distribution is decomposed into: Based on the update rule of the variational Bayesian method, there are the following update steps for each variational distribution: In the formula, is the complete distribution of the system identification problem: The update expression for each variational distribution is: Update expression of q(Z): Update expression of q(Q): where Ψ is given by: Ψ=∫(from t -g(z t-1 )) T (With t -g(z t-1 ))q(z t ,With t-1 )dz t d t-1 The expected value of Q is calculated as: The update expression of q(R) is: Where Φ is given by: The expected value of R is: q(ω t )’s update expression: Among them, B is given as Β=∫(y t -h(z t )) T (y t -h(z t ))q(z t )dz t S5. Use the smoothing algorithm to solve the estimated augmented state, that is, Then we can get the estimated values ​​of the randomized unknown parameters at different times. use The estimated value of the unknown parameter is finally obtained by the average method, that is,

Images