A State Estimation Method for ATR Engine Control System Based on Adaptive Robust UKF
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2025-07-31
- Publication Date
- 2026-08-14
AI Technical Summary
[0004]本发明的目的是提供一种基于自适应鲁棒UKF的ATR发动机控制系统状态估计方法,以解决现有非线性系统状态估计方法在复杂噪声环境中鲁棒性不足、对异常测量值敏感、难以应对非高斯噪声的影响等问题
[0055]与现有技术相比,本发明具有以下技术特点:
Smart Images

Figure CN120798544B_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of ATR engine control, specifically relating to a state estimation method for an ATR engine control system based on adaptive robust UKF. Background Technology
[0002] The ATR (Air Turbo-Rocket) engine is a hybrid propulsion system with a high thrust-to-weight ratio and strong adaptability, widely used in aerospace and military fields. Due to its unique structure and performance, the ATR engine often operates in complex environments, requiring precise real-time monitoring and control of its state to ensure system safety and stability. However, the operating state of the ATR engine is highly nonlinear, its operating environment is variable, and it is also affected by high noise interference and measurement errors. Traditional filtering methods struggle to guarantee high-precision state estimation in complex noise environments. Furthermore, existing filters typically use fixed parameters and cannot adaptively adjust according to real-time noise conditions, resulting in insufficient robustness and adaptability in practical applications. Current Kalman filters and their derivatives, when dealing with state estimation of complex nonlinear systems, make overly idealistic noise assumptions (usually assuming a Gaussian distribution). In real-world non-Gaussian noise and anomalous data conditions, filter performance degrades significantly, failing to meet the high-precision state estimation requirements of the ATR engine.
[0003] To address these challenges, researchers have proposed several improvements, such as adaptive filters and robust filters, aiming to enhance the robustness of the filters and their adaptability to complex noise. However, these methods still suffer from inflexible parameter tuning, high computational cost, and inability to handle state estimation of high-dimensional complex systems. Therefore, there is an urgent need for an ATR engine state estimation method that can adaptively cope with the influence of complex noise environments and sensor error values, while maintaining high accuracy and good robustness. Summary of the Invention
[0004] The purpose of this invention is to provide a state estimation method for an ATR engine control system based on adaptive robust UKF, in order to solve the problems of insufficient robustness of existing nonlinear system state estimation methods in complex noise environments, sensitivity to abnormal measurement values, and difficulty in dealing with the influence of non-Gaussian noise.
[0005] To achieve the above objectives, the present invention employs the following technical solution:
[0006] A state estimation method for an ATR engine control system based on adaptive robust UKF includes:
[0007] Establish a nonlinear discrete model for the ATR engine control system, including establishing a state-space model for the ATR engine control system, determining the specific parameters of engine state variables, engine input variables and engine output variables, discretizing the state-space model and calculating the noise covariance matrix.
[0008] The time update process includes initializing the engine state variables and the state error covariance matrix, selecting sampling points and constructing the state variables of the sampling points, and predicting the estimated values of the engine state variables and the error covariance matrix.
[0009] The measurement update process includes updating the state variables of the sampling points, transforming the updated state variables of the sampling points using measurement equation functions, calculating the estimated value of engine output and covariance matrix through unscented transformation, and calculating the state-measurement cross-covariance matrix.
[0010] The process of calculating the Kalman gain matrix using the maximum correlation entropy criterion includes defining a pseudo-measurement matrix, using the pseudo-measurement matrix to convert the nonlinear engine output into a linear form, defining a cost function, and solving for the Kalman filter gain matrix based on the cost function.
[0011] The process of updating the state estimate and covariance matrix includes updating the engine state quantity estimates and updating the state error covariance matrix.
[0012] Furthermore, the method also includes:
[0013] The robust state estimator design process involves: first, designing the batch regression form of the robust state estimator, and then rewriting it into a compact form to obtain the observation vector, regression matrix, error vector, and error covariance matrix; minimizing the robust scale of the residual vector, including defining a robust scale estimator and constructing a weight function to solve for the robust scale estimator.
[0014] Furthermore, construct with σ y The Gaussian kernel function with kernel bandwidth σ is used to define the cost function; kernel bandwidth σ y The expression is as follows:
[0015]
[0016] Where y(k) is the engine output. R(k) represents the estimated engine output at time k-1 based on the prediction at time k; R(k) is the measurement noise covariance matrix at time k, and the superscript -1 indicates matrix inversion; ||t|| M =(t T Mt) 1 / 2 .
[0017] Furthermore, selecting sampling points and constructing the state variables for those sampling points includes:
[0018]
[0019] In the formula, χ i (k-1|k-1) represents the state variable of the i-th sampling point at time k-1, where i = 0, 1, ..., 2n; (·) i This indicates the i-th column or i-th row of the matrix enclosed in parentheses; P(k-1) represents the state estimate and state error covariance matrix at time k-1; n is the system state dimension, and λ is the composite scaling factor, λ = α. 2 (n+ε)-n; 0<α<1, ε is the scaling factor;
[0020] Solving the problem using singular value decomposition.
[0021] Furthermore, the cost function is as follows:
[0022]
[0023] In the formula, α′ and β′ are adjustable weighting coefficients. Therefore, σ y is a Gaussian kernel function with a kernel bandwidth; x(k) represents the engine state variables. Let represent the state estimate at time k-1; P(k|k-1) is the error covariance matrix of the prediction at time k-1 to time k; and y(k) is the engine output. H(k) represents the estimated engine output at time k-1, and H(k) is the pseudo-measurement matrix.
[0024] Furthermore, the Kalman filter gain matrix is solved based on the cost function, as follows:
[0025]
[0026] In the formula, This represents the partial derivative, where ξ is an intermediate variable, and the expression is as follows:
[0027]
[0028] Let α′=1, Substituting into the above equation, we get:
[0029]
[0030] Summarized as follows:
[0031]
[0032] In the formula, The Kalman filter gain matrix is expressed as follows:
[0033]
[0034] Furthermore, the robust state estimator is expressed as:
[0035]
[0036] In the formula, x(k) represents the engine state quantity. This represents the state estimate at time k-1; For measuring the noise sequence, I is the identity matrix; y(k) is the engine output. H(k) represents the estimated engine output at time k-1 for time k; H(k) is the pseudo-measurement matrix.
[0037] The robust state estimator is rewritten in the following compact form:
[0038]
[0039] Among them, the observation vector Regression Matrix Error vector And the error covariance matrix is:
[0040]
[0041] Where E(·) represents the expectation operation. P(k|k-1) is the measurement noise matrix after unscented transformation, and P(k|k-1) is the error covariance matrix of the prediction at time k-1 to time k.
[0042] If the robust scale of the residual vector can be minimized, then the state estimation result is robust, i.e.:
[0043]
[0044] In the formula, Let k be the estimated value of the engine state variables at time k. For robust scaling estimators; Let r be the residual vector. i (x(k)) represents the i-th residual component in the residual vector, i = 1, 2, ..., m; m is the observation vector. The dimension of.
[0045] Furthermore, the robust scaling estimator is constructed as follows:
[0046]
[0047] In the formula, δ∈[0,1], and the function ρ(·) satisfies A bounded even function;
[0048] Robust scaling estimator Solve iteratively:
[0049]
[0050] In the formula, For relative standard residuals, the weighting function is:
[0051]
[0052] in, observation vector The i-th component.
[0053] A terminal device includes a processor, a memory, and a computer program stored in the memory; when the processor executes the computer program, it implements the state estimation method for the ATR engine control system based on adaptive robust UKF.
[0054] A computer-readable storage medium storing a computer program; when executed by a processor, the computer program implements the state estimation method for an ATR engine control system based on adaptive robust UKF.
[0055] Compared with the prior art, the present invention has the following technical features:
[0056] 1. This invention employs an adaptive filter kernel bandwidth adjustment mechanism, which can adaptively adjust according to changes in noise in the environment, making it suitable for complex environments with a mixture of Gaussian and non-Gaussian noise.
[0057] 2. This invention combines the maximum correlation entropy criterion to improve the filter's ability to handle high-order moment characteristics of nonlinear systems and enhance the accuracy of state estimation, especially in noisy environments.
[0058] 3. By introducing a robust S-estimator, this invention effectively isolates the influence of abnormal measurement data, ensuring that even when the sensor malfunctions or measurement is abnormal, the filter can still provide stable state estimation results, thus guaranteeing the safe operation of the engine under extreme conditions.
[0059] 4. This invention reduces the computational burden while maintaining high estimation accuracy and fast response performance by employing batch regression and optimized iterative algorithms.
[0060] 5. This invention is particularly applicable to the atmospheric and external operating modes of ATR engines, and can perform effective state estimation in turbojet mode and rocket mode, adapting to the state monitoring needs of engines in various operating environments. Attached Figure Description
[0061] Figure 1 This is a schematic diagram of the ATR engine structure;
[0062] Figure 2 This is a flowchart illustrating the method of the present invention. Detailed Implementation
[0063] 1. Working principle of ATR engine.
[0064] The structural schematic diagram of a monocomponent liquid propellant ATR engine is shown below. Figure 1 As shown. Unlike the combustion chamber of a turbojet engine, the ATR engine uses a gas generator component from a rocket engine, thus it can be seen as a combination of a turbojet engine and a rocket engine. The working principle of the ATR engine is as follows: In the gas generator, the single-component liquid propellant is catalytically decomposed into high-temperature, high-pressure fuel-rich gas. The fuel-rich gas flows through the turbine, thereby driving the turbine components, which in turn drive the compressor to do work on the outside air entering through the intake. The pressurized high-pressure air enters the mixing chamber through the bypass duct, where it is fully mixed with the fuel-rich gas at the turbine outlet and undergoes chemical combustion in the combustion chamber to generate high-temperature, high-pressure gas. Finally, the high-temperature, high-pressure gas is fully expanded in the tailpipe and discharged at high speed, thereby generating thrust.
[0065] Compared to traditional turbojet engines, the ATR engine separates the airflow path through the compressor from the high-temperature, high-pressure gas flow path through the turbine, so that the turbine inlet parameters are not affected by the free flow. This achieves decoupling of the compressor airflow path and the turbine gas flow path in terms of thermodynamic parameters. This special structure and working method gives the ATR engine better hypersonic performance.
[0066] Figure 1 The definitions of ATR engine component names and section numbers are given, and the names of each section number are shown in Table 1.
[0067] Table 1. ATR Engine Section Numbers and Names
[0068]
[0069] 2. Adaptive Robust UKF Design.
[0070] like Figure 2 As shown, the state estimation method for an ATR engine control system based on adaptive robust UKF provided by this invention comprises the following steps:
[0071] Step 1: Establish a nonlinear discrete model of the ATR engine control system, including the following sub-steps:
[0072] Step 1.1: Establish a state-space model for the ATR engine control system:
[0073]
[0074] In the formula, the superscript dots of the parameters indicate the differential of the parameters, and the number of dots is the order of the differential, the same below; x is the engine state variable, u is the engine input variable, y is the engine output variable, n, m, and p are the number of state variables, internal input variables, external input variables, and output variables of the state-space model, respectively, and n < m < p, f(·) represents the state nonlinear function, and g(·) represents the measurement nonlinear function.
[0075] Step 1.2: Determine the specific parameters of engine state variables, engine input variables, and engine output variables.
[0076]
[0077] In the formula, the superscript T in (·) indicates transpose, the same below; N1 is the low-pressure rotor speed, N2 is the high-pressure rotor speed, W F Main fuel flow rate, A8 is the nozzle throat area, T t7 SM represents the total temperature at the combustion chamber outlet. c This refers to the compressor surge margin.
[0078] Step 1.3: Discretize the state-space model.
[0079] x(k)=f k-1 (x(k-1),u(k))+w(k-1)
[0080] y(k)=g k (x(k),u(k))+v(k)
[0081] In the formula, x(k), y(k), and u(k) are the discretized representations of the engine state variables x, engine output variables y, and engine input variables u at time k, and f k-1 (·), g k (·) represents the discretized representation of the state nonlinear function and measurement nonlinear function at time k-1 and time k, x(k-1) is the engine state quantity at time k-1, w(k-1) is the system process noise at time k-1, and v(k) is the system measurement noise at time k.
[0082] Step 1.4, calculate the noise covariance matrix:
[0083] E[w(k-1)wT [(k-1)]=Q(k-1)
[0084] E[v(k)v T [(k)]=R(k)
[0085] In the formula, Q(k-1) is the system process noise covariance matrix at time k-1, R(k) is the measurement noise covariance matrix at time k, E[·] represents the expected value of the calculated variable, and the superscript T of the parameter indicates transpose, the same below.
[0086] Step 2, the time update process, includes the following sub-steps:
[0087] Step 2.1, Variable initialization:
[0088]
[0089] In the formula, x0 is the engine state variable at the initial moment. P0 represents the estimated value of the engine state variables at the initial moment, and P0 is the state error covariance matrix at the initial moment. In this scheme, the superscript ∧ of the parameter indicates the estimated value of the parameter.
[0090] Step 2.2: Select 2n+1 sampling points and construct the state variables for the sampling points:
[0091]
[0092] In the formula, χ i (k-1|k-1) represents the state variable of the i-th sampling point at time k-1, i = 0, 1, ..., 2n; (·)i represents the i-th column or i-th row of the matrix in parentheses, the same below; for example yes The i-th column or the i-th row; P(k-1) represents the state estimate and state error covariance matrix at time k-1; n is the system state dimension, and λ is the composite scaling factor, given by the following equation:
[0093] λ=α 2 (n+ε)-n
[0094] In the formula, the parameter α determines the distribution of the sampling points, and is usually chosen to be 0 < α < 1; ε is the scaling factor, which is usually set to 3-n.
[0095] Solve The Cholesky decomposition method is commonly used, but it requires the covariance matrix to be positive definite. Due to uncertainties, noise, and truncation errors in numerical solutions during engine state estimation, the positive definiteness of the covariance matrix is often not satisfied, leading to decomposition failure and a rapid decline in filter performance. This invention solves for the square root of the covariance matrix using Singular Value Decomposition (SVD), avoiding the failure problem caused by positive definiteness and thus ensuring the reliability of state estimation.
[0096] Step 2.3: Predict the estimated engine state variables and error covariance matrix at time k:
[0097]
[0098]
[0099] In the formula, P(k|k-1) represents the estimated engine state variables and error covariance matrix predicted at time k-1 for time k, and χ² represents the value of the engine state variables predicted at time k. i (k|k-1) represents the state quantity predicted by the i-th sampling point at time k-1 for time k. Let S and W be the weighting coefficients of the state estimate and error covariance corresponding to the i-th sampling point, respectively, and their expressions are as follows:
[0100]
[0101] In the formula, β is a parameter related to the prior knowledge distribution of the engine state variable x(k) at time k, which is generally set to 2 under a Gaussian distribution.
[0102] Step 3, the measurement update process, includes the following sub-steps:
[0103] Step 3.1: Update the state variables of the sampling points using the state estimate and error covariance predicted at time k-1.
[0104]
[0105] Step 3.2: Use the measurement equation function to transform and update the state variables of the sampled points:
[0106] γ i (k|k-1)=h k (χ i (k∣k-1)),i=0...2n
[0107] In the formula, h k (·) represents the measurement equation function, γ i (k|k-1) represents the state variables of the sampled points after the transformation.
[0108] Step 3.3: Calculate the estimated engine output and covariance matrix using unscented transformation.
[0109]
[0110] in, P represents the estimated engine output at time k-1 for prediction at time k. yy R(k) represents the covariance matrix, and R(k) is the measurement noise covariance matrix at time k.
[0111] Step 3.4, calculate the state-measurement cross-covariance matrix P. xy (k):
[0112]
[0113] Step 4, the process of calculating the Kalman gain matrix using the maximum correlation entropy criterion, includes the following sub-steps:
[0114] Step 4.1, define the pseudo-measurement matrix as follows:
[0115] H(k)=(P -1 (k∣k-1)P xy (k)) T
[0116] Where the superscript -1 of the parameter indicates inversion, such as P -1 (k|k-1) denotes the inverse matrix of P(k|k-1), and the same applies below; P(k|k-1) is the error covariance matrix of the prediction at time k-1 to time k.
[0117] Theorem 1 (Statistical Linear Regression Theorem) states that for a given nonlinear function that is approximately linearized as y = g(x) ≈ Ax + c, where A is the coefficient matrix and c is the constant term, the coefficient solution of the minimum weighted sum of squared errors is:
[0118]
[0119] in, Let P be the mean of x and y. xy Let P be the cross covariance matrix of x and y. xx Let x be the covariance matrix.
[0120] And the mean of the error term e = y - Ax - c and variance P ee satisfy:
[0121]
[0122] Where P yy Let y be the covariance matrix. Let c be the mean.
[0123] Step 4.2, according to Theorem 1, use a pseudo-measurement matrix to convert the nonlinear engine output into a linear form:
[0124]
[0125] In the formula, To measure the noise sequence, N(·) represents a Gaussian distribution. The measurement noise matrix after unscented transformation is expressed as follows:
[0126]
[0127] Step 4.3, define the cost function:
[0128] Correlation entropy is a measure of the local similarity between two random variables, encompassing their higher-order moments. Currently, this method has been successfully applied to linear and nonlinear filtering, demonstrating its advantages in processing non-Gaussian signals. Consider two joint density functions F... X,Y(x,y) The correlation entropy of random variables X and Y is:
[0129] V(X,Y)=E[κ(X,Y)]=∫κ(x,y)dF X,Y (x,y)
[0130] In the formula, E represents the expectation operator, and κ(·,·) represents a positive definite bounded Mercer kernel function, given by the following formula:
[0131]
[0132] In the formula, e = xy is the random variable error, and σ > 0 is the kernel function bandwidth.
[0133] In many practical situations, the joint density function of random variables cannot be accurately obtained. This invention approximates the correlation entropy function using a finite set of data samples and a sample mean estimator.
[0134]
[0135] In the formula, N is the number of samples, and e(i) = x(i) - y(i). It is taken from F X,Y Using the sample data (x, y), this formula reveals the positive definite boundedness of the Gaussian correlation entropy, and it reaches its maximum value only when X = Y. When performing a Taylor expansion on the Gaussian kernel function, we obtain:
[0136]
[0137] Analysis shows that the correlation entropy is the weighted sum of all even moments of a random variable XY. When a suitable kernel bandwidth is chosen, the correlation entropy captures the higher-order moment characteristics of the variable. For Gaussian signals, the mean and covariance are sufficient to describe the entire distribution characteristics of the signal. Therefore, traditional Kalman filtering and its derivatives that satisfy the minimum mean square error criterion can achieve optimal estimation of Gaussian systems. However, when non-Gaussian signals exist in the system state or measurement equations, the maximum correlation entropy criterion is more suitable for handling higher-order information, thereby improving the accuracy of state estimation.
[0138] Given an error data sequence The cost function for maximum relevance entropy is defined as follows:
[0139]
[0140] Suppose x is a vector parameter to be estimated in an adaptive system. The parameter estimation problem based on the maximum correlation entropy criterion can be transformed into the following optimization problem:
[0141]
[0142] In the formula, Let Ω represent the optimal estimate of the state, and let Ω represent the feasible set of parameters.
[0143] Inspired by the maximum correlation entropy criterion, this invention defines a weighted combined cost function that includes state and measurement output errors. It uses a weighted least squares method to handle system noise and the maximum correlation entropy criterion to handle non-Gaussian sensor measurement noise. The aim is to capture the higher-order moment characteristics of sensor noise as much as possible through the maximum correlation entropy criterion, thereby improving the estimation accuracy under non-Gaussian noise.
[0144] The cost function is defined as follows:
[0145]
[0146] In the formula, α′ and β′ are adjustable weighting coefficients, ||t|| M =(t T Mt) 1 / 2 ,For example:
[0147] Therefore, σ y The Gaussian kernel function with kernel bandwidth.
[0148] Gaussian kernel bandwidth σ yThe kernel bandwidth has a significant impact on filter performance. Currently, kernel bandwidth selection is still based on experience and trial-and-error methods, lacking adaptability and robustness in the face of large variations in the operating environment or noise characteristics. This invention proposes an adaptive adjustment method for the kernel function bandwidth based on the error between actual measurements and estimates. This adaptive law enables robust adaptive adjustment of the filter, ensuring good filtering performance even with large variations in environmental noise characteristics. Its expression is as follows:
[0149]
[0150] Based on the maximum correlation entropy criterion, the engine state variables estimate at time k can be transformed into the following optimization problem:
[0151]
[0152] When the partial derivative of the cost function is zero, the state is the optimal solution.
[0153] Step 4.4, Solve for the Kalman filter gain matrix based on the cost function:
[0154]
[0155] In the formula, This represents the partial derivative, where ξ is an intermediate variable, and the expression is as follows:
[0156]
[0157] Let α′=1, Substituting into the above equation, we get:
[0158]
[0159] Summarized as follows:
[0160]
[0161] In the formula, The Kalman filter gain matrix is expressed as follows:
[0162]
[0163] Step 5, the robust state estimator design process, includes the following sub-steps:
[0164] Step 5.1, let the robust state estimator be in batch regression form:
[0165]
[0166] In the formula, For measuring the noise sequence, I is the identity matrix.
[0167] The robust state estimator is rewritten in the following compact form:
[0168]
[0169] Among them, the observation vector Regression Matrix Error vector And the error covariance matrix is:
[0170]
[0171] Where E(·) represents the expectation operation. This is the measurement noise matrix after unscented transformation.
[0172] If the robust scale of the residual vector can be minimized, then the state estimation result is robust, i.e.:
[0173]
[0174] In the formula, Let k be the estimated value of the engine state variables at time k. As a robust scaling estimator, minimizing the robust scaling estimator can ensure isolation of outlier measurements (engine output y); for example, the least squares median (LMS) estimator will minimize the median of the absolute residuals; the weighted least squares estimator (WLS) will minimize the standard deviation of the residuals; Let r be the residual vector. i (x(k)) represents the i-th residual component in the residual vector, i = 1, 2, ..., m; m is the observation vector. The dimension of.
[0175] Step 5.2: Define a robust scaling estimator for the robust state estimator to minimize the robust scaling of the residual vector.
[0176]
[0177] In the formula, δ∈[0,1], and the function ρ(·) satisfies A bounded even function.
[0178] Robust scaling estimator Generally, iterative solutions are used:
[0179]
[0180] In the formula, For relative standard residuals, the weighting function is:
[0181]
[0182] in, observation vector The i-th component.
[0183] In one embodiment of the present invention, the selected ρ(·) function is a double square function, defined as follows:
[0184]
[0185] In the formula, c is the generalized relative residual threshold.
[0186] Step 5.3, calculate the weight function as follows:
[0187]
[0188] This yields a robust state estimator, which is applied to the front end of an adaptive UKF filter. The reweighting function of the robust scaling estimator is then used to perform robust scaling correction on the residual vector, thereby isolating the influence of outlier measurements.
[0189] The robust scaling in the form of traditional standard deviation corresponds to an unbounded ρ function estimator, i.e. And since δ = 1, this method will use an equal weight for all observations. Weighting reduces robustness to outliers. Our method reweights all observations, assigning a weight of 1 to normal observations and a weight less than 1 to outliers. The further an outlier is from the normal value, the smaller its weight, down to 0, thus isolating it. Minimizing the residual robust scale of the robust state estimator has proven to provide high-precision robust estimation results.
[0190] Step 6, updating the state estimate and covariance matrix, includes the following sub-steps:
[0191] Step 6.1, update the engine state quantity estimates:
[0192]
[0193] In the formula, ε is an intermediate variable, and the parameter subscript t represents the result of the t-th iteration with respect to that parameter, for example... Represents the result obtained in the t-th iteration. Solving the state variable matrix is essentially a fixed-point iterative equation, which can generally be solved iteratively. Once the required accuracy is met, the iterative solution approximates the state estimate. To satisfy the algorithm's recursive form and reduce computational burden, this invention adopts a single-iteration approach. If x(k) is approximated, then the intermediate variable ξ can be simplified to:
[0194]
[0195] at this time It can be approximated by the following formula:
[0196]
[0197] in, There exists another form of expression; the proof is given in Theorem 1:
[0198]
[0199] Step 6.2, Update the state error covariance matrix:
[0200]
[0201] Theorem 1 When the Gaussian kernel bandwidth σ y As the value approaches infinity, the filter converges to the standard UKF filter.
[0202] Proof: The time update steps in the above filter derivation are completely consistent with the standard UKF; only proof is required.
[0203] When σ y As the value approaches infinity, x(k) and P(k) in the measurement update step can be the same as the update algorithm in the standard UKF.
[0204] Lemma 2: If square matrix A is invertible, square matrix C is invertible, and matrix A+BCD is invertible, then
[0205] (A+BCD) -1 =A -1 -A -1 B(C -1 +DA -1 B)DA -1
[0206] Consideration
[0207] Let A = P -1 (k|k-1),B=ξH T (k), Since the invertibility condition of Lemma 2 is satisfied, we can obtain the following based on Lemma 2:
[0208]
[0209] Assume A=I, B=I, D=H(k)P(k|k-1)ξH T (k), applying Lemma 2 in reverse, we get:
[0210]
[0211] This formula is also a formula. Find the Kalman gain matrix Another form
[0212] When σ y As ξ approaches infinity, ξ approaches 1. Then H(k) and... Substituting into the above equation, we get:
[0213]
[0214] Therefore, it can be seen that the Kalman gain matrix and the standard UKF gain matrix in this formula are different. The calculations are the same. Consider the state covariance matrix P(k):
[0215]
[0216] Will and Substituting into the above equation, we get:
[0217]
[0218] Based on It can be known Substituting into the above equation, we get:
[0219]
[0220] The covariance matrix and the standard UKF covariance matrix The calculations are the same. Therefore, combining the above proof process, it can be seen that when σ y As the threshold of infinity approaches infinity, the adaptive UKF filter proposed in this paper converges to the standard UKF filter.
[0221] Theorem 1 is now proven.
[0222] When the system measurement suddenly has a large outlier or shot noise, i.e., ||y(k)|| ∞ As ξ approaches infinity, ξ approaches zero. And P(k) = P(k|k-1). The filter can maintain the state of the previous moment in a timely manner and has a certain degree of robustness to the noise.
[0223] The above embodiments are only used to illustrate the technical solutions of this application, and are not intended to limit them. Although this application has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features. Such modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of this application, and should all be included within the protection scope of this application.
Claims
1. A state estimation method for an ATR engine control system based on adaptive robust UKF, characterized in that, include: Establish a nonlinear discrete model for the ATR engine control system, including establishing a state-space model for the ATR engine control system, determining the specific parameters of engine state variables, engine input variables and engine output variables, discretizing the state-space model and calculating the noise covariance matrix. The time update process includes initializing the engine state variables and the state error covariance matrix, selecting sampling points and constructing the state variables of the sampling points, and predicting the estimated values of the engine state variables and the error covariance matrix. The measurement update process includes updating the state variables of the sampling points, transforming the updated state variables of the sampling points using measurement equation functions, calculating the estimated value of engine output and covariance matrix through unscented transformation, and calculating the state-measurement cross-covariance matrix. The process of calculating the Kalman gain matrix using the maximum correlation entropy criterion includes defining a pseudo-measurement matrix, using the pseudo-measurement matrix to convert the nonlinear engine output into a linear form, defining a cost function, and solving for the Kalman filter gain matrix based on the cost function. The process of updating the state estimate and covariance matrix includes updating the engine state quantity estimates and updating the state error covariance matrix.
2. The state estimation method for an ATR engine control system based on adaptive robust UKF as described in claim 1, characterized in that, The method further includes: The robust state estimator design process involves: first, designing the batch regression form of the robust state estimator, and then rewriting it into a compact form to obtain the observation vector, regression matrix, error vector, and error covariance matrix; minimizing the robust scale of the residual vector, including defining a robust scale estimator and constructing a weight function to solve for the robust scale estimator.
3. The state estimation method for an ATR engine control system based on adaptive robust UKF as described in claim 1, characterized in that, Construct with σ y The Gaussian kernel function with kernel bandwidth σ is used to define the cost function; kernel bandwidth σ y The expression is as follows: Where y(k) is the engine output. R(k) represents the estimated engine output at time k-1 based on the prediction at time k; R(k) is the measurement noise covariance matrix at time k, and the superscript -1 indicates matrix inversion; ||t|| M =(t T Mt) 1 / 2 .
4. The state estimation method for an ATR engine control system based on adaptive robust UKF as described in claim 1, characterized in that, Selecting sampling points and constructing their state variables includes: In the formula, χ i (k-1|k-1) represents the state variable of the i-th sampling point at time k-1, i = 0, 1, ..., 2n; (·)i represents the i-th column or i-th row of the matrix in parentheses; P(k-1) represents the state estimate and state error covariance matrix at time k-1; n is the system state dimension, and λ is the composite scaling factor, λ = α. 2 (n+ε)-n; 0<α<1, ε is the scaling factor; Solving using the singular value decomposition method 5. The state estimation method for an ATR engine control system based on adaptive robust UKF according to claim 1, characterized in that, The cost function is as follows: In the formula, α′ and β′ are adjustable weighting coefficients. Therefore, σ y is a Gaussian kernel function with a kernel bandwidth; x(k) represents the engine state variables. Let represent the state estimate at time k-1; P(k|k-1) is the error covariance matrix of the prediction at time k-1 to time k; and y(k) is the engine output. H(k) represents the estimated engine output at time k-1, and H(k) is the pseudo-measurement matrix.
6. The state estimation method for an ATR engine control system based on adaptive robust UKF according to claim 1, characterized in that, The Kalman filter gain matrix is solved based on the cost function, as follows: In the formula, This represents the partial derivative, where ξ is an intermediate variable, and the expression is as follows: Let α′=1, Substituting into the above equation, we get: Summarized as follows: In the formula, The Kalman filter gain matrix is expressed as follows:
7. The state estimation method for an ATR engine control system based on adaptive robust UKF according to claim 2, characterized in that, The robust state estimator is represented as: In the formula, x(k) represents the engine state quantity. This represents the state estimate at time k-1; For measuring the noise sequence, I is the identity matrix; y(k) is the engine output. H(k) represents the estimated engine output at time k-1 for time k; H(k) is the pseudo-measurement matrix. The robust state estimator is rewritten in the following compact form: Among them, the observation vector Regression Matrix Error vector And the error covariance matrix is: Where E(·) denotes the expectation operation. P(k|k-1) is the measurement noise matrix after unscented transformation, and P(k|k-1) is the error covariance matrix of the prediction at time k-1 to time k. If the robust scale of the residual vector can be minimized, then the state estimation result is robust, i.e.: In the formula, Let k be the estimated value of the engine state variables at time k. For robust scaling estimators; Let r be the residual vector. i (x(k)) represents the i-th residual component in the residual vector, i = 1, 2, ..., m; m is the observation vector. The dimension of.
8. The state estimation method for an ATR engine control system based on adaptive robust UKF according to claim 7, characterized in that, The robust scaling estimator is constructed as follows: In the formula, δ∈[0,1], and the function ρ(·) satisfies A bounded even function; Robust scaling estimator Solve iteratively: In the formula, For relative standard residuals, the weighting function is: in, observation vector The i-th component.
9. A terminal device, comprising a processor, a memory, and a computer program stored in the memory; characterized in that, When the processor executes the computer program, it implements the state estimation method for the ATR engine control system based on adaptive robust UKF as described in any one of claims 1-8.
10. A computer-readable storage medium storing a computer program; characterized in that, When the computer program is executed by the processor, it implements the state estimation method for the ATR engine control system based on adaptive robust UKF as described in any one of claims 1-8.
Citation Information
Patent Citations
Gas pipeline parameter estimation method based on improved ARUKF
CN110532517A
Self-adaptive unscented Kalman filter state estimation method with noise estimator
CN111985093A