基于moesp和n4sid改进的子空间模型辨识方法

基于moesp和n4sid改进的子空间模型辨识方法

0模型的子空间辨识方法传递函数模型和状态空间模型是描述工农业对象的动态特性的两种常见模型。较之传递函数模型,状态空间模型具有如下优点:1)可以方便地描述多输入多输出系统;2)可以采用统一表达方式描述不同维数系统;3)作为一种时域分析方法,可以直接在时域中进行数值计算。其次,与外源自回归滑动平均模型(auto-regressivemovingaverageexogenous,ARMAX)相比,状态空间模型以其特有的模型结构,而无需求解复杂的丢番图方程,在各种预测控制算法中得到广泛应用。因此,如何根据现场操作信号直接求得热工对象的状态空间模型,一直是控制领域研究的重要课题。传统的系统辨识方法有预测误差法和辅助变量法。这2种方法都以优化为基础,需要进行迭代计算,辨识结果对优化算法的初始条件比较敏感,而且对多变量系统辨识的参数化也比较困难。20世纪90年代提出的子空间辨识方法可以克服这些缺点。近20年来,子空间辨识方法在理论和应用上都得到了巨大的发展。子空间辨识方法综合了系统理论,线性代数和统计学3方面的思想,大部分子空间辨识方法如:MOESP、N4SID和CAV方法,可以用1种统一的形式表示。子空间辨识方法的实质,就是基于这个统一的表达形式,对特定加权矩阵进行分解从而获取系统矩阵的方法。本文结合MOESP与N4SID算法在计算系统矩阵时各自的优点,提出了1种改进的子空间辨识方法(subspaceidentificationmethod,SIM),该方法先采用MOESP计算系统矩阵A和输出矩阵C;再根据矩阵A和C,采用N4SID得到输入矩阵B和前馈矩阵D。然后将改进算法应用于某天然气电站和Alstom气化炉模型的辨识中,并将本文提出的改进SIM算法与MOESP和N4SID算法的计算结果进行了比较。1模型描述和分析1.1模型的生成和编码一个线性时不变系统的离散状态空间模型可以描述为式中:uk∈Rm×1为第k时刻m维模型输入向量;yk∈Rl×1为第k时刻l维模型输出向量;m、l分别代表模型的输入和输出变量个数;xk∈Rn×1为第k时刻模型的状态向量,n为模型的阶数;wk∈Rn×1和vk∈Rl×1分别代表模型的不可以观测的过程噪声和测量噪声,是与模型输入uk无关的零均值白噪声序列;A∈Rn×n为模型的系统矩阵;B∈Rn×m为模型的输入矩阵;C∈Rl×n为模型的输出矩阵;D∈Rl×n为模型的前馈矩阵。1.2[bq]噪声范围为了保证辨识结果具有良好的统计特性,做如下假设:1)系统是渐进稳定的,A的特征值严格在单位圆内;2)(A,C)可观测,(A,[BQ1/2])可控;3)输入是确定性的拟平稳序列,且和过程噪声和测量噪声无关;4)过程噪声和测量噪声都是平稳的零均值白噪声序列,且式中:Q∈Rn×n、S∈Rn×l和R∈Rl×l为噪声wk和vk的协方差矩阵;δpq为克罗内克函数(即δpq=0,当p≠q;δpq=1,当p=q)。2矩阵和投影的定义2.1角toeplitch矩阵文献中所给出的输入输出矩阵等式在子空间辨识方法的发展中起到了非常重要的作用。该公式如下所示:其中增广观测矩阵确定性下三角Toeplitz矩阵Hid:随机性下三角Toeplitz矩阵His:输入和输出的Hankel矩阵:本文中假设j→∞。为了书写方便,作如下定义:式中p和f分别表示过去和将来。同样可以得到过程噪声wk和测量噪声vk的Hankel矩阵,分别记作M0|i-1、Mi|2i-1、N0|i-1和N0|2i-1,简记为Mp、Mf、Np和Nf。最后定义状态序列Xi为2.2矩阵补空间的正交投影首先定义p行j列矩阵A∈Rp×j和q行j列矩阵B∈Rq×j。定义1:矩阵A的行空间在矩阵B的行空间上的正交投影记作:定义2:矩阵A的行空间在矩阵B的行空间的正交补空间上的正交投影记作:式中:⋅†表示Moore-Penrose伪逆;B⊥表示B的行空间的正交补空间,即(B⊥)TB=0。3改进的子空间状态的空间识别方法3.1改进的moesp算法所有的子空间辨识算法都包含2个重要步骤,如图1所示。第一步,对上文中定义的Hankel矩阵进行加权投影,然后对投影进行奇异值分解(singularvaluedecomposition,SVD),得到系统的观测矩阵Γi或者系统状态序列Xi的估计Xi;第二步,由投影分解得到的观测矩阵Γi或者由系统状态序列Xi的估计Xi确定系统矩阵(A、B、C、D)和系统噪声的协方差矩阵Q、S和R。MOESP和N4SID算法的主要区别在于矩阵A和C的求解。N4SID方法首先需要通过投影获得系统状态估计,然后根据系统的状态估计构造一个最小二乘问题,求解得到矩阵A和C。而MOESP算法只需先确定系统的观测矩阵Γi,由观测矩阵直接得到矩阵A和C。MOESP和N4SID的另一个区别是矩阵B和D的计算。MOESP算法在求解B和D时涉及到1个计算量较大的最小二乘求解问题,需要构造一些大的矩阵,这会导致算法计算时间的显著增加。而N4SID算法可以在计算矩阵A和C的同时得到矩阵B和D。考虑2种算法的优缺点,本文提出了1种改进的辨识算法。该方法的改进之处在于:利用了MOESP算法可以根据系统的观测矩阵直接计算出系统矩阵A和输出矩阵C的优点,先计算出矩阵A和C,然后采用N4SID算法计算输入矩阵B和前馈矩阵D。该方法既能避免MOESP算法在计算矩阵B和D时需要构建大矩阵的缺点,又能避免N4SID算法在计算矩阵A和C时需要求解线性最小二乘的问题,降低了算法的复杂性。3.2加权矩阵oi的svd分解利用上文中的正交投影定义,将Yf的行空间正交投影到Uf的行空间的正交补空间Uf⊥上:由于噪声和输出是不相关的,因此有Mf/Uf⊥=Mf,Nf/Uf⊥=Nf,可以得到:选择恰当的加权矩阵W1和W2,并将W1和W2分别左乘和右乘于公式(14)的两端,得到:W1和W2需满足如下条件:1)W1Γi的秩等于Γi的秩;2)Xi/Uf⊥W2的秩等于Xi秩;如果上述条件都能得到满足,则定义矩阵Oi:对矩阵Oi进行SVD分解,得到:于是可以得到W1Γi=U1S11/2,Xi/Uf⊥W2=S11/2V2T。由文献可知,选取恰当的权矩阵W2,可得到系统状态序列的卡尔曼滤波估计:Xide=fXi/Uf⊥W2。3.3i和oi+1:2矩阵的形成有关MOESP和N4SID算法详细的推导和证明可参考文献[9-12],本文基于上述文献中的方法,给出了改进的SIM算法的计算步骤如下:1)计算矩阵Oi和Oi+1:2)对Oi进行SVD分解得到3)根据S1中非零奇异值的个数确定系统的阶数;4)计算Γi和Γi-1,Γi=W1†U1S11/2,代表由Γi的前(i-1)×l行数据组成的矩阵;5)由Γi直接计算求得A和C,,C=Γi(1:l,:)(Matlab表示形式),是由Γi的后(i-1)×l行数据组成的矩阵;6)确定系统的状态变量序列,7)求解下列线性等式得到B和D:8)由式(21)的残差确定系统噪声协方差矩阵Q,S和R。4模拟计算示例4.1模型辨识和算法比较为了验证本文提出的改进SIM的优越性,本文以某天然气电站和Alstom气化炉模型为例,分别采用MOESP、N4SID和改进SIM对其进行模型辨识,并采用文中4.2中的性能指标对3种算法的计算结果进行比较。由于本文所提出的改进SIM结合了MOESP和N4SID的优点,为了方便描述,在下文的某些图中用MOESN4表示本文提出的改进SIM。4.2浮点数运算次数1)模型总的相对误差(relativeerror,RE):式中:ye为对象辨识模型的输出估计;y为对象实际输出;norm(⋅)表示⋅的1-范数。2)模型各输出的均方根误差(rootmeansquareerror,RMSE):式中:为对象辨识模型的第i个样本的第j个输出估计;yi,j为对象的第i个样本的第j个实际输出。3)算法计算时所需的浮点数运算次数。4)算法计算时所需的CPU运算时间。4.3模型建立与辨识图2、3分别为法国某120MW天然气电站的输入输出采样数据曲线,采样时间为1228.8s,样本数为200个,该数据可由文献获得。本文以该电站的天然气流量变化量dW1、汽轮机阀门开度变化量dµ1、过热器喷水流量变化量dW2、烟气档板开度变化量dµ2和空气流量变化量dW3为输入,以主蒸汽压力变化量dp、主蒸汽温度变化量dT1和再热蒸汽温度变化量dT2为输出,建立1个5输入3输出的天然气电站机炉协调动态模型,并采用上述3种子空间辨识算法对该组数据进行辨识。通过计算可知,3种辨识算法所得模型均为5阶,计算结果如图3所示,3种模型均能够较好地逼近天然气电站输出。表1中给出了3种辨识模型的估计误差,与MOESP算法相比较,MOESN4和N4SID算法的辨识精度更高。3种辨识算法在计算天然气电站机炉协调动态模型的系统矩阵时所做的浮点数运算次数和CPU运算时间如图4所示。MOESN4算法所做的浮点数运算最少,计算速度也最快。但由于本算例中的样本数较少,模型的阶次也较低,MOESN4算法较其他2种算法在计算效率方面的优越性没有得到充分的体现。为此,本文在4.4节中采用2000组输入输出数据对Alstom气化炉进行辨识,较好地体现了MOESN4算法优越性。4.4模型建立与辨识Alstom气化炉模型是由Alstom能源技术中心开发的,以一台87MW空气气化炉为对象建立的全工况动态模型。通过与现场实测数据的校核,该模型能够有效地描述了气化炉的各种特性,模型具有很高的准确度。根据这个气化炉模型,Alstom公司发起了基准测试问题。Alstom气化炉有5个输入变量(煤量,石灰石量,空气量,水蒸气量和排渣量)和4个输出变量(合成气压力,合成气温度,床料量和合成气热值)。石灰石用于脱去煤中的硫,它的质量流量被设定为给煤量的1/10。于是,模型变成为了一个4输入、4输出的对象,输入和输出变量如表2所示。该模型包含了干燥、脱硫、气化、质量守恒和能量守恒等许多重要的过程,是一个高耦合和强非线性的多变量模型。为了充分激励出该模型的特性,选取伪随机噪声为激励信号,其幅值服从以100%稳定工况点输入量为中心的高斯(正态)分布,对4个输入同时进行激励。本文共采集了3000组输入输出信号,采样周期为1s,其中前2000组输入输出信号用于计算系统矩阵,如图5、6所示,剩余的1000组输入输出信号用于测试所得辨识模型的辨识精度,如图7、8所示。通过计算可以得到,Alstom气化炉模型的3种辨识模型均为14阶。本文在图8中给出3种辨识模型和Alstom气化炉模型的输出比较,3种模型均能够较好地辨识Alstom气化炉。3种模型的各种辨识误差如表3所示,较之其他2种算法,MOESN4算法同样具有较好的辨识精度。但是,通过比较3种算法的浮点数运算次数和CPU运算时间,如图9所示,可以发现,算例2中MOESN4的浮点数运算次数和CPU运算时间较其他2种算法有明显的降低。算例2充分说明了MOESN4能够在保证较好辨识精度的前提下,提高了原有算法的计算效率,特别是在大容量数据样本条件下,能够有效降低浮点数运算次数和CPU运算时间。4.5算例2:表1子空间辨识模型算法是在有无穷多输入输出样本的假设下发展起来,因此样本容量对子空间辨识算法的辨识精度有着重要的影响。在算例1中,由于文献中只给出200组输入输出样本,从表1中可以看出3种算法的相对误差和各输出的均方根误差均较大。相比之下,算例2中采用了2000组输入输出数据对模型进行辨识,3种算法的辨识精度较算例1均有较大提高。但并不是样本容量越大越好,样本容量过大容易导致较大的计算负荷,样本容量的大小要根据被辨识对象的特性(如模型阶次等)和辨识精度选取。5模型辨识精度对比本文提出了一种改进的子空间辨识方法(MOESN4算法),该方法利用了MOESP算法可以根据系统观测矩阵直接计算出系统矩阵A和输出矩阵C的优点,先计算出矩阵A和C,然后采用N4SID算法计算输入矩阵B和前馈矩阵D。并将该算法应用于某天然气电站和Alstom气化炉模型的