
你在搜李雅普诺夫一阶系统/二阶系统稳定性分析或者自适应自抗扰控制Simulink实现的时候大概率已经被理论书里的V函数、正定矩阵、半负定这些概念劝退过一轮了。这很正常我自己当年啃这些内容的时候也差不多是同样的状态——公式能看懂但一回到Simulink里就不知道从哪下手更不知道调出来的参数到底在干什么。这篇东西就是想把这中间的沟填上一点我会从一阶系统的李雅普诺夫稳定性分析讲起再做二阶系统的V函数构造然后过渡到自适应自抗扰控制ADRC的设计思路最后把整套东西在Simulink里面跑通附上我实际调试时踩过的坑和总结出来的参数整定经验。适合正在做控制理论课程设计、研究生课题里需要用到稳定性证明或者工作中要用ADRC做扰动抑制但不想只看PPT的读者。1. 为什么李雅普诺夫方法值得从一阶二阶系统开始练手很多初学者有个误区觉得李雅普诺夫稳定性分析是纯理论课内容跟实际仿真关系不大。但我做过几个工程项目之后发现恰恰相反——李雅普诺夫方法是在仿真发散时最快帮你定位问题根源的工具之一。你辛辛苦苦搭好一个自适应控制器仿真跑了0.3秒就数值爆炸了这时候如果只会盯着波形看大概率调一晚上都找不出原因。但如果你有一阶、二阶系统的李雅普诺夫分析功底你会习惯性地先检查V函数是否正定、V̇是否满足负定条件然后很快定位到是自适应律里的某个增益符号搞反了还是ESO带宽给得太高导致观测误差项在V̇里贡献了正项。1.1 从一阶系统开始的核心逻辑一阶系统之所以是入门的最小可行案例是因为它的状态量只有一个V函数的选择几乎是唯一的不需要面对多变量系统里V函数构造方式不唯一这个让人头疼的问题。考虑最简单的标量系统ẋ -a·x d(t)其中a 0是系统固有衰减系数d(t)是外部扰动。如果d(t) 0这个系统的解析解是x(t) x(0)·e^(-at)所有人都知道它稳定。但问题在于你凭什么严格证明它稳定用李雅普诺夫方法取V ½x²显然V 0当x ≠ 0然后求导V̇ x·ẋ -a·x² ≤ 0当a 0时V̇是负定的所以系统在李雅普诺夫意义下渐近稳定。这个推导简单到让人觉得这不是废话吗。但它的价值在于建立了一种思维模式把稳定性问题转化为能不能找到一个能量函数让它沿着系统轨迹单调下降。1.2 为什么实际工程里要重新用李雅普诺夫方法纯粹的一阶线性系统当然不需要这么绕。但一旦加入不确定性、非线性项、时变参数甚至控制器本身带有自适应调节机制解析解就不存在了这时候李雅普诺夫方法是少数几条能给出严格保证的路径。比如说一阶系统带不确定参数θ和未知扰动ẋ θ·u d(t)其中u是控制输入θ的真实值未知但范围可以估计。你设计自适应控制器u -k·x - θ̂·x这里的θ̂是θ的估计值这时候整个闭环系统变成了(x, θ̃)的二维系统θ̃ θ̂ - θ是参数估计误差。怎么证明整体稳定还是用李雅普诺夫方法但V函数要扩展为V ½x² ½·θ̃²/γ其中γ是自适应增益。对这个V求导然后把自适应律设计成θ̂̇ γ·x²即θ̃̇ γ·x²假设θ是常值就能消掉交叉项得到V̇ -k·x² ≤ 0。这就是自适应控制里最经典的李雅普诺夫导出自适应律的思路——不是先拍脑袋定自适应律再证明稳定而是从V̇要满足负半定的需求里反推自适应律应该长什么样。这个方法在Simulink里实现起来非常直观下面第三节我会给完整的模型拆解。2. 一阶系统V函数构造的最小可行案例与Simulink验证先搭一个最简单的一阶系统Simulink模型来验证上面的理论推导。这个模型虽然简单但它包含了后续做ADRC时要用到的所有基础要素被控对象、控制器、扰动注入、信号观测。2.1 一阶系统V函数在Simulink里的直观呈现模型结构大致是这几块一个积分器表示状态x一个Gain模块表示-a·x一个Step或正弦信号作为扰动d(t)注入。控制器部分先用最简单的比例控制u -k·x然后我们把V ½x²和V̇ x·ẋ同时在模型里算出来用Scope观察。具体操作步骤新建Simulink模型从Continuous库拖一个Integrator模块命名x。Integrator的输出接一个Gain模块增益设为-1对应a 1再接回Integrator的输入构成负反馈闭环。在Integrator输入端用Add模块叠加一个Step扰动Step Time设为5Final Value设为0.5。用Fcn模块或者直接搭乘法器计算V ½·x²。用Derivative模块对V求导得到V̇——注意这里只是仿真验证用实际控制器里不要用Derivative模块噪声会放大到你怀疑人生。仿真参数设置Solver用ode45相对误差1e-6仿真时长20秒初始条件x(0) 1。跑完仿真你会看到三组曲线x指数衰减到0附近受到Step扰动后有一个小凸起然后恢复V单调下降并趋于0V̇始终在0以下、只在扰动注入瞬间有一个短暂波动。这个结果跟理论推导完全吻合。这里有个值得注意的细节理论上的V̇恒小于等于0但在仿真里扰动注入的瞬间V̇可能会短暂为正因为Step信号是不连续的等效于给系统注入了一个瞬时的能量冲击。这是连续系统理想模型和数值仿真之间的固有差异不表示控制器有问题。2.2 一阶系统仿真的V函数可视化技巧如果你想让自己对V函数的理解更直观可以做一个更有意思的可视化用XY Graph模块X轴接xY轴接V你会看到一条开口向上的抛物线V ½x²。然后运行仿真观察点沿着抛物线滑向原点。这个可视化在课程答辩里非常加分。不过我要提醒一下XY Graph在新版本MATLAB里被归类为legacy不影响使用但如果你的MATLAB版本比较新R2020a之后也可以用Simulink Dashboard里的Scope替代或者把数据导出到MATLAB工作区后用plot命令画图这样排版更可控。我个人的习惯是Simulink里跑数据然后导出到工作区用MATLAB脚本统一出图。因为Simulink原生的Scope在论文/报告里截图的清晰度不够而且不方便统一坐标轴范围。具体做法是在模型里加To Workspace模块变量名设为x_out和V_outSave format选Timeseries跑完后在工作区执行figure; subplot(2,1,1); plot(x_out.Time, x_out.Data); grid on; xlabel(Time (s)); ylabel(x); subplot(2,1,2); plot(V_out.Time, V_out.Data); grid on; xlabel(Time (s)); ylabel(V 0.5*x^2);这样出来的图直接能放进报告里。2.3 带自适应律的一阶系统模型搭法接下来把问题升级被控对象改成ẋ θ·u d(t)其中θ 1.5是未知常值我们要设计自适应控制器让x收敛到0同时在线估计θ。控制器结构u -k·x - θ̂·x自适应律θ̂̇ γ·x²这里的γ是自适应增益k是比例增益。在Simulink里实现对象部分Integrator的输入是θ·u d(t)用一个Gain模块值为1.5在乘法后面模拟真实对象参数θ。控制器部分用MATLAB Function模块写u -kx - theta_hatx输入x和theta_hat。自适应律部分θ̂̇ γ·x²用Fcn模块算γ·x²然后接Integrator得到θ̂。注意真实对象里的θ和控制器估计的θ̂是两个完全独立的信号路径千万不要在模型中把它们连在一起。仿真初始条件x(0) 1θ̂(0) 0。参数取k 2γ 0.5。跑完之后你会看到两个关键现象一是x收敛到0二是θ̂从0逐渐上升到1.5附近但上升过程伴随着震荡——这是自适应控制里非常典型的参数估计与状态收敛互相耦合的现象。振荡的幅度和γ直接相关γ越大收敛越快但振荡越剧烈γ太小则收敛很慢。这个矛盾就是自适应控制调参的核心痛点后面做ADRC时还会再遇到。3. 二阶系统的V函数进阶能量视角与相轨迹验证二阶系统的李雅普诺夫分析比一阶系统上了一个台阶因为状态变量有两个V函数的构造有了自由度——你可以选不同的V只要满足正定且V̇负定都能证明稳定性但不同V给出的物理直觉完全不同。这里我推荐从机械系统的能量视角切入因为这是最自然、也最容易在Simulink里验证的路径。3.1 质量-弹簧-阻尼系统的V函数构造考虑经典的二阶系统m·ẍ c·ẋ k·x 0其中m、c、k分别是质量、阻尼系数、弹簧刚度都为正。这是控制理论里最经典的被对象也是你后面做位置伺服控制、振动抑制、乃至ADRC仿真的常见对象。把它写成状态空间形式ẋ₁ x₂ẋ₂ -(k/m)·x₁ - (c/m)·x₂这里x₁ x是位移x₂ ẋ是速度。系统的总机械能是V ½k·x₁² ½m·x₂²第一项是弹簧的弹性势能第二项是质量的动能。显然V ≥ 0且只有x₁ 0且x₂ 0时V 0所以V是正定函数。对V求导V̇ k·x₁·x₁̇ m·x₂·x₂̇ k·x₁·x₂ m·x₂·(-(k/m)x₁ - (c/m)x₂) k·x₁·x₂ - k·x₂·x₁ - c·x₂² -c·x₂² ≤ 0妙就妙在弹性势能项和动能项的交叉项正好抵消了最后只剩下阻尼耗散项-c·x₂²。这就是能量从机械能转换为热能系统的总能量单调下降最终停在平衡点这个物理事实的数学表达。用Simulink搭这个模型非常顺手两个Integrator串联构成二阶对象然后自己算V函数看它是否单调下降。3.2 相轨迹与V函数等值线联合分析二阶系统的分析要比一阶系统多一个工具相平面。把x₂作为纵轴、x₁作为横轴画相轨迹你会看到系统从任意初始状态出发相轨迹螺旋收敛到原点欠阻尼情况或直接滑向原点过阻尼情况。而V ½k·x₁² ½m·x₂² E这个方程在相平面上是一族椭圆等值线每个椭圆代表一个固定的能量水平。把相轨迹和V的等值线叠在一张图上看你能直观理解李雅普诺夫定理的几何含义相轨迹从外层椭圆穿向内层椭圆能量在不断减小。我在实际做课程展示时很喜欢用这个图因为它能让没学过控制理论的人也能一眼看出系统在往能量低的方向走。Simulink里实现方法用XY Graph模块X轴接x₁Y轴接x₂跑完就能看到相轨迹。如果想看到等值线把仿真数据导出后加一行MATLAB画图代码[X1, X2] meshgrid(-1.5:0.1:1.5, -1.5:0.1:1.5); E 0.5*100*X1.^2 0.5*1*X2.^2; % 以k100, m1为例 contour(X1, X2, E, 20); hold on; plot(x1_out, x2_out, r, LineWidth, 1.5);你会看到红色相轨迹从外圈椭圆一圈一圈地切向原点每一圈穿过的椭圆等值线对应的能量都在下降。这个画面比我在这里写一百句话都有用。3.3 二阶系统V函数构造的两个常见误区第二个坑是V̇是负半定的V̇ -c·x₂²只在x₂ 0时为0这种情况只能推出系统稳定Lyapunov stability要推出渐近稳定还需要用到LaSalle不变集原理——即从集合{x₂ 0}出发的最大不变集只有原点。很多教程草草一句所以系统渐近稳定就带过了但你要知道这里缺失的论证是什么。在Simulink仿真里你可以直接验证如果初始条件里的x₁ ≠ 0但x₂ 0比如拉开的弹簧静止释放系统不会停在x₂ 0那条线上因为弹簧力会推动质量运动所以系统的状态不会永远停留在V̇ 0的集合里最终还是会回到原点。第二个常见的坑是强行套用V ½x₁² ½x₂²这种正方形求和形式的二次型而不顾系统的物理结构。对于质量-弹簧-阻尼系统如果取V ½x₁² ½x₂²无加权你会得到V̇ x₁·x₂ x₂·(-(k/m)x₁ - (c/m)x₂) (1 - k/m)·x₁·x₂ - (c/m)·x₂²当k/m ≠ 1时交叉项(1 - k/m)·x₁·x₂无法保证消掉V̇的符号就不确定了。这说明V函数不是随便选的它必须跟系统的能量结构匹配。这也是为什么我一直强调从物理能量视角出发选V函数而不是从凑一个二次型出发——前者有物理直觉兜底后者纯靠运气。4. 从稳定性分析到自适应自抗扰控制ESO与控制律的设计逻辑一阶、二阶系统的李雅普诺夫分析打底之后就可以进入本文的重头戏自适应自抗扰控制ADRC。很多人第一次接触ADRC是被它的抗扰名头吸引的觉得它可以不管对象模型就干活。这句话对了一半——ADRC确实不需要精确模型但它对哪些东西需要未知、哪些东西必须已知是有明确界定的如果你搞不清楚这个边界仿真模型就是一锅粥。4.1 扩张状态观测器ESO到底在观测什么ADRC的灵魂是扩张状态观测器Extended State Observer, ESO。以二阶系统为例假设被控对象可以写成ẍ f(x, ẋ, d, t) b·u其中f(x, ẋ, d, t)包含了内部非线性动态、模型不确定性和外部扰动统称总扰动b是控制增益可以理解为输入通道的不确定系数但通常用估计值b₀代替。ADRC的核心思路是把总扰动f扩张成一个新的状态变量x₃ f然后构造一个三阶观测器去同时估计x₁、x₂、x₃ẑ₁ z₂ β₁·(y - z₁)ẑ₂ z₃ β₂·(y - z₁) b₀·uẑ₃ β₃·(y - z₁)其中β₁、β₂、β₃是观测器增益通常按带宽参数化β₁ 3ω₀β₂ 3ω₀²β₃ ω₀³ω₀被称为观测器带宽。z₁跟踪输出yz₂跟踪速度ẋz₃跟踪总扰动f。控制律设计为u (u₀ - z₃) / b₀其中u₀是某个基础控制器比如PD控制器的输出。这样做的意图很清晰把z₃总扰动估计值从控制量里减掉相当于用一个负的估计来抵消真实扰动——如果z₃ ≈ f那么被控对象就变成了ẍ ≈ u₀一个完美的二阶积分器。于是你只需要对u₀设计一个简单的PD控制器就能控制它。这就是ADRC的动态线性化或者说反馈线性化思想。在Simulink里实现ESO最直接的办法是用MATLAB Function模块写一个三阶积分器表达式或者拆成三个Integrator模块——前者代码短、易维护后者纯模块化、适合演示底层原理看你个人习惯。4.2 自适应机制加在哪里从固定增益到在线调节理论中的ADRC往往假设b₀已知且固定但工程实际中b₀经常是不精确的——比如你控制电机的力矩系数、控制无人机时的气动增益都只是粗略估计。当b₀与实际b偏差较大时ESO的补偿效果会打折扣甚至出现极限环振荡。自适应ADRC的思路就是把这个增益失配也当作一种不确定性用在线估计的方式去处理。具体自适应机制的选择有几种自适应方案核心思想优点缺点增益b₀在线估计估计真实b并实时更新b₀直接解决增益失配需要持续激励条件否则参数不可辨识控制增益自适应用自适应律补偿残余扰动结构简单易与李雅普诺夫分析结合自适应律带宽受系统动态限制模糊/神经网络自适应用万能逼近器估计非线性函数适应复杂非线性工程落地难度大参数多对于Simulink仿真入门我个人推荐第二种在PD控制器输出的基础上叠加一个自适应补偿项用李雅普诺夫方法导出自适应律——就是把第一节里那个反推自适应律的思路搬到二阶系统上。具体做法是取滑模面s ė λ·e其中e r - x是跟踪误差然后设计自适应律让参数估计误差带来的项在V̇里被消掉。这跟我在第一节一阶系统里做的推导是同构的只是状态从一维变成了二维。4.3 一个具体的设计实例二阶位置跟踪自适应ADRC拿一个具体的仿真案例来说被控对象是质量-弹簧-阻尼系统加上一个非线性扰动项目标是从初始位置0跟踪一个方波信号。数学描述m·ẍ -c·ẋ - k·x f_nl(x) b·u d(t)其中f_nl(x) sin(2x)是非线性弹簧力d(t) 0.5·sin(3t)是外部扰动m 1, c 1, k 5, b 2。设计以下控制结构ESO估计总扰动包括f_nl、参数偏差、d(t)带宽ω₀ 30。PD控制器u₀ kp·e kd·ėkp 100kd 20。自适应参数θ_filet用于补偿残余的常值扰动自适应律θ̂̇ γ·ss是滑模面γ 5。这个结构在Simulink里搭起来大约需要20分钟但跑通之后你会非常清楚地看到三个现象一是ESO的z₃波形会跟踪f_nl d(t)的总和二是有了补偿之后跟踪误差比纯PD控制小一个数量级三是自适应项能进一步把稳态误差压到接近零。5. Simulink仿真模型搭建架构、参数与调参全流程理论讲完了下面进入纯干货阶段——模型怎么搭、参数怎么设、调参顺序是什么。我尽量按我实际动手时的操作顺序来写你照着做就能复现。5.1 仿真模型总体架构与模块清单整个模型分成四层被控对象层、ESO层、控制律层、自适应层。推荐用Bus信号把各层之间的接口整清楚别用Goto/From满天飞的方式——短期看省事长期看维护成本极高。核心模块清单如下层级模块配置要点被控对象Integrator ×2x1初值0x2初值0被控对象MATLAB FunctionPlant写微分方程右侧x2; (-cx2 - kx1 sin(2x1) bu d(t))/m被控对象Sine Wave注入扰动d(t)Amplitude 0.5Frequency 3ESOMATLAB FunctionESO输入y, u输出z1, z2, z3控制律MATLAB FunctionControl输入r, z1, z2, z3, theta_hat输出u自适应Fcn Integrator Gain计算γ·s再积分得到theta_hat信号记录To Workspace ×4分别记录y, r, u, z3, theta_hat一个关键的工程建议把b₀、ω₀、kp、kd、γ这些参数全部定义成MATLAB工作区的变量Simulink模块里的Gain参数直接填变量名。这样你就不用每次调参都去双击模块改数字直接在工作区改完跑仿真就行。如果你用的是R2020a之后版本推荐用Simulink的Parameter Configuration方式统一管理或者干脆写一个初始化脚本set_params.m把参数集中放在文件头部。5.2 关键参数取值与整定顺序ADRC参数整定的顺序我建议严格遵循从内到外、从后到前的原则。这个顺序我踩了不少坑才总结出来先整定ESO的ω₀。把控制律暂时退化成u -z3/b0只做扰动补偿观察z1对y的跟踪速度。ω₀太小z1滞后严重后面全白搭ω₀太大观测器对噪声敏感Simulink里表现为z3波形高频抖动。经验值ω₀从对象响应速度的35倍起步逐步加大到z1能快速跟踪y且z3不抖为止。再整定PD的kp、kd。整定ESO时可以把PD设为纯P控制把kp从小往大加看到系统开始振荡就回退20%然后加入kd抑制超调。这里有个经验法则kp决定闭环刚度kd决定阻尼先刚度后阻尼。最后加自适应项。γ从很小的值比如0.1开始往上加每加一次就观察θ̂的收敛速度和s的大小。γ太大θ̂本身开始振荡并且会跟ESO的观测动态互相激发整个系统出现低频振荡。以我上面那个例子实测下来一组稳定的参数是ω₀ 30kp 100kd 20γ 5b₀ 2仿真采样步长0.001秒定步长求解器。这套参数跑出来方波跟踪的上升时间约为0.15秒超调小于5%稳态误差在0.01以内。5.3 仿真求解器设置定步长还是变步长这个问题在ADRC仿真里特别关键。ESO本质是高频观测器它需要足够小的仿真步长才能稳定。我的建议是ADRC仿真一律用定步长求解器比如ode4或者ode5步长按控制系统最高频率的1/501/100来选。如果你用变步长求解器ode45仿真会自动在ESO动态剧烈的地方加密步长结果是你无法复现实际嵌入式系统里的固定采样行为。况且最终你要把控制器部署到真实的MCU或者DSP上控制周期是固定的比如1kHz定步长仿真才贴近工程现实。定步长具体设置方法Simulink模型窗口 - Solver - Solver options - Type选Fixed-stepSolver选ode4四阶龙格库塔Fixed-step size填0.001。如果你用了MATLAB Function模块而且里面写了连续状态比如ESO的三个积分状态用ode4是够用的如果模型规模大、非线性强可以考虑ode5提高精度但仿真时间会明显变长。6. 仿真踩坑记录代数环、参数发散与离散化问题最后这部分是我实际调试中遇到的三个最典型的问题每个都花了不止一个晚上才解决。我把排查过程写出来希望能帮你省下这些时间。6.1 MATLAB Function模块引起的代数环用MATLAB Function模块写ESO时很多人第一次跑会弹出一个warning存在Algebraic Loop代数环。原因是ESO的输入里有y和u而u又由ESO的输出z3计算而来形成了一个u → ESO → z3 → u的环。Simulink默认会尝试用迭代法解这个代数环但迭代会显著拖慢仿真速度严重时直接导致仿真失败或者结果振荡。我的处理办法有两种。第一种最治本把ESO模型拆成基于积分器的状态方程实现用Integrator模块的输出作为状态输入输出信号分开这样Simulink就能以积分环节断开代数环。第二种做法是给ESO的输入u加一个单位延迟Unit Delay模块人为打破代数环。这相当于把控制量延迟一个仿真步长对于采样时间在毫秒级的系统影响非常小。工程上很多嵌入式实现本来就是上一周期的u参与本周期ESO计算反而更接近实际代码行为。6.2 自适应增益引发的参数发散自适应项刚加上去的时候系统100%会先发散一段时间——至少在我做的十几个不同工程里还没有一次是第一次跑就稳定的。最典型的症状是前0.1秒一切正常然后θ̂突然往一个方向狂奔比如涨到几百同时x也开始发散。原因通常是γ太大加上持续激励条件不满足或者符号错误。排查思路的建议先用纯PDESO不加自适应项跑通确认基础层没问题然后加入自适应项但把γ设到极小0.01再逐步增大γ并用Scope同时观察θ̂和s滑模面如果θ̂出现高频振荡说明γ和ESO带宽之间产生了耦合需要降低ω₀或γ两者之一。我还习惯在θ̂的积分器上加一个饱和限幅Saturation模块比如限制在[-10, 10]这样即使调试过程中参数出问题也不至于一次仿真就把数据搞到天文数字方便观察局面。6.3 ESO带宽与采样步长的匹配问题最后一个高频坑ESO的ω₀设得很大但仿真步长不够小结果就是离散化误差导致的数值发散。ESO作为高频观测器它的动态时间尺度是1/ω₀仿真步长至少要小于1/(10·ω₀)才能保证数值稳定。比如ω₀ 30步长至少要小于0.0033秒而我上面推荐0.001秒是留了安全余量的。你可以做个简单的实验对比同样的参数变步长ode45默认容差跑出来可能看起来很正常但一旦切成0.01秒的定步长系统立刻发散。这正是离散化误差的威力。我以前给一个课程设计做辅导的时候一个学生折腾了三天没解决的问题就是这个——他把ESO带宽设成了1000而仿真步长是0.001理论上1/(10·1000) 0.0001秒的步长要求他差了10倍。把ω₀降到300或者把步长改成0.0001秒问题立刻消失。6.4 一个工程建议把仿真当实验台把调参当实验设计最后分享一点个人体会。我刚接触ADRC的时候总觉得调参是一个试错的过程调到一个稳定点就算完事。后来做多了才意识到这就是一个典型的实验设计问题——每一个参数组合都是一个实验点你要做的不是碰运气找到能用的那个点而是理解每个参数在系统中所扮演的角色知道往哪个方向调会产生什么效果。比如ω₀在ESO里的角色是观测带宽它的物理意义是观测器对输出误差的放大程度对应z3对总扰动的跟踪速度。kp和kd的角色是闭环刚度和阻尼对应跟踪误差在PD层面的响应速度。γ的角色是参数自适应速度对应未知参数搜索的快慢。三个参数分别作用于三个时间尺度ω₀管快变量kp/kd管中速变量γ管慢变量。调参时只要按这个时间尺度的顺序来先快后慢、先观测后控制、先固定后自适应大部分问题都能有条理地解决而不是靠瞎试。这套从李雅普诺夫稳定性分析到自适应自抗扰控制的仿真链路我前后带过不少学生和同事复现过。每次有人卡住十有八九不是理论不懂而是仿真实现里的某个细节没处理好——代数环没破、步长没匹配、参数整定顺序反了。希望这篇东西能帮你把这几个坑提前填平。如果你在复现过程中遇到了这篇没覆盖到的问题欢迎带着你的模型参数和波形来交流。