简介一阶倒立摆系统作为控制理论中的经典力学模型其建模、仿真与控制分析资源面向自动控制、机械工程及机器人领域学习者帮助理解动态系统稳定性与反馈控制核心方法。内容先介绍由可移动支点和质量点构成的基本结构及直立稳定目标再阐述基于重力、惯性、摩擦力等因素的一阶非线性动力学方程推导并演示如何用Simulink构建仿真模型设置输入输出与系统参数通过Scope观察摆杆角度变化设定不同初始条件与外部扰动进行仿真验证。之后重点说明PID控制参数整定并拓展滑模控制、自适应控制、LQR等高级策略便于按需选择优化。资源包仅14KB已有3841人学习下载适合快速获取从系统建模到控制器设计的完整思路与实现技巧为机器人稳定控制等应用奠定基础。1. 一阶倒立摆建模仿真先明白它为什么开环必倒一阶倒立摆是控制领域最典型的“看着简单、上手翻车”的对象一辆小车加一根摆杆自由度只有两个动力学上却同时具备开环不稳定、欠驱动、强非线性三个特征。做建模与仿真分析时很多人第一步就把牛顿方程抄错——铰链反力、符号约定、摆角以竖直向上还是向下为基准任何一处错了后续仿真要么直接发散要么控制器怎么调都只能撑两秒。这篇文章按一条完整可复现的路径展开从拉格朗日方程推出运动方程在平衡点线性化得到状态空间模型再用 Simulink 搭出保留非线性的仿真模型最后用级联 PID 和 LQR 把系统压住。控制部分的意义在于验证模型重心始终落在建模与仿真分析这条主线上。适合正在做课程设计、数学建模竞赛控制题或者从零接触机器人平衡控制的工程师对照落地所有代码和参数都按可直接复现的标准给出只要不换符号约定跑通不需要额外调试。2. 一阶倒立摆数学建模从拉格朗日方程到状态空间2.1 建模方法选型为什么放弃牛顿法一阶倒立摆的推导几乎都会从牛顿第二定律开场但对着既有平动又有转动的摆杆牛顿法必须先把铰链处的约束力当作未知量列出来水平、垂直各写一个方程再联立消元。这个过程里只要某个力的方向画反后边全盘皆错而且中间变量一多别人也很难复核。我一般直接用拉格朗日方程选小车位移 x 和摆角 θ 作为广义坐标写出动能 T 和势能 V代入 L T − V 求偏导约束力根本不会出现在方程里。做数学建模竞赛或课程报告时从能量出发推导也更好审——评审能顺着公式一步步核对中间没有来路不明的消元。建模参数先统一放在下表后面所有代码共用这一套符号和取值。符号物理含义仿真取值M小车质量1.0 kgm摆杆质量按集中质量处理0.1 kgl摆杆质心到转轴的距离0.5 mθ摆角竖直向上为 0向右偏为正初始 0.1 radx小车位移向右为正0 mF作用在小车上的外力由控制器给出注意 θ 的方向定义直接决定后面控制增益的符号。很多人栽在“摆角相对向上还是向下、顺时针还是逆时针”上本文统一取竖直向上为 0、向右偏为正后面每一个公式、每一段代码都跟这个约定走中途不要切换。2.2 拉格朗日方程推导与非线性运动方程摆杆质心坐标是 (x l·sinθ, l·cosθ)对时间求导得到质心速度平方 ẋ² 2l·ẋ·θ̇·cosθ l²·θ̇²。系统动能 T 由小车平动、摆杆平动和摆杆绕自身质心转动三项组成势能 V mgl·cosθ竖直向上时势能最大拉格朗日量为L ½(Mm)ẋ² ml·ẋ·θ̇·cosθ ½(ml²J)·θ̇² − mgl·cosθ把 L 分别代入广义坐标 x 和 θ 的拉格朗日方程整理后得到两个非线性运动方程(Mm)ẍ ml·θ̈·cosθ − ml·θ̇²·sinθ F 式 1ml·ẍ·cosθ (ml²J)·θ̈ − mgl·sinθ 0 式 2J 是摆杆绕自身质心的转动惯量。仿真分析里最常用 J 0 的集中质量近似相当于把整根杆的质量压到质心趋势与实物一致且推导省事如果摆杆明显不是细杆把 (ml²J) 整体记为 I 代入即可公式结构不变。以下统一取 I ml²。2.3 在平衡点线性化得到状态空间模型在竖直向上的平衡点附近取 θ ≈ 0令 cosθ ≈ 1、sinθ ≈ θ忽略 θ̇² 等高阶小量式 1 和式 2 退化为(Mm)ẍ ml·θ̈ Fml·ẍ ml²·θ̈ − mgl·θ 0从第二个方程解出 θ̈ gθ/l − ẍ/l代回第一个方程得到 ẍ F/M − mg·θ/M。取状态向量 [x; ẋ; θ; θ̇]标准状态空间表达式如下这段 MATLAB 代码可以直接运行验证M 1.0; m 0.1; l 0.5; g 9.81; A [0 1 0 0; 0 0 -m*g/M 0; 0 0 0 1; 0 0 (Mm)*g/(M*l) 0]; B [0; 1/M; 0; -1/(M*l)]; C [1 0 0 0; 0 0 1 0]; D zeros(2,1); sys ss(A, B, C, D); fprintf(开环极点: \n); disp(eig(A)); fprintf(能控性矩阵的秩: %d\n, rank(ctrb(A, B)));A 矩阵第二行第三列的 −mg/M 是摆杆偏角通过重力分量对小车加速度的反作用第四行第三列的 (Mm)g/(Ml) 是重力项构成的正反馈这一项正是一阶倒立摆开环不稳定的根源。C 矩阵取前两行让输出同时观测小车位移和摆角这样后续做全状态反馈时不需要额外设计观测器。2.4 开环极点与能控性先看清系统有多糟运行上面的代码eig(A) 得到 0、0、±4.65。一对零极点来自小车的纯积分特性——没有外力时小车匀速漂移±4.65 是摆杆的失稳模态对应时间常数约 0.22 秒意味着一个微小扰动之后摆杆在零点几秒内就会明显倒下。所有控制方案的本质就是把这对正负极点中位于右半平面的那个拉回左半平面。能控性矩阵的秩为 4说明四个状态都可以由外力 F 控制到任意目标这是后面 LQR 全状态反馈的前提。反过来如果实物上只测摆角 θ 一个量需要单独验算能观性通常小车编码器给 x、摆杆编码器给 θ两个量都测全状态反馈才能直接成立。3. 一阶倒立摆仿真Simulink 两种搭法与仿真发散排查3.1 最快闭环State-Space 模块五分钟跑通线性模型拿到 2.3 节的 A、B、C、D 后最快的验证方式是在 Simulink 里放一个 State-Space 模块。双击填入四个矩阵Initial conditions 填 [0;0;0.1;0] 表示给摆杆一个初始偏角控制量用一个 Gain 模块实现 u −K·[x;ẋ;θ;θ̇]把四个状态量全部引出接到矩阵增益 K 上做负反馈Scope 看 x 和 θ 两条曲线。这个方案只对线性模型成立如果矩阵是在工作区里用 ss 命令建好的也可以直接用 LTI System 模块加载对象省去手抄矩阵。求解器先用默认的 ode45 就能跑但建议把 Max step size 手动改成 1e-3。原因在于模型里有一个 4.65 rad/s 的失稳极点仿真步长过大时每一步的数值误差都会被这个极点指数放大最终结果看起来像控制器失效其实是求解器精度不够。3.2 保留非线性的 S-Function 模型线性模型适合快速验证控制器结构但初始摆角到 0.5 rad 左右时线性化误差已经不可忽略往实物移植前必须在非线性模型上再测一轮。常见做法是写一个 Level-1 S-Function把式 1、式 2 联立解成显式的一阶微分方程组由 Simulink 的积分器完成数值求解。下面这段代码直接可用模块参数 M、m、l、g 在 S-Function 模块的参数列表里按顺序传入function [sys,x0,str,ts] pend_sfun(t,x,u,flag,M,m,l,g) % 状态: x1小车位移 x2车速 x3摆角 theta x4角速度 switch flag case 0 [sys,x0,str,ts] mdlInitializeSizes(); case 1 sys mdlDerivatives(t,x,u,M,m,l,g); case 3 sys x; % 全状态输出 case {2,4,9} sys []; otherwise error([unhandled flag ,num2str(flag)]); end function [sys,x0,str,ts] mdlInitializeSizes() sizes simsizes; sizes.NumContStates 4; sizes.NumDiscStates 0; sizes.NumOutputs 4; sizes.NumInputs 1; sizes.DirFeedthrough 0; % 输出不经过输入避免代数环 sizes.NumSampleTimes 1; sys simsizes(sizes); x0 [0; 0; 0.1; 0]; % 初始摆角 0.1 rad str []; ts [0 0]; % 连续系统 function dx mdlDerivatives(~,x,u,M,m,l,g) F u(1); th x(3); w x(4); I m*l^2; % 集中质量近似 D (Mm)*I - (m*l*cos(th))^2; dx zeros(4,1); dx(1) x(2); dx(2) (I*F - m^2*g*l^2*sin(th)*cos(th) I*m*l*w^2*sin(th)) / D; dx(3) x(4); dx(4) (m*g*l*sin(th) - m*l*dx(2)*cos(th)) / I;mdlInitializeSizes 里的 NumContStates 必须为 4对应四个连续状态ts [0 0] 表示连续系统而不是离散采样DirFeedthrough 置 0 是因为输出直接取状态 x不经输入 u这个标志位写错会出现代数环警告。mdlDerivatives 里第四个式子用到了 dx(2)这是 θ̈ 依赖 ẍ 的动力学耦合必须先算 ẍ 再算 θ̈两行顺序不能颠倒。3.3 仿真发散按这张表逐项排查仿真发散在倒立摆调试里几乎人人会遇到而且大半不是控制器问题。下表是排错优先级最高的四类情况现象常见原因处理方式输出瞬间变成 Inf/NaN步长过大或初值越界Max step size 降到 1e-4~1e-3检查初始摆角是否在 ±π/2 内振荡幅度指数级增长反馈符号接反核对 θ 定义方向把对应增益取反而不是加大增益曲线呈高频锯齿求解器步长与模型动态不匹配换 ode4 固定步长步长不高于 1e-3反复报代数环错误控制器输出直接参与自身计算在回路中插入 Memory 或 Unit Delay 模块断开代数量提示仿真发散时先跑开环。把控制器增益全部置零给一个 0.1 rad 初始摆角如果模型输出有界摆下去后在最低点附近来回摆动说明模型本身没问题问题在反馈回路里。求解器选择上纯线性模型 ode45 足够一旦模型里加入摩擦、饱和、齿隙等刚性环节换 ode15s 这类变步长隐式求解器发散概率会明显下降。4. 一阶倒立摆控制律设计级联 PID 与 LQR 的整定路径4.1 结构先行只控摆角还是一并控位置只让摆杆不倒一个 PD 就够了控制律写成 F Kp_θ·θ Kd_θ·ω。关键在符号按本文 θ 向右偏为正的约定这两个增益必须为正。直觉解释是摆杆向右倒时小车要先向右加速去“接住”它这和我们习惯的位置回路方向相反所以角度项看起来像正反馈。如果换成 θ 逆时针为正的定义两个增益全部取反。课程设计和竞赛里几乎都要同时压住小车位置这时把位置回路叠在外面就成了最典型的级联 PID 控制结构内环 PD 管摆角外环位置回路把整个角度回路当成执行机构控制律为F Kp_θ·θ Kd_θ·ω − Kp_x·(x − x_ref) − Kd_x·ẋ位置项是常规负反馈角度项保持正的 PD。先调内环再调外环是两个回路参数整定的基本顺序。4.2 整定顺序与完整闭环仿真下面这段 MATLAB 脚本直接把控制器和 2.2 节非线性模型封进一个函数用 ode45 做闭环仿真比 Simulink 调起来更直观适合先扫参数% 一阶倒立摆闭环仿真内环 PD 控角 外环 PD 控位 M 1.0; m 0.1; l 0.5; g 9.81; Kp_t 40; Kd_t 6; % 内环摆角、角速度 Kp_x 2; Kd_x 1.5; % 外环位置、速度 x_ref 0; x0 [0; 0; 0.1; 0]; [t, z] ode45((t,z) closedPendulum(t,z,M,m,l,g, ... Kp_t,Kd_t,Kp_x,Kd_x,x_ref), [0 10], x0, ... odeset(MaxStep, 1e-3)); subplot(2,1,1); plot(t, z(:,3)*180/pi); ylabel(theta (deg)); subplot(2,1,2); plot(t, z(:,1)); ylabel(x (m)); xlabel(t (s)); function dz closedPendulum(~,z,M,m,l,g,Kp_t,Kd_t,Kp_x,Kd_x,x_ref) x z(1); vx z(2); th z(3); w z(4); F Kp_t*th Kd_t*w - Kp_x*(x - x_ref) - Kd_x*vx; I m*l^2; D (Mm)*I - (m*l*cos(th))^2; ax (I*F - m^2*g*l^2*sin(th)*cos(th) I*m*l*w^2*sin(th)) / D; ath (m*g*l*sin(th) - m*l*ax*cos(th)) / I; dz [vx; ax; w; ath]; end初始参数按 Kp_θ40、Kd_θ6 起调内环闭环自然频率约 7.6 rad/s阻尼比约 0.8再叠加 Kp_x2、Kd_x1.5 的位置回路。整定顺序建议如下只保留角度项给定 0.1 rad 初始偏角Kp_θ 从 20 往上加回零太慢就加大出现振荡就加 Kd_θ目标是 3 秒内摆角收敛到 0.5° 以内。加入位置项Kp_x 从 1 起调、Kd_x 取 0.8 左右观察 x 是否在 5 秒内回到 ±1 cm。外环带宽保持在内环的三分之一到五分之一。如果小车还没到位摆角就开始抖先降 Kp_x而不是继续加内环增益。一旦出现“越控越倒”优先查符号而不是扫参数。把 F 里任一单向反向整定就永远不收敛。注意以上增益依赖 θ 的符号约定和集中质量近似。换模型参数后先跑一遍开环仿真确认失稳模态的数值再按上述顺序重调。4.3 LQR 参数设计Q、R 怎么给不踩坑全状态可测时LQR 比手调 PID 省事得多。基于 2.3 节的 A、B直接调用 lqrQ diag([100 1 200 10]); % x, vx, theta, omega 的权重 R 1; K lqr(A, B, Q, R); eig(A - B*K) % 验证闭环极点Q 的对角元素对应四个状态的重要程度摆角权重最高50~300位置其次10~100角速度是阻尼项5~20小车速度权重最小。R 从 1 起调R 越小控制越猛但控制量高频抖动的风险也越大实物上通常不敢比 0.1 更小。有一个容易吓到人的现象按本文 θ 约定lqr 返回的 K 第三、四列是负数这恰恰是对的。u −K·x 展开后负的 θ 系数会把角度项变成与直觉一致的正反馈方向。看到负号不要“好心”去修正直接放进 u −K·x 即可。PID 与 LQR 的取舍可以看下表对比项级联 PIDLQR对模型依赖低方向对即可调稳高A、B 不准则增益失真整定成本逐项试凑依赖经验Q、R 两个矩阵一次成型鲁棒性参数拉偏时衰减慢固定增益下对摆长变化敏感实现成本两个 PD 叠加MCU 上极简需要全状态反馈和矩阵乘适用阶段实物联调、快速验证性能优化、竞赛指标冲刺当负载 m 或摆长 l 变化时失稳模态的固有频率 ω_p sqrt((Mm)g/(Ml)) 会跟着变固定增益不再最优。常见的工程做法是在线辨识 ω_p按频率查表切换 PID 或 LQR 增益这就是自适应频率控制在一阶倒立摆系统上的落地形态可以先在 Simulink 里用 Lookup Table 做增益调度验证稳定后再搬到实物。5. 一阶倒立摆移植实物前的离散化与验证技巧5.1 控制器离散化与采样周期选择实物控制器运行在 MCU 上所有控制律必须写成差分方程。常见做法是用 c2d 把连续模型离散化再用 lqrd 直接求解离散域的最优增益Ts 0.005; sysd c2d(ss(A,B,eye(4),zeros(4,1)), Ts, zoh); Kd lqrd(A, B, Q, R, Ts); % Q/R 沿用连续域设置zoh 即零阶保持器对应实物中 DAC 或 PWM 在每个采样周期内保持输出不变。采样周期从失稳极点倒推4.65 rad/s 对应时间常数约 0.22 秒闭环带宽通常设计在 10~20 rad/s采样频率低于带宽 10 倍时离散化误差就会吃掉稳定裕度。电机控制里常见的 1 kHz 已接近下限一般取 1~5 ms摆杆惯量很大的场景才敢放到 10 ms 以上。5.2 执行器饱和与抗积分饱和实物电机的 PWM 输出有上限仿真里要在 u 后面加 Saturation 模块饱和值取电机连续出力的 ±80%。外环若引入积分项饱和期间积分仍在累积撤销饱和后控制器会大幅过冲这是积分饱和的典型表现。简单的 clamping 抗饱和写法如下% 外环 PI 积分项Ts 为控制周期 integral_x integral_x Ki_x*(x_ref - x)*Ts; F Kp_t*th Kd_t*w - Kp_x*(x - x_ref) - Kd_x*vx - integral_x; if abs(F) Fmax integral_x integral_x - Ki_x*(x_ref - x)*Ts; % 冻结积分 end饱和判断必须用叠加后的合力 F而不是只看积分项本身否则冻结逻辑会在控制器正常输出时误触发。这个细节在 Simulink 里经常被忽略实物上板却发现小车反复冲过目标点。5.3 上电前必跑的五项验证仿真收敛只是起点往实物移植前建议把以下五项跑完一是大扰动测试初始摆角 0.3 rad 时控制器仍能在 3 秒内回稳二是参数拉偏把 M 在 0.6~1.4 kg、m 在 0.05~0.2 kg 范围内扫一遍记录稳定边界三是噪声注入给位置量测加 ±1 mm、角度量测加 ±0.5° 的白噪声看控制量抖动幅度是否超过执行器承受范围四是饱和检查max(|F|) 必须低于电机连续出力的 80%五是离散裕度离散闭环极点的模值全部小于 0.95并留 10% 以上的余量。实物上的摩擦力、齿隙和传感器延迟会把仿真里留的稳定裕度吃掉一大半这五项里的余量至少按两倍准备不要卡着稳定边界上电。本文还有配套的精品资源点击获取