当前位置：首页 > news >正文

【MPC】实战：基于MPC的车辆自适应巡航控制 (ACC) 系统设计

news 2025/6/26 15:37:41

前言

ACC 简介：车辆自适应巡航控制系统（Adaptive Cruise Control, ACC）是一种智能驾驶辅助技术，通过雷达或摄像头实时监测前方路况，自动调整车速以保持安全行驶。当前方无车辆时，ACC 会按照驾驶员设定的目标车速进行定速巡航，减少驾驶疲劳；当检测到前方有车辆时，系统会自动调整车速，保持预设的安全车距，并在必要时减速甚至制动，确保行车安全。

一、定速巡航控制器

1.1 系统模型

设 $v$ 为车辆的纵向速度， $v_{res}$ 为目标速度，设 $e=v-v_{res}$ ， $\dot{v}_{res} = 0$ ，有：
$\begin{align*} \dot{e}(t) &= \dot{v}(t) - \dot{v}_{res}(t) \\ &=[F(t) - F_{aero}(t) - F_{r}(t)]/m \end{align*}$
其中 $m = 1230 k g$ 为车辆总质量， $F$ 为总牵引力，
$F_{aero} = \frac{1}{2} \rho C_dA_Fv^2$ 为空气阻力（设风速为 0 ）， $\rho = 1.206 kg/m^3$ 为空气密度， $C_d \approx 0.3$ 为空气阻力系数， $A_F = 1.6 m^2$ 为迎风面积，
$F_r = fmg$ 为滚动阻力， $2.5\times10^{-5}v$ 为滚阻系数， $g=9.8 m/s^2$ 为重力加速度。

设控制系统采样时间为 $t_s$ ，在 $t_0,~t_0+t_s]$ 时域内对两端做定积分，有：
$\begin{align*} \int_{t_0}^{t_0+t_s} \dot{e}(t)~dt &= \int_{t_0}^{t_0+t_s} [F(t) - F_{aero}(t) - F_{r}(t)]/m~dt \end{align*} \tag{A.1}$
由积分中值定理可知：
$\int_{t_0}^{t_0+t_s} [F(t) - F_{aero}(t) - F_{r}(t)]/m~dt = \frac{F(\xi_1) - F_{aero}(\xi_2) - F_{r}(\xi_3)}{m}t_s$
使用前向欧拉近似，有 $F(\xi_1) = F(t_0)$ ， $F_{aero}(\xi_2) = F_{aero}(t_0)$ ， $F_{r}(\xi_3) = F_{r}(t_0)$ ，式 (A.1) 变为：
$e(t_0+t_s) - e(t_0) = \frac{F(t_0) - F_{aero}(t_0) - F_{r}(t_0)}{m}t_s$
设 $t_0$ 时刻为状态 $k$ ，离散系统可写为：
$e_{k+1} = e_k + \frac{t_s}{m}F_k - \frac{ F_{aero_k} + F_{r_k}}{m}t_s \tag{A.2}$

1.2 MPC 控制器设计

将系统 (A.2) 记为：
$x_{k+1} = Ax_k + Bu_k + D \tag{A.3}$
其中 $x_k = e_k$ ， $u_k = F_k$ ， $A = 1$ ， $B = t_s/m$ ， $D = - (F_{aero_k} + F_{r_k})t_s/m$ ，取 $t_s = 0.01s$ .
为了加速平缓，且考虑路面附着极限，设总控制输入约束为：
$\le u_k \le 0.4mg \approx 4800$
针对以上系统，显然我们有两种选择：（1）忽略非线性项，设计鲁棒MPC；（2）设计非线性MPC。

1.2.1 鲁棒MPC

（详细参考【MPC】模型预测控制笔记 (5)：抗干扰鲁棒MPC）
定义名义系统为：
$z_{k+1} = Az_k + Bv_k \tag{A.4}$
定义误差为 $r_k = x_k-z_k$ ，设最终控制输入 $u_k = v_k - Kr_k$ ，误差状态方程为：
$e_{k+1} = (A-BK)e_k + D$
取 $Q = 500$ ， $R = 0.0001$ ，使用LQR计算得最优增益 $K = 2216$ .
计算得终端代价权重 $P = 27755$ .
设控制器在误差 $r_k|<0.5 m/s$ 下工作， $∣ Ke ∣ < 1108$ .
故设MPC控制约束为：
$-3692\le v_k \le 3692$
取控制时域 $N = 50$ ，在无约束下，无约束MPC计算得到的最优反馈增益 $K_{mpc} = [1~0~0~\cdots~0](\mathcal{H}^T \mathcal{Q} \mathcal{H} + \mathcal{R})^{-1}\mathcal{H}^T \mathcal{Q} \mathcal{G} \approx 235$ .
显然，在名义系统中，整个数域都是不变集，取不变集 $\Omega = [-10,10]$ ，以满足 $-K_{mpc}\Omega$ 在控制约束范围内。
根据以上约束设计终端不等式MPC.
通过Simulink-Carsim联合仿真：
在这里插入图片描述

function [vk, rk] = accRobustMPC(xCur, H, G, Qp, Rp, N, A, B)persistent lb ub tmpReshape tmpAin tmpbin options zCur V;if isempty(zCur)lb = -3692*ones(N,1);ub = -lb;tmpReshape = kron(ones(1, N-1), zeros(1));tmpReshape = [tmpReshape, eye(1)];tmpAin = [1; -1];tmpbin = [10; 10];options = optimoptions('quadprog','Algorithm', 'active-set', 'MaxIterations', 200, 'Display','none');zCur = xCur;V = zeros(N, 1);endHp = 2 * (H' * Qp * H + Rp);fp = 2 * zCur' * G' * Qp * H;Hp = 0.5 * (Hp + Hp');Ain = tmpAin * tmpReshape * H;bin = tmpbin - tmpAin * tmpReshape * G * zCur;V = quadprog(Hp, fp, Ain, bin, [], [], lb, ub, V, options);vk = V(1);zCur = A*zCur + B*vk;rk = xCur - zCur;end

设初始车速为 $25 m / s$ ，目标跟踪车速为 $30 m / s$ ，有以下结果：
在这里插入图片描述
名义系统 (A.4) 中状态 $z_k$ 趋于 0，但我们把模型中非线性项都视为了扰动，真实系统可稳定但误差不趋于0，最终有界.。

1.2.2 非线性MPC

详细参考【MPC】模型预测控制笔记 (7)：非线性MPC 设计非线性MPC：
在这里插入图片描述

function uk  = nonlinearMPC(xCur, vRes, N)persistent lb ub UX;if isempty(ub)lb = -4800*ones(2*N,1);ub = -lb;UX = zeros(2*N, 1);endoptions = optimoptions('fmincon', 'Algorithm', 'sqp', 'Display', 'none');nonlconFun = @(UX) nonlcon(UX, xCur, vRes);UX = fmincon(@objfun, UX, [], [], [], [], lb, ub, nonlconFun, options);uk = UX(1);
end%% 代价函数
function f = objfun(UX) % XU = [u0, u1, ..., u(N-1), x1, x2, .., xN]N = length(UX)/2;U = UX(1:N);X = UX(N+1:end);% 定义权重矩阵Q = 500;R = 0.0001;Qp = kron(eye(N), Q);Rp = kron(eye(N), R);f = X'*Qp*X + U'*Rp*U;
end
%% 非线性约束
function [c,ceq] = nonlcon(UX, xCur, vRes)N = length(UX)/2;n = 1; % x维数p = 1; % u维度U = UX(1:N);X = UX(N+1:end);c = []; % 非线性不等式约束ceq = zeros(n*(N), 1); % 非线性等式约束ceq(1:n) = X(1:n) - model(xCur, U(1:p), vRes); % x1k = f(x0k, u0k);for i = 2:Nceq( (i-1)*n + 1 : i*n ) = X( (i-1)*n+1 : i*n) - model( X((i-2)*n+1 : (i-1)*n), U((i-1)*p + 1:i*p), vRes);end
endfunction xPlus = model(xCur, uCur, vRes)v = xCur + vRes;D = -(0.5*1.206*0.3*1.6*v^2 + (0.004 + 2.5e-5*v)*1230*9.8)*0.01/1230;xPlus = xCur + 0.01/1230*uCur + D;
end