动静态数据混合驱动的Hammerstein非线性工业系统简约灰箱子空间辨识方法

动静态数据混合驱动的Hammerstein非线性工业系统简约灰箱子空间辨识方法

　　技术领域

　　本发明属于控制理论与控制工程非线性系统辨识领域，涉及动态/静态数据混合驱动的Hammerstein非线性工业系统简约灰箱子空间辨识方法。

　　背景技术

　　Hammerstein非线性系统是一种模块化的非线性系统，由静态非线性模块和线性动态模块组成，能够有效的反应大多数实际输入非线性工业过程的系统特性，已经被大量应用于对实际工业过程进行描述，例如高炉炼铁系统、电池管理系统、无线电力传输系统等。在过去的几十年里，业内开展了许多关于Hammerstein非线性系统辨识的研究工作，多种成熟的辨识方法被相继提出，如最大似然估计法、最小二乘法、随机逼近法和子空间辨识法(SIMs)。

　　由于状态空间能够有效描述系统内部动态特性，方便控制器的设计。具有状态空间线性动态模块的Hammerstein非线性系统辨识成为控制届的研究热点，相应的子空间辨识方法得到了控制界的广泛研究。在过去的二十年里，许多Hammerstein非线性系统子空间辨识方法被相继提出，现有方法包括两大类，即非迭代方法和迭代方法。(1)非迭代方法包括基于多变量输出误差状态空间的子空间辨识方法(MOESP)、基于数值子空间状态空间系统的子空间方法(N4SID)、以及基于奇偶空间的子空间方法。(2)迭代方法包括基于固定点迭代的子空间辨识方法、基于高斯-牛顿迭代的子空间辨识方法。以上方法在机械系统、无线电力传输系统、流程工业系统得到很好的应用。但是现有方法全部是动态数据驱动的，也即只采用动态数据来进行Hammerstein非线性系统模型辨识，属于黑箱子空间辨识方法。由于动态数据仅含有高频系统信息，不能提供系统的低频信息，导致现有黑箱系统辨识方法难以建立高精度的系统模型。

　　发明内容

　　有鉴于此，本发明的目的在于提供一种动态/静态数据混合驱动的新型灰箱简约子空间辨识方法对Hammerstein非线性系统进行精确辨识，以克服现有黑箱子空间辨识方法不能提供系统低频信息，导致模型精度低的问题。

　　为达到上述目的，本发明提供如下技术方案：

　　一种动静态数据混合驱动的Hammerstein非线性工业系统简约灰箱子空间辨识方法，具体包括以下步骤：

　　S1：采集系统动态数据和静态数据；

　　S2：选择和处理Hammerstein非线性系统模型，得到预测模型；

　　S3：构建Hammerstein非线性系统动态简约模型；

　　S4：构建Hammerstein非线性系统静态简约模型；

　　S5：利用递阶拉格朗日最优加权方法求解系统参数融合辨识。

　　进一步，所述步骤S1具体包括：对实际工业过程进行数据采集，获得系统动态数据和静态数据其中，u(t)，y(t)分别为动态输入输出数据，N为动态输入输出数据的数量，分别为静态输入输出数据，M为静态输入输出数据的数量。

　　进一步，所述步骤S2具体包括以下步骤：

　　S21：对非线性特性的实际工业过程，采用具有状态空间模块的Hammerstein非线性模型进行描述，具体形式如下：

　　

　　f(t)＝ωz(t)(2)

　　ω＝[ω1，ω2，…，ωr](3)

　　z(t)＝[f1(u(t))，f2(u(t))，…，fr(u(t))]T(4)

　　其中，u(t)∈R，y(t)∈R，e(t)∈R分别表示状态向量、采集输入数据、采集输出数据、新息向量。(A，B，C，K)表示系统矩阵。nx表示模型阶次，新息为零均值高斯噪声，非线性输入函数f(t)∈R为已知基函数fi(u(t))∈R的线性组合，ωi∈R表示未知系数。

　　S22：对上述选择模型进行处理可以得到预测模型，以方便后续模型参数辨识，具体模型形式如下：

　　

　　其中：

　　同时对系统的可辨识性做以下一般假设：

　　A1线性系统是可观测与可实现的。

　　A2的特征值位于单位圆内。

　　A3ω的第一个非零元素为正，||ω||2＝1(||.||2定义欧氏范数)。

　　本发明从动态和稳态数据中估计未知参数(A，B，C，ω)，可以提高仅由动态数据估计模型参数的精度。

　　进一步，步骤S3中，构建Hammerstein非线性系统动态简约模型，具体包括以下步骤：

　　S31：对新息模型进行p次迭代处理，得到新的状态方程：

　　

　　其中：

　　

　　

　　zp(t)＝[zT(t-p)，zT(t-p+1)，…，zT(t-1)]T(9)

　　yp(t)＝[yT(t-p)，yT(t-p+1)，…，yT(t-1)]T(10)

　　需要说明的是p可以根据实际系统特性任意选择，一般在10-60范围内进行选择。

　　S32：将新的状态方程带入输出，根据的特性对系统的输出y(t)进行重写得到新的系统输出函数：

　　(t)＝P0(t)θ0+P1(t)θ1+e(t)(11)

　　其中：

　　P0(t)＝[zT(t-1)，zT(t-2)，…，zT(t-p)] (12)

　　P1(t)＝[yT(t-1)，yT(t-2)，…，yT(t-p)](13)

　　

　　

　　由于θ0中未知参数的数目是pr，而系统参数和非线性参数向量ω的数目是p+r。对于包括高度非线性模型在内的大多数系统而言，过参数化模型(11)与系统参数和ω相比较，它包含了冗余参数，这会导致不必要的高方差与模型失真，特别是当系统使用小而嘈杂的数据集时。

　　S33：对θ0进行分解处理，重写过参数化模型(11)，得到两个动态简约模型：

　　y(t)＝P2(t)θ2+P1(t)θ1+e(t)(16)

　　＝P3(t)ωT+P1(t)θ1+e(t) (17)

　　其中：

　　P2(t)＝[ωz(t-1)，ωz(t-2)，…，ωz(t-p)](18)

　　

　　

　　且θ1和θ2包括所有扩展的马尔科夫参数(i＝1，…，p)，如下所示：

　　

　　

　　相应的系统马尔科夫参数定义为：

　　hi＝CAi-1K(23)

　　gi＝CAi-1B(24)

　　需要说明的是这样的分解处理能够避免对冗余参数的估计，从而提高估计精度。

　　S34：对动态数据和动态简约模型进行矩阵化处理：

　　Y＝[y(1)，y(2)，…，y(N)]T(25)

　　E＝[e(1)，e(2)，…，e(N)] (26)

　　

　　

　　

　　

　　在矩阵化的基础上，对(11)、(16)和(17)进行迭代运算，得到新的矩阵化动态简约模型如下：

　　Y＝φ0θ0+φ1θ1+E(31)

　　Y＝φ2θ2+φ1θ1+E(32)

　　Y＝φ3ωT+φ1θ1+E(33)

　　由此，采用动态简约模型(32)估计线性扩展马尔可夫参数向量θ2；用动态简约模型(33)估计非线性参数ω；用辅助模型(31)估计扩展(32)和(33)中的冗余参数，即马尔可夫参数θ1。

　　进一步，步骤S4中，构建Hammerstein非线性系统静态简约模型，具体包括以下步骤：

　　S41：对应所得静态输入输出代入式(11)，整理后得到稳态输出的另一种表示形式：

　　

　　其中：

　　

　　

　　

　　S42：对静态数据进行矩阵化处理：

　　

　　

　　得到稳态参数G和Q的一致最小二乘解：

　　

　　S43：由于参数G的函数表达式(36)和Q的函数表达式(37)包含过参数化的模型参数，将过参数化模型(37)重写为两个静态简约模型：

　　Q＝K1ωT(41)

　　

　　其中：

　　

　　K2＝ωi[1，…，1](44)

　　其中，Qi是Q的第i个参数，可以从Q中任意选择；

　　采用静态简约模型(41)估计非线性参数向量ω；采用静态简约模型(42)估计线性扩展马尔可夫参数向量θ2。

　　进一步，步骤S5中，利用递阶拉格朗日最优加权方法求解系统参数融合辨识，具体包括以下步骤：

　　S51：采用投影方法估计辅助模型(31)中的θ1；

　　将动态简约模型(31)中的输出Y正交投影到φ0的正交补上，进一步进行最小二乘估计，得到θ1的估计值如下：

　　

　　其中：

　　

　　

　　

　　其中，IN表示维度为N的单位矩阵；

　　S52：构建基于加权矩阵简约模型的多正则化框架：

　　

　　

　　其中，τ1和τ2是与稳态数据的等式约束相关的未知拉格朗日乘子；

　　S53：求导计算：取对和τ1的一阶导数，取对和τ2的一阶导数，得到：

　　

　　

　　

　　

　　其中，和的定义相似；

　　S54：采用递阶拉格朗日迭代方法估计ω和θ2；

　　S55：估计线性系统参数。

　　进一步，步骤S54中，采用递阶拉格朗日迭代方法估计ω和θ2，具体包括：在第(k+1)次迭代时，用第k次迭代时的估计值替换第(k+1)次迭代时φ2和K2中的未知参数向量ω，用于估计θ2，用代替φ3和K1中的未知参数向量θ2来估计ω；对执行规范化操作：

　　

　　初始值ω(0)和θ2(0)在0-1之间任意生成。

　　进一步，步骤S55中，估计线性系统参数，具体包括：

　　S551：从估计马尔科夫参数提取系统马尔科夫参数：

　　

　　

　　S552：执行以下奇异值分解(SVD)：

　　

　　其中，和分别是S的最大nx特征向量、U的左特征向量和VT的右特征向量；由于系统是可观察和可实现的，显然SS＝0。

　　S553：从和中提取可观测矩阵和可控性矩阵

　　

　　

　　S554：从和中直接检索系统矩阵：

　　

　　

　　

　　本发明的有益效果在于：本发明采用基于分解的简约模型，避免对额外中间参数的估计，降低估计模型的方差、提高模型精度；另外本发明采用最优加权的多重正则化框架对静态数据和动态数据进行融合，实现对模型参数的辨识求解。相比现有黑箱子空间辨识，本发明大大提高了Hammerstein非线性系统的模型精度、降低了估计参数方差。

　　本发明的其他优点、目标和特征在某种程度上将在随后的说明书中进行阐述，并且在某种程度上，基于对下文的考察研究对本领域技术人员而言将是显而易见的，或者可以从本发明的实践中得到教导。本发明的目标和其他优点可以通过下面的说明书来实现和获得。

　　附图说明

　　为了使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明作优选的详细描述，其中：

　　图1为本发明方法的流程框图；

　　图2为本实施例采用的WPT电路拓扑图；

　　图3为WTP系统动态输入输出数据和稳态输入输出数据仿真图；

　　图4真实WTP系统训练结果于测试结果对比图。

　　具体实施方式

　　以下通过特定的具体实例说明本发明的实施方式，本领域技术人员可由本说明书所揭露的内容轻易地了解本发明的其他优点与功效。本发明还可以通过另外不同的具体实施方式加以实施或应用，本说明书中的各项细节也可以基于不同观点与应用，在没有背离本发明的精神下进行各种修饰或改变。需要说明的是，以下实施例中所提供的图示仅以示意方式说明本发明的基本构想，在不冲突的情况下，以下实施例及实施例中的特征可以相互组合。

　　请参阅图1～图4，将本发明方法应用到无线电能传输(wireless powertransfer，WPT)系统中研究，以说明本发明的有效性和优越性。本实施例研究的WPT系统由相控全桥逆变器和接收机组成，电路拓扑图和实物结构图分别如图2，其中具体参数如表1所示。

　　表1 WPT实验装置主要参数

　　具体实施过程如下：

　　S1：对无线电能传输系统进行数据采集。

　　输入为WPT系统逆变器相移，输出为负载电压。采样周期选择为T＝0.1ms。为了保证系统持续激励，在[0，0.4]范围内任意生成输入变量u(t)，以获取动态输入输出数据在[0，0.4]范围内以0.01为间隔并且持续3.8s，作为变量以获取静态输入输出数据为了验证本发明方法的有效性，数据被分为训练数据和测试数据2部分。训练数据包括5000个动态数据和11780个稳态数据组成的训练数据，测试数据由1000个动态数据和3800个采样稳态数据组成。训练数据和测试数据如图3所示。

　　S2：对具有非线性特性的无线电能传输系统，采用具有状态空间模块的Hammerstein非线性模型进行描述，并且进一步得到预测模型形式，以方便后续辨识。需要说明的是Hammerstein模型的基函数可以根据系统特性选取不同的非线性函数，如径向基函数、铰链函数、多项式函数等。由于未知系数提供了自由度，能够保证非线性函数能够精确逼近系统特性，因此基函数可以简单地选取为多项式函数。另外r可在[2，11]范围内选择。可以通过对比实验获取最优值，也即在不同r取值情况下，以系统静态输出和模型静态输出两者的误差作为选择依据，误差最小对应的r即为最优值。本例中非线性基函数个数最优值为7。

　　S3：对新息模型进行p次迭代处理，能够得到Hammerstein非线性系统动态简约模型，该模型能够有效处理动态数据。其中p可以在应的[10，60]的范围任意选择。可以通过对比实验获取最优值，也即在不同p取值情况下，以系统静态输出和模型静态输出两者的误差作为选择依据，误差最小对应的p即为最优值。本实施例中非线性基函数个数最优值为27。

　　S4：将对应所得静态输入输出数据进行处理，整理后可以得到矩阵化静态数据，进一步可以得到稳态参数G和Q的一致最小二乘解。从而可以得到Hammerstein非线性系统静态简约模型。

　　S5：首先，通过投影估计得到θ1，接着将动态简约模型中的输出Y正交投影到φ0的正交补，进一步进行最小二乘估计得到θ1的估计值。再采用具有随机初值的递阶拉格朗日迭代方法估计ω和θ2。其中ω的估计值如表2所示。

　　表2 WPT系统参数

　　

　　然后，从估计所得θ1和θ2提取系统马尔可夫参数，并且构造相应矩阵，执行奇异值分解(SVD)。最后从而从和中提取得到系统矩阵参数A、B、C。所得结果如图4所示。

　　基于所得模型得到WPT系统的模型训练输出数据和模型测试输出数据如图4所示。模型训练输出数据和模型测试输出数据的拟合度如表3所示。拟合度定义为：

　　

　　其中y(T)和分别代表采样实际输出和模型输出，是y(T)的平均值。可以看出，使用本发表所提出的方法对WPT系统进行参数估计取得了较好的效果。所得模型模型精度和理想值基本重合、拟合度较高，建立了高精度的WPT模型。

　　表3 WPT系统拟合度

　　使用本发明提出的动态/静态数据混合驱动的新型灰箱简约子空间辨识方法对WPT系统进行精确辨识，采用基于分解的简约模型，避免对额外中间参数的估计，降低估计模型的方差、提高模型精度；另外该方法采用最优加权的多重正则化框架对静态数据和动态数据进行融合，实现对模型参数的辨识求解，克服现有黑箱子空间辨识方法不能提供系统低频信息，导致模型精度低的问题。得到精确的WPT模型。可以在其它复杂系统进行推广应用。

　　最后说明的是，以上实施例仅用以说明本发明的技术方案而非限制，尽管参照较佳实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，可以对本发明的技术方案进行修改或者等同替换，而不脱离本技术方案的宗旨和范围，其均应涵盖在本发明的权利要求范围当中。