简介:本资源是一套面向控制工程与无人机方向初学者及进阶学习者的四旋翼滑模控制MATLAB仿真完整实现,聚焦非线性系统鲁棒控制设计痛点,适用于课程设计、毕业设计及科研原型验证。压缩包共7个文件(5个.m主程序脚本、1个Simulink模型.mdl、1个备份.asv),总大小仅13KB,轻量但结构完整:包含动力学plant模型、双层滑模控制器(ctrl1/ctrl2)、积分器int模块及结果可视化plot脚本,Simulink模型直观呈现闭环控制架构,各m文件分工明确、注释清晰,便于理解滑模函数设计、边界层抖振抑制及姿态跟踪实现逻辑。目前已有1127人学习下载,读者可直接运行复现四旋翼俯仰/横滚/偏航角响应曲线、控制输入时序图及状态收敛过程,快速掌握滑模控制在实际飞控中的建模—设计—仿真—分析全流程。
1. 四旋翼滑模控制:从理论到仿真的完整闭环
如果你正在做四旋翼飞行器的控制算法研究,尤其是想尝试一些比经典PID更“硬核”的现代控制方法,滑模控制(Sliding Mode Control, SMC)大概率会出现在你的备选清单里。它以其对模型不确定性和外部干扰的强鲁棒性而闻名,听起来就像是专为四旋翼这种强耦合、非线性、易受扰的系统量身定做的。但理论上的美好,往往在第一步仿真时就可能碰壁:状态方程怎么建?滑模面参数怎么调?抖振怎么抑制?Simulink模型怎么搭才能既清晰又高效?最后那一堆Plot图,又该怎么解读才能证明你的控制器真的有效,而不是在自欺欺人?
我花了相当一段时间,才把这一整套流程从零到一跑通。网上能找到的代码和模型往往支离破碎,要么只有理论推导没有程序,要么程序跑出来结果诡异却无从调试。今天,我就把自己搭建的这套四旋翼滑模控制Matlab/Simulink仿真程序的完整思路、核心代码、模型架构以及结果分析,毫无保留地分享出来。这不是一个简单的“Hello World”式Demo,而是一个包含了姿态与位置控制、参数可调、带有抗抖振设计、并能输出完整分析图表的实战项目。无论你是刚开始接触滑模控制的研究生,还是需要在工程中验证算法可行性的工程师,这套东西都能让你少走很多弯路。
2. 四旋翼建模与滑模控制律设计:核心原理拆解
在动笔写代码和拖拽Simulink模块之前,我们必须把控制对象和控制器的数学模型搞清楚。这一步是地基,地基不稳,后面的仿真再漂亮也是空中楼阁。
2.1 四旋翼的动力学模型简化与状态定义
四旋翼是一个典型的欠驱动系统(6个自由度,只有4个输入),其完整的动力学模型相当复杂。为了专注于控制算法本身,我们通常采用基于牛顿-欧拉方程推导的简化模型,并做如下合理假设:
- 机体是刚体且结构对称。
- 重心与几何中心重合。
- 桨叶升力与电机转速的平方成正比。
- 忽略空气阻力矩的高阶项。
在这样的假设下,我们可以将系统分为**姿态环(内环)和位置环(外环)**进行解耦控制。这是最经典也最实用的级联控制结构。
我们定义以下状态变量:
- 位置与速度:
P = [x, y, z]^T,V = [v_x, v_y, v_z]^T,在地面惯性坐标系下。 - 姿态与角速度:
Θ = [φ, θ, ψ]^T(滚转、俯仰、偏航角),Ω = [p, q, r]^T(机体坐标系下的角速度)。
系统的控制输入是四个电机的转速,但更直接的控制量是它们组合产生的总升力 U1和三个力矩 U2, U3, U4(分别对应滚转、俯仰、偏航)。
经过推导,可以得到如下形式的动力学方程:
位置通道:m * ddot(P) = R * [0; 0; U1] - [0; 0; m*g]其中,m是质量,g是重力加速度,R是由姿态角Θ构成的旋转矩阵(从机体坐标系到惯性坐标系)。这个方程告诉我们,位置的变化最终是通过调整姿态Θ来改变拉力矢量的方向,从而产生水平方向的分力来实现的。
姿态通道:I * dot(Ω) + Ω × (I * Ω) = [U2; U3; U4]其中,I是机体的惯性张量矩阵,×表示叉乘。这个方程描述了角加速度与施加力矩之间的关系。
在实际仿真编程时,我们会使用这些方程的一阶形式(状态空间方程),并处理好旋转矩阵R和科里奥利项Ω × (I * Ω)的计算。这是后续所有控制器设计的起点。
2.2 滑模控制器的设计思路与关键公式
滑模控制的精髓在于“滑模面”的设计。一旦系统状态被“拉”到这个预设的曲面上,就会沿着它滑向平衡点,并且对外部干扰和参数变化具有不变性。
我们以高度通道(z方向)的控制为例,来具体说明设计过程。假设我们的控制目标是让四旋翼稳定在期望高度z_d。
- 定义误差:
e_z = z - z_d。 - 设计滑模面:我们通常选择线性滑模面,其形式为
s = dot(e) + λ * e。对于高度通道,我们定义:s_z = dot(e_z) + λ_z * e_z = (v_z - dot(z_d)) + λ_z * (z - z_d)。 其中λ_z > 0是一个设计参数,它决定了系统状态在滑模面上的收敛速度。λ越大,收敛越快,但可能带来更大的控制输入。 - 设计趋近律:为了让系统状态从任意初始点到达滑模面
s=0,我们需要设计控制律。最常用的是指数趋近律:dot(s) = -ε * sign(s) - k * s。 其中,sign()是符号函数,ε > 0和k > 0是控制器参数。-k*s项保证状态指数趋近滑模面,-ε*sign(s)项用于克服干扰,保证到达条件。 - 推导控制量:将
s_z的导数dot(s_z)与动力学方程联立。dot(s_z) = dot(v_z) - ddot(z_d) + λ_z * (v_z - dot(z_d))。 从位置通道动力学方程中,我们可以提取出z方向的方程:m * dot(v_z) = U1 * cosφ * cosθ - m*g。 令dot(s_z) = -ε_z * sign(s_z) - k_z * s_z,并将上面的dot(v_z)代入,即可反解出总升力U1的表达式:U1 = (m / (cosφ * cosθ)) * [ g + ddot(z_d) - λ_z*(v_z - dot(z_d)) - ε_z*sign(s_z) - k_z*s_z ]。 这里有一个细节:cosφ * cosθ在分母上,当姿态角接近90度时会导致奇异,这也是为什么四旋翼通常不做大机动翻滚的原因之一。在实际程序中,我们需要对这个值做一个安全限幅。
对于姿态角(φ, θ, ψ)的控制,设计流程完全类似。只是动力学方程变成了姿态通道的方程,控制量变成了力矩U2, U3, U4。最终,我们会得到四个控制量的计算公式,它们都是关于状态误差、滑模面s及其参数的函数。
注意:符号函数
sign(s)是抖振的主要来源。为了减轻抖振,工程上普遍采用饱和函数sat(s/Φ)或连续近似函数(如s/(|s|+δ))来代替它。其中Φ是边界层厚度,δ是一个很小的正数。这属于“准滑模控制”,在牺牲一点点鲁棒性的前提下,极大平滑了控制信号。在我的仿真程序中,就采用了饱和函数sat()来处理。
3. Matlab仿真程序架构:核心脚本与函数详解
有了理论公式,我们就可以开始用Matlab代码将它们实现。我的程序通常由几个核心的.m脚本和函数文件组成,结构清晰,便于修改和调试。
3.1 主仿真脚本main_smc_quadrotor.m
这个脚本是仿真的总指挥,负责设置参数、调用求解器、运行仿真和绘制结果。它的典型结构如下:
%% 四旋翼滑模控制仿真主程序 clear; close all; clc; %% 1. 初始化参数 % 物理参数 params.m = 1.2; % 质量 (kg) params.g = 9.81; % 重力加速度 params.Ixx = 0.023; params.Iyy = 0.023; params.Izz = 0.046; % 转动惯量 (kg*m^2) params.arm_length = 0.225; % 机臂长度 (m) % 滑模控制器参数 % 高度通道 smc_params.z.lambda = 1.5; smc_params.z.epsilon = 0.8; smc_params.z.k = 2.0; smc_params.z.phi = 0.1; % 边界层厚度 % 姿态通道 (以滚转角为例) smc_params.phi.lambda = 8.0; smc_params.phi.epsilon = 5.0; smc_params.phi.k = 10.0; smc_params.phi.phi = 0.05; % 期望轨迹 t_sim = 20; % 仿真时间 % 期望位置:前5秒悬停于(0,0,1),5-15秒沿x轴移动到(2,0,1),最后5秒悬停 t_vec = linspace(0, t_sim, 1000); z_d = 1 * ones(size(t_vec)); x_d = zeros(size(t_vec)); x_d(t_vec>=5 & t_vec<15) = 2 * (t_vec(t_vec>=5 & t_vec<15)-5)/10; % 匀速运动 x_d(t_vec>=15) = 2; % 期望偏航角:保持0度 psi_d = 0 * ones(size(t_vec)); % 初始状态 init_state = [0; 0; 0; ... % x, y, z 0; 0; 0; ... % vx, vy, vz 0; 0; 0; ... % phi, theta, psi (rad) 0; 0; 0]; % p, q, r (rad/s) %% 2. 调用ODE求解器进行仿真 % 将参数打包 sim_data.params = params; sim_data.smc_params = smc_params; sim_data.ref_t = t_vec; sim_data.ref_x = x_d; sim_data.ref_z = z_d; sim_data.ref_psi = psi_d; % 使用ode45求解微分方程 options = odeset('RelTol', 1e-6, 'AbsTol', 1e-9); [t, state_history] = ode45(@(t, x) quadrotor_dynamics_smc(t, x, sim_data), ... [0, t_sim], init_state, options); %% 3. 提取并处理仿真结果 % state_history的每一列对应init_state的一个状态 x_pos = state_history(:,1); y_pos = state_history(:,2); z_pos = state_history(:,3); vx = state_history(:,4); vy = state_history(:,5); vz = state_history(:,6); phi = state_history(:,7); theta = state_history(:,8); psi = state_history(:,9); p = state_history(:,10); q = state_history(:,11); r = state_history(:,12); % 插值得到对应时间点的期望值 x_d_interp = interp1(t_vec, x_d, t); z_d_interp = interp1(t_vec, z_d, t); psi_d_interp = interp1(t_vec, psi_d, t); % 计算控制输入历史(需要在动力学函数中记录) % 假设动力学函数通过全局变量或额外输出返回了控制量U % 这里需要根据你的函数具体实现来获取,例如: global U_history; % 或者将动力学函数修改为返回额外参数 %% 4. 绘制结果图形 (Plot部分将在第5章详细展开) % plot_results(t, state_history, U_history, ref_data);这个脚本的关键在于调用了ode45求解器,而微分方程的具体内容则封装在quadrotor_dynamics_smc这个函数中。
3.2 核心动力学函数quadrotor_dynamics_smc.m
这个函数是仿真的心脏,它根据当前状态和时间,计算状态的导数(dxdt)。它必须严格按照状态方程和滑模控制律来编写。
function dxdt = quadrotor_dynamics_smc(t, x, sim_data) % 输入:t - 当前时间, x - 当前状态向量(12维), sim_data - 包含参数和参考轨迹的结构体 % 输出:dxdt - 状态导数向量 % 解包参数 params = sim_data.params; smc_params = sim_data.smc_params; m = params.m; g = params.g; I = diag([params.Ixx, params.Iyy, params.Izz]); % 解包当前状态 pos = x(1:3); % [x; y; z] vel = x(4:6); % [vx; vy; vz] angles = x(7:9); % [phi; theta; psi] omega = x(10:12); % [p; q; r] phi = angles(1); theta = angles(2); psi = angles(3); p = omega(1); q = omega(2); r = omega(3); % 获取当前时刻的期望值(通过插值) ref_t = sim_data.ref_t; x_d = interp1(ref_t, sim_data.ref_x, t, 'linear', 'extrap'); z_d = interp1(ref_t, sim_data.ref_z, t, 'linear', 'extrap'); psi_d = interp1(ref_t, sim_data.ref_psi, t, 'linear', 'extrap'); % 计算期望值的一阶、二阶导数(这里假设轨迹已知且可微,简单示例中常设为零或根据轨迹公式计算) dot_x_d = 0; ddot_x_d = 0; % 根据你的轨迹设计修改 dot_z_d = 0; ddot_z_d = 0; dot_psi_d = 0; ddot_psi_d = 0; %% 滑模控制器计算 % --- 位置通道(外环,生成期望姿态角)--- % 高度z控制,计算总升力U1 e_z = pos(3) - z_d; dot_e_z = vel(3) - dot_z_d; s_z = dot_e_z + smc_params.z.lambda * e_z; % 使用饱和函数替代符号函数 sat_s_z = sat(s_z, smc_params.z.phi); U1 = (m / (cos(phi)*cos(theta))) * (g + ddot_z_d ... - smc_params.z.lambda * dot_e_z ... - smc_params.z.epsilon * sat_s_z ... - smc_params.z.k * s_z); % 安全限制,防止除零或U1为负 U1 = max(0.1*m*g, U1); % 至少保持10%的升力 % 水平位置x,y控制,生成期望的滚转和俯仰角指令 % 以x方向为例,设计虚拟控制量ux e_x = pos(1) - x_d; dot_e_x = vel(1) - dot_x_d; s_x = dot_e_x + smc_params.x.lambda * e_x; sat_s_x = sat(s_x, smc_params.x.phi); ux = (m / U1) * (ddot_x_d ... - smc_params.x.lambda * dot_e_x ... - smc_params.x.epsilon * sat_s_x ... - smc_params.x.k * s_x); % 同理计算y方向的uy % ... % 根据ux, uy和期望偏航角psi_d,解算期望的滚转角phi_d和俯仰角theta_d % 这是一个代数变换: [ux; uy] = [cos(psi_d), sin(psi_d); -sin(psi_d), cos(psi_d)] * [theta_d; -phi_d] % 因此: R_psi = [cos(psi_d), sin(psi_d); -sin(psi_d), cos(psi_d)]; angle_cmd = R_psi \ [ux; uy]; % 求解线性方程组 theta_d = angle_cmd(1); phi_d = -angle_cmd(2); % 注意符号 %% --- 姿态通道(内环,计算力矩U2, U3, U4)--- % 滚转角phi控制 e_phi = phi - phi_d; dot_e_phi = p - 0; % 假设期望角速度dot(phi_d)=0 s_phi = dot_e_phi + smc_params.phi.lambda * e_phi; sat_s_phi = sat(s_phi, smc_params.phi.phi); % 根据姿态动力学方程推导控制力矩U2 % Ixx * dot(p) = U2 + (Iyy - Izz)*q*r % 令 dot(s_phi) = -epsilon*sat - k*s, 且 dot(s_phi) = dot(p) + lambda*dot(e_phi) % 可解得: U2 = I(1,1) * ( -smc_params.phi.lambda * dot_e_phi ... - smc_params.phi.epsilon * sat_s_phi ... - smc_params.phi.k * s_phi ) ... - (I(2,2)-I(3,3)) * q * r; % 俯仰角theta控制 (计算U3) % ... 类似phi的计算 % 偏航角psi控制 (计算U4) % ... 类似phi的计算,但动力学方程不同 %% 计算状态导数 % 1. 线速度导数 (来自位置动力学) R = rotation_matrix(phi, theta, psi); % 需要实现一个旋转矩阵函数 thrust_vector = R * [0; 0; U1]; gravity_vector = [0; 0; -m*g]; accel = (thrust_vector + gravity_vector) / m; % 2. 角速度导数 (来自姿态动力学) I_inv = inv(I); moments = [U2; U3; U4]; coriolis = cross(omega, I*omega); % 科里奥利项 omega_dot = I_inv * (moments - coriolis); % 3. 位置导数就是速度 pos_dot = vel; % 4. 姿态角导数 (欧拉角微分方程,注意奇异性) % 对于小角度或非极端机动,常用以下近似: angles_dot = [1, sin(phi)*tan(theta), cos(phi)*tan(theta); 0, cos(phi), -sin(phi); 0, sin(phi)/cos(theta), cos(phi)/cos(theta)] * omega; % 注意:当theta接近±90度时,此矩阵奇异。在实际飞行或高机动仿真中,应使用四元数。 % 组装状态导数向量 dxdt = [pos_dot; accel; angles_dot; omega_dot]; % 可选:记录控制输入U1-U4,用于后续绘图 persistent U_log; if isempty(U_log) U_log = []; end U_log = [U_log; t, U1, U2, U3, U4]; assignin('base', 'U_history', U_log); % 存入基础工作区 end % 饱和函数定义 function y = sat(s, phi) if abs(s) <= phi y = s / phi; else y = sign(s); end end % 旋转矩阵函数 function R = rotation_matrix(phi, theta, psi) Rz = [cos(psi), -sin(psi), 0; sin(psi), cos(psi), 0; 0, 0, 1]; Ry = [cos(theta), 0, sin(theta); 0, 1, 0; -sin(theta), 0, cos(theta)]; Rx = [1, 0, 0; 0, cos(phi), -sin(phi); 0, sin(phi), cos(phi)]; R = Rz * Ry * Rx; % 通常旋转顺序为Z-Y-X (偏航-俯仰-滚转) end这个函数内容非常密集,是算法实现的核心。它清晰地展示了如何将滑模控制律嵌入到四旋翼的动力学模型中。其中,姿态环的期望指令phi_d,theta_d是由位置环的滑模控制器产生的,这体现了级联控制的思想。
实操心得:在编写这个函数时,最大的坑是单位一致性和坐标系转换。确保所有角度都用弧度制,惯性矩单位正确,并且旋转矩阵
R的乘法顺序(是Rz*Ry*Rx还是Rx*Ry*Rz)与你定义的欧拉角旋转顺序完全匹配。一个符号错误或顺序错误就足以让仿真结果变得毫无意义。我建议在写完函数后,先用一个简单的开环恒定输入测试一下,看看四旋翼是否在重力作用下正确下坠,或者给一个力矩看姿态角是否按预期旋转,这是验证模型正确性的第一步。
4. Simulink模型搭建:可视化与模块化实现
虽然纯Matlab脚本灵活,但对于复杂系统、需要快速调整参数或直观查看信号流图的场景,Simulink更具优势。我的Simulink模型遵循模块化设计原则,将四旋翼模型、控制器、参考轨迹生成器等分开。
4.1 整体模型架构
模型顶层通常包含以下几个主要子系统:
Reference Trajectory Generator:根据时间t生成期望的位置[x_d, y_d, z_d]和偏航角psi_d,及其一阶、二阶导数。可以用Matlab Function模块或Lookup Table实现。Sliding Mode Controller:这是核心控制器子系统。输入是当前状态State和期望值Ref,输出是四个控制量[U1, U2, U3, U4]。内部就是实现了第3.2节中quadrotor_dynamics_smc函数里的控制律计算部分。Quadrotor Plant Model:四旋翼被控对象模型。输入是四个控制量,输出是12个状态。内部使用Matlab Function模块或S-Function来实现动力学方程dxdt = f(x, U),并使用积分器模块(如Integrator)进行数值积分。强烈建议使用Matlab Function模块,因为它可以直接调用.m文件中的函数,便于和脚本代码保持同步。Scope和To Workspace:用于实时观察和记录仿真数据。将关键信号(如位置误差、姿态角、控制输入)连接到Scope,同时用To Workspace模块将数据保存到Matlab工作区,供后续精细绘图和分析。
4.2 控制器子系统内部细节
在Sliding Mode Controller子系统中,你需要用Simulink模块搭建出滑模面的计算、饱和函数的实现以及最终控制量的计算。
- 滑模面计算:使用
Add、Gain(增益,即λ参数)、Derivative(求导,慎用,易引入噪声)或直接使用误差的微分信号(如果从状态中可得)来构建s = dot(e) + λ*e。对于轨迹跟踪,dot(e)通常等于(当前速度 - 期望速度)。 - 饱和函数实现:可以用
Matlab Function模块写几行代码实现sat(s/Φ),也可以用基础的Relational Operator(比较)、Switch(切换)和Gain模块来搭建。例如: |s| <= Φ ? 是 -> 输出 s/Φ 否 -> 输出 sign(s) - 控制量计算:按照推导出的公式,使用
Product、Divide、Sum、Trigonometric Function(三角函数)等模块组合计算。特别注意分母可能为零的情况,用MinMax或Saturation模块进行保护。
4.3 被控对象模型子系统
在Quadrotor Plant Model子系统中,最简洁高效的方式是使用一个Matlab Function模块。该模块的输入是控制向量U和当前状态x(来自积分器的反馈),输出是状态导数dxdt。模块内部的代码就是quadrotor_dynamics_smc函数中只包含动力学计算、不包含控制律的那部分。然后外接一组Integrator模块,对dxdt积分得到x,形成闭环。
Simulink建模避坑指南:
- 求解器选择:对于这种非线性模型,建议使用变步长求解器,如
ode45(Dormand-Prince) 或ode15s(刚性系统)。在Model Configuration Parameters中设置合适的相对容差(如1e-6)和绝对容差(如1e-9)。- 代数环:如果你的控制律计算中,
U1依赖于当前的phi,theta,而phi,theta又由包含U1的动力学方程积分得到,这可能会形成代数环。Simulink会报错。解决方法是在反馈回路中插入一个Memory模块或Unit Delay模块,打破代数环。这相当于引入了一个计算步长的延迟,在仿真步长足够小时是可接受的。- 信号维度:务必用
Mux和Demux模块清晰地管理向量信号,并用Signal Dimensions显示功能检查维度是否正确,避免维度不匹配的错误。- 初始状态:双击
Integrator模块,设置正确的初始值(如高度z0=0,其他为0),这对应了四旋翼的起飞初始状态。
5. 结果分析与Plot图解读:如何判断控制器好坏
仿真跑完了,生成了海量的数据。怎么画图才能最有说服力地展示你的滑模控制器性能?一堆乱七八糟的曲线堆在一起是毫无意义的。我通常从以下几个维度来绘制和分析Plot图,它们构成了评价控制器的“证据链”。
5.1 轨迹跟踪性能图
这是最直观的图,展示四旋翼是否能跟上期望的轨迹。
- 三维轨迹对比图:使用
plot3函数,在同一张图上绘制期望轨迹(x_d, y_d, z_d)和实际轨迹(x, y, z)。用不同颜色和线型区分。这张图能宏观展示跟踪效果。figure(‘Position‘, [100, 100, 800, 600]); plot3(x_ref, y_ref, z_ref, ‘b--‘, ‘LineWidth‘, 2, ‘DisplayName‘, ‘期望轨迹‘); hold on; plot3(x_actual, y_actual, z_actual, ‘r-‘, ‘LineWidth‘, 1.5, ‘DisplayName‘, ‘实际轨迹‘); xlabel(‘X (m)‘); ylabel(‘Y (m)‘); zlabel(‘Z (m)‘); legend; grid on; axis equal; title(‘四旋翼三维轨迹跟踪效果‘); - 各通道位置/姿态跟踪曲线:将时间
t作为横轴,分别绘制x, y, z, φ, θ, ψ的期望值与实际值随时间的变化曲线。通常用子图(subplot)排列。这张图能清晰显示稳态误差、超调量和响应速度。figure; subplot(3,2,1); plot(t, x_ref, ‘b--‘, t, x_actual, ‘r-‘); legend(‘期望x‘, ‘实际x‘); ylabel(‘x (m)‘); grid on; subplot(3,2,2); plot(t, y_ref, ‘b--‘, t, y_actual, ‘r-‘); legend(‘期望y‘, ‘实际y‘); ylabel(‘y (m)‘); grid on; % ... 依次绘制z, phi, theta, psi
5.2 误差与滑模面演化图
这部分是滑模控制特有的分析,用于验证控制器是否按理论工作。
- 跟踪误差图:绘制
e_x, e_y, e_z, e_phi, e_theta, e_psi随时间的变化。理想情况下,所有误差应渐近收敛到零或一个很小的边界层内。这是控制器性能的直接量化。 - 滑模面变量图:绘制各个通道的滑模面
s随时间的变化。这是最关键的一张图之一。理论证明,当系统进入滑模运动后,s应保持在零附近一个很小的范围内(边界层厚度Φ内)高频抖振。你的图应该显示出,在初始阶段后,s迅速收敛并稳定在边界层内。如果s一直发散或振幅很大,说明控制器参数(ε, k, λ, Φ)设置不当。figure; subplot(2,1,1); plot(t, s_z, ‘b-‘); ylabel(‘s_z‘); grid on; title(‘高度通道滑模面‘); hold on; yline(smc_params.z.phi, ‘r--‘); yline(-smc_params.z.phi, ‘r--‘); % 画出边界层 legend(‘s_z‘, ‘边界层‘); subplot(2,1,2); plot(t, s_phi, ‘b-‘); ylabel(‘s_\phi‘); grid on; title(‘滚转通道滑模面‘); hold on; yline(smc_params.phi.phi, ‘r--‘); yline(-smc_params.phi.phi, ‘r--‘);
5.3 控制输入分析图
控制输入U1, U2, U3, U4的曲线同样重要。
- 控制输入时间序列图:展示四个控制量随时间的变化。观察它们是否平滑、是否有物理可实现(如U1不能为负),以及抖振现象是否明显。使用饱和函数后,抖振应大幅减弱,表现为在稳态值附近的高频小幅度波动。
- 控制输入频谱分析(可选):如果对抖振研究深入,可以对稳态后的
U2或U3信号做快速傅里叶变换(FFT),查看其高频分量。滑模控制带来的抖振通常分布在较宽的频带上。
5.4 参数敏感性测试与对比图(进阶)
为了体现滑模控制的鲁棒性,你可以设计对比实验。
- 参数摄动:在仿真中,故意让控制器使用的模型参数(如质量
m、转动惯量I)与实际被控对象模型的参数有10%-20%的差异。然后运行相同的轨迹跟踪任务。 - 加入外部干扰:在动力学方程中,加入持续或脉冲形式的风扰力矩或力。
- 绘制对比图:将标称情况(无摄动无干扰)下的误差曲线,与存在参数摄动和干扰下的误差曲线画在同一张图上。如果滑模控制器设计良好,两条曲线应该非常接近,这表明控制器对模型不确定性和干扰不敏感。你可以同时对比一下PID控制器在相同扰动下的表现,通常PID的误差会大很多。
结果解读经验:不要只看轨迹跟踪“像不像”。一个看起来跟踪不错的曲线,可能隐藏着巨大的、高频率的控制输入(抖振),这在实物电机上是无法实现的,会迅速烧毁电机。因此,必须将跟踪误差图、滑模面图和控制输入图结合起来看。理想的仿真结果是:跟踪误差小且收敛快,滑模面迅速进入并保持在边界层内,控制输入连续、平滑且幅值合理。如果控制输入出现剧烈的锯齿波,即使跟踪效果好,也意味着需要调整边界层厚度
Φ或趋近律参数ε, k,在鲁棒性和平滑性之间取得更好的平衡。
本文还有配套的精品资源,点击获取