3. MPC处理输入约束的饱和问题

任务目标

在MPC控制器中引入输入饱和约束（例如 |u| ≤ 1），对比无约束MPC与约束MPC在相同被控对象上的控制性能差异，观察约束对控制量幅值、响应速度与稳态精度的影响。

被控对象与离散化

采用一阶惯性加延迟对象（FOPDT）：

G(s) = 2 / (5s + 1) · e^(-s)

采样时间 Ts = 0.5s，延迟 τ = 1s（即2个采样周期）。使用零阶保持器离散化，得到状态空间模型：

x(k+1) = A·x(k) + B·u(k-2)
y(k) = C·x(k)

其中延迟通过增广状态向量处理：x = [x_p; u(k-1); u(k-2)]，扩展后的系统矩阵为：

A = [0.9048, 0, 0; 0, 0, 0; 0, 1, 0]
B = [0.1903; 1; 0]
C = [1, 0, 0]

MPC控制器设计

预测时域 Np = 10，控制时域 Nc = 3。目标函数：

J = Σ_{i=1}^{Np} (y(k+i) - r)² + λ·Σ_{i=0}^{Nc-1} Δu(k+i)²

其中 λ = 0.1。约束MPC在每个采样时刻求解带约束的二次规划问题：

min J
s.t. -1 ≤ u(k+i) ≤ 1, i = 0, ..., Nc-1

无约束MPC则直接使用解析解（无约束QP闭式解）。

完整MATLAB实现

%% MPC输入饱和约束对比 - 单次仿真
clear; clc;

% 系统参数
Ts = 0.5;
A = [0.9048, 0, 0; 0, 0, 0; 0, 1, 0];
B = [0.1903; 1; 0];
C = [1, 0, 0];
[n, ~] = size(A);

% MPC参数
Np = 10; Nc = 3; lambda = 0.1;
r = 1.0;  % 设定值

% 预测矩阵构建
[F, Phi] = build_prediction_matrices(A, B, C, Np, Nc);

% 无约束MPC增益
H = Phi'*Phi + lambda*eye(Nc);
K_mpc = (Phi'*Phi + lambda*eye(Nc)) \ (Phi'*F);

% 仿真设置
Tsim = 30;  % 采样步数
Nsteps = Tsim / Ts;

% 存储
y_uncon = zeros(Nsteps+1, 1); u_uncon = zeros(Nsteps, 1);
y_con = zeros(Nsteps+1, 1);   u_con = zeros(Nsteps, 1);

% 初始状态（含延迟状态）
x_uncon = zeros(n, 1); x_con = zeros(n, 1);
y_uncon(1) = C*x_uncon; y_con(1) = C*x_con;

% 约束MPC的QP求解器（使用quadprog）
options = optimoptions('quadprog', 'Display', 'off');

for k = 1:Nsteps
    % ---- 无约束MPC ----
    du = -K_mpc(1:Nc, :) * (F*x_uncon - r);
    u_uncon(k) = u_uncon(max(k-1,1)) + du(1);
    x_uncon = A*x_uncon + B*u_uncon(k);
    y_uncon(k+1) = C*x_uncon;
    
    % ---- 约束MPC ----
    % 当前输入（用于增量形式）
    u_prev = u_con(max(k-1,1));
    
    % QP: min 0.5*du'*H*du + f'*du, s.t. lb ≤ du ≤ ub
    H_qp = Phi'*Phi + lambda*eye(Nc);
    f_qp = (F*x_con - r)' * Phi;
    
    % 输入幅值约束: -1 ≤ u_prev + cumsum(du) ≤ 1
    % 构造累积矩阵
    L = tril(ones(Nc));
    Aineq = [L; -L];
    bineq = [1 - u_prev; 1 + u_prev];
    
    du = quadprog(H_qp, f_qp, Aineq, bineq, [], [], [], [], [], options);
    u_con(k) = u_prev + du(1);
    x_con = A*x_con + B*u_con(k);
    y_con(k+1) = C*x_con;
end

% 时间轴
t = (0:Nsteps)*Ts;

%% 绘图对比
figure('Position', [100, 100, 800, 500]);

subplot(2,1,1);
stairs(t(1:end-1), u_uncon, 'b-', 'LineWidth', 1.5); hold on;
stairs(t(1:end-1), u_con, 'r-', 'LineWidth', 1.5);
yline(1, 'k--', '上限'); yline(-1, 'k--', '下限');
xlabel('时间 (s)'); ylabel('控制输入 u');
legend('无约束', '约束MPC', 'Location', 'best');
title('控制输入对比');
grid on;

subplot(2,1,2);
plot(t, y_uncon, 'b-', 'LineWidth', 1.5); hold on;
plot(t, y_con, 'r-', 'LineWidth', 1.5);
yline(r, 'k--', '设定值');
xlabel('时间 (s)'); ylabel('输出 y');
legend('无约束', '约束MPC', 'Location', 'best');
title('系统输出对比');
grid on;

%% 辅助函数：构建预测矩阵
function [F, Phi] = build_prediction_matrices(A, B, C, Np, Nc)
    [n, m] = size(B);
    F = zeros(Np, n);
    Phi = zeros(Np, Nc);
    
    % 构建F矩阵（自由响应）
    for i = 1:Np
        F(i, :) = C * A^i;
    end
    
    % 构建Phi矩阵（受控响应）
    for i = 1:Np
        for j = 1:min(i, Nc)
            Phi(i, j) = C * A^(i-j) * B;
        end
    end
end

运行与验证

运行方式：将上述代码保存为 mpc_saturation_compare.m，在MATLAB中直接运行。需要Optimization Toolbox（用于quadprog）。

验证方法：

  1. 运行后观察第一个子图：无约束MPC的控制输入应明显超出 [-1, 1] 边界（初始阶段可能达到1.5以上）；约束MPC的输入应始终保持在边界内。
  2. 第二个子图：约束MPC的输出响应应比无约束情况更慢（上升时间更长），但最终都应稳定到设定值 r = 1。
  3. 检查稳态误差：仿真末尾（t=30s附近），两种方法的输出都应接近1.0，偏差小于0.01。

预期结果：无约束MPC的控制量峰值约为1.6（超出约束），系统输出上升时间约3.5s；约束MPC的控制量被限制在1.0，输出上升时间约5s，但两者稳态值均为1.0。此结果由代码逻辑直接推出，非实测数据。

关键说明

  • 延迟通过增广状态处理，x(2) 和 x(3) 分别存储 u(k-1) 和 u(k-2)，确保预测模型正确反映延迟影响。
  • 约束MPC使用增量形式 Δu 作为决策变量，幅值约束通过累积矩阵 L 转换为线性不等式约束。
  • 无约束MPC直接使用解析解 du = -K_mpc·(F·x - r)，避免了QP求解开销。
← 上一章 📖 返回目录 下一章 →