CN102867129B - 基于可变数据长度最大信息量-可信度准则的飞行器建模方法 - Google Patents

Abstract

本发明公开了一种基于可变数据长度最大信息量-可信度准则的飞行器建模方法，用于解决现有的最大信息量准则不考虑数据长度而导致飞行试验给出的气动模型和参数验证正确性差的技术问题。技术方案是通过在最大信息量准则考虑数据长度并且引入可信度参数，可以根据飞行器的不同飞行试验因素修正建模准则；对测量方差估计R_j和R_j+1的U-D分解，得到了标量模型选择和验证判别式。便于直接根据飞行试验数据建立飞行器气动力、力矩模型，避免了最大信息量准则未直接考虑数据长度导致用不同飞行试验数据建立和验证气动模型不正确的技术问题。

Description

基于可变数据长度最大信息量-可信度准则的飞行器建模方法

技术领域

本发明涉及一种飞行器建模方法，特别是涉及一种基于可变数据长度最大信息量-可信度准则的飞行器建模方法。

背景技术

根据飞机气动模型和参数不仅可以确定飞机的操纵稳定性，还可为地面和空中仿真器提供正确的数学模型；验证飞机气动参数的风洞实验和理论计算结果；为飞机控制系统的设计和改进提供基本数据；鉴定定型飞机的飞行品质；研究高性能飞机的飞行品质；进行飞机失事的事故分析等等；准确地建立飞机数学模型问题与通过基本定律、定理等机理建模的理论方法截然不同，主要根据实验所得的输入和输出数据建立模型，其基本理论依据为非线性系统辨识学和非线性飞行动力学；当飞机作小迎角小扰动飞行时，气动力和力矩模可以用台劳级数展开取一次项，即Bryan模型表示。当马赫数、高度一定时，这一模型是线性定常模型，此模型因为形式简单而一直沿用至今，成为气动数学模型的基石；采用这种模型，飞行器系统辨识就成了对已知数学模型的系统参数估计了；现代战斗机、战术导弹在作战时需要较大机动、过失速甚至尾旋，其迎角可以从十几度、几十度直至一百多度，已不能采用线性模型；飞机大迎角形成的脱体涡、分离涡所引起的非定常下洗流场、使得定常模型也不能再适用了。研究在大迎角下飞行器的非定常、非线性气动模型已成为当前飞机研制的迫切需要的问题。然而，非线性气动力的辨识异常复杂，它是一般的非线性系统辨识问题，输入量与状态之间的函数关系很难确定，需要对模型进行辨识；模型辨识的关键是建模判据和优选算法，对于给定的结构形式，应用建模判据来确定模型的最优阶数并从侯选模型中选出最优模型；由于实测数据含有噪声，建模判据不能仅仅考察对现有数据的拟合误差大小，而且综合考虑其它因素，否则将会使模型不正确；通常，建模判据应能使优选出的模型具有以下特点：1．模型很好地拟合现有飞行数据；2．模型各项有明显的物理意义；3．模型能预测类似条件下的实测数据；4．在性能相当的条件下阶次最低；最常用的模型辨识方法是逐步回归法，其原理是逐项将影响显著性的预报因子选入，并将影响小的因子剔除，建立回归方程的方法；这一方法计算简单、比较实用；但这一方法有两个明显的缺点：一是选择标准由人而定，而且没有给出结果的可信度；二是误差积累大，容易漏选和误选；为此，人们对要求较高的飞行器模型辨识问题常常采用最大信息量准则AIC方法，但是该方法处理速度慢，信噪比较小时模型辨识精度差；由于在非线性情况下，只能对飞机非线性方程进行数值积分，进行灵敏度矩阵计算和迭代计算，从而使计算的复杂程度和计算量比线性估计高得多，同时也使模型输出与实验数据之间的拟合变得更加困难，特别是当飞行数据长度不同时，现有AIC准则没有直接考虑不同的数据长度，常常会导致飞行试验给出的气动模型和参数验证不正确。

发明内容

为了克服现有最大信息量准则不考虑数据长度而导致飞行试验给出的气动模型和参数验证正确性差的不足，本发明提供一种基于可变数据长度最大信息量-可信度准则的飞行器建模方法。该方法通过分析数据长度的影响，并且引入可信度参数对最大信息量准则进行了修正，得到了新的模型辨识判据，由新判据建立了指数建模，直接可以用于飞行器的飞行试验建模和模型验证，可以避免根据飞行试验建立和验证飞行器大迎角模型存在的技术问题。

本发明解决其技术问题所采用的技术方案是：一种基于可变数据长度最大信息量-可信度准则的飞行器建模方法，其特点是包括以下步骤：

步骤一、飞行试验待确定的飞行器候选模型的状态方程为

$\stackrel{\·}{x}\left(t\right)=f\left\{f_0\right[x\left(t\right),{\mathrm{\Ω}}_0],f_1[x\left(t\right),{\mathrm{\θ}}_1],...,f_q[{\mathrm{\θ}}_q,x\left(t\right)],t\}---\left(1\right)$

观测方程为

$\left\{\begin{array}{c}y\left(t\right)=g[x\left(t\right),\mathrm{\Ω},t]=g\left\{g_0\right[x\left(t\right),{\mathrm{\Ω}}_0],g_1[x\left(t\right),{\mathrm{\θ}}_1],...,g_q[{\mathrm{\θ}}_q,x\left(t\right)],t\}\\ z\left(t_k\right)=y\left(t_k\right)+v\left(k\right)\end{array}\right.---\left(2\right)$

(1)、(2)式中，x(t)为n维状态向量；y(t)为m维观测向量；f{f₀[x(t)，Ω₀]，f₁[x(t)，θ₁]，...，f_q[θ_q，x(t)]，t}、g{g₀[x(t)，Ω₀]，g₁[x(t)，θ₁]，...，g_q[θ_q，x(t)]，t}为表达式已知的待确定模型结构函数，f₀[x(t)，Ω₀]、g₀[x(t)，Ω₀]为根据物理概念必须选入的模型，f_i[x(t)，θ_i]、g_i[x(t)，θ_i](i＝1，2，…，q)为候选模型，z(t_k)为在t_k时刻对y(t_k)的测量值；Ω为未知维数的参数向量，Ω₀为已知维数的参数向量；v(k)为测量噪声，假定方差为R_k的零均值高斯白噪声；f_i[x(t)，θ_i]、g_i[x(t)，θ_i](i＝1，2，…，q)是否在模型中出现及Ω₀、θ_i(i＝1，2，…，q)的取值需要辨识，q为已知的候选模型个数；

由于对飞行器的模型结构准确度要求较高，给出以下最大信息量-可信度准则：

$\left|\frac{\ln L_2[x\left(t\right),{\mathrm{\Ω}}_2,t]-\ln L_1[x\left(t\right),{\mathrm{\Ω}}_1,t]}{\ln L_1[x\left(t\right),{\mathrm{\Ω}}_1,t]}\right|>\mathrm{\δ},---\left(3\right)$

式中，L₁[x(t)，Ω₁，t]、L₂[x(t)，Ω₂，t]为取不同的参数向量Ω₁和Ω₂的极大似然函数，δ为给定正实数、表示可信度参数，

$\left\{\begin{array}{c}\ln L_1[x\left(t\right),{\mathrm{\Ω}}_1,t]=-\frac{\mathrm{mN}}{2}\ln \left(\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\stackrel{\‾}{\mathrm{\η}}}_1^T\left(k\right)R_1^{-1}\left(k\right){\stackrel{\‾}{\mathrm{\η}}}_1\left(k\right)\right)-\frac{\mathrm{mN}}{2}[\ln \left(\frac{2\mathrm{\π}}{\mathrm{mN}}\right)+1]-\frac{1}{2}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}\ln \left|R_1\left(k\right)\right|\\ \ln L_2[x\left(t\right),{\mathrm{\Ω}}_2,t]=-\frac{\mathrm{mN}}{2}\ln \left(\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\stackrel{\‾}{\mathrm{\η}}}_2^T\left(k\right)R_2^{-1}\left(k\right){\stackrel{\‾}{\mathrm{\η}}}_2\left(k\right)\right)-\frac{\mathrm{mN}}{2}[\ln \left(\frac{2\mathrm{\π}}{\mathrm{mN}}\right)+1]-\frac{1}{2}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}\ln \left|R_2\left(k\right)\right|\end{array}\right.---\left(4\right)$ ， ${\stackrel{\‾}{\mathrm{\η}}}_1\left(k\right)=z\left(t_k\right)-g[x\left(t_x\right),{\mathrm{\Ω}}_1,t_k],$ ${\stackrel{\‾}{\mathrm{\η}}}_2\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\Ω}}_2,t_k],$ $R_1=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\stackrel{\‾}{\mathrm{\η}}}_1\left(k\right){\stackrel{\‾}{\mathrm{\η}}}_1^T\left(k\right),$ $R_2=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\stackrel{\‾}{\mathrm{\η}}}_2\left(k\right){\stackrel{\‾}{\mathrm{\η}}}_2^T\left(k\right),$ Ω₁和Ω₂为不同的参数向量，N为数据长度，ln为自然对数符号；

步骤二、根据假定f₀[x(t)，Ω₀]、g₀[x(t)，Ω₀]、Ω₀=Ω₀已经通过优选算法选入模型，并由以下算法迭代计算得到：

令j=0，1，2，…，q，假定f_j[x(t)，θ_j]、g_j[x(t)，θ_j]、Ω_j已经选入模型，按照以下方式选择其它候选模型：

求（4）式极大值，迭代计算：

$\mathrm{\Δ}{\mathrm{\Ω}}_j=A_j^{-1}{b}_j---\left(5\right)$

以及

$R_j=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\stackrel{\‾}{\mathrm{\η}}}_j\left(k\right){\stackrel{\‾}{\mathrm{\η}}}_j^T\left(k\right),$ ${\stackrel{\‾}{\mathrm{\η}}}_j\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\Ω}}_j,t_k]---\left(6\right)$

(5）、（6）式中： $\mathrm{\Δ}{\mathrm{\Ω}}_j={\mathrm{\Ω}}_j-{\widehat{\mathrm{\Ω}}}_j,$ $b_j=\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\left(\frac{\∂y}{\∂{\mathrm{\Ω}}_j^T}\right)}^T{R}_j^{-1}[z\left(t_k\right)-y\left(t_k\right)],$

$A_j=\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\left(\frac{\∂y}{\∂{\mathrm{\Ω}}_j^T}\right)}^T{R}_j^{-1}\frac{\∂y}{\∂{\mathrm{\Ω}}_j^T}=B_j^T{P}_j^{-1}{B}_j,$ $B_j^T=[{\left(\frac{\∂y\left(t_1\right)}{\∂{\mathrm{\Ω}}_j^T}\right)}^T,{\left(\frac{\∂y\left(t_2\right)}{\∂{\mathrm{\Ω}}_j^T}\right)}^T,\·\·\·,{\left(\frac{\∂y\left(t_N\right)}{\∂{\mathrm{\Ω}}_j^T}\right)}^T]$

$P_j^{-1}=\mathrm{diag}[R_j^{-1},R_j^{-1},\·\·\·R_j^{-1}],$

当两个飞行试验数据长度分别为N、M、设 ${\mathrm{\Ω}}_{j+1}=\left[\begin{array}{c}{\mathrm{\Ω}}_j\\ {\mathrm{\θ}}_{j+1}\end{array}\right],$ θ_j+1的选入或剔除及验证条件为：θ_j+1的选入或剔除模型验证条件为：当

$\frac{\ln \left|R_{\mathrm{Nj}}\left(k\right)\right|-\ln \left|R_{N(j+1)}\left(k\right)\right|}{|m(\ln 2\mathrm{\&pi;}+1)+\ln |R_{\mathrm{Nj}}\left(k\right)\left|\right|}>\mathrm{\&delta;}$ 及 $\left|\frac{m}{2}(N-M)\right[1+\ln \left(2\mathrm{\&pi;}\right)]+\ln \frac{{\left|R_N\left(k\right)\right|}^N}{{\left|R_M\left(k\right)\right|}^M}|<\mathrm{\&epsiv;}---\left(7\right)$

时，θ_j+1、f_j+1[x(t)，θ_j+1]、g_j+1[x(t)，θ_j+1]选入模型，且 ${\mathrm{\&Omega;}}_{j+1}=\left[\begin{array}{c}{\mathrm{\&Omega;}}_j\\ {\mathrm{\&theta;}}_{j+1}\end{array}\right];$ 否则剔除f_j+1[x(t)，θ_j+1]、g_j+1[x(t)，θ_j+1]候选项，且Ω_j+1=Ω_j；

(7)式中：ε为给定正数， $R_{\mathrm{Nj}}=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j^T\left(k\right),$ $R_{N(j+1)}=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}^T\left(k\right),$

$R_{\mathrm{Mj}}=\frac{1}{M}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j^T\left(k\right),$ $R_{M(j+1)}=\frac{1}{M}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}^T\left(k\right),$

${\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_j,t_k],$ ${\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_{j+1},t_k];$

步骤三、由于飞行器测量向量y的维数m较大，采用Gram-Schmidt正交化方法对R_Nj、R_Mj、R_N(j+1)和R_M(j+1)进行U-D分解，R_Nj、R_Mj、R_N(j+1)和R_M(j+1)的U-D分解分别为： $R_{\mathrm{Nj}}=U_{\mathrm{RNj}}{D}_{\mathrm{RNj}}{U}_{\mathrm{RNj}}^T,$ $R_{\mathrm{Mj}}=U_{\mathrm{RMj}}{D}_{\mathrm{RMj}}{U}_{\mathrm{RMj}}^T,$ $R_{N(j+1)}=U_{\mathrm{RN}(j+1)}{D}_{\mathrm{RN}(j+1)}{U}_{\mathrm{RN}(j+1)}^T,$ $R_{M(j+1)}=U_{\mathrm{RM}(j+1)}{D}_{\mathrm{RM}(j+1)}{U}_{\mathrm{RM}(j+1)}^T,$

式中，U_RNj、U_RMj、U_RN(j+1)、U_RM(j+1)为单位上三角阵；

D_RNj=diag[d_RNj(1)，d_RNj(2)，…，d_RNj(m)]，D_RN(j+1)=diag[d_RN(j+1)(1)，d_RN(j+1)(2)，…，d_RN(j+1)(m)]，

D_RMj=diag[d_RMj(1)，d_RMj(2)，…，d_RMj(m)]，D_RM(j+1)=diag[d_RM(j+1)(1)，d_RM(j+1)(2)，…，d_RM(j+1)(m)]；

diag为对角符号；

模型验证的最大信息量准则写成：当

$\left\{\begin{array}{c}\frac{\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}[\ln d_{\mathrm{RNj}}\left(i\right)-\ln d_{\mathrm{RN}(j+1)}\left(i\right)]}{|m(\ln 2\mathrm{\&pi;}+1)+\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\ln d_{\mathrm{RNj}}\left(i\right)|}>\mathrm{\&delta;}\\ \left|\frac{m}{2}(N-M)\right[1+\ln \left(2\mathrm{\&pi;}\right)]+N\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\ln d_{\mathrm{RNj}}\left(i\right)-M\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\ln d_{\mathrm{RMj}}\left(i\right)|<\mathrm{\&epsiv;}\end{array}\right.---\left(8\right)$

成立时，θ_j+1、f_j+1[x(t)，θ_j+1]、g_j+1[x(t)，θ_j+1]选入模型正确，且 ${\mathrm{\&Omega;}}_{j+1}=\left[\begin{array}{c}{\mathrm{\&Omega;}}_j\\ {\mathrm{\&theta;}}_{j+1}\end{array}\right];$ 否则剔除f_j+1[x(t)，θ_j+1]、g_j+1[x(t)，θ_j+1]候选项，且Ω_j+1=Ω_j。

本发明的有益效果是：由于通过在最大信息量准则考虑数据长度并且引入可信度参数，可以根据飞行器的不同飞行试验因素修正建模准则；对测量方差估计R_j和R_j+1的U-D分解，得到了标量模型选择和验证判别式，便于直接根据飞行试验数据建立飞行器气动力、力矩模型，避免了最大信息量准则未直接考虑数据长度导致用不同飞行试验数据建立和验证气动模型不正确的技术问题。

下面结合具体实施方式对本发明作详细说明。

具体实施方式

本发明基于可变数据长度最大信息量-可信度准则的飞行器建模方法具体步骤如下：

1、许多飞行器在迎角小于60度时常用候选模型形式为：

$\stackrel{\&CenterDot;}{x}\left(t\right)=\mathrm{\&Phi;}\left({\mathrm{\&Omega;}}_0\right)f_0\left[x\left(t\right)\right]+{\mathrm{\&theta;}}_1{f}_1\left[x\left(t\right)\right]+...+{\mathrm{\&theta;}}_q{f}_q\left[x\left(t\right)\right]---\left(1\right)$

$\left\{\begin{array}{c}y\left(t\right)=g[x\left(t\right),\mathrm{\&Omega;},]=\mathrm{\&Psi;}\left({\mathrm{\&Omega;}}_0\right)g_0\left[x\left(t\right)\right]+{\mathrm{\&theta;}}_1{g}_1\left[x\left(t\right)\right]+...+g_q[{\mathrm{\&theta;}}_q,x\left(t\right)]\\ z\left(t_k\right)=y\left(t_k\right)+v\left(k\right)\end{array}\right.---\left(2\right)$

(1)、(2)式中，，x(t)为n维状态向量；y(t)为m维观测向量，Φ(Ω₀)f₀[x(t)]、Ψ(Ω₀)g₀[x(t)]为根据物理概念必须选入的模型，θ_if_i[x(t)]、θ_ig_i[x(t)](i＝1，2，…，q)为候选模型，z(t_k)为在t_k时刻对y(t_k)的测量值；Ω为未知维数的参数向量，Ω₀为已知维数的参数向量；v(k)为测量噪声，假定方差为R_k的零均值高斯白噪声；θ_if_i[x(t)]、θ_ig_i[x(t)](i＝1，2，…，q)是否在模型中出现及Ω₀、θ_i(i＝1，2，…，q)的取值需要辨识·，q为已知的候选模型个数；

通常对飞行器的模型结构准确度要求较高，本发明给出以下最大信息量-可信度准则：

$\left|\frac{\ln L_2[x\left(t\right),{\mathrm{\&Omega;}}_2,t]-\ln L_1[x\left(t\right),{\mathrm{\&Omega;}}_1,t]}{\ln L_1[x\left(t\right),{\mathrm{\&Omega;}}_1,t]}\right|>\mathrm{\&delta;},---\left(3\right)$

式中，L₁[x(t)，Ω₁，t]、L₂[x(t)，Ω₂，t]为取不同的参数向量Ω₁和Ω₂的极大似然函数，δ为给定正实数、表示可信度参数，

$\left\{\begin{array}{c}\ln L_1[x\left(t\right),{\mathrm{\&Omega;}}_1,t]=-\frac{\mathrm{mN}}{2}\ln \left(\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_1^T\left(k\right)R_1^{-1}\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_1\left(k\right)\right)-\frac{\mathrm{mN}}{2}[\ln \left(\frac{2\mathrm{\&pi;}}{\mathrm{mN}}\right)+1]-\frac{1}{2}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}\ln \left|R_1\left(k\right)\right|\\ \ln L_2[x\left(t\right),{\mathrm{\&Omega;}}_2,t]=-\frac{\mathrm{mN}}{2}\ln \left(\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_2^T\left(k\right)R_2^{-1}\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_2\left(k\right)\right)-\frac{\mathrm{mN}}{2}[\ln \left(\frac{2\mathrm{\&pi;}}{\mathrm{mN}}\right)+1]-\frac{1}{2}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}\ln \left|R_2\left(k\right)\right|\end{array}\right.---\left(4\right)$ ， ${\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_1\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_1,t_k],$ ${\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_2\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_2,t_k],$ $R_1=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_1\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_1^T\left(k\right),$ Ω₁和Ω₂为不同的参数向量，N为数据长度，ln为自然对数符号；

2、根据假定f₀[x(t)，Ω₀]、g₀[x(t)，Ω₀]、Ω₀=Ω₀已经通过优选算法选入模型，并由以下算法迭代计算得到：

令j=0，1，2，…，q，假定f_j[x(t)，θ_j]、g_j[x(t)，θ_j]、Ω_j已经选入模型，按照以下方式选择其它候选模型：

求（4）式极大值，迭代计算：

$\mathrm{\&Delta;}{\mathrm{\&Omega;}}_j=A_j^{-1}{b}_j---\left(5\right)$

以及

$R_j=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j^T\left(k\right),$ ${\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_j,t_k]---\left(6\right)$

（5）、（6）式中： $\mathrm{\&Delta;}{\mathrm{\&Omega;}}_j={\mathrm{\&Omega;}}_j-{\widehat{\mathrm{\&Omega;}}}_j,$ $b_j=\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\left(\frac{\&PartialD;y}{\&PartialD;{\mathrm{\&Omega;}}_j^T}\right)}^T{R}_j^{-1}[z\left(t_k\right)-y\left(t_k\right)],$

$A_j=\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\left(\frac{\&PartialD;y}{\&PartialD;{\mathrm{\&Omega;}}_j^T}\right)}^T{R}_j^{-1}\frac{\&PartialD;y}{\&PartialD;{\mathrm{\&Omega;}}_j^T}=B_j^T{P}_j^{-1}{B}_j,$ $B_j^T=[{\left(\frac{\&PartialD;y\left(t_1\right)}{\&PartialD;{\mathrm{\&Omega;}}_j^T}\right)}^T,{\left(\frac{\&PartialD;y\left(t_2\right)}{\&PartialD;{\mathrm{\&Omega;}}_j^T}\right)}^T,\&CenterDot;\&CenterDot;\&CenterDot;,{\left(\frac{\&PartialD;y\left(t_N\right)}{\&PartialD;{\mathrm{\&Omega;}}_j^T}\right)}^T]$

$P_j^{-1}=\mathrm{diag}[R_j^{-1},R_j^{-1},\&CenterDot;\&CenterDot;\&CenterDot;R_j^{-1}],$

当两个飞行试验数据长度分别为N、M、设 ${\mathrm{\&Omega;}}_{j+1}=\left[\begin{array}{c}{\mathrm{\&Omega;}}_j\\ {\mathrm{\&theta;}}_{j+1}\end{array}\right],$ θ_j+1的选入或剔除及验证条件为：θ_j+1的选入或剔除模型验证条件为：当

$\frac{\ln \left|R_{\mathrm{Nj}}\left(k\right)\right|-\ln \left|R_{N(j+1)}\left(k\right)\right|}{|m(\ln 2\mathrm{\&pi;}+1)+\ln |R_{\mathrm{Nj}}\left(k\right)\left|\right|}>\mathrm{\&delta;}$ 及 $\left|\frac{m}{2}(N-M)\right[1+\ln \left(2\mathrm{\&pi;}\right)]+\ln \frac{{\left|R_N\left(k\right)\right|}^N}{{\left|R_M\left(k\right)\right|}^M}|<\mathrm{\&epsiv;}---\left(7\right)$

时，θ_j+1、θ_jf_j[x(t)]、θ_jg_j[x(t)]选入模型，且 ${\mathrm{\&Omega;}}_{j+1}=\left[\begin{array}{c}{\mathrm{\&Omega;}}_j\\ {\mathrm{\&theta;}}_{j+1}\end{array}\right];$ 否则剔除θ_jf_j[x(t)]、θ_jg_j[x(t)]候选项，且Ω_j+1=Ω_j；

(7)式中：ε为给定正数， $R_{\mathrm{Nj}}=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j^T\left(k\right),$ $R_{N(j+1)}=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}^T\left(k\right),$

$R_{\mathrm{Mj}}=\frac{1}{M}\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j^T\left(k\right),$ $R_{M(j+1)}=\frac{1}{M}\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}^T\left(k\right),$

${\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_j,t_k],$ ${\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_{j+1},t_k];$

3、通常飞行器测量向量y的维数m较大，采用Gram-Schmidt正交化方法对R_Nj、R_Mj、R_N(j+1)和R_M(j+1)进行U-D分解，R_Nj、R_Mj、R_N(j+1)和R_M(j+1)的U-D分解分别为： $R_{\mathrm{Nj}}=U_{\mathrm{RNj}}{D}_{\mathrm{RNj}}{U}_{\mathrm{RNj}}^T,$ $R_{\mathrm{Mj}}=U_{\mathrm{RMj}}{D}_{\mathrm{RMj}}{U}_{\mathrm{RMj}}^T,$ $R_{N(j+1)}=U_{\mathrm{RN}(j+1)}{D}_{\mathrm{RN}(j+1)}{U}_{\mathrm{RN}(j+1)}^T,$ $R_{M(j+1)}=U_{\mathrm{RM}(j+1)}{D}_{\mathrm{RM}(j+1)}{U}_{\mathrm{RM}(j+1)}^T,$

式中，U_RNj、U_RMj、U_RN(j+1)、U_RM(j+1)为单位上三角阵；

D_RNj=diag[d_RNj(1)，d_RNj(2)，…，d_RNj(m)]，D_RN(j+1)=diag[d_RN(j+1)(1)，d_RN(j+1)(2)，…，d_RN(j+1)(m)]，

D_RMj=diag[d_RMj(1)，dR_Mj(2)，…，d_RMj(m)]，D_RM(j+1)=diag[d_RM(j+1)(1)，d_RM(j+1)(2)，…，d_RM(j+1)(m)]；

diag为对角符号；

模型验证的最大信息量准则可写成：当

$\left\{\begin{array}{c}\frac{\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}[\ln d_{\mathrm{RNj}}\left(i\right)-\ln d_{\mathrm{RN}(j+1)}\left(i\right)]}{|m(\ln 2\mathrm{\&pi;}+1)+\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\ln d_{\mathrm{RNj}}\left(i\right)|}>\mathrm{\&delta;}\\ \left|\frac{m}{2}(N-M)\right[1+\ln \left(2\mathrm{\&pi;}\right)]+N\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\ln d_{\mathrm{RNj}}\left(i\right)-M\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\ln d_{\mathrm{RMj}}\left(i\right)|<\mathrm{\&epsiv;}\end{array}\right.---\left(8\right)$

成立时，θ_j+1、θ_jf_j[x(t)]、θ_jg_j[x(t)]选入模型正确，且 ${\mathrm{\&Omega;}}_{j+1}=\left[\begin{array}{c}{\mathrm{\&Omega;}}_j\\ {\mathrm{\&theta;}}_{j+1}\end{array}\right];$ 否则剔除θ_jf_j[x(t)]、θ_jg_j[x(t)]候选项，且Ω_j+1=Ω_j。

Claims (1)

1.一种基于可变数据长度最大信息量-可信度准则的飞行器建模方法，其特征在于包括以下步骤：

步骤一、飞行试验待确定的飞行器候选模型的状态方程为

$\stackrel{\&CenterDot;}{x}\left(t\right)=f\left\{f_0\right[x\left(t\right),{\mathrm{\&Omega;}}_0],f_1[x\left(t\right),{\mathrm{\&theta;}}_1],...,f_q[{\mathrm{\&theta;}}_q,x\left(t\right)],t\}---\left(1\right)$

观测方程为

$\left\{\begin{array}{c}y\left(t\right)=g[x\left(t\right),\mathrm{\&Omega;},t]=g\left\{g_0\right[x\left(t\right),{\mathrm{\&Omega;}}_0],g_1[x\left(t\right),{\mathrm{\&theta;}}_1],...,g_q[{\mathrm{\&theta;}}_q,x\left(t\right)],t\}---\left(2\right)\\ z\left(t_k\right)=y\left(t_k\right)+v\left(k\right)\end{array}\right.---\left(2\right)$

(1)、(2)式中，x(t)为n维状态向量；y(t)为m维观测向量；f{f₀[x(t),Ω₀],f₁[x(t),θ₁],...,f_q[θ_q,x(t)],t}、g{g₀[x(t),Ω₀],g₁[x(t),θ₁],...,g_q[θ_q,x(t)],t}为表达式已知的待确定模型结构函数，f₀[x(t),Ω₀]、g₀[x(t),Ω₀]为根据物理概念必须选入的模型，f_i[x(t),θ_i]、g_i[x(t),θ_i]为候选模型，其中，i＝1,2,…,q；z(t_k)为在t_k时刻对y(t_k)的测量值；Ω为未知维数的参数向量，Ω₀为已知维数的参数向量；v(k)为测量噪声，假定方差为R_k的零均值高斯白噪声；f_i[x(t),θ_i]、g_i[x(t),θ_i]是否在模型中出现及Ω₀、θ_i的取值需要辨识，q为已知的候选模型个数；

由于对飞行器的模型结构准确度要求较高，给出以下最大信息量-可信度准则：

$\left|\frac{\ln L_2[x\left(t\right),{\mathrm{\&Omega;}}_2,t]-\ln L_1[x\left(t\right),{\mathrm{\&Omega;}}_1,t}{\ln L_1[x\left(t\right),{\mathrm{\&Omega;}}_1,t]}\right|>\mathrm{\&delta;}---\left(3\right)$

式中，L₁[x(t),Ω₁,t]、L₂[x(t),Ω₂,t]为取不同的参数向量Ω₁和Ω₂的极大似然函数，δ为给定正实数、表示可信度参数，

$\left\{\begin{array}{c}\ln L_1[x\left(t\right),{\mathrm{\&Omega;}}_1,t]=-\frac{\mathrm{mN}}{2}\ln \left(\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_1^T\left(k\right)R_1^{-1}\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_1\left(k\right)\right)-\frac{\mathrm{mN}}{2}[\ln \left(\frac{2\mathrm{\&pi;}}{\mathrm{mN}}\right)+1]-\frac{1}{2}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}\ln \left|R_1\left(k\right)\right|\\ \ln L_2[x\left(t\right),{\mathrm{\&Omega;}}_2,t]=-\frac{\mathrm{mN}}{2}\ln \left(\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_2^T\left(k\right)R_2^{-1}\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_2\left(k\right)\right)-\frac{\mathrm{mN}}{2}[\ln \left(\frac{2\mathrm{\&pi;}}{\mathrm{mN}}\right)+1]-\frac{1}{2}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}\ln \left|R_2\left(k\right)\right|\end{array}\right.---\left(4\right)$ ， ${\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_1\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_1,t_k],{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_2\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_2,t_k],R_1=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_1\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_1^T\left(k\right),$ Ω₁和Ω₂为不同的参数向量，N为数据长度，ln为自然对数符号；

步骤二、根据假定f₀[x(t),Ω₀]、g₀[x(t),Ω₀]、Ω₀＝Ω₀已经通过优选算法选入模型，并由以下算法迭代计算得到：

令j＝0,1,2,…,q，假定f_j[x(t),θ_j]、g_j[x(t),θ_j]、Ω_j已经选入模型，按照以下方式选择其它候选模型：

求(4)式极大值,迭代计算：

${\mathrm{\&Delta;\&Omega;}}_j=A_j^{-1}{b}_j---\left(5\right)$

以及

$R_j=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j^T\left(k\right),{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_j,t_k]---\left(6\right)$

(5)、(6)式中： ${\mathrm{\&Delta;\&Omega;}}_j={\mathrm{\&Omega;}}_j-{\widehat{\mathrm{\&Omega;}}}_j,b_j=\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\left(\frac{\&PartialD;y}{{\&PartialD;\mathrm{\&Omega;}}_j^T}\right)}^T{R}_j^{-1}[z\left(t_k\right)-y\left(t_k\right)],$

$A_j=\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\left(\frac{\&PartialD;y}{{\&PartialD;\mathrm{\&Omega;}}_j^T}\right)}^T{R}_j^{-1}\frac{\&PartialD;y}{{\&PartialD;\mathrm{\&Omega;}}_j^T}=B_j^T{P}_j^{-1}{B}_j,B_j^T=[{\left(\frac{\&PartialD;y\left(t_1\right)}{{\&PartialD;\mathrm{\&Omega;}}_j^T}\right)}^T,{\left(\frac{\&PartialD;y\left(t_2\right)}{{\&PartialD;\mathrm{\&Omega;}}_j^T}\right)}^T,...,{\left(\frac{\&PartialD;y\left(t_N\right)}{{\&PartialD;\mathrm{\&Omega;}}_j^T}\right)}^T]$

$P_j^{-1}=\mathrm{diag}\left[\begin{array}{cccc}{R}_j^{-1},& R_j^{-1},& ...& R_j^{-1}\end{array}\right],$

当两个飞行试验数据长度分别为N、M、设 ${\mathrm{\&Omega;}}_{j+1}=\left[\begin{array}{c}{\mathrm{\&Omega;}}_j\\ {\mathrm{\&theta;}}_{j+1}\end{array}\right],$ θ_j+1的选入或剔除及验证条件为：当

$\frac{\ln \left|R_{\mathrm{Nj}}\left(k\right)\right|-\ln \left|R_{N(j+1)}\left(k\right)\right|}{|m(\ln 2\mathrm{\&pi;}+1)+\ln |R_{\mathrm{Nj}}\left(k\right)\left|\right|}>\mathrm{\&delta;}$ 及 $\left|\frac{m}{2}(N-M)\right[1+\ln \left(2\mathrm{\&pi;}\right)]+\ln \frac{{\left|R_N\left(k\right)\right|}^N}{{\left|R_M\left(k\right)\right|}^M}|<\mathrm{\&epsiv;}---\left(7\right)$ 时，θ_j+1、f_j+1[x(t),θ_j+1]、g_j+1[x(t),θ_j+1]选入模型，且 ${\mathrm{\&Omega;}}_{j+1}=\left[\begin{array}{c}{\mathrm{\&Omega;}}_j\\ {\mathrm{\&theta;}}_{j+1}\end{array}\right];$ 否则剔除f_j+1[x(t),θ_j+1]、g_j+1[x(t),θ_j+1]候选项，且Ω_j+1＝Ω_j；

(7)式中：ε为给定正数， $R_{\mathrm{Nj}}=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j^T\left(k\right),R_{N(j+1)}=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}^T\left(k\right),$

$R_{\mathrm{Mj}}=\frac{1}{M}\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j^T\left(k\right),R_{M(j+1)}=\frac{1}{M}\underset{k=1}{\overset{M}{\mathrm{\&Sigma;}}}{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}\left(k\right){\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}^T\left(k\right),$

${\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_j\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_j,t_k],{\stackrel{\&OverBar;}{\mathrm{\&eta;}}}_{j+1}\left(k\right)=z\left(t_k\right)-g[x\left(t_k\right),{\mathrm{\&Omega;}}_{j+1},t_k];$

步骤三、由于飞行器测量向量y的维数m较大，采用Gram-Schmidt正交化方法对R_Nj、R_Mj、R_N(j+1)和R_M(j+1)进行U-D分解，R_Nj、R_Mj、R_N(j+1)和R_M(j+1)的U-D分解分别为： $R_{\mathrm{Nj}}=U_{\mathrm{RNj}}{D}_{\mathrm{RNj}}{U}_{\mathrm{RNj}}^T,R_{\mathrm{Mj}}=U_{\mathrm{RMj}}{D}_{\mathrm{RMj}}{U}_{\mathrm{RMj}}^T,R_{N(j+1)}=U_{\mathrm{RN}(j+1)}{D}_{\mathrm{RN}(j+1)}{U}_{\mathrm{RN}(j+1)}^T,$ $R_{M(j+1)}=U_{\mathrm{RM}(j+1)}{D}_{\mathrm{RM}(j+1)}{U}_{\mathrm{RM}(j+1)}^T,$

式中，U_RNj、U_RMj、U_RN(j+1)、U_RM(j+1)为单位上三角阵；

D_RNj＝diag[d_RNj(1),d_RNj(2),…,d_RNj(m)]，D_RN(j+1)＝diag[d_RN(j+1)(1),d_RN(j+1)(2),…,d_RN(j+1)(m)]，D_RMj＝diag[d_RMj(1),d_RMj(2),…,d_RMj(m)]，D_RM(j+1)＝diag[d_RM(j+1)(1),d_RM(j+1)(2),…,d_RM(j+1)(m)]；diag为对角符号；

模型验证的最大信息量准则写成：当

$\begin{array}{c}\frac{\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}[\ln d_{\mathrm{RNj}}\left(i\right)-\ln d_{\mathrm{RN}(j+1)}\left(i\right)]}{|m(\ln 2\mathrm{\&pi;}+1)+\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\ln d_{\mathrm{RNj}}\left(i\right)}>\mathrm{\&delta;}\\ \left|\frac{m}{2}(N-M)\right[1+\ln \left(2\mathrm{\&pi;}\right)]+N\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\ln d_{\mathrm{RNj}}\left(i\right)-M\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\ln d_{\mathrm{RMj}}\left(i\right)|<\mathrm{\&epsiv;}\end{array}---\left(8\right)$

成立时，θ_j+1、f_j+1[x(t),θ_j+1]、g_j+1[x(t),θ_j+1]选入模型正确，且 ${\mathrm{\&Omega;}}_{j+1}=\left[\begin{array}{c}{\mathrm{\&Omega;}}_j\\ {\mathrm{\&theta;}}_{j+1}\end{array}\right];$ 否则剔除f_j+1[x(t),θ_j+1]、g_j+1[x(t),θ_j+1]候选项，且Ω_j+1＝Ω_j。