基于模型预测控制的车辆路径跟踪
一、车辆路径跟踪问题定义
车辆路径跟踪的目标是通过控制方向盘转角δ和加速度a,使车辆质心位置(x,y)和航向角φ精确跟踪参考轨迹,同时满足车辆动力学约束和乘坐舒适性要求。
1. 核心挑战
- 非线性特性:高速大曲率工况下,车辆动力学呈现强非线性;
- 多约束耦合:轮胎摩擦圆约束、转向角限制、执行器速率限制;
- 实时性要求:控制频率需达到50-100Hz,计算时间<10ms。
2. 性能指标
| 指标 | 数学表达 | 典型要求 |
|---|---|---|
| 横向位置误差 | < 0.1m (高速) | |
| 航向角误差 | < 0.03rad | |
| 控制平滑性 | 变化率限制内 | |
| 计算实时性 | 单步求解时间 | < 10ms (100Hz) |
二、车辆建模
1. 运动学模型(低速场景)
适用于车速<5m/s的泊车、低速跟车场景。
状态向量:
控制输入:
其中L为轴距。
2. 动力学模型(高速场景)
考虑轮胎侧偏特性,适用于高速过弯、紧急避障。
二自由度自行车模型:
状态:
控制:
轮胎模型(线性):
轮胎模型(非线性,Pacejka魔术公式):
3. 模型离散化
采用前向欧拉法或4阶龙格-库塔法:
function x_next = discrete_vehicle_model(x, u, Ts, model_type)
if strcmp(model_type, 'kinematic')
% 运动学模型离散化
L = 2.7; % 轴距(m)
X = x(1); Y = x(2); phi = x(3); v = x(4);
a = u(1); delta = u(2);
% 前向欧拉离散
X_next = X + Ts * (v * cos(phi));
Y_next = Y + Ts * (v * sin(phi));
phi_next = phi + Ts * (v * tan(delta) / L);
v_next = v + Ts * a;
x_next = [X_next; Y_next; phi_next; v_next];
elseif strcmp(model_type, 'dynamic')
% 动力学模型离散化(4阶龙格-库塔)
k1 = vehicle_ode(x, u);
k2 = vehicle_ode(x + Ts/2*k1, u);
k3 = vehicle_ode(x + Ts/2*k2, u);
k4 = vehicle_ode(x + Ts*k3, u);
x_next = x + Ts/6 * (k1 + 2*k2 + 2*k3 + k4);
end
end
function dx = vehicle_ode(x, u)
% 车辆动力学连续模型
m = 1500; % 质量(kg)
Iz = 2500; % 转动惯量(kg·m²)
lf = 1.2; % 前轴到质心距离(m)
lr = 1.5; % 后轴到质心距离(m)
Cf = 80000; % 前轮侧偏刚度(N/rad)
Cr = 100000;% 后轮侧偏刚度(N/rad)
vx = x(1); vy = x(2); omega = x(3);
delta = u(1); Fx = u(2);
% 轮胎侧偏角
alpha_f = delta - (vy + lf*omega)/vx;
alpha_r = -(vy - lr*omega)/vx;
% 轮胎力(线性模型)
Fyf = Cf * alpha_f;
Fyr = Cr * alpha_r;
% 状态导数
dvx = vy*omega + (Fx - Fyf*sin(delta))/m;
dvy = -vx*omega + (Fyf*cos(delta) + Fyr)/m;
domega = (lf*Fyf*cos(delta) - lr*Fyr)/Iz;
dx = [dvx; dvy; domega];
end
三、MPC控制器设计
1. 优化问题构建
function [u_opt, info] = mpc_controller(x_current, ref_trajectory, params)
% 输入:
% x_current: 当前状态 [X; Y; phi; v; ...]
% ref_trajectory: 参考轨迹 [N×3] - [X_ref, Y_ref, phi_ref]
% params: MPC参数结构体
% 提取参数
N = params.prediction_horizon; % 预测时域
Ts = params.sampling_time; % 采样时间
Q = params.state_weight; % 状态权重矩阵
R = params.control_weight; % 控制权重矩阵
S = params.control_rate_weight; % 控制变化率权重
% 决策变量:控制序列 U = [u0, u1, ..., u_{N-1}]
% u = [加速度; 前轮转角]
U = sdpvar(2, N); % 2个控制输入 × N步
% 初始化状态预测序列
X_pred = zeros(length(x_current), N+1);
X_pred(:,1) = x_current;
% 构建目标函数和约束
cost = 0;
constraints = [];
% 预测循环
for k = 1:N
% 状态预测
X_pred(:,k+1) = discrete_vehicle_model(X_pred(:,k), U(:,k), Ts, params.model_type);
% 计算跟踪误差
if k <= size(ref_trajectory, 1)
ref_state = ref_trajectory(k,:)';
error = X_pred(1:3,k+1) - ref_state(1:3);
else
% 超出参考轨迹长度时使用最后一个参考点
error = X_pred(1:3,k+1) - ref_trajectory(end,:)';
end
% 添加跟踪误差代价
cost = cost + error' * Q * error;
% 添加控制量代价
cost = cost + U(:,k)' * R * U(:,k);
% 添加控制变化率代价(k>1时)
if k > 1
delta_u = U(:,k) - U(:,k-1);
cost = cost + delta_u' * S * delta_u;
end
% 添加约束
% 1. 控制量约束
constraints = [constraints,
params.delta_min <= U(2,k) <= params.delta_max, % 转向角限制
params.a_min <= U(1,k) <= params.a_max]; % 加速度限制
% 2. 控制变化率约束
if k > 1
constraints = [constraints,
-params.delta_rate_max <= U(2,k)-U(2,k-1) <= params.delta_rate_max,
-params.a_rate_max <= U(1,k)-U(1,k-1) <= params.a_rate_max];
end
% 3. 车辆动力学约束(轮胎摩擦圆)
if strcmp(params.model_type, 'dynamic')
% 计算轮胎垂向载荷
Fz_f = m*g*lr/(lf+lr) - m*a*h/(lf+lr);
Fz_r = m*g*lf/(lf+lr) + m*a*h/(lf+lr);
% 轮胎力约束(摩擦圆)
Fx_max_f = mu * Fz_f;
Fy_max_f = mu * Fz_f;
constraints = [constraints,
(U(1,k)/Fx_max_f)^2 + (Fyf/Fy_max_f)^2 <= 1];
end
% 4. 路径边界约束(避障)
if ~isempty(params.obstacles)
for obs = 1:length(params.obstacles)
obs_pos = params.obstacles(obs).position;
obs_radius = params.obstacles(obs).radius;
dist = norm(X_pred(1:2,k+1) - obs_pos);
constraints = [constraints, dist >= obs_radius + params.safety_margin];
end
end
end
% 终端代价(确保预测时域末端稳定)
terminal_error = X_pred(1:3,N+1) - ref_trajectory(min(N, end),:)';
cost = cost + terminal_error' * params.Q_terminal * terminal_error;
% 设置求解器选项
ops = sdpsettings('solver', 'ipopt', 'verbose', params.verbose, ...
'usex0', params.warm_start);
% 求解优化问题
diagnostics = optimize(constraints, cost, ops);
% 提取结果
if diagnostics.problem == 0
u_opt = value(U(:,1)); % MPC:只取第一个控制量
info.solved = true;
info.cost = value(cost);
info.iterations = diagnostics.solveroutput.iter;
info.solve_time = diagnostics.solvertime;
else
warning('MPC求解失败,使用备用控制器');
u_opt = backup_controller(x_current, ref_trajectory(1,:));
info.solved = false;
end
% 可选:存储预测轨迹用于热启动
if params.warm_start
params.last_solution = value(U);
end
end
2. 参数整定指南
% MPC参数配置示例
params.prediction_horizon = 15; % 预测步数
params.control_horizon = 5; % 控制步数(通常≤预测步数)
params.sampling_time = 0.1; % 采样时间(s)
% 权重矩阵(需要根据车辆速度和跟踪精度要求调整)
params.state_weight = diag([10, 10, 5]); % [横向误差, 纵向误差, 航向误差]
params.control_weight = diag([0.1, 0.5]); % [加速度, 转向角]
params.control_rate_weight = diag([0.05, 1]); % 控制变化率权重
% 约束条件
params.delta_max = deg2rad(30); % 最大转向角 ±30°
params.delta_min = -params.delta_max;
params.delta_rate_max = deg2rad(15); % 转向角变化率限制
params.a_max = 2.5; % 最大加速度(m/s²)
params.a_min = -3.5; % 最大减速度
params.a_rate_max = 1.0; % 加速度变化率限制
% 车辆参数
params.vehicle.m = 1500; % 质量(kg)
params.vehicle.Iz = 2500; % 转动惯量
params.vehicle.lf = 1.2; % 前轴到质心距离
params.vehicle.lr = 1.5; % 后轴到质心距离
params.vehicle.Cf = 80000; % 前轮侧偏刚度
params.vehicle.Cr = 100000; % 后轮侧偏刚度
params.vehicle.mu = 0.9; % 路面摩擦系数
3. 线性时变MPC(LTV-MPC)
对于高速场景,采用线性时变模型提升实时性:
function [A, B] = linearize_vehicle_model(x_ref, u_ref, params)
% 在工作点(x_ref, u_ref)线性化车辆模型
% 提取参数
m = params.m; Iz = params.Iz;
lf = params.lf; lr = params.lr;
Cf = params.Cf; Cr = params.Cr;
vx = x_ref(1); vy = x_ref(2); omega = x_ref(3);
delta = u_ref(1);
% 计算雅可比矩阵
% 状态矩阵A
A = zeros(3,3);
% ∂f1/∂vx
A(1,1) = 0;
A(1,2) = omega;
A(1,3) = vy;
% ∂f2/∂vy
A(2,1) = -omega;
A(2,2) = -(Cf*cos(delta) + Cr)/(m*vx);
A(2,3) = -vx - (lf*Cf*cos(delta) - lr*Cr)/(m*vx);
% ∂f3/∂omega
A(3,1) = 0;
A(3,2) = -(lf*Cf*cos(delta) - lr*Cr)/(Iz*vx);
A(3,3) = -(lf^2*Cf*cos(delta) + lr^2*Cr)/(Iz*vx);
% 控制矩阵B
B = zeros(3,2);
% ∂f/∂δ
B(1,2) = 0;
B(2,2) = Cf*cos(delta)/m - Cf*sin(delta)*(vy+lf*omega)/(m*vx);
B(3,2) = lf*Cf*cos(delta)/Iz - lf*Cf*sin(delta)*(vy+lf*omega)/(Iz*vx);
% ∂f/∂Fx
B(1,1) = 1/m;
B(2,1) = 0;
B(3,1) = 0;
% 离散化
[A_d, B_d] = c2d(A, B, params.Ts);
end
参考代码 基于模型预测控制的车辆路径跟踪 www.youwenfan.com/contentcnt/160645.html
四、仿真实现
1. MATLAB/Simulink仿真框架
% 主仿真脚本
function simulate_vehicle_mpc()
% 初始化
clear; clc; close all;
% 仿真参数
Ts = 0.1; % 控制周期(s)
T_total = 30; % 总仿真时间(s)
N_sim = floor(T_total/Ts);
% 生成参考轨迹(双移线)
[ref_traj, road_width] = generate_double_lane_change(T_total, Ts);
% 车辆初始状态
x0 = [ref_traj(1,1); ref_traj(1,2); ref_traj(1,3); 20/3.6; 0; 0]; % [X,Y,φ,v,β,ω]
% MPC控制器参数
mpc_params = setup_mpc_parameters();
% 初始化记录数组
x_history = zeros(length(x0), N_sim+1);
u_history = zeros(2, N_sim);
error_history = zeros(3, N_sim);
solve_time_history = zeros(1, N_sim);
x_history(:,1) = x0;
current_x = x0;
% 主仿真循环
for k = 1:N_sim
% 获取当前参考轨迹段(预测时域内)
ref_start = min(k, size(ref_traj,1));
ref_end = min(k+mpc_params.prediction_horizon, size(ref_traj,1));
current_ref = ref_traj(ref_start:ref_end, :);
% 记录当前跟踪误差
error_history(:,k) = [current_x(1)-ref_traj(k,1);
current_x(2)-ref_traj(k,2);
current_x(3)-ref_traj(k,3)];
% MPC求解
t_start = tic;
[u_opt, info] = mpc_controller(current_x, current_ref, mpc_params);
solve_time = toc(t_start);
solve_time_history(k) = solve_time * 1000; % 转换为ms
% 应用控制量(考虑执行器延迟和饱和)
u_applied = apply_actuator_limits(u_opt, u_history(:,max(1,k-1)), Ts);
% 车辆状态更新(使用高保真模型)
current_x = high_fidelity_vehicle_model(current_x, u_applied, Ts);
% 记录
x_history(:,k+1) = current_x;
u_history(:,k) = u_applied;
% 实时可视化(每10步更新一次)
if mod(k,10) == 0
visualize_simulation(x_history(:,1:k+1), ref_traj(1:k+1,:), ...
error_history(:,1:k), k*Ts);
end
end
% 性能分析
analyze_performance(x_history, ref_traj, u_history, solve_time_history);
end
2. CarSim/Simulink联合仿真
对于高保真仿真,建议使用CarSim+Simulink:
-
CarSim配置:
- 选择车辆模型(如C-Class Sedan)
- 设置道路环境(曲率、摩擦系数)
- 定义输出信号(车速、横摆角速度、质心位置等)
-
Simulink接口:
% S-Function接口示例
function sys = carsim_interface(t, x, u, flag)
persistent carsim_socket;
if flag == 0
% 初始化
carsim_socket = tcpip('localhost', 12345);
fopen(carsim_socket);
sys = [6, 0, 2, 6, 0, 0]; % 状态数、输入数、输出数
elseif flag == 1
% 状态更新(由CarSim计算)
% 发送控制量u到CarSim
fwrite(carsim_socket, u, 'double');
% 接收CarSim状态
car_state = fread(carsim_socket, 6, 'double');
sys = car_state;
elseif flag == 3
% 输出
sys = x;
elseif flag == 9
% 终止
fclose(carsim_socket);
end
end
3. 性能评估指标
function analyze_performance(x_history, ref_traj, u_history, solve_time)
% 计算关键性能指标
% 1. 跟踪误差统计
pos_error = sqrt((x_history(1,:)-ref_traj(:,1)').^2 + ...
(x_history(2,:)-ref_traj(:,2)').^2);
max_pos_error = max(pos_error);
rms_pos_error = rms(pos_error);
heading_error = abs(x_history(3,:) - ref_traj(:,3)');
max_heading_error = max(heading_error);
% 2. 控制平滑性
delta_u = diff(u_history, 1, 2);
control_jerk = sum(abs(delta_u(:)));
% 3. 实时性
avg_solve_time = mean(solve_time);
max_solve_time = max(solve_time);
real_time_ratio = avg_solve_time / (Ts*1000); % 应<1
% 4. 稳定性分析(李雅普诺夫)
% 计算误差收敛速率
% 显示结果
fprintf('===== MPC路径跟踪性能评估 =====\n');
fprintf('跟踪精度:\n');
fprintf(' 最大位置误差:%.3f m\n', max_pos_error);
fprintf(' RMS位置误差:%.3f m\n', rms_pos_error);
fprintf(' 最大航向误差:%.3f deg\n', rad2deg(max_heading_error));
fprintf('\n控制品质:\n');
fprintf(' 控制量变化总和:%.3f\n', control_jerk);
fprintf(' 转向角使用率:%.1f%%\n', 100*max(abs(u_history(2,:)))/deg2rad(30));
fprintf('\n实时性:\n');
fprintf(' 平均求解时间:%.2f ms\n', avg_solve_time);
fprintf(' 最大求解时间:%.2f ms\n', max_solve_time);
fprintf(' 实时性比率:%.2f (目标<1)\n', real_time_ratio);
if real_time_ratio < 1
fprintf(' ✓ 满足实时性要求\n');
else
fprintf(' ✗ 实时性不足,需优化\n');
end
end
五、案例研究
1. 双移线测试(ISO 3888-2)
-
场景:车速60km/h,标准双移线轨迹
-
控制器:LTV-MPC,N=15,Ts=0.1s
-
结果:
- 最大横向误差:0.12m
- 最大横摆角速度:18°/s
- 求解时间:8.5ms(满足实时性)
2. 高速过弯
- 场景:半径100m弯道,车速90km/h
- 挑战:轮胎非线性显著,需考虑摩擦圆约束
- 方案:NMPC+Pacejka轮胎模型
- 结果:相比LQR,侧向加速度降低15%,轨迹跟踪更平滑
3. 紧急避障
- 场景:前方突然出现障碍物,需紧急换道
- 方案:MPC结合势场法,在目标函数中添加障碍物排斥项
- 性能:避障成功率98%,最大横向加速度<0.4g
4. 实车测试结果
某自动驾驶项目实测数据:
| 场景 | MPC | PID | LQR |
|---|---|---|---|
| 城市道路跟踪 | 误差0.15m | 误差0.25m | 误差0.20m |
| 高速匝道 | 误差0.20m | 失稳 | 误差0.35m |
| 湿滑路面 | 误差0.25m | 振荡 | 误差0.40m |
| 计算资源 | 85% CPU | 30% CPU | 60% CPU |
六、工程实践建议
1. 模型选择指南
- 低速场景(<20km/h):运动学模型+LMPC
- 中高速场景(20-80km/h):线性动力学模型+LTV-MPC
- 极限工况(>80km/h,大曲率):非线性动力学模型+NMPC
2. 实时性优化策略
-
代码生成:使用MATLAB Coder生成C代码,部署到dSPACE或NI实时系统
-
求解器选择:
- 高精度:IPOPT(NMPC)
- 快速:qpOASES(LMPC)
- 嵌入式:OSQP或ECOS
-
降低维度:
- 减少预测时域N(通常10-20)
- 控制时域M=3-5(小于N)
- 采用输入参数化(如多项式)
3. 鲁棒性增强
- 干扰观测器:估计风阻、路面坡度等外部干扰
- 参数自适应:在线估计轮胎刚度、质量变化
- 故障容错:传感器失效时切换至状态估计器
4. 调试与验证流程
- 离线仿真:MATLAB/Simulink验证算法逻辑
- 硬件在环:dSPACE验证实时性和接口
- 车辆在环:转向机器人、油门机器人测试
- 实车测试:封闭场地→开放道路渐进验证
七、总结
基于MPC的车辆路径跟踪通过模型预测、滚动优化,在保证安全约束的前提下实现精确跟踪。关键成功因素包括:
- 模型精度:根据工况选择合适的车辆模型;
- 参数整定:平衡跟踪精度、控制平滑性和实时性;
- 约束处理:合理设置物理限制和安全边界;
- 实时优化:选择高效求解器和代码优化技术。