面向阿克曼移动机器人（自行车模型）的LQR（最优二次型调节器）的路径跟踪方法

news2024/11/17 3:49:12

        线性二次调节器（Linear Quadratic Regulator，LQR）是针对线性系统的最优控制方法。LQR 方法标准的求解体系是在考虑到损耗尽可能小的情况下, 以尽量小的代价平衡其他状态分量。一般情况下，线性系统在LQR 控制方法中用状态空间方程描述，性能能指标函数由二次型函数描述。

LQR 方法存在以下优点：

  1. 最小能量消耗和最高路径跟踪精度。
  2. 求解时能够考虑多状态情况。
  3. 鲁棒性较强。

缺点：

  1. 控制效果和模型精确程度有很大相关性。
  2. 实时计算状态反馈矩阵和控制增益。

一、系统模型

1.1 车辆模型

        一般来说阿克曼移动机器人可以简化为自行车模型，是一个非线性时变系统，工程上一般通过在工作平衡点附近差分线性化转化为线性系统来来分析和控制，具体就不推导了，我直接给出模型。

1.2 线性系统状态反馈控制示意图

状态反馈是线性能控线性系统镇定的一个有效方法，主要是通过极点配置方法寻找一组非正的闭环极点使得闭环系统大范围渐进稳定。

 A，B，C分别代表系统矩阵、输入矩阵和输出矩阵，K是待设计的状态反馈增益。

二、控制器设计

2.1 代价函数泛函设计

最优控制里，代价函数一般设计为控制性能和控制代价的范数加权和，形式如下

其中，期望和实际的误差系统定义为

2.2 最优状态反馈控制律推导

当想要状态与期望状态之间的误差越差越小，同时控制消耗更少的能量。求解极小值点时，新定义拉格朗日函数如下 

在拉格朗日函数基础上对各个优化变量的一阶导为零 ，得

 当\Delta_{x(t)}=0时候，

推导LQR控制律时候，设\lambda(t)=P_tx(t) ,求偏导后可得

由 得

 由于状态量初始不为零，只能是

又由于当上述方程成立时候，P收敛到了期望的范围 ，\dot{P}为零，因此得到立卡提方程形式的矩阵微分方程

综上，通过迭代或者近似方法求解上述立卡提方程后，最优的控制律为

 2.3 连续立卡提方程求解流程

For~~ iteration=1: iteration_{max}\\ ~~~~~~~\bar{P}=Q+A^{\top}PB(R+B^TPB)^{-1}B^{\top}PA\\ ~~~~~~~criteria=\overline\lambda(\vert\bar{P}-P\vert)\\ ~~~~~~~IF~~criteria<\epsilon\\ ~~~~~~~~~~Break\\ ~~~~~~P=\bar{P}

三、具体实现代码

3.1 main函数

close all
clear;
clc; 
cx = [];
cy= [];
y0 = @(t_step) 10*sin(2 * t_step + 1);
x0_dot= @(t_step) 5 * 2 * cos(2 * t_step + 1);
x0 = @(t_step) -40*cos(t_step + 0.5);
for theta=0:pi/200:2*pi
    cx(end + 1) = x0(theta);
    cy(end + 1) = y0(theta);
end
refer_path_primary= [cx', cy'];
x = refer_path_primary(:, 1)';
y = refer_path_primary(:, 2)';
points = [x; y]';
ds = 0.1 ;%等距插值处理的间隔
distance = [0, cumsum(hypot(diff(x, 1), diff(y, 1)))]';
distance_specific = 0:ds:distance(end);
hypot(diff(x, 1), diff(y, 1));
diff(x, 1);
diff(y, 1);
s = 0:ds:distance(end);
refer_path= interp1(distance, points, distance_specific, 'spline');
figure(1)
% 绘制拟合曲线
plot(refer_path(:, 1), refer_path(:, 2),'b-',x, y,'r.');
hold on
refer_path_x = refer_path(:,1);  % x
refer_path_y = refer_path(:,2); % y
for i=1:length(refer_path)
    if i==1
         dx = refer_path(i + 1, 1) - refer_path(i, 1);
         dy = refer_path(i + 1, 2) - refer_path(i, 2);
         ddx = refer_path(3, 1) + refer_path(1, 1) - 2 * refer_path(2, 1);
         ddy = refer_path(3, 2) + refer_path(1, 2) - 2 * refer_path(2, 2);
    elseif  i==length(refer_path)
         dx = refer_path(i, 1) - refer_path(i - 1, 1);
         dy = refer_path(i, 2) - refer_path(i - 1, 2);
         ddx = refer_path(i, 1) + refer_path(i - 2, 1) - 2 * refer_path(i - 1, 1);
         ddy = refer_path(i, 2) + refer_path(i - 2, 2) - 2 * refer_path(i - 1, 2);
    else
          dx = refer_path(i + 1, 1) - refer_path(i, 1);
          dy = refer_path(i + 1, 2) - refer_path(i, 2);
          ddx = refer_path(i + 1, 1) + refer_path(i - 1, 1) - 2 * refer_path(i, 1);
          ddy = refer_path(i + 1, 2) + refer_path(i - 1, 2) - 2 * refer_path(i, 2);
    end
    refer_path(i,3)=atan2(dy, dx);%yaw
    refer_path(i,4)=(ddy * dx - ddx * dy) / ((dx ^ 2 + dy ^ 2) ^ (3 / 2));
end
figure(2)
plot(refer_path(:, 3),'b-');
figure(3)
plot(refer_path(:, 4),'b-')
% 
%%目标及初始状态
L=2;%车辆轴距
v=2;%初始速度
dt=0.05;%时间间隔
goal=refer_path(end,1:2);
x_0=refer_path_x(1);
y_0=refer_path_y(1);
psi_0 = refer_path(1, 3);
% %运动学模型
ugv = KinematicModel(x_0, y_0, psi_0, v, dt, L);
Q = eye(3) * 3.0;
R = eye(2) * 2.0;
robot_state = zeros(4, 1);
step_points = length(refer_path(:, 1));
for i=1:1:step_points
    robot_state(1)=ugv.x;
    robot_state(2)=ugv.y;
    robot_state(3)=ugv.psi;
    robot_state(4)=ugv.v;
    [e, k, ref_yaw, min_idx] = calc_track_error(robot_state(1), robot_state(2), refer_path);
    ref_delta = atan2(L * k, 1);
    [A, B] = state_space( robot_state(4), ref_delta, ref_yaw, dt, L);
    delta = lqr(robot_state, refer_path, min_idx, A, B, Q, R);
    delta = delta + ref_delta;
    [ugv.x, ugv.y, ugv.psi, ugv.v] = update(0, delta, dt, L, robot_state(1), robot_state(2),robot_state(3), robot_state(4));
    ugv.record_x(end + 1 ) = ugv.x;
    ugv.record_y(end + 1 ) = ugv.y;
    ugv.record_psi(end + 1 ) = ugv.psi;
     ugv.record_phy(end + 1 ) = ref_delta;
end
figure(4)
% 绘制拟合曲线
%  scr_size = get(0,'screensize');
% set(gcf,'outerposition', [1 1 scr_size(4) scr_size(4)]);
plot(ugv.record_x , ugv.record_y, Color='m',LineStyle='--',LineWidth=2);
axis([-40,40,-40,40])
grid on
hold on
% 绘制车辆曲线
axis equal
for ii = 1:1:length(ugv.record_x)
       
        h =  PlotCarbox(ugv.record_x(ii), ugv.record_y(ii), ugv.record_psi(ii), 'Color', 'r',LineWidth=2);
        h1 = plot(ugv.record_x(1:ii), ugv.record_y(1:ii),'Color', 'b');
        th1 = text(ugv.record_x(ii), ugv.record_y(ii)+10, ['#car', num2str(1)], 'Color', 'm');
        set(th1, {'HorizontalAlignment'},{'center'});
        h2 = PlotCarWheels(ugv.record_x(ii), ugv.record_y(ii), ugv.record_psi(ii),ugv.record_phy(ii),'k',LineWidth=2);
        h3 = plot(ugv.record_x(1:ii) , ugv.record_y(1:ii), Color='b',LineStyle='-',LineWidth=4);
        drawnow
        delete(h); delete(h1);delete(th1);delete(h3);
        for jj = 1:1:size(h2)
            delete(h2{jj});
        end
end
% 
function [P] =  cal_Ricatti(A, B, Q, R)
 Qf = Q;
 P = Qf;
  iter_max = 100;
  Eps = 1e-4;
 for step = 1:1:iter_max
     P_bar = Q + A' * P * A - A' * P * B * pinv(R + B' * P *B) * B' * P * A;
     criteria = max(abs(P_bar - P));
     if  criteria < Eps
         break;
     end
     P = P_bar;
 end

end
%%LQR控制器
function[u_star]=lqr(robot_state, refer_path, s0, A, B, Q, R)
    x = robot_state(1:3) - refer_path(s0,1:3)';
    P =  cal_Ricatti(A, B, Q, R);
    K= -pinv(R + B' * P * B) * B' * P * A;
    u = K * x;%状态反馈
    u_star = u(2);
end

function [e, k, yaw, min_idx]=calc_track_error(x, y, refer_path)
    p_num = length(refer_path);
    d_x = zeros(p_num, 1);
    d_y = zeros(p_num, 1);
    d = zeros(p_num, 1);
    for i=1:1:p_num 
        d_x(i) = refer_path(i, 1) - x; 
        d_y(i) = refer_path(i, 2) - y;
    end
    for i=1:1:p_num 
        d(i) = sqrt(d_x(i) ^2 + d_y(i) ^ 2) ;
    end
    [~, min_idx] = min(d); 
    yaw = refer_path(min_idx, 3);
    k = refer_path(min_idx, 4);
    angle= normalize_angle(yaw - atan2(d_y(min_idx), d_x(min_idx)));
    e = d(min_idx);
    if angle < 0
        e = e*(-1);
    end
end
%%将角度取值范围限定为[-pi,pi]
function [angle]=normalize_angle(angle)
    while angle > pi
        angle = angle - 2*pi;
    end
    while angle < pi
        angle = angle + 2*pi;
    end
end
function [x_next, y_next, psi_next, v_next] = update(a, delta_f, dt, L, x, y, psi, v)
     x_next = x + v * cos(psi) * dt;
     y_next = y + v * sin(psi) * dt;
     psi_next = psi + v / L * tan(delta_f) * dt;
     v_next = v + a * dt;
end
function [A, B]=state_space(v, ref_delta, ref_yaw, dt, L)
    A=[ 1.0, 0.0,  -v * dt * sin(ref_yaw);
           0.0, 1.0,   v *  dt * cos(ref_yaw);
           0.0, 0.0,   1.0  ];
    B =[ dt * cos(ref_yaw), 0;
           dt * sin(ref_yaw), 0;
           dt * tan(ref_delta) / L, v * dt / (L * cos(ref_delta) * cos(ref_delta))];
end
function h = PlotCarbox(x, y, theta, varargin)
Params = GetVehicleParams();

carbox = [-Params.Lr -Params.Lb/2; Params.Lw+Params.Lf -Params.Lb/2; Params.Lw+Params.Lf Params.Lb/2; -Params.Lr Params.Lb/2];
carbox = [carbox; carbox(1, :)];

transformed_carbox = [carbox ones(5, 1)] * GetTransformMatrix(x, y, theta)';
h = plot(transformed_carbox(:, 1), transformed_carbox(:, 2), varargin{:});

end
function hs = PlotCarWheels(x, y, theta, phy, varargin)
Params = GetVehicleParams();

wheel_box = [-Params.wheel_radius -Params.wheel_width / 2;
    +Params.wheel_radius -Params.wheel_width / 2;
    +Params.wheel_radius +Params.wheel_width / 2;
    -Params.wheel_radius +Params.wheel_width / 2];

front_x = x + Params.Lw * cos(theta);
front_y = y + Params.Lw * sin(theta);
track_width_2 = (Params.Lb - Params.wheel_width / 2) / 2;

boxes = {
    TransformBox(wheel_box, x - track_width_2 * sin(theta), y + track_width_2 * cos(theta), theta);
    TransformBox(wheel_box, x + track_width_2 * sin(theta), y - track_width_2 * cos(theta), theta);
    TransformBox(wheel_box, front_x - track_width_2 * sin(theta), front_y + track_width_2 * cos(theta), theta+phy);
    TransformBox(wheel_box, front_x + track_width_2 * sin(theta), front_y - track_width_2 * cos(theta), theta+phy);
};

hs = cell(4, 1);
for ii = 1:4
    hs{ii} = fill(boxes{ii}(:, 1), boxes{ii}(:, 2), varargin{:});
end

end

function transformed = TransformBox(box, x, y, theta)
transformed = [box; box(1, :)];
transformed = [transformed ones(5, 1)] * GetTransformMatrix(x, y, theta)';
transformed = transformed(:, 1:2);
end
function mat = GetTransformMatrix(x, y, theta)
mat = [ ...
    cos(theta) -sin(theta), x; ...
    sin(theta) cos(theta), y; ...
    0 0 1];
end

3.2 运动学结构体： 

classdef  KinematicModel<handle
    properties
        x;
        y;
        psi;
        v;
        dt;
        L;
        record_x;
        record_y;
        record_psi;
        record_phy;
    end
    methods
        function self=KinematicModel(x, y, psi, v, dt, L)
            self.x=x;
            self.y=y;
            self.psi=psi;
            self.v = v;
            self.L = L;
            % 实现是离散的模型
            self.dt = dt;
            self.record_x = [];
            self.record_y= [];
            self.record_psi= [];
            self.record_phy= [];
        end
end
end

四、仿真参数和效果

 4.1 参数配置

%期望轨迹
y0 = @(t_step) 10*sin(2 * t_step + 1);
x0_dot= @(t_step) 5 * 2 * cos(2 * t_step + 1);
L=2;%车辆轴距
v=2;%初始速度
dt=0.05;%时间间隔
Q = eye(3) * 3.0;
R = eye(2) * 2.0;
robot_state = zeros(4, 1);
VehicleParams.Lw = 2.8 * 2; % wheel base
VehicleParams.Lf = 0.96 * 2; % front hang
VehicleParams.Lr = 0.929 * 2; % rear hang
VehicleParams.Lb = 1.942 * 2; % width
VehicleParams.Ll = VehicleParams.Lw + VehicleParams.Lf + VehicleParams.Lr; % length
VehicleParams.f2x = 1/4 * (3*VehicleParams.Lw + 3*VehicleParams.Lf - VehicleParams.Lr);
VehicleParams.r2x = 1/4 * (VehicleParams.Lw + VehicleParams.Lf - 3*VehicleParams.Lr);
VehicleParams.radius = 1/2 * sqrt((VehicleParams.Lw + VehicleParams.Lf + VehicleParams.Lr) ^ 2 / 4 + VehicleParams.Lb ^ 2);
VehicleParams.a_max = 0.5;
VehicleParams.v_max = 2.5;
VehicleParams.phi_max = 0.7;
VehicleParams.omega_max = 0.5;
% for wheel visualization
VehicleParams.wheel_radius = 0.32*2;
VehicleParams.wheel_width = 0.22*2;
iter_max = 100;
Eps = 1e-4;

  4.1 仿真效果

本文来自互联网用户投稿，该文观点仅代表作者本人，不代表本站立场。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如若转载，请注明出处：http://www.coloradmin.cn/o/1883032.html

如若内容造成侵权/违法违规/事实不符，请联系多彩编程网进行投诉反馈，一经查实，立即删除！

相关文章

Docker 一篇到位

Docker 一篇到位

目录 01. Docker使用导航 02. Build Share Run 样例 03. 理解容器 04. 安装 Docker 05. Docker 样例&#xff08;常见命令使用&#xff09; 下载镜像 启动容器 修改页面 保存镜像 docker commit docker save docker load 分享社区 docker login docker tag do…
阅读更多...
浅聊JavaScript中的栈（stack）

浅聊JavaScript中的栈（stack）

前言 这篇文章结合leetcode题目讲解一下栈这种结构 第20题&#xff1a;20. 有效的括号 - 力扣&#xff08;LeetCode&#xff09; 第LCR-147题LCR 147. 最小栈 - 力扣&#xff08;LeetCode&#xff09; 栈&#xff08;stack&#xff09; 在讲解题目之前&#xff0c;我们先来…
阅读更多...
天池大赛Higress插件官方demo详细部署+调试

天池大赛Higress插件官方demo详细部署+调试

天池大赛Higress插件官方demo详细部署调试 契机 ⚙ 使用Higress AI网关优化AI调用成本。就是基于向量召回相似问题的缓存&#xff0c;降低LLM API调用成本。就是开发一个网关插件做QA缓存嘛。前文已经成功复现了hello-world插件&#xff0c;这次结合官方提供的AI-Cache插件自…
阅读更多...
二叉树遍历练习题

二叉树遍历练习题

2.已知某二叉树的前序遍历序列为5 7 4 9 6 2 1&#xff0c;中序遍历序列为4 7 5 6 9 1 2&#xff0c;则其后序遍历序列为&#xff08; &#xff09; A.4 2 5 7 6 9 1 B.4 2 7 5 6 9 1 C.4 7 6 1 2 9 5 D.4 7 2 9 5 6 1 答案&#xff1a;C 解析&#xff1a; 通过前序遍历找…
阅读更多...
失眠焦虑植物神经紊乱应该怎么改善饮食？

失眠焦虑植物神经紊乱应该怎么改善饮食？

在快节奏的现代社会中&#xff0c;越来越多的人受到植物神经紊乱的困扰&#xff0c;尤其是失眠、焦虑、胸闷气短等症状频发。这些症状不仅影响日常生活和工作效率&#xff0c;还可能引发一系列健康问题。今天&#xff0c;我们就来谈谈如何通过调整饮食来改善这些症状。 饮食调整…
阅读更多...
Java [ 基础 ] 方法引用 ✨

Java [ 基础 ] 方法引用 ✨

✨探索Java基础✨ Java基础&#xff1a;方法引用 方法引用是Java 8中引入的一种新特性&#xff0c;它使得代码更加简洁和易读。方法引用提供了一种可以直接引用已有方法作为Lambda表达式的替代方案。本文将深入介绍方法引用的基本概念、使用方法、具体实例及其在实际开发中的…
阅读更多...
HCIA4.26-5.10

HCIA4.26-5.10

OSPF ——开放式最短路径优先协议 无类别链路状态IGP动态路由协议 距离矢量协议 运行距离矢量协议的路由器会周期性的泛洪自己的路由表&#xff0c;通过路由之间的交互&#xff0c;每台路由器都从相邻的路由器学习到路由条目&#xff0c;随后加载进自己的路由表中。对于网络…
阅读更多...
简单分享下python封装

简单分享下python封装

目录&#xff1a; 一、简介&#xff0c;什么是封装 二、封装的优点与好处 三、封装的示例 四、可封装的场景 一、简介&#xff0c;什么是封装 通俗理解&#xff1a;封装&#xff0c;简而言之&#xff0c;就是把数据&#xff08;变量&#xff09;和操作这些数据的方法&#xff0…
阅读更多...
全球AI新闻速递7.1

全球AI新闻速递7.1

全球AI新闻速递 1.科大讯飞发布讯飞星火 V4.0。 2.成都人形机器人创新中心&#xff1a;基于视觉扩散架构的人形机器人任务生成式模型 R-DDPRM。 3.安徽省人形机器人产业创新中心获批&#xff0c;将打造国内首创、世界领先研究基地。 4.亳州牵手华为打造华佗中医药大模型。 …
阅读更多...
系统安全及应用（命令）

系统安全及应用（命令）

目录 一、账号安全控制 1.1 系统账号清理 1.2 密码安全控制 1.3 历史记录控制 1.4 终端自动注销 二、系统引导和登陆控制 2.1 限制su命令用户 2.2 PAM安全认证 示例一&#xff1a;通过pam 模块来防止暴力破解ssh 2.3 sudo机制提升权限 2.3.1 sudo命令&#xff08;ro…
阅读更多...
SRS流媒体服务器安装与推拉流测试

SRS流媒体服务器安装与推拉流测试

一、首先打开SRS Github https://github.com/ossrs/srs二、在Linux系统下安装SRS 1. git clone https://github.com/ossrs/srs2、克隆完后进入trunk文件夹 cd srs/trunk3.进行环境配置 ./configure 检测当前环境是否满足 4.进行编译 make 5.编译完成后运行 ./objs/srs …
阅读更多...
统计学期末名词解释说明

统计学期末名词解释说明

基本名词解释P值、显著性、显著水平、样本量、三大分布等 20基本名词详细解释&#xff1a; √P值&#xff1a; 相关名词&#xff1a;显著性、显著性差异、0.01水平显著、0.05水平显著。 P值&#xff0c;也称显著性值或者Sig.值&#xff0c;用于描述某件事情发生的概率情况&a…
阅读更多...
器件频频更换为哪桩

器件频频更换为哪桩

曾想象&#xff0c;在一家大型研发型企业里有如下案例&#xff1a; 硬件工程师设计电路选择了器件库中的某器件&#xff0c;在批量试产产品时&#xff0c;却发现没有库存&#xff0c;即时申请采购&#xff0c;却发现货期相当长&#xff0c;一时难以采购&#xff0c;甚至根本不…
阅读更多...
【ElementPlus源码】Scrollbar 滚动条

【ElementPlus源码】Scrollbar 滚动条

文章目录 thumbclickThumbHandlerstartDragmouseMoveDocumentHandlermouseUpDocumentHandlerclickTrackHandler其他 barScrollbar导出的方法noresize更新滚动条相关属性 utilsruntime.tsbuildProps 看源码时候做的笔记。若有问题&#xff0c;请指出&#xff01; 路径相关格式请…
阅读更多...
什么是协程？协程和线程的区别

什么是协程？协程和线程的区别

文章目录 前置知识应用程序和内核阻塞和非阻塞同步和异步并发和并行IO 发展历史同步编程异步多线程/进程异步消息 回调函数&#xff08;响应式编程&#xff09; 协程协程基本概念go 示例代码协程和线程的区别 个人简介 前置知识 在了解协程前&#xff0c;我们先理解一些相关的…
阅读更多...
VBA数据库解决方案第十二讲：如何判断数据库中数据表是否存在

VBA数据库解决方案第十二讲：如何判断数据库中数据表是否存在

《VBA数据库解决方案》教程&#xff08;版权10090845&#xff09;是我推出的第二套教程&#xff0c;目前已经是第二版修订了。这套教程定位于中级&#xff0c;是学完字典后的另一个专题讲解。数据库是数据处理的利器&#xff0c;教程中详细介绍了利用ADO连接ACCDB和EXCEL的方法…
阅读更多...
平安养老险陕西分公司参加2024上半年省级单位驻富平帮扶团联席会

平安养老险陕西分公司参加2024上半年省级单位驻富平帮扶团联席会

6月28日&#xff0c;平安养老险陕西分公司工会副主席武媛携驻村工作队赴富平县庄里镇永安村参加2024上半年度省级单位驻富平帮扶团联席会议。 会议由省委金融办副主任、省委金融工委委员李嘉辉及省委金融办选派挂职干部、富平县副县长席玮共同主持。 会上&#xff0c;席玮县长带…
阅读更多...
全球AI新闻速递6.28

全球AI新闻速递6.28

全球AI新闻速递 1.首款 Transformer 专用 AI 芯片 Sohu 登场。 2.钉钉&#xff1a;宣布对所有AI大模型厂商开放&#xff0c;首批7家接入。 3.华为联合清华大学发布《AI 终端白皮书》。 4.国家卫生健康委&#xff1a;推动AI技术在制定个性化营养、运动干预方案中的应用。 …
阅读更多...
地下水电站3D虚拟仿真展示平台

地下水电站3D虚拟仿真展示平台

借助先进的VR技术&#xff0c;我们将水电站的每一个角落、每一处细节都以三维全景的形式真实呈现。您可以自由穿梭于水电站的各个区域&#xff0c;无论是发电机组、巍峨的水坝&#xff0c;还是错综复杂的输水管道&#xff0c;都近在咫尺。感受水流的澎湃力量&#xff0c;聆听机…
阅读更多...
中文TeX，各种数学符号和表格

中文TeX，各种数学符号和表格

\documentclass{article} \usepackage{amsmath,amssymb,amsfonts} \usepackage{CJKutf8} \begin{document}\begin{CJK}{UTF8}{gkai}%正文放在此行下与\end{CJK}之间就行你好, LaTeX!平方根 $\sqrt{x}$立方根 $\sqrt[3]{x}$分数的代码是 $\frac{a}{b}$求和的代码是 $\sum_{i1}^{…
阅读更多...
推荐文章
最新文章

Copyright @ 2022~2023