自编Matlab代码实现MPC基本算法

自编Matlab代码实现MPC的基本算法

本文属于入门级的教程，适合刚入坑的萌新看。在写之前我在网上找了一下相关的文章，发现没有一个写得很直观的答案，刚好最近在研究MPC，于是决定自己写一个，针对最简单的线性约束系统，推导其优化的过程原理，揭示其简单的结构。另外关于MPC的理论，可以自行参考相关教材，这里不再赘述。

首先要介绍一下一个基于Matlab的工具MPT3，这是一个优化工具箱，同样也能调用MPC，感兴趣的同学可以自己点进去看一下。今天要用的模型也是来自这里面的一个demo，可以参考这个网址https://www.mpt3.org/UI/RegulationProblem

一、模型准备

这里用最简单的线性模型来进行解释，用差分方程进行描述如下：
$$x^{+}=Ax+Bu$$
其中，$A=\left[\begin{array}{cc}1& 1\\ 0& 1\end{array}\right]$， $B=\left[\begin{array}{c}1\\ 0.5\end{array}\right]$. 初始状态$x_0={\left[\begin{array}{cc}4& 0\end{array}\right]}^{\prime }$. 然后是状态量以及控制量的约束限制：
$$\left[\begin{array}{c}-5\\ -5\end{array}\right]\le x\le \left[\begin{array}{c}5\\ 5\end{array}\right],-1\le u\le 1$$
一般的线性MPC中用的是二次优化，分别包含了针对状态量的代价函数$x_k^{\prime }Q{x}_k$项以及针对控制量的代价函数$u_k^{\prime }R{u}_k$:
$$J_{(k,N_p)}=\mathrm{m}\mathrm{i}\mathrm{n}\sum _{t=k}^{k+N_p}(x_t^{\prime }Q{x}_t+u_t^{\prime }R{u}_t)$$
表示的是在时刻$t$，滚动时域步长为$N_p$的优化目标函数，目的是要让状态量$x$以及控制量$u$尽量的小，即使用最少的能量使系统状态量尽快稳定到零点（或者是稳定到某一个状态，这也是一般优化问题的目的）。其中代价系数矩阵的取值分别为：
$$Q=\left[\begin{array}{cc}1& 0\\ 0& 1\end{array}\right],R=1$$
其取值可以根据设计要求在状态量与控制量的大小之间平衡。预测时域的长度选为$N_p=5$，模型及参数准备完毕开始进行下一步。

二、两个关键问题

2.1 未来状态量的变量替换

为了使用Matlab中的quadprog函数进行在线二次优化，求解Matlab中的优化问题，需要对优化目标进行转换，变换为关于控制量$u_t$的优化函数，回顾优化目标函数，里面包含了在$N_p$步后的未来系统状态量$x_t$，这里需要将其描述为关于优化目标$u_t$的函数。由于是线性系统，通过不停的迭代可以很容易得到$x_t$关于$u_t$的表达式：
x k + 1 ∣ k = A x k + B u k x k + 2 ∣ k = A x k + 1 ∣ k + B u k + 1 = A 2 x k + A B u k + B u k + 1 … … x k + N p ∣ k = A x k + N p − 1 ∣ k + B u k + N p − 1 = A N p