基于数据驱动的子空间模型辨识预测控制方法与流程

1.本发明涉及工业控制技术领域，尤其涉及一种基于数据驱动的子空间模型辨识预测控制方法，可广泛应用与过程控制，如石油、冶炼等过程工业，以及电力电子等行业。


背景技术：

2.模型预测控制是起源于20世纪70年代的一类应用于过程控制领域的算法。该方法在学术界和各种世纪工业应用中得到了广泛的研究，它具有三个基本特征：预测模型、滚动优化和反馈校正。然而，随着科学技术的进步，现代工业装备的规模越来越大，建模工作也越来越复杂，如何有效地建立控制系统的动态模型变得非常有必要。传统的工业预测控制都是采用输入输出模型，包括参数模型和非参数模型，但是为了进一步提高控制性能和控制精度，学术界和工业界普遍认为应该采用状态空间模型，这样近些年所发展起来的现代滤波理论和控制器设计方法就可以发挥作用，而且实际工业过程中往往模型结构未知，只有可测输入输出数据已知，很难得到较复杂的状态空间模型。


技术实现要素：

3.为此，本发明在于克服现有技术中很难得到较复杂的状态空间模型的技术问题，提供一种基于数据驱动的子空间模型辨识预测控制方法，在未知系统模型的情况下，通过子空间模型辨识方法求得状态空间系统矩阵的系数，在此基础上，对被控对象进行控制并且更新模型，最终使得被控对象模型辨识准确且稳定运行，达到良好的运行状态。
4.为解决上述技术问题，本发明提供了一种基于数据驱动的子空间模型辨识预测控制方法，包括以下步骤，
5.步骤一：以离线状态下被控对象的输入输出数据为对象，通过子空间模型辨识方法对被控对象进行离线辨识，获取被控对象的状态空间方程的系统参数，根据该系统参数确定出所述被控对象的离线状态空间模型；
6.步骤二：将获得的所述离线状态空间模型作为预测模型，并基于模型预测控制策略控制系统，在控制过程中改变所述被控对象的特征值，并采集所述被控对象特征值改变后的在线输入输出数据；
7.步骤三：以步骤二中采集到的所述在线输入输出数据为对象，通过子空间模型辨识方法对被控对象进行在线辨识，获取所述被控对象的在线状态空间模型；
8.步骤四：不断更新所述在线状态空间模型，直至根据所述在线状态空间模型控制的系统输出值与所述被控对象的真实输出值的差值小于预设值时，确定出该时刻下的所述在线状态空间模型为所述被控对象的最终预测控制模型。
9.在本发明的一个实施例中，步骤一中获取所述被控对象的状态空间方程的系统参数包括，
10.s11、实验环境下给定所述被控对象在离线状态下的输入输出数据，根据该输入输
出数据建立hankel矩阵k为hankel矩阵的行数；0k-1为hankel矩阵第一列的第一个到最后一个参数；
11.s12、基于多变量输出误差状态空间算法moesp处理hankel矩阵以计算出所述被控对象的状态空间方程的系统参数，该系统参数包括状态空间方程的常数矩阵a、b、c、d；
12.所述状态空间方程为：
13.其中，为状态矢量，为控制输入，为输出矢量，为常数矩阵；n为状态矢量的维数；m为控制输入的维数；p为输出矢量的维数。
14.在本发明的一个实施例中，步骤s12中，建立hankel矩阵后,根据所述状态空间方程迭代出hankel矩阵的以下输入输出方程：
15.y
0|k-1
＝οkx0 ψku
0|k-1
ꢀꢀꢀꢀ
(1)
16.其中，为扩展可观测矩阵；
17.为toeplitz矩阵；
18.为初始状态矩阵；
19.n为hankel矩阵的列数。
20.在本发明的一个实施例中，步骤s12中，计算所述系统参数包括，
21.对hankel矩阵进行lq分解，得到右上方为零的下三角矩阵，
[0022][0023]
对l
22
进行奇异值分解得到系统状态矢量的维数n；
[0024]
根据奇异值分解得到的结果计算出矩阵a、b、c和d；
[0025]
其中，且l
11
,l
22
为下三角矩阵；
[0026]
互为正交矩阵。
[0027]
在本发明的一个实施例中，l
22
的奇异值分解为：
[0028][0029]
σ1为一个对角矩阵；v1和v2为奇异值分解后得到的矩阵；
[0030]
根据以下等式确定出矩阵c：
[0031]
c＝οk(1:p,1:n)
ꢀꢀꢀꢀ
(4)；
[0032]
通过求解以下线性方程得到矩阵a：
[0033]
οk(1:p(k-1),1:n)a＝οk(p 1:kp,1:n)
ꢀꢀꢀꢀ
(5)。
[0034]
在本发明的一个实施例中，根据以下等式确定矩阵b和d的线性方程，
[0035][0036]
通过最小二乘法求解以上式(6)，求解出矩阵b和d。
[0037]
在本发明的一个实施例中，步骤二中，模型预测控制策略包括，
[0038]
s21、构建以下优化模型函数：
[0039][0040]
其中，w＞0和r＞0为权重矩阵；
·
(s i|s)表示系统状态的矢量x(t)或者控制输入u(t)在s时刻时t＝s i的值；i表示步长；
[0041]
s22、在线求解所述优化模型函数的最优值；
[0042]
s23、将步骤s22中求解得到的最优值作为模型预测控制器的输出，该输出用于控制被控对象。
[0043]
在本发明的一个实施例中，基于最小-最大方法求解所述优化函数模型的最优值，其包括，
[0044]
求解所述优化模型函数的最大界v(x(s|s))，其满足以下关系：
[0045][0046]
求解所述最大界v(x(s|s))的最小值，该最小值用作所述优化模型函数的最优值。
[0047]
在本发明的一个实施例中，求解所述最大界v(x(s|s))的最小值时通过舒尔补定理将不等式的问题优化转化为矩阵不等式的优化。
[0048]
在本发明的一个实施例中，求解所述优化模型函数的最优值时考虑函数的以下约束：
[0049][0050]
其中，表示约束条件中输入的最大值。
[0051]
本发明的上述技术方案相比现有技术具有以下优点：
[0052]
本发明所述的基于数据驱动的子空间模型辨识预测控制方法，采用子空间辨识和
模型预测控制相结合，以及离线子空间辨识与在线子空间辨识相结合，通过子空间辨识方法辨识离线状态下的输入输出数据即可获得离线状态空间模型，随后使用该离线状态空间模型作为预测模型，采集在线输入输出数据，并通过子空间辨识方法辨识该在线输入输出数据获得在线状态空间模型，通过不断更新在线状态空间模型，直至根据在线状态空间模型控制的系统输出值与真实输出值的差值小于预设值时，确定出该时刻下的在线状态空间模型为预测控制模型，最后使用该预测控制模型控制系统。该方法既能辨识更新后的状态空间模型，又能实现闭环控制系统的稳定，最终使得被控对象模型辨识准确且稳定运行，达到良好的运行状态，简化了建模过程，计算量适中适用于获取复杂模型的控制系统。
附图说明
[0053]
为了使本发明的内容更容易被清楚的理解，下面根据本发明的具体实施例并结合附图，对本发明作进一步详细的说明，其中
[0054]
图1为本发明基于数据驱动的子空间模型辨识预测控制方法的流程图；
[0055]
图2为本发明实施例中具体应用场景下辨识模型和实际模型的频率响应幅值误差图；
[0056]
图3为本发明实施例中具体应用场景下辨识模型控制输入示意图；
[0057]
图4为本发明实施例中具体应用场景下辨识模型下电感电流值示意图；
[0058]
图5为本发明实施例中具体应用场景下辨识模型下输出电压值示意图。
具体实施方式
[0059]
下面结合附图和具体实施例对本发明作进一步说明，以使本领域的技术人员可以更好地理解本发明并能予以实施，但所举实施例不作为对本发明的限定。
[0060]
参照图1所示，本发明公开一种基于数据驱动的子空间模型辨识预测控制方法，包括以下步骤，步骤一：以离线状态下被控对象的输入输出数据为对象，通过子空间模型辨识方法对被控对象进行离线辨识，获取被控对象的状态空间方程的系统参数，根据该系统参数确定出所述被控对象的离线状态空间模型；步骤二：将获得的所述离线状态空间模型作为预测模型，并基于模型预测控制策略控制系统，在控制过程中改变所述被控对象的特征值，并采集所述被控对象特征值改变后的在线输入输出数据；步骤三：以步骤二中采集到的所述在线输入输出数据为对象，通过子空间模型辨识方法对被控对象进行在线辨识，获取所述被控对象的在线状态空间模型；步骤四：不断更新所述在线状态空间模型，直至根据所述在线状态空间模型控制的系统输出值与所述被控对象的真实输出值的差值小于预设值时，确定出该时刻下的所述在线状态空间模型为所述被控对象的最终预测控制模型。
[0061]
本发明所述的基于数据驱动的子空间模型辨识预测控制方法，采用子空间辨识和模型预测控制相结合，以及离线子空间辨识与在线子空间辨识相结合，通过子空间辨识方法辨识离线状态下的输入输出数据即可获得离线状态空间模型，随后使用该离线状态空间模型作为预测模型，采集在线输入输出数据，并通过子空间辨识方法辨识该在线输入输出数据获得在线状态空间模型，通过不断更新在线状态空间模型，直至根据在线状态空间模型控制的系统输出值与真实输出值的差值小于预设值时，确定出该时刻下的在线状态空间模型为预测控制模型，最后使用该预测控制模型控制系统。该方法既能辨识更新后的状态
空间模型，又能实现闭环控制系统的稳定，最终使得被控对象模型辨识准确且稳定运行，达到良好的运行状态，简化了建模过程，计算量适中适用于获取复杂模型的控制系统。
[0062]
具体的，步骤一中获取所述被控对象的状态空间方程的系统参数包括，
[0063]
s11、实验环境下给定所述被控对象在离线状态下的输入输出数据，根据该输入输出数据建立hankel矩阵k为hankel矩阵的行数；0k-1为hankel矩阵第一列的第一个到最后一个参数；
[0064]
s12、基于多变量输出误差状态空间算法moesp处理hankel矩阵以计算出所述被控对象的状态空间方程的系统参数，该系统参数包括状态空间方程的常数矩阵a、b、c、d；
[0065]
所述状态空间方程为：
[0066]
其中，为状态矢量，为控制输入，为输出矢量，为常数矩阵；n为状态矢量的维数；m为控制输入的维数；p为输出矢量的维数。
[0067]
例如，实验环境下给定以下输入输出数据(u(0)u(1)
…
u(k n-2))，(y(0)y(1)
…
y(k n-2))；其中，为控制输入，为输出矢量；m为控制输入的维数；p为输出矢量的维数；n为hankel矩阵的列数。根据以上输入输出数据建立如下所示的hankel矩阵
[0068][0069][0070]
建立hankel矩阵后,根据所述状态空间方程迭代出：
[0071][0072]
为了简化以上公式1，定义：
[0073][0074]
扩展可观测矩阵οk和toeplitz矩阵：
[0075][0076]
由此公式1可表达为如下：
[0077]
yk(t)＝οkx(t) ψkuk(t),t＝0,1,
…ꢀꢀꢀꢀ
(公式2)；
[0078]
hankel矩阵u
0|k-1
和y
0|k-1
可以分别表示为由uk(t)和yk(t)组成的形式：u
0|k-1
＝[uk(0) uk(1)
ꢀ…ꢀ
uk(n-1)]，y
0|k-1
＝[yk(0) yk(1)
ꢀ…ꢀ
yk(n-1)]，将公式2进一步表达为如下形式：
[0079]y0|k-1
＝οkx0 ψku
0|k-1
ꢀꢀꢀꢀ
(1)
[0080]
为初始状态矩阵。式(1)是具有初始状态x0的矩阵输入输出方程，hankel矩阵u
0|k-1
和y
0|k-1
表示过去的输入输出。
[0081]
被控对象状态空间方程的系统参数的求取，采用子空间辨识方法的多变量输出误差状态空间算法，具体实施步骤和结果如下：
[0082]
首先对hankel矩阵u
0|k-1
和y
0|k-1
进行lq分解得到右上方为零的下三角矩阵，
[0083][0084]
其中，且l
11
,l
22
为下三角矩阵；互为正交矩阵。
[0085]
由式(2)可以分别得到hankel矩阵u
0|k-1
和y
0|k-1
的表达式：
[0086][0087][0088]
由于l
11
是非奇异矩阵，公式3a变换得到将代入公式3b可以写成：
[0089][0090]
通过比较式(1)和公式3b，都是有关y
0|k-1
的表达式，由此可知：
[0091][0092]
需要注意的是，虽然公式4右边是正交和，但其左边是直和，所以其中的两个量不一定是正交的，这意味着
[0093]
将公式4右乘以q2得到οkx0q2＝l
22
，其中已知rank(οk)＝n，等于rank(l
22
)，接下来通过奇异值分解(svd)对进行分解得到系统的状态矢量的维数n(即系统的阶数)。
[0094]
的奇异值分解为：
[0095][0096]
σ1为一个对角矩阵；已知οkx0q2＝l
22
，然后我们有οkx0q2＝u1σ
1v1t
。定义οk为并且n＝dimσ1，随后就容易得到矩阵c：
[0097]
c＝οk(1:p,1:n)
ꢀꢀꢀꢀ
(4)；
[0098]
通过求解如下的线性方程可以得到矩阵a
[0099]
οk(1:p(k-1),1:n)a＝οk(p 1:kp,1:n)
ꢀꢀꢀꢀ
(5)。
[0100]
接下来考虑矩阵b和d的估计，由于公式4左乘得到对前面公式中右乘q1：
[0101][0102]
式(6)给出了关于矩阵b和d的线性方程，然后采用最小二乘法求解出矩阵b和d。
[0103]
综上，通过子空间辨识方法的多变量输出误差状态空间算法计算出状态空间方程的系统参数，即计算出状态空间方程的矩阵a、b、c和d。
[0104]
作为本发明实施例的进一步改进，基于以上方法在模型结构未知的情况下确定出预测模型，随后基于该预测模型采用模型预测控制策略进行系统控制，其包括，
[0105]
s21、构建以下优化模型函数：
[0106][0107]
其中，w＞0和r＞0为权重矩阵；
·
(s i|s)表示系统状态的矢量x(t)或者控制输入u(t)在s时刻时t＝s i的值；i表示步；
[0108]
s22、在线求解所述优化模型函数的最优值；
[0109]
s23、将步骤s22中求解得到的最优值作为模型预测控制器的输出，该输出用于控制被控对象。
[0110]
具体的，基于最小-最大方法求解所述优化函数模型的最优值，其包括，
[0111]
首先，求解所述优化模型函数的最大界v(x(s|s))，其满足以下关系：
[0112][0113]
再，求解所述最大界v(x(s|s))的最小值，该最小值用作所述优化模型函数的最优值。
[0114]
其中，求解所述优化模型函数的最优值时考虑函数的以下约束：
[0115][0116]
其中，表示约束条件中输入的最大值。
[0117]
即在线求解的优化问题为：
[0118][0119]
为了求解所述最大界v(x(s|s))的最小值，采用状态反馈控制律，定义u(s i|s)＝fx(s i|s),i≥0，其中f为控制增益，定义一个函数v(x)＝x
t
px,p＞0，引入一个标量变量γ，并且考虑最小化γ，该γ满足v(x(s|s))≤γ，则问题转化为求解γ的最小值，通过舒尔补定理将不等式的优化问题转化为矩阵不等式的优化。
[0120]
定义矩阵q＝γp-1
，由舒尔补定理可知v(x(s|s))≤γ等价于如下的矩阵不等式：
[0121][0122]
将定义的反馈控制律等式代入式(8)，进一步得到：
[0123]
x(s i|s)
t
[(a bf)
t
p(a bf)-p f
t
rf w]x(s i|s)≤0
ꢀꢀꢀꢀ
(12)；
[0124]
由以下式(13)可以保证式(12)满足i≥0所有要求：
[0125]
(a bf)
t
p(a bf)-p f
t
rf w≤0
ꢀꢀꢀꢀ
(13)；
[0126]
定义f＝yq-1
，式(13)代入p＝γq-1
和f＝yq-1
，在式(13)两边乘以q，通过舒尔补定理，得到以下线性矩阵不等式：
[0127][0128]
式(14)中符号"*"表示对称位置的分块。
[0129]
基于此，在线求解的优化问题为：
[0130][0131]
转换成求解以下优化问题：
[0132][0133]
下面介绍一个具体的应用场景：
[0134]
考虑升压dc-dc变换器电路模型，选择状态空间向量x
t
＝(i
l v)，状态空间动态模型：
[0135][0136]
其中，rs＝(1-d)2r r
l
d(1-d)(rc||r)；vg为直流输入电压，为电压变化值；x为状态向量的直流值，叠加的交流摄动值。同样的，假设占空比随着周期的变化，也就是d为稳态直流占空比，为叠加的变化值。
[0137]
如图2所示，基于本发明的子空间辨识方法更新后的模型与系统实际模型一致，波德幅值图差异小，由此可知，本发明的子空间辨识方法的多变量输出误差状态空间算法moesp可以更准确地识别被控对象模型。由图3可以看出，当使用升压变压器模型时，仿真在247ms，旧模型的负载电阻由5ω变为50ω。模型更改后，收集更新后的输入输出数据，新数据集的长度为30组，在子空间辨识方法mpc下，更新后的模型可以继续控制和稳定运行。
[0138]
本领域内的技术人员应明白，本技术的实施例可提供为方法、系统、或计算机程序产品。因此，本技术可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本技术可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、cd-rom、光学存储器等)上实施的计算机程序产品的形式。
[0139]
本技术是参照根据本技术实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
[0140]
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
[0141]
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
[0142]
显然，上述实施例仅仅是为清楚地说明所作的举例，并非对实施方式的限定。对于所属领域的普通技术人员来说，在上述说明的基础上还可以做出其它不同形式变化或变动。这里无需也无法对所有的实施方式予以穷举。而由此所引伸出的显而易见的变化或变动仍处于本发明创造的保护范围之中。