模型预测控制快速求解算法

模型预测控制快速求解算法

模型预测控制（Model Predictive Control，MPC）是一种基于在线计算的控制优化算法，能够统一处理带约束的多参数优化控制问题。当被控对象结构和环境相对复杂时，模型预测控制需选择较大的预测时域和控制时域，因此大大增加了在线求解的计算时间，同时降低了控制效果。从现有的算法来看，模型预测控制通常只适用于采样时间较大、动态过程变化较慢的系统中。因此，研究快速模型预测控制算法具有一定的理论意义和应用价值。

虽然MPC方法为适应当今复杂的工业环境已经发展出各种智能预测控制方法，在工业领域中也得到了一定应用，但是算法的理论分析和实际应用之间仍然存在着一定差距，尤其在多输入多输出系统、非线性特性及参数时变的系统和结果不确定的系统中。预测控制方法发展至今，仍然存在一些问题，具体如下：

①模型难以建立。模型是预测控制方法的基础，因此建立的模型越精确，预测控制效果越好。尽管模型辨识技术已经在预测控制方法的建模过程中得以应用，但是仍无法建立非常精确的系统模型。

②在线计算过程不够优化。预测控制方法的一大特征是在线优化，即根据系统当前状态、性能指标和约束条件进行在线计算得到当前状态的控制律。在在线优化过程中，当前的优化算法主要有线性规划、二次规划和非线性规划等。在线性系统中，预测控制的在线计算过程大多数采用二次规划方法进行求解，但若被控对象的输入输出个数较多或预测时域较大时，该优化方法的在线计算效率也会无法满足系统快速性需求。而在非线性系统中，在线优化过程通常采用序列二次优化算法，但该方法的在线计算成本相对较高且不能完全保证系统稳定，因此也需要不断改进。

③误差问题。由于系统建模往往不够精确，且被控系统中往往存在各种干扰，预测控制方法的预测值和实际值之间一定会产生误差。虽然建模误差可以通过补偿进行校正，干扰误差可以通过反馈进行校正，但是当系统更复杂时，上述两种校正结合起来也无法将误差控制在一定范围内。

模型预测控制区别于其它算法的最大特征是处理多变量多约束线性系统的

能力，但随着被控对象的输入输出个数的增多，预测控制方法为保证控制输出的精确性，往往会选取较大的预测步长和控制步长，但这样会大大增加在线优化过程的计算量，从而需要更多的计算时间。因此，预测控制方法只能适用于采样周期较大且动态变化过程较慢的系统中。为使预测控制方法能在更多场合中应用，快速模型预测控制算法成为了一个新的研究方向。

国内外研究现状

近年来人们对预测控制算法的不足有了越来越清晰的认识，为了将该算法应用到更多领域，越来越多的学者对其进行不断研究和改进。阅读近些年国内外核心期刊的文献可知，人们对预测控制方法能够在更大更复杂的系统中应用寄予了很高期望，同时也在其不足方面做了很多探究和尝试，发展出了多种智能预测控制方法。快速模型预测控制算法作为目前智能预测控制方法之一，其研究方向主要有以下几个方面：

（1）显式模型预测控制（Explicit Model Predictive Control，EMPC）

2002年Alberto Bemporad等学者[1]提出了显式模型预测控制，该方法在预测控制基础之上，在线性时不变系统优化求解过程中引入多参数二次规划理论，对系统的状态区域进行凸划分，根据最优控制问题计算得到状态分区上相应的控制律。EMPC将模型预测控制的在线计算过程转化为离线和在线计算相结合的过程，大大减少了算法的在线计算时间，弥补了MPC方法反复在线计算的不足，EMPC 也在电力电子[2]、电机控制[3]等领域得到了很好应用。

但随着被控对象问题规模（如输入、状态维数、约束等）的增大，EMPC算法在离线计算过程中所求的状态分区数会呈指数倍增加，而状态分区数的增加不仅会导致存储状态分区所需的存储空间增加，还会导致EMPC算法在线查找最优解所需的计算时间增加，因此该算法很难适用于状态约束较多（状态变量往往不超过5）、预测步长较大的复杂系统。

基于以上原因，许多学者提出用近似的显式模型预测控制方法[4,5]来代替精确的EMPC算法，即通过牺牲一定的控制精度来降低计算过程中的复杂度，从而简化整个求解过程。如文献[6]提出了显式模型预测控制的多尺度近似方法，通

过引入分段线性插入法和自适应分层函数近似法，运用重心插值理论对EMPC离线计算出的状态空间进行网格划分，得到近似的状态分区和近似控制律。基于同样的思想，文献[7]提出了显式模型预测控制多胞体近似方法，主要利用双描述法对最优控制问题进行近似处理，再通过重心插值得到近似控制律。文献[8]基于小波的多分辨率分析提出了近似EMPC，通过二次插值和网格划分得到低复杂度且可保证系统稳定的近似控制律。上述文献提出的方法均能在误差允许的范围内保证系统的控制性能，在一定程度上解决了EMPC随着问题规模增大而带来的复杂度和存储空间增大的问题。

（2）模型预测控制的简化算法

MPC算法采用的是在线滚动优化的控制策略，但随着工业模型和环境越来越复杂，其在线计算量越来越大，所以限制了MPC算法在动态变化较快系统中的应用。为减少MPC在线优化求解的计算量，有学者考虑对参数进行优化，提出了预测控制的简化算法，如将参数分块化的blocking技术[9]，其思想是将越远离当前时刻的控制输入越粗略计算，从而减少在线计算量。在此基础上，又有学者提出了移动blocking方法，其核心是限制系统优化变量个数同时增大系统的有效输入步长[10]。相较于blocking技术，文献[11]提出的简化方法主要把约束分为起作用约束集、不起作用约束集和不确定约束集，对不起作用的约束集进行忽略操作，对起作用约束集和不确定约束集进行优化计算，从而降低在线计算量，加快控制进程。

（3）改进的在线优化算法

模型预测控制的核心是采用反复的迭代优化进行在线求解，选取适当的在线优化方法可以提高在线计算速度。近年来，有学者试图对模型预测控制的标准形式做适当变形或者近似处理，继而降低预测控制方法的在线计算量。如文献[12]提出了一种用扩展的牛顿拉夫逊（Newton-Raphson）方法来代替传统的二次规划方法，当问题规模增大时，不仅能够保证优化问题总是收敛，还可以有效解决MPC在线计算量过大的问题。文献[13]提出了将表存储和在线优化相结合的部分枚举法（Partial Enumeration，PE），对于规模较大的线性模型不仅有很好的控制

效果，而且求解速度是传统MPC 算法的5倍以上。文献[14]提出了基于降精度求

解准则的序列二次规划法，主要通过牺牲一定的精确度来降低在线优化所需的迭

代次数，从而提高在线求解速度。文献[15]提出了可用于预测控制在线优化过程

的有效内点法，该内点法不仅可以减少代码量，而且可以使在线计算时间提高

2-5倍。

理论基础

凸集与凸函数

首先给出凸优化理论中凸集、凸函数和仿射的定义。

凸集：假设C 为n 维实数空间n R 中的集合，若C 中任意两点之间的线段仍

然在集合C 中，即对于任意12,x x C ∈且对任意实数[]0,1λ∈，都满足：

12(1)x x C λλ+-∈

（3-1）则称集合C 为凸集，反之则为非凸集。

满足凸集的集合具有以下性质：

（1）若集合1C 、2C 都是凸集则集合{}12121122|,,C C x x y y y C y C +==+∈∈也

是凸集。

（2）两个或多个凸集的交集仍为凸集。

（3）空集也是凸集。

凸函数：设()f x 是定义在非空凸集n C R ∈上的函数，对于凸集C 中任意的1

x 和2x 以及任意实数[]0,1λ∈均满足

1212((1))()(1)()f x x f x f x λλλλ+-≤+-

（3-2）则称()f x 为凸函数，反之则为非凸函数（或称凹函数）。若对12x x ≠和任意

实数[]0,1λ∈均满足

1212((1))()(1)()f x x f x f x λλλλ+-<+-

（

3-3）则称()f x 为严格凸函数。

同样，根据凸函数的定义可以推导出凸函数的几个性质，具体如下：

（1）若1(),,()m f x f x L 是定义在凸集C 上的凸函数，则11()()

m m f x f x αα++L 也是凸集C 上的凸函数，其中α为任意实数。

（2）若函数()f x 是定义在凸集C 上的凸函数，则对于任意实数α，函数

()f x α也是凸集C 上的凸函数。

（3）若()g x 是单调递增的凸函数，()f x 也是凸函数，则复合函数[]()g f x 也

是凸函数。

（4）若()f x 和()g x 都是定义在凸集C 中非单调递增的凸函数，则()()

f x

g x 也是定义在凸集C 中的凸函数。

仿射：对于集合C 中任意两点1x 和2x 且任意实数R λ∈，若满足

12(1)x x C λλ+-∈

（3-4）则称集合C 是仿射的。也就是说，通过集合n C R ?中任意两点的直线仍在

集合C 中。

凸优化问题

凸优化问题是指性能指标函数和不等式约束函数均为凸函数且由约束条件

所得到的集合是凸集的最优化问题。由于凸优化问题的局部最优解就是全局最优

解，因此对于一个实际问题，如果能表示成凸优化问题，则意味着该问题可以得

到彻底解决，对实际应用有着十分重要的意义。凸优化问题有多种形式，其基本

形式可表示为：

min ()

.. ()0 1,, ()0 1,,i i f x s t g x i m h x i p

≤===L L

（

3-5）其中()f x 为性能指标函数，()0i g x ≤为不等式约束条件，()0i h x =为等式约束条件，()0 1,,()0 1,,i i g x i m C x h x i p ≤=??

=??≤=??L L 为可行域（也称可行集或约束集）。

如果上述问题中性能指标函数()f x 和不等式约束()i g x 都是仿射函数，则该

问题又称线性规划问题（Linear Programming ，LP ），其形式可表示为：

min .. T c x d

s t Gx h Ax b

+≤= （3-6）

其中,m n p n G R A R ??∈∈，对于线性规划问题而言，通常将性能指标函数中常

数项d 忽略，因为其不影响最优解的集合。该问题的几何解释可由图表示，图中

阴影多面体为可行域P ，若性能指标函数中T c x 为线性的，则其等位曲线是和c

正交的超平面（图中的虚线），点*x 对即为最优点。

图 线性规划的几何解释

如果性能指标函数()f x 是凸二次型且约束条件是仿射函数，则该问题又称

二次规划问题（Quadratic Programming ，QP ）。对于二次规划问题，其标准形式

可表示为： 1

min 2.. T T x Qx q x r

s t Ax b Gx h

++=≤

（3-7）

其中,,,,,n n m n

m p n p Q S q R A R b R G R h R ??+∈∈∈∈∈∈。对于该问题，最优解的

求解即是在多面体上极小化一个凸二次函数，其几何解释如图所示，该问题的多

面体可行域P 如图阴影所示，性能指标函数的等位线如图虚线所示，则点*x 即

为最优点。

图 二次规划问题的几何解释

对于凸优化问题的求解，理论上很容易得到其最优解，因为其局部极小点就

是全局极小点，如果求解的凸优化问题中性能指标函数是严格凸函数，则求得的

极小点是唯一极小点。

基于凸优化的快速模型预测控制

原始对偶内点法求解QP 问题

对于标准二次规划问题（3-7），本算法采用原始对偶内点法进行求解。相较

于有效集等方法，原始对偶内点法无需热启动，且通过适当优化即可保证算法在

5-25步内迭代便可精确求解。

首先引入松弛变量p s R ∈，将标准二次规划问题转化为： 1

min 2.. 0

T T x Qx q x

s t Ax b Gx s h

s +=+=≥

（3-8）对于式（3-8），定义双变量m y R ∈且满足等式约束条件，p z R ∈满足不等式

约束，求导后得到其最优化条件（KKT 条件）：

0,00

0,1,,T T i i Ax b

Gx s h

s z Qx q G z A y z s i p =+=≥≥+++===L

（

3-9）

其中0i i z s =为互补松弛条件。

模型预测控制问题

模型预测控制算法需要依靠计算机实现求解，其处理的问题形式一般是离散

控制问题而非连续问题。对任意线性时不变随机控制问题，其动态方程可表示为：

(1)()()()x t Ax t Bu t t ω+=++ （3-10）

其中0,1,t =L 表示时间，()n x t R ∈为状态变量，()m u t R ∈为控制输入，

()n t R ω∈为扰动项，n n A R ?∈和n m B R ?∈分别为状态和输入矩阵。

预测控制方法的核心任务往往是求解系统的性能指标函数，在优化求解过程

中，当该函数的值达到极大或极小时，表明该情况下的系统处于最优状态，而此

时的控制输入即为系统的最优控制律，因此系统最优控制律的选择往往取决于系

统的性能指标函数。可以定义上述问题中的性能指标函数为：

()()()()

()1

()()T T T T final t J x Qx u Ru x t T Q x t T ττττττ+-==++++∑

（3-11）其中{}12,,,n Q diag q q q =L 为状态权重矩阵，{}12,,,m R diag r r r =L 为输入权重

矩阵，final Q 为状态终端权重矩阵，且均为对角阵;T 为预测步长。由于实际被控

系统的输入输出往往存在局限性，需考虑状态变量和控制变量的约束条件，该问

题的模型预测控制形式可表示为：

()()()()()()()min max

min max

min .. (1) ,,1 J

s t x A x B u t t T u u u x x x τττττωττττ+=++=+-≤≤≤≤L

（3-12）

其中(1),,()x t x t T ++L 为状态变量，(1),,()u t u t T ++L 为控制变量。对于该

二次规划问题，求解后的最优解可表示为**(),,(1)u t u t T +-L 和

**(),,(1)x t x t T +-L ，预测控制方法仅取最优控制序列的第一项*()u τ作为控制输

入。

障碍内点法

为更快地求解上式（3-12）中的最优控制问题，利用其特殊结构，采用改进

后的原障碍内点法可以降低问题求解的复杂度，从而达到快速求解的目的。首先

定义全局优化变量[](),(1),,(1),()z u t x t u t T x t T =++-+L 对式子（3-12）进行重

写，代入化简后可转化为：

min .. T T z Hz g z

s t Cz b Pz h

+=≤ （3-13

）

矩阵H 是由各权重矩阵Q 、R 和final Q 组成的块对角矩阵，C 是由状态矩阵A 、

B 和单位矩阵I 组成的稀疏矩阵，矩阵h 由状态和输入的最大值和最小值组成。

利用原障碍内点法对重组后的QP 问题进行求解。原障碍内点法的实质是通

过引入障碍参数，将含不等式约束的最优控制问题转化为等式约束的凸优化问

题。因此，令障碍参数κ>0，替换（3-13）中的不等式约束后得到其近似问题为：

()

min .. T T z Hz g z z s t Cz b κφ++=

（3-14）

其中()()1log lT k

T

i i i z h p z φ++=--∑，i p 表示矩阵P 的行。此时，式（3-13）已经

转化为由确定目标和线性等式约束所组成的凸优化问题。同时由式（3-14）易知，障碍参数κ越接近于零，近似性能指标函数越接近于精确的性能指标函数，求得的最优解也越精确。

对于式（3-14），采用不可行初始点牛顿法进行求解。定义双变量Tn v R ∈且满足上述等式约束条件，求导后得到最优化条件： 0

20T T Cz b J

Hz g P d C v z κ-=?=++++?

（3-15）其中1()T i i i

d h p z =-，

T

i p 表示矩阵P 的行，T P d κ是()z κφ的梯度。此时，

优化问题（3-14）的求解可等价为对优化条件（3-15）的求解。对于当前优化条件，假设当前时刻全局优化变量z 不可行（z 满足不等式约束，但不满足等式约束条件），则需找到一个方向z ?和v ?分别使得z z +?和v v +?满足（至少近似满足）最优化条件，即*z z z +?≈、*v v v +?≈。然后在优化条件中用z z +?代替*z 、v v +?代替*v ，并利用泰勒公式的一阶近似，将上述问题（3-15）进行整理后可得关于梯度z ?和v ?的一组线性方程：

22()20T T T T z H P diag d P C Hz g P d C v v C Cz b κκ?????

++++??=-???????-??????

（

3-16）其中2()T P diag d P κ为()z κφ的Hessian 阵，同时将上述方程等式右边变量分别记为对偶残差2T T d r Hz g P d C v κ=+++，原残差p r Cz b =-，残差向量(,)((,),(,))d p r z v r x v r x v =。对式（3-16）进行计算后可得到梯度z ?和v ?的值，再利用回溯直线搜索法得到步长(]0,1s ∈，同时可保证更新后的优化点满足不等式Pz h <。最后利用步长s 更新原始和对偶变量:z z s z =+?，:v v s v =+?。该迭代修正过程不断重复进行直到残余向量的二范数小于设定的误差阈值。上述问题若严格可行，则需要有限步数的计算便能获得一定的原始可行性，而一旦=0p r ，那么在以后的迭代过程中残余向量会一直等于0，且变量z 和v 将趋于优化点。

改进后的原障碍内点法

K来为加快求解，本算法采用固定的障碍参数κ和固定的最大牛顿迭代步

max

简化计算。在求解优化问题（3-16）时，原障碍内点法分别对一系列不同障碍参数x的优化问题进行依次求解，再分别比较不同障碍参数下的控制性能得出最优解，因此会花费较多的求解时间。大量仿真表明只需使得性能指标函数达到最小就能得到很好的控制律，而非精确求解QP问题，因此对于优化问题（3-16），障碍参数κ越接近于零，近似性能指标函数越接近于精确的性能指标函数J，求得的最优解也越精确。固定障碍参数x能够减少预测控制问题次迭代求解所需的牛顿迭代次数（最大牛顿迭代步），从而加快凸优化问题的求解。

在标准的牛顿法中，只有当残差达到足够小或者迭代达到上限时则停止迭代寻优。若问题不是发散的，往往最大迭代步

K在3到10之间就可以求出最优

max

解。因此，固定的迭代步

K可能会使残差未能达到足够小，优化变量没有达到

max

最优，但这并不影响最终控制效果，还能提高在线求解速度。

[1] Bemporad A, Morari M, Dua V, Pistikopoulos E N. The Explicit Linear Quadratic Regulator for Constrained Systems[J]. Automatica,2002,38(1):3-20.

[2] Beccuti A GPapafotiou G, Roberto Frasca, Marari M. Explicit Hybrid Model Predictive Control of the Boost DC-DC Converter[A]. IEEE Power Power Electronics Specialists Conference[C], Orlando, Florida, USA,2007:2503-2509.

[3] Mariethoz S, Morari M. Explicit Model-Predictive Control of a PWM Imverter With an LCL Filter[J].IEEE Transactions on Industrial Electronics,2009,56(2):389-399.

[4] Nguyen H, Olaru s, Hovd patchy approximation of explicit model predictive control[J].Intermational Jourmal of Control,2012,85(12):1929-1941.

[5] Bemporad A, Filippi C. An algorithm for approximate multiparametric convex programming[J].Computational Optimization and Applications,2006,35(1):87-108. [6] Summers S, Jones C N, Lygeros J, Morari multiscale approximation scheme for explicit model predictive control with stability, feasibility, and performance guarantees[C].48th EEE Conference on Decision and Control and 28th Chinese Control Conference[C],2009:6327-6332.

[7] Jones CN, Morari M. Polytopic Approximation of Explicit Model Predictive Controllers[J]. IEEE Transactions on Automatic Control,2010,55(11):2542-2553. [8] Summers S, Jones CN, Morari multiresolution approximation method for fast explicit model predictive control[J]. IEEE Transactions on Automatic Control,2011,56(11):2530-2541.

[9] Ricker N L. Use of quadratic programming for constrained intermal model control[J]. Industrial and Engineering Chemistry Process Design and Development,1985,24(4):925-938.

[10] Valencia-Palomo GPelegrinis M, Rossiter J A, Gondhalekar move-blocking strategy to improve tracking in predictive control[A]. American Control Conference[C], Baltimore, USA,2010:

[11]李德伟，席裕庚，秦辉.预测控制等效集结优化策略的研究[J].自动化学报，2007，33（3）：302-308.

[12]Kouvaritakis B，Cannon M，Rossiter J needs QP for linear MPC anyway[J].Automatica，2002，38（5）：879-884.

[13]Pannocchia J B，Wright S ，large scale model predictive control by partial enumeration[J].Automatica，2007，43（5）：852-860.

[14]万娇娜，邵之江，王可心，王志强，钱积新.基于降精度求解准则的快速最优控制问题求解算法

[15] Domahidi A, Zgraggen A U, Zeilinger M N, Morari M, Jones C N. Efficient Interior Point Methods for Multistage Problems Arising in Receding Horizon Control[A]. IEEE Conference on Decision and Control[C], Maui, USA,2012:668-674.