【无人机控制】基于MATLAB的欠驱动无人机控制算法仿真，针对四旋翼飞行器。该框架比较了激进轨迹、测量噪声和外部干扰下的几何控制、微分平坦度的控制、INDI 和 NMPC
✅作者简介：热爱科研的Matlab仿真开发者，擅长毕业设计辅导、数学建模、数据处理、建模仿真、程序设计、完整代码获取、论文复现及科研仿真。
🍎 往期回顾关注个人主页：Matlab科研工作室
👇 关注我领取海量matlab电子书和数学建模资料
🍊个人信条：做科研，博学之、审问之、慎思之、明辨之、笃行之，是为：博学慎思，明辨笃行。
🔥 内容介绍
一、研究背景与仿真框架设计
欠驱动四旋翼无人机是典型的四输入六自由度欠驱动系统，凭借结构简单、机动灵活的优势，在巡检、测绘、竞速等场景得到广泛应用。在激进轨迹跟踪、传感器测量噪声、外部强风扰等复杂工况下，不同控制算法的性能差异极大，传统基于欧拉角的PID控制在高动态场景下极易出现姿态失稳、轨迹跟踪发散的问题，完全无法适配高复杂度作业需求。
本研究搭建统一的MATLAB/Simulink仿真框架，针对欠驱动四旋翼的全非线性动力学模型，在完全一致的硬件参数、仿真步长、测试工况下，系统性对比几何控制、微分平坦度控制、INDI（增量非线性动态逆）、NMPC（非线性模型预测控制）四类主流先进控制算法的性能差异，覆盖激进轨迹跟踪、测量噪声干扰、外部风扰三大核心测试场景，量化不同算法的适用边界，为工程场景下的控制器选型提供精准的量化依据。
仿真框架采用统一的四旋翼物理参数配置：整机质量1.2kg，旋翼臂长0.25m，转动惯量对角矩阵为diag(0.0023, 0.0023, 0.0042) kg·m²，所有算法的仿真步长统一设置为1ms，传感器采样频率统一为200Hz，从根源上消除仿真条件差异带来的性能对比偏差，保证测试结果的公平性与可复现性。
二、四类先进控制算法核心原理实现
1. 几何控制算法
几何控制直接在特殊欧氏群SE(3)上推导控制律，完全不使用欧拉角或四元数的局部参数化表示，从根源上规避了传统欧拉角控制的奇异点问题，非常适合大角度机动场景。控制律直接基于四旋翼的旋转矩阵、位置、速度误差设计，将位置环输出的期望合力映射到机体坐标系，通过姿态误差的动态补偿直接生成期望力矩，实现全姿态的稳定跟踪。该算法无需对动力学模型做线性化近似，在大角度机动场景下依然可以保证全局渐近稳定性，代码实现简洁，计算量极低，非常适合机载嵌入式平台部署。
2. 微分平坦度控制
四旋翼系统已被严格证明是微分平坦系统，选择位置坐标与偏航角作为平坦输出，系统所有的状态量与输入量都可以表示为平坦输出及其各阶导数的代数函数，无需求解复杂的微分方程即可完成轨迹规划与控制律设计。在轨迹跟踪场景下，直接通过期望轨迹的位置、速度、加速度、加加速度、加加加速度的高阶微分，即可直接解算出四旋翼所需的总推力与期望姿态角，大幅简化轨迹跟踪的控制逻辑，对平滑连续的轨迹具备天然的优秀跟踪性能。
3. INDI增量非线性动态逆控制
INDI算法是在传统非线性动态逆的基础上改进得到的增量式控制方法，不直接依赖精确的系统动力学模型，而是利用传感器实时测量得到的当前状态量，计算控制量的增量，仅需要非常粗略的模型先验信息即可实现高精度控制。该算法对模型参数不确定性的鲁棒性极强，即使四旋翼的升力系数、转动惯量等参数出现较大偏差，依然可以维持稳定的控制性能，非常适合模型参数难以精准标定的工程场景。
4. NMPC非线性模型预测控制
NMPC算法在每个控制周期内，基于四旋翼的全非线性动力学模型，在满足所有物理约束的前提下，在线求解有限时间域内的最优控制序列，直接将无人机的最大角速度、最大加速度、电机转速饱和等物理约束显式嵌入优化问题，从根源上避免控制量超出执行机构的饱和边界，在多约束场景下可以生成最优的控制输出，是目前理论上性能上限最高的先进控制算法之一。
⛳️ 运行结果
📣 部分代码
%% plot_paper_figures.m
% Paper-style figure entry point. This file reads saved benchmark results and
% exports editable PDF plus high-resolution PNG without rerunning simulations.
%
% Edit the settings in Sections 1-5 for paper layout, colors, line widths,
% labels, legends, and output size.
clear; clc; close all;
%% ========================================================================
%% 1. What to export
figuresToMake = ["disturbance_monte_carlo", "trajectory_sweep_3d"];
outputDir = fullfile(pwd, "results", "paper_figures");
savePdf = true;
savePng = true;
%% ========================================================================
%% 2. Global paper style
P.fontName = "Times New Roman";
P.fontSize = 8;
P.axisLineWidth = 0.5;
P.gridAlpha = 0.18;
P.titleFontSize = 9;
P.labelFontSize = 8;
P.legendFontSize = 7;
P.export.pdfResolution = 1200;
P.export.pngResolution = 600;
P.export.pdfContentType = "vector"; % set "image" if vector PDF is slow
P.export.pngContentType = "image";
P.export.padding = "tight";
%% ========================================================================
%% 3. Controller display style
P.controllerOrder = ["lee", "johnson", "lu", ...
"sun_dfbc", "sun_dfbc_indi", "sun_nmpc", "sun_nmpc_indi", ...
"tal", "geometric_indi"];
P.controllerLabel.lee = "Geometric Control";
P.controllerLabel.johnson = "Log-Geometric Control";
P.controllerLabel.lu = "On-Manifold MPC";
P.controllerLabel.sun_dfbc = "DFBC";
P.controllerLabel.sun_dfbc_indi = "DFBC+INDI";
P.controllerLabel.sun_nmpc = "NMPC";
P.controllerLabel.sun_nmpc_indi = "NMPC+INDI";
P.controllerLabel.tal = "DF-based INDI";
P.controllerLabel.geometric_indi = "Geometric INDI";
P.controllerColor.lee = [0.000, 0.447, 0.741];
P.controllerColor.johnson = [0.850, 0.325, 0.098];
P.controllerColor.lu = [0.929, 0.694, 0.125];
P.controllerColor.sun_dfbc = [0.494, 0.184, 0.556];
P.controllerColor.sun_dfbc_indi = [0.494, 0.184, 0.556];
P.controllerColor.sun_nmpc = [0.466, 0.674, 0.188];
P.controllerColor.sun_nmpc_indi = [0.466, 0.674, 0.188];
P.controllerColor.tal = [0.301, 0.745, 0.933];
P.controllerColor.geometric_indi = [0.635, 0.078, 0.184];
P.controllerLineStyle.lee = "-";
P.controllerLineStyle.johnson = "-";
P.controllerLineStyle.lu = "-";
P.controllerLineStyle.sun_dfbc = "--";
P.controllerLineStyle.sun_dfbc_indi = "-";
P.controllerLineStyle.sun_nmpc = "--";
P.controllerLineStyle.sun_nmpc_indi = "-";
P.controllerLineStyle.tal = "-";
P.controllerLineStyle.geometric_indi = "-";
%% ========================================================================
%% 4. Disturbance Monte Carlo boxchart figure
P.mc.resultDir = fullfile(pwd, "results", "disturbance_monte_carlo", ...
"random_gust");
P.mc.trajNames = ["figure8_horizontal", "helix_flip"];
P.mc.controllerNames = strings(0,1); % empty means all saved controllers
P.mc.levelNames = ["low", "medium", "high"];
P.mc.levelLabel.low = "Low disturbance";
P.mc.levelLabel.medium = "Medium disturbance";
P.mc.levelLabel.high = "High disturbance";
P.mc.failureRmseThreshold = 5; % [m], set inf to keep all finite runs
P.mc.figureWidthCm = 16;
P.mc.figureHeightCm = 9;
P.mc.boxWidth = 0.10;
P.mc.boxLineWidth = 0.55;
P.mc.boxFaceAlpha = 0.25;
P.mc.medianColor = [0.80, 0.10, 0.10];
P.mc.outlierMarker = "o";
P.mc.outlierMarkerSize = 3.5;
P.mc.legendLocation = "northwest";
P.mc.legendNumColumns = 1;
P.mc.showLegend = true;
P.mc.yLabel = "RMS position tracking error (m)";
P.mc.xLabel = "";
P.mc.title.figure8_horizontal = "Horizontal figure-eight";
P.mc.title.helix_flip = "Helix flip";
P.mc.yLimit.figure8_horizontal = [];
P.mc.yLimit.helix_flip = [];
P.mc.disturbanceLabelFontSize = 6.5;
%% ========================================================================
%% 5. Trajectory sweep 3-D figure
P.sweep.resultFile = "latest"; % "latest" or a .mat file path
P.sweep.trajNames = ["figure8_horizontal", "helix_flip"];
P.sweep.figureWidthCm = 16;
P.sweep.figureHeightCm = 7;
P.sweep.tileSpacing = "compact";
P.sweep.padding = "compact";
P.sweep.actualColor = [0.000, 0.447, 0.741];
P.sweep.refColor = [0.850, 0.325, 0.098];
P.sweep.actualLineStyle = "-";
P.sweep.refLineStyle = "--";
P.sweep.actualLineWidth = 0.9;
P.sweep.refLineWidth = 0.8;
P.sweep.viewAzEl = [35, 25];
P.sweep.axisEqual = true;
P.sweep.legendLocation = "northeast";
P.sweep.showLegend = true;
P.sweep.xLabel = "x (m)";
P.sweep.yLabel = "y (m)";
P.sweep.zLabel = "z (m)";
P.sweep.title.figure8_horizontal = "Horizontal figure-eight";
P.sweep.title.helix_flip = "Helix flip";
%% ========================================================================
%% 6. Run selected paper plots
P.outputDir = outputDir;
P.savePdf = savePdf;
P.savePng = savePng;
applyPaperDefaults(P);
if any(figuresToMake == "disturbance_monte_carlo")
makePaperDisturbanceBoxFigures(P);
end
if any(figuresToMake == "trajectory_sweep_3d")
makePaperTrajectorySweepFigure(P);
end
fprintf('Paper figures saved to:\n %s\n', outputDir);
%% ========================================================================
%% Local plotting functions
function applyPaperDefaults(P)
set(groot, ...
'defaultAxesFontSize', P.fontSize, ...
'defaultAxesFontName', char(P.fontName), ...
'defaultAxesLineWidth', P.axisLineWidth, ...
'defaultAxesLabelFontSizeMultiplier', 1, ...
'defaultAxesTitleFontSizeMultiplier', 1, ...
'defaultTextFontName', char(P.fontName), ...
'defaultLegendFontName', char(P.fontName));
end
function makePaperDisturbanceBoxFigures(P)
resultFile = fullfile(P.mc.resultDir, "disturbance_benchmark_results.mat");
if ~isfile(resultFile)
error("Monte Carlo result not found: %s.", resultFile);
end
S = load(resultFile, 'results');
results = S.results;
trajNames = selectedNames(P.mc.trajNames, unique(results.Trajectory, ...
'stable'));
controllerNames = selectedNames(P.mc.controllerNames, ...
unique(results.Controller, 'stable'));
controllerOrder = appendMissing(P.controllerOrder, controllerNames);
levelNames = selectedNames(P.mc.levelNames, ...
unique(results.DisturbanceLevel, 'stable'));
levelLabels = disturbanceLevelLabels(levelNames, P.mc.levelLabel);
outDir = fullfile(P.outputDir, "disturbance_monte_carlo");
ensureDir(outDir);
for iTraj = 1:numel(trajNames)
trajName = trajNames(iTraj);
mask = string(results.Trajectory) == trajName ...
& ismember(string(results.Controller), controllerNames) ...
& ismember(string(results.DisturbanceLevel), levelNames) ...
& results.IsFinite;
if isfinite(P.mc.failureRmseThreshold)
mask = mask & isfinite(results.RMSE) ...
& results.RMSE <= P.mc.failureRmseThreshold;
end
if ~any(mask)
warning("No Monte Carlo rows selected for trajectory %s.", trajName);
continue;
end
fig = paperFigure(P.mc.figureWidthCm, P.mc.figureHeightCm);
ax = axes(fig);
hold(ax, 'on');
grid(ax, 'on');
ax.GridAlpha = P.gridAlpha;
legendHandles = gobjects(0);
legendLabels = strings(0,1);
nLevels = numel(levelNames);
nCtrl = numel(controllerOrder);
offsets = groupedOffsets(nCtrl, 0.70);
for iLevel = 1:nLevels
for iCtrl = 1:nCtrl
ctrlName = controllerOrder(iCtrl);
rowMask = mask ...
& string(results.DisturbanceLevel) == levelNames(iLevel) ...
& string(results.Controller) == ctrlName;
y = results.RMSE(rowMask);
y = y(isfinite(y));
if isempty(y)
continue;
end
x = ones(numel(y),1)*(iLevel + offsets(iCtrl));
color = controllerColor(P, ctrlName, iCtrl);
h = boxchart(ax, x, y, ...
'BoxFaceColor', color, ...
'BoxFaceAlpha', P.mc.boxFaceAlpha, ...
'BoxEdgeColor', color, ...
'BoxMedianLineColor', P.mc.medianColor, ...
'WhiskerLineColor', [0.15 0.15 0.15], ...
'MarkerStyle', char(P.mc.outlierMarker), ...
'MarkerColor', color, ...
'MarkerSize', P.mc.outlierMarkerSize, ...
'BoxWidth', P.mc.boxWidth, ...
'LineWidth', P.mc.boxLineWidth);
if iLevel == 1
legendHandles(end+1,1) = h; %#ok<AGROW>
legendLabels(end+1,1) = controllerLabel(P, ctrlName); %#ok<AGROW>
end
end
end
xlim(ax, [0.5, nLevels + 0.5]);
xticks(ax, 1:nLevels);
ylabel(ax, P.mc.yLabel, 'FontSize', P.labelFontSize);
xlabel(ax, P.mc.xLabel, 'FontSize', P.labelFontSize);
title(ax, titleFor(P.mc.title, trajName), ...
'FontSize', P.titleFontSize, 'Interpreter', 'none');
formatAxes(ax, P);
applyOptionalYLimit(ax, P.mc.yLimit, trajName);
addDisturbanceLevelLabels(ax, levelLabels, P.mc.disturbanceLabelFontSize);
xticklabels(ax, repmat({''}, 1, nLevels));
if P.mc.showLegend && ~isempty(legendHandles)
legend(ax, legendHandles, cellstr(legendLabels), ...
'Location', char(P.mc.legendLocation), ...
'NumColumns', P.mc.legendNumColumns, ...
'FontSize', P.legendFontSize, ...
'Interpreter', 'none');
end
savePaperFigure(fig, fullfile(outDir, "mc_" + trajName), ...
P.mc.figureWidthCm, P.mc.figureHeightCm, P);
end
end
function makePaperTrajectorySweepFigure(P)
resultFile = string(P.sweep.resultFile);
if resultFile == "latest"
resultFile = latestSweepResultsFile();
end
if ~isfile(resultFile)
error("Trajectory sweep result not found: %s.", resultFile);
end
S = load(resultFile, 'results');
results = S.results;
trajNames = selectedNames(P.sweep.trajNames, unique(results.Trajectory, ...
'stable'));
data = loadSweepData(results, trajNames);
if isempty(data)
error("No saved trajectory run data found for the selected sweep.");
end
outDir = fullfile(P.outputDir, "trajectory_sweep");
ensureDir(outDir);
fig = paperFigure(P.sweep.figureWidthCm, P.sweep.figureHeightCm);
tl = tiledlayout(fig, 1, numel(data), ...
'TileSpacing', char(P.sweep.tileSpacing), ...
'Padding', char(P.sweep.padding));
for i = 1:numel(data)
ax = nexttile(tl);
p = data(i).log.p;
pd = data(i).log.pd;
plot3(ax, p(1,:), p(2,:), p(3,:), ...
'LineStyle', char(P.sweep.actualLineStyle), ...
'Color', P.sweep.actualColor, ...
'LineWidth', P.sweep.actualLineWidth, ...
'DisplayName', 'actual');
hold(ax, 'on');
plot3(ax, pd(1,:), pd(2,:), pd(3,:), ...
'LineStyle', char(P.sweep.refLineStyle), ...
'Color', P.sweep.refColor, ...
'LineWidth', P.sweep.refLineWidth, ...
'DisplayName', 'reference');
grid(ax, 'on');
ax.GridAlpha = P.gridAlpha;
if P.sweep.axisEqual
axis(ax, 'equal');
end
view(ax, P.sweep.viewAzEl(1), P.sweep.viewAzEl(2));
set(ax, 'ZDir', 'reverse');
xlabel(ax, P.sweep.xLabel, 'FontSize', P.labelFontSize);
ylabel(ax, P.sweep.yLabel, 'FontSize', P.labelFontSize);
zlabel(ax, P.sweep.zLabel, 'FontSize', P.labelFontSize);
title(ax, titleFor(P.sweep.title, data(i).name), ...
'FontSize', P.titleFontSize, 'Interpreter', 'none');
setTrajectoryLimits(ax, [p, pd]);
formatAxes(ax, P);
if P.sweep.showLegend && i == numel(data)
legend(ax, 'Location', char(P.sweep.legendLocation), ...
'FontSize', P.legendFontSize, 'Interpreter', 'none');
end
end
savePaperFigure(fig, fullfile(outDir, "trajectory_sweep_3d"), ...
P.sweep.figureWidthCm, P.sweep.figureHeightCm, P);
end
function fig = paperFigure(widthCm, heightCm)
fig = figure('Color', 'w', 'Units', 'centimeters', ...
'Position', [0, 0, widthCm, heightCm], ...
'PaperPositionMode', 'auto', ...
'InvertHardcopy', 'off');
end
function savePaperFigure(fig, stem, widthCm, heightCm, P)
[folder, ~, ~] = fileparts(stem);
ensureDir(folder);
set(fig, 'Units', 'centimeters');
set(fig, 'Position', [0, 0, widthCm, heightCm]);
drawnow;
if P.savePdf
exportgraphics(fig, stem + ".pdf", ...
'ContentType', P.export.pdfContentType, ...
'BackgroundColor', 'white', ...
'Resolution', P.export.pdfResolution, ...
'Width', widthCm, ...
'Height', heightCm, ...
'Padding', P.export.padding, ...
'Units', 'centimeters');
end
if P.savePng
exportgraphics(fig, stem + ".png", ...
'ContentType', P.export.pngContentType, ...
'BackgroundColor', 'white', ...
'Resolution', P.export.pngResolution, ...
'Width', widthCm, ...
'Height', heightCm, ...
'Padding', P.export.padding, ...
'Units', 'centimeters');
end
end
function data = loadSweepData(results, trajNames)
data = struct('name', {}, 'time', {}, 'log', {}, 'par', {});
for i = 1:numel(trajNames)
idx = find(string(results.Trajectory) == trajNames(i), 1);
if isempty(idx) || strlength(string(results.ErrorMessage(idx))) > 0
continue;
end
matFile = string(results.MatFile(idx));
if ~isfile(matFile)
continue;
end
S = load(matFile, 'time', 'log', 'par');
S.name = trajNames(i);
data(end+1) = S; %#ok<AGROW>
end
end
function resultFile = latestSweepResultsFile()
rootDir = fullfile(pwd, "results", "main_trajectory_sweep");
files = dir(fullfile(rootDir, "*", "main_trajectory_sweep_results.mat"));
if isempty(files)
error("No trajectory sweep result found under %s.", rootDir);
end
[~, idx] = max([files.datenum]);
resultFile = string(fullfile(files(idx).folder, files(idx).name));
end
function names = selectedNames(value, defaultNames)
names = string(value);
names = names(strlength(names) > 0);
if isempty(names)
names = string(defaultNames);
end
names = names(:).';
end
function values = appendMissing(values, observed)
values = string(values(:)).';
observed = string(observed(:)).';
values = [values, observed(~ismember(observed, values))];
values = values(ismember(values, observed));
end
function offsets = groupedOffsets(n, width)
if n <= 1
offsets = 0;
else
offsets = linspace(-width/2, width/2, n);
end
end
function color = controllerColor(P, ctrlName, idx)
field = matlab.lang.makeValidName(char(ctrlName));
if isfield(P.controllerColor, field)
color = P.controllerColor.(field);
else
color = lines(max(idx, 1));
color = color(end,:);
end
end
function label = controllerLabel(P, ctrlName)
field = matlab.lang.makeValidName(char(ctrlName));
if isfield(P.controllerLabel, field)
label = string(P.controllerLabel.(field));
else
label = string(ctrlName);
end
end
function txt = titleFor(titleStruct, name)
field = matlab.lang.makeValidName(char(name));
if isfield(titleStruct, field)
txt = string(titleStruct.(field));
else
txt = string(name);
end
end
function applyOptionalYLimit(ax, yLimits, trajName)
field = matlab.lang.makeValidName(char(trajName));
if isstruct(yLimits) && isfield(yLimits, field) ...
&& ~isempty(yLimits.(field))
ylim(ax, yLimits.(field));
end
end
function labels = disturbanceLevelLabels(levelNames, labelStruct)
labels = string(levelNames);
for i = 1:numel(labels)
field = matlab.lang.makeValidName(char(labels(i)));
if isstruct(labelStruct) && isfield(labelStruct, field)
labels(i) = string(labelStruct.(field));
end
end
end
function addDisturbanceLevelLabels(ax, labels, fontSize)
fig = ancestor(ax, 'figure');
n = numel(labels);
if n == 0
return;
end
ax.Units = 'normalized';
pos = ax.Position;
bottomPad = 0.10;
ax.Position = [pos(1), pos(2) + bottomPad, ...
pos(3), max(0.1, pos(4) - bottomPad)];
boxW = 0.92*pos(3)/n;
boxH = 0.05;
boxY = max(0.01, pos(2) + 0.025);
for i = 1:n
xCenter = pos(1) + pos(3)*(i - 0.5)/n;
annotation(fig, 'textbox', ...
[xCenter - boxW/2, boxY, boxW, boxH], ...
'String', char(labels(i)), ...
'HorizontalAlignment', 'center', ...
'VerticalAlignment', 'top', ...
'Interpreter', 'none', ...
'EdgeColor', 'none', ...
'FitBoxToText', 'off', ...
'FontSize', fontSize);
end
end
function setTrajectoryLimits(ax, pAll)
span = max(max(pAll, [], 2) - min(pAll, [], 2));
margin = 0.1*max(span, 1);
xlim(ax, [min(pAll(1,:))-margin, max(pAll(1,:))+margin]);
ylim(ax, [min(pAll(2,:))-margin, max(pAll(2,:))+margin]);
zlim(ax, [min(pAll(3,:))-margin, max(pAll(3,:))+margin]);
end
function formatAxes(ax, P)
ax.FontName = char(P.fontName);
ax.FontSize = P.fontSize;
ax.LineWidth = P.axisLineWidth;
ax.TickLabelInterpreter = 'none';
box(ax, 'on');
end
function ensureDir(folder)
if strlength(string(folder)) == 0
return;
end
if ~exist(folder, 'dir')
mkdir(folder);
end
end
🔗 参考文献
[1] T. Lee, M. Leok, and N. H. McClamroch, “Geometric Tracking Control of a Quadrotor UAV on SE(3),” Mar. 10, 2010, arXiv: arXiv:1003.2005. doi: 10.48550/arXiv.1003.2005.
[2] T. Lee, M. Leok, and N. H. McClamroch, “Control of Complex Maneuvers for a Quadrotor UAV using Geometric Methods on SE(3),” Nov. 12, 2010, arXiv: arXiv:1003.2005. doi: 10.48550/arXiv.1003.2005.
[3] F. Bullo and R. M. Murray, “Proportional derivative (PD) control on the Euclidean group,” 1995.
[4] Y. Yu, S. Yang, M. Wang, C. Li, and Z. Li, “High performance full attitude control of a quadrotor on SO (3),” in 2015 IEEE International Conference on Robotics and Automation (ICRA), Seattle, WA, USA: IEEE, 2015, pp. 1698–1703. doi: 10.1109/icra.2015.7139416.
[5] D. Brescianini, M. Hehn, and R. D’Andrea, “Nonlinear Quadrocopter Attitude Control: Technical Report,” ETH Zurich, 2013. doi: 10.3929/ETHZ-A-009970340.
[6] M. Faessler, A. Franchi, and D. Scaramuzza, “Differential Flatness of Quadrotor Dynamics Subject to Rotor Drag for Accurate Tracking of High-Speed Trajectories,” IEEE Robotics and Automation Letters, vol. 3, no. 2, pp. 620–626, Apr. 2018, doi: 10.1109/LRA.2017.2776353.
[7] J. Sola, “Quaternion kinematics for the error-state Kalman filter,” arXiv:1711.02508 [cs], Nov. 2017, Accessed: Sep. 26, 2020. [Online]. Available: http://arxiv.org/abs/1711.02508
[8] D. Brescianini and R. D’Andrea, “Tilt-Prioritized Quadrocopter Attitude Control,” IEEE Transactions on Control Systems Technology, vol. 28, no. 2, pp. 376–387, Mar. 2020, doi: 10.1109/TCST.2018.2873224.
[9] J. Johnson and R. Beard, “Globally-Attractive Logarithmic Geometric Control of a Quadrotor for Aggressive Trajectory Tracking,” Dec. 01, 2021, arXiv: arXiv:2109.07025. doi: 10.48550/arXiv.2109.07025.
[10] J. Sola, J. Deray, and D. Atchuthan, “A micro Lie theory for state estimation in robotics,” Dec. 08, 2021, arXiv: arXiv:1812.01537. doi: 10.48550/arXiv.1812.01537.
[11] E. Tal and S. Karaman, “Accurate Tracking of Aggressive Quadrotor Trajectories Using Incremental Nonlinear Dynamic Inversion and Differential Flatness,” IEEE Trans. Contr. Syst. Technol., vol. 29, no. 3, pp. 1203–1218, May 2021, doi: 10.1109/tcst.2020.3001117.
[12] S. Sun, A. Romero, P. Foehn, E. Kaufmann, and D. Scaramuzza, “A Comparative Study of Nonlinear MPC and Differential-Flatness-Based Control for Quadrotor Agile Flight,” Feb. 23, 2022, arXiv: arXiv:2109.01365. Accessed: May 27, 2022. [Online]. Available: http://arxiv.org/abs/2109.01365
[13] G. Lu, W. Xu, and F. Zhang, “On-Manifold Model Predictive Control for Trajectory Tracking on Robotic Systems,” IEEE Transactions on Industrial Electronics, vol. 70, no. 9, pp. 9192–9202, Sep. 2023, doi: 10.1109/TIE.2022.3212397.
🍅更多创新智能优化算法模型和应用场景可扫描关注
🌟机器学习/深度学习类：BP、SVM、RVM、DBN、LSSVM、ELM、KELM、HKELM、DELM、RELM、DHKELM、RF、SAE、LSTM、BiLSTM、GRU、BiGRU、PNN、CNN、XGBoost、LightGBM、TCN、BiTCN、ESN、Transformer、模糊小波神经网络、宽度学习等等均可~
方向涵盖风电预测、光伏预测、电池寿命预测、辐射源识别、交通流预测、负荷预测、股价预测、PM2.5浓度预测、电池健康状态预测、用电量预测、水体光学参数反演、NLOS信号识别、地铁停车精准预测、变压器故障诊断
🌟组合预测类：CNN/TCN/BiTCN/DBN/Transformer/Adaboost结合SVM、RVM、ELM、LSTM、BiLSTM、GRU、BiGRU、Attention机制类等均可（可任意搭配非常新颖）~
🌟分解类：EMD、EEMD、VMD、REMD、FEEMD、TVFEMD、CEEMDAN、ICEEMDAN、SVMD、FMD、JMD等分解模型均可~
🌟路径规划类：旅行商问题（TSP）、车辆路径问题（VRP、MVRP、CVRP、VRPTW等）、无人机三维路径规划、无人机协同、无人机编队、机器人路径规划、栅格地图路径规划、多式联运运输问题、 充电车辆路径规划（EVRP）、 双层车辆路径规划（2E-VRP）、 油电混合车辆路径规划、 船舶航迹规划、 全路径规划规划、 仓储巡逻、公交车时间调度、水库调度优化、多式联运优化等等~
🌟小众优化类：生产调度、经济调度、装配线调度、充电优化、车间调度、发车优化、水库调度、三维装箱、物流选址、货位优化、公交排班优化、充电桩布局优化、车间布局优化、集装箱船配载优化、水泵组合优化、解医疗资源分配优化、设施布局优化、可视域基站和无人机选址优化、背包问题、 风电场布局、时隙分配优化、 最佳分布式发电单元分配、多阶段管道维修、 工厂-中心-需求点三级选址问题、 应急生活物质配送中心选址、 基站选址、 道路灯柱布置、 枢纽节点部署、 输电线路台风监测装置、 集装箱调度、 机组优化、 投资优化组合、云服务器组合优化、 天线线性阵列分布优化、CVRP问题、VRPPD问题、多中心VRP问题、多层网络的VRP问题、多中心多车型的VRP问题、 动态VRP问题、双层车辆路径规划（2E-VRP）、充电车辆路径规划（EVRP）、油电混合车辆路径规划、混合流水车间问题、 订单拆分调度问题、 公交车的调度排班优化问题、航班摆渡车辆调度问题、选址路径规划问题、港口调度、港口岸桥调度、停机位分配、机场航班调度、泄漏源定位、冷链、时间窗、多车场等、选址优化、港口岸桥调度优化、交通阻抗、重分配、停机位分配、机场航班调度、通信上传下载分配优化、微电网优化、无功优化、配电网重构、储能配置、有序充电、MPPT优化、家庭用电、电/冷/热负荷预测、电力设备故障诊断、电池管理系统（BMS）SOC/SOH估算（粒子滤波/卡尔曼滤波）、 多目标优化在电力系统调度中的应用、光伏MPPT控制算法改进（扰动观察法/电导增量法）、电动汽车充放电优化、微电网日前日内优化、储能优化、家庭用电优化、供应链优化\智能电网分布式能源经济优化调度，虚拟电厂，能源消纳，风光出力，控制策略，多目标优化，博弈能源调度，鲁棒优化等等均可~
🌟 无人机应用方面：无人机路径规划、无人机控制、无人机编队、无人机协同、无人机任务分配、无人机安全通信轨迹在线优化、车辆协同无人机路径规划
🌟通信方面：传感器部署优化、通信协议优化、路由优化、目标定位优化、Dv-Hop定位优化、Leach协议优化、WSN覆盖优化、组播优化、RSSI定位优化、水声通信、通信上传下载分配
🌟信号处理方面：信号识别、信号加密、信号去噪、信号增强、雷达信号处理、信号水印嵌入提取、肌电信号、脑电信号、信号配时优化、心电信号、DOA估计、编码译码、变分模态分解、管道泄漏、滤波器、数字信号处理+传输+分析+去噪、数字信号调制、误码率、信号估计、DTMF、信号检测
🌟电力系统方面： 微电网优化、无功优化、配电网重构、储能配置、有序充电、MPPT优化、家庭用电、电/冷/热负荷预测、电力设备故障诊断、电池管理系统（BMS）SOC/SOH估算（粒子滤波/卡尔曼滤波）、 多目标优化在电力系统调度中的应用、光伏MPPT控制算法改进（扰动观察法/电导增量法）、电动汽车充放电优化、微电网日前日内优化、储能优化、家庭用电优化、供应链优化\智能电网分布式能源经济优化调度，虚拟电厂，能源消纳，风光出力，控制策略，多目标优化，博弈能源调度，鲁棒优化
🌟原创改进优化算法（适合需要创新的同学）：原创改进2025年的波动光学优化算法WOO以及三国优化算法TKOA、白鲸优化算法BWO等任意优化算法均可，保证测试函数效果，一般可直接核心