模型预测控制：设定点跟踪（Set Point Tracking）

模型预测控制（Model Predictive Control, MPC）不仅可以用于系统稳定性问题，还可以用于设定点跟踪问题（Set Point Tracking），即系统需要跟踪一个动态变化的参考输入。本文将介绍如何将设定点跟踪问题转化为稳定性问题，并使用MPC解决。

系统描述

首先，我们假设系统模型为一个离散时间线性时不变（LTI）系统：

$x_{k+1} = A x_k + B u_k$

其中：

$x_k$ 是时间步 $k$ 的系统状态向量
$u_k$ 是时间步 $k$ 的控制输入向量
$A$ 是状态转移矩阵
$B$ 是输入矩阵

对于设定点跟踪问题，目标是让系统状态 $x_k$ 跟踪一个参考状态 $x_r$ ，同时控制输入 $u_k$ 跟踪一个参考输入 $u_r$ 。

转化为稳定性问题

设定误差状态为：

$e_k = x_k - x_r$

误差控制输入为：

$v_k = u_k - u_r$

将其代入系统状态方程：

$e_{k+1} = x_{k+1} - x_r = A x_k + B u_k - x_r = A (x_k - x_r) + B (u_k - u_r) + A x_r + B u_r - x_r$

考虑到在稳态时 $A x_r + B u_r = x_r$ ，则有：

$e_{k+1} = A e_k + B v_k$

这样，我们将设定点跟踪问题转化为误差状态的稳定性问题。

优化问题

优化问题的目标是最小化误差状态和控制输入的代价函数，同时满足系统的约束。代价函数可以表示为：

$\sum_{i=0}^{N-1} \left( e_{k+i|k}^T Q e_{k+i|k} + v_{k+i|k}^T R v_{k+i|k} \right) + e_{k+N|k}^T P e_{k+N|k}$

其中：

$N$ 是预测时域长度
$Q$ 是误差状态权重矩阵
$R$ 是误差控制输入权重矩阵
$P$ 是终端误差状态权重矩阵

约束条件为：

$\begin{align*} e_{k+i+1|k} &= A e_{k+i|k} + B v_{k+i|k} \quad \forall i = 0, \ldots, N-1 \\ x_{min} &\leq e_{k+i|k} + x_r \leq x_{max} \quad \forall i = 1, \ldots, N \\ u_{min} &\leq v_{k+i|k} + u_r \leq u_{max} \quad \forall i = 0, \ldots, N-1 \end{align*}$

Python代码示例

以下是一个简单的Python代码示例，展示如何实现设定点跟踪的线性MPC：

import numpy as np
from scipy.linalg import solve_discrete_are
from scipy.optimize import minimize
import matplotlib.pyplot as plt# Define the system model
A = np.array([[1.1]])
B = np.array([[1.0]])
Q = np.array([[1.0]])
R = np.array([[0.1]])
N = 10  # Prediction horizon# Define constraints
x_min = -10.0
x_max = 10.0
u_min = -1.0
u_max = 1.0# Terminal weight matrix P
P = solve_discrete_are(A, B, Q, R)# Initial state and reference state
x0 = np.array([0.0])
x_r = np.array([5.0])
u_r = np.array([0.5])def cost_function(v, x0, x_r, u_r):e = x0 - x_rcost = 0.0for i in range(N):e = A @ e + B * v[i]cost += e.T @ Q @ e + v[i] * R * v[i]cost += e.T @ P @ ereturn cost# Optimization
v0 = np.zeros(N)
constraints = [{'type': 'ineq', 'fun': lambda v: v + u_r - u_min},{'type': 'ineq', 'fun': lambda v: u_max - (v + u_r)}]result = minimize(cost_function, v0, args=(x0, x_r, u_r), constraints=constraints)
v_opt = result.x# Apply the optimal control input
x = x0
x_history = [x]
u_history = []
for v in v_opt:u = v + u_ru_history.append(u)x = A @ x + B * ux_history.append(x)# Convert x_history to a 1D array of scalars
x_history = [x.item() for x in x_history]# Plot results
plt.plot(range(N+1), x_history, label='State')
plt.axhline(y=x_r, color='r', linestyle='--', label='Reference state')
plt.xlabel('Time step')
plt.ylabel('State')
plt.title('Linear MPC for Set Point Tracking')
plt.legend()
plt.grid(True)
plt.show()plt.plot(range(N), u_history, label='Control input')
plt.axhline(y=u_r, color='r', linestyle='--', label='Reference input')
plt.xlabel('Time step')
plt.ylabel('Control input')
plt.title('Control Input for Set Point Tracking')
plt.legend()
plt.grid(True)
plt.show()

在这里插入图片描述

代码解释

系统模型和参数定义：
- 定义系统的状态矩阵 $A$ 和输入矩阵 $B$ 。
- 定义代价函数中的权重矩阵 $Q$ 和 $R$ 。
- 设置预测时域长度 $N$ 。
- 定义状态和控制输入的约束。
终端权重矩阵计算：
- 使用 solve_discrete_are 计算离散时间代数Riccati方程的解，得到终端权重矩阵 $P$ 。
初始状态和参考状态：
- 定义系统的初始状态 $x_0$ 和参考状态 $x_r$ 。
- 定义参考控制输入 $u_r$ 。
代价函数：
- 定义代价函数 cost_function，用于计算给定误差控制输入序列 $v$ 的代价。
优化：
- 使用 minimize 函数进行优化，得到最优误差控制输入序列 $v_{opt}$ 。
应用控制输入：
- 计算实际控制输入 $u = v + u_r$ ，并更新系统状态 $x$ 。
- 记录系统状态和控制输入的历史数据。
绘制结果：
- 绘制系统状态和参考状态的变化曲线。
- 绘制控制输入和参考输入的变化曲线。