首页 / 资讯中心 / 文章详情

自动控制原理实验全攻略:MATLAB与MWORKS从建模到校正

自动控制原理实验全攻略:MATLAB与MWORKS从建模到校正 ★ FEATURED ARTICLE
简介这份文档面向具备基础控制理论与编程能力的高年级本科生、研究生及工程技术人员围绕自动控制原理实验教学系统讲解控制系统的建模、分析与校正方法。内容涵盖控制系统数学模型的MWORKS与MATLAB实现、时域与频域分析、典型环节动态性能与稳定性分析、频率特性测量、系统串联校正以及EF-Servo集成与滤波六个实验每个实验均配有目的、内容、步骤与思考题帮助读者在动手实践中掌握MATLAB与MWORKS工具的使用。资源包为1个docx文档大小约1.43MB结构完整、条理清晰便于按实验模块逐步学习。目前已有676人学习下载。读者可借助其中的传递函数、状态空间、Simulink建模及频域校正等具体案例对照操作并完成思考题从而加深对控制系统工作原理的理解提升系统设计与分析能力。1. 从传递函数到仿真波形这套自动控制原理实验资源到底能干什么很多控制理论的学习者都有过这样的体验课堂上把劳斯判据、奈奎斯特判据背得滚瓜烂熟一到要用 MATLAB 或 MWORKS 把传递函数跑出阶跃响应曲线就不知道从哪下手。这份自动控制原理实验资源就是冲着这个断层来的。它把六个实验串成一条完整的动手链路从传递函数、零极点增益、状态空间三种模型的建立与互相转换到二阶系统的时域指标测量再到频率特性测量和串联校正设计最后落到 EF-Servo 集成与滤波。每个实验都配了明确的实验目的、操作步骤和思考题工具链覆盖 MATLAB 和 MWORKS 两套环境。适合正在上控制理论实验课的高年级本科生和研究生也适合需要快速把控制理论捡起来做系统分析的工程技术人员。如果你之前只会在 Simulink 里拖模块、不会用脚本批量分析或者反过来只会敲命令、不会搭仿真模型这份资源正好能补上你缺的那一半。2. 数学模型建立与转换tf、zpk、ss 三种模型怎么选怎么转2.1 三种 LTI 模型各自的适用场景线性定常系统在 MATLAB 和 MWORKS 里有四种描述方式传递函数模型、零极点增益模型、状态空间模型以及 Simulink 结构图。前三种是数学表达式统称 LTI 模型都有连续和离散两种类型。选哪种模型取决于你接下来要做什么分析。传递函数模型最直观适合做频域分析和串联校正设计。零极点增益模型在分析稳定性时最方便因为极点的位置直接决定系统是否稳定一眼就能看出来。状态空间模型适合处理多输入多输出系统也适合做状态反馈和观测器设计。Simulink 结构图则适合验证非线性环节或复杂系统的时域响应。常见做法是先用 tf 建立传递函数需要判断稳定性时转成 zpk 看极点分布需要做状态反馈时转成 ss。三种模型之间可以通过 tf()、zpk()、ss() 三个函数互相转换。需要注意的是传递函数模型和零极点增益模型是唯一的但状态空间模型不唯一同一个系统可以有无穷多种状态空间实现。2.2 建立模型与转换的完整操作先看传递函数模型的建立。假设系统传递函数为 G(s) 6 / (s³ 6s² 11s 6)在 MWORKS 的 Julia 环境中这样写using TyControlSystems using TyBase # 分子系数按 s 的降幂排列 num 6 # 分母系数按 s 的降幂排列 den [1 6 11 6] # 建立传递函数模型 G tf(num, den)这里 num 和 den 都是按 s 降幂排列的系数向量。如果分子是多项式而不是常数比如 s1就写成 [1 1]。分母的系数必须完整缺项要补零比如 s³ 6s² 11s 6 对应 [1 6 11 6]不能省略任何一项。再看零极点增益模型。同一个系统写成零极点形式是 G(s) 6 / [(s1)(s2)(s3)]建立方式如下# 零点为空因为分子是常数 z [] # 极点位置 p [-1, -2, -3] # 增益 k 6 # 建立零极点增益模型 G zpk(z, p, k)注意无零点时 z 要设成空数组不能省略。极点向量里每个元素对应一个极点重极点就重复写比如两个 -1 极点写成 [-1, -1]。模型转换是实验一的核心操作。下面这段代码演示了从传递函数出发转成零极点增益模型和状态空间模型再把状态空间模型转回传递函数using TyControlSystems using TyBase clear() clc() # 原始传递函数 G(s) 5 / (s^3 4s^2 5s 2) num 5 den [1 4 5 2] Gtf tf(num, den) # 转换为零极点增益模型 Gzpk zpk(Gtf) print(Gzpk) # 转换为状态空间模型 Gss ss(Gtf) print(Gss) # 从状态空间模型转回传递函数 Gtf1 tf(Gss) print(Gtf1)运行后会看到转回传递函数的系数和原始系数在数值上几乎一致但会有极小的浮点误差比如 8.88e-16 这种量级。这是正常的数值计算误差不影响分析结果。如果你发现转换后的传递函数和原始传递函数差异很大那就要检查是不是在转换过程中把某个参数搞错了。2.3 子系统连接串联、并联、反馈实验一还要求对两个子系统做串联、并联和反馈连接。假设两个子系统的传递函数分别是 G₁(s) (s1) / (s² 4s 4) 和 G₂(s) 1 / (s1)操作如下using TyControlSystems # 子系统1 num1 [1 1] den1 [1 4 4] G1 tf(num1, den1) # 子系统2 num2 1 den2 [1 1] G2 tf(num2, den2) # 串联连接 G_series G1 * G2 # 并联连接 G_parallel G1 G2 # 反馈连接默认负反馈 G_feedback feedback(G1, G2)串联用乘法并联用加法反馈用 feedback() 函数。feedback(G1, G2) 表示前向通道是 G1、反馈通道是 G2 的负反馈连接。如果是正反馈需要加参数指定具体写法查一下当前版本的文档确认。提示串联和并联的顺序不影响结果但反馈连接中前向通道和反馈通道的位置不能互换否则传递函数完全不同。3. 时域与频域分析从阶跃响应曲线里读出性能指标3.1 二阶系统时域分析的参数扫描方法实验二的核心是研究阻尼比 ζ 和无阻尼振荡频率 ωₙ 对二阶系统动态性能的影响。典型二阶系统的闭环传递函数是G(s) ωₙ² / (s² 2ζωₙs ωₙ²)分析思路是固定一个参数、扫描另一个参数观察阶跃响应曲线的变化。先固定 ωₙ 5让 ζ 分别取 0、0.5、0.707、1、1.5using TyControlSystems using TyPlot w 5 zita_list [0, 0.5, 0.707, 1, 1.5] for zita in zita_list num w^2 den [1, 2*zita*w, w^2] G tf(num, den) step(G) hold(on) end这段代码把五条阶跃响应曲线画在同一张图上。ζ 0 时系统等幅振荡ζ 0.5 和 0.707 时系统欠阻尼、有超调ζ 1 时临界阻尼、无超调ζ 1.5 时过阻尼、响应缓慢。从曲线上可以直接读出上升时间、峰值时间、超调量和调节时间。再固定 ζ 0.5让 ωₙ 分别取 0.5、1、3、5zita 0.5 w_list [0.5, 1, 3, 5] for w in w_list num w^2 den [1, 2*zita*w, w^2] G tf(num, den) step(G) hold(on) endωₙ 越大响应越快上升时间和峰值时间越短但超调量基本不变。这就是为什么在设计控制系统时阻尼比决定超调量无阻尼振荡频率决定响应速度。3.2 状态空间模型的时域响应求解实验二还要求对状态空间模型求解状态响应和输出响应。给定系统ẋ [0 1 0; 0 0 1; -6 -10 -4] x [0; 0; 1] u y [6 0 0] x初始状态为零输入为单位阶跃时用 step() 求解using TyControlSystems using TyPlot A [0 1 0; 0 0 1; -6 -10 -4] B [0; 0; 1] C [6 0 0] D 0 G ss(A, B, C, D) t 0:0.01:10 # 单位阶跃响应 y, t_out, x step(G, t) plot(t_out, x, :) plot(t_out, y, -)step() 返回三个值输出响应 y、时间向量 t_out、状态向量 x。如果 G 是状态空间模型x 就是状态轨迹如果是传递函数模型x 返回空矩阵。输入为单位斜坡时用 lsim()u t y, t_out, x lsim(G, u, t) plot(t_out, x, :) plot(t_out, y, -)lsim() 的第二个参数是输入向量长度必须和时间向量一致。这里 u t 就是单位斜坡信号。如果要分析初始状态响应用 initial() 函数传入初始状态向量 x0 即可。3.3 频域分析奈奎斯特图与伯德图实验二还涉及频域稳定性分析。给定开环传递函数 G(s) 2 / (s³ 9s² 14s)绘制奈奎斯特图using TyControlSystems using TyPlot num 2 den [1 9 14 0] G tf(num, den) # 绘制奈奎斯特图 nyquist(G) # 绘制伯德图并显示稳定裕度 bode(G) margin(G)奈奎斯特判据的核心是看开环频率特性曲线是否包围 (-1, j0) 点。如果包围闭环系统不稳定如果不包围闭环系统稳定。用 pole() 函数可以直接求出闭环系统的极点验证奈奎斯特判据的结论G_cl feedback(G, 1) p pole(G_cl)如果所有极点的实部都小于零闭环系统稳定。两种方法的结果应该一致如果不一致检查一下开环传递函数是否写错或者反馈极性是否搞反。注意nyquist() 画出的曲线在低频段和高频段的行为需要特别关注。低频段对应 s → 0高频段对应 s → ∞。如果曲线穿过 (-1, j0) 点附近需要放大看细节。4. 典型环节电路模拟与频率特性测量从电路到传递函数4.1 典型环节的模拟电路与传递函数对应关系实验三和实验四的核心是把模拟电路和传递函数对应起来。比例环节的传递函数是 K积分环节是 1/(Ts)比例积分环节是 K(1 1/(Ts))比例微分环节是 K(1 Ts)惯性环节是 K/(Ts1)比例积分微分环节是 K(1 1/(Ts) Td·s)。以积分环节为例模拟电路中 R₀ 100kΩC 1μFR 10kΩ时间常数 T R₀C 0.1s。在 MWORKS 中建立传递函数并仿真阶跃响应using TyControlSystems using TyPlot # 积分环节 G(s) 1 / (0.1s) num 1 den [0.1 0] G tf(num, den) step(G)惯性环节的电路参数是 R₁ 200kΩR₀ 200kΩC 1μF传递函数为 G(s) 1 / (0.2s 1)num 1 den [0.2 1] G tf(num, den) step(G)比例积分环节的传递函数是 G(s) K(1 1/(Ts))其中 K R₁/R₀T R₀C。如果 R₁ 100kΩR₀ 100kΩC 1μF则 K 1T 0.1sK 1 T 0.1 num K * [T 1] den [T 0] G tf(num, den) step(G)比例微分环节的实际电路是近似的传递函数为 G(s) K(1 Ts)其中 K (R₁R₂)/R₀T R₁R₂C/(R₁R₂)。如果 R₁ 10kΩR₂ 10kΩR₀ 10kΩC 10μF则 K 2T 0.05sK 2 T 0.05 num K * [T 1] den [1] G tf(num, den) step(G)4.2 频率特性测量与传递函数求取实验四要求用实验方法测量一阶惯性环节和典型二阶系统的频率特性曲线然后根据曲线求取传递函数。一阶惯性环节 G(s) K/(Ts1) 的幅相频率特性曲线是一个半圆。从半圆的直径可以得到放大倍数 KK r(0)。在曲线上取一点 ωₖ可以确定时间常数 T tan φ(ωₖ) / ωₖ。用 MWORKS 绘制一阶惯性环节的频率特性using TyControlSystems using TyPlot K 1 T 0.2 num K den [T 1] G tf(num, den) # 绘制奈奎斯特图 nyquist(G) # 绘制伯德图 bode(G)从奈奎斯特图上读出半圆的直径应该等于 K 1。在某个频率点 ωₖ 处读出相角 φ(ωₖ)代入 T tan φ(ωₖ) / ωₖ 计算时间常数应该接近 0.2。对于典型二阶系统开环传递函数为 G(s) K / [s(Ts1)]频率特性测量需要记录不同频率下的幅值和相角然后画出伯德图从图上读出穿越频率和相位裕度。提示实验测量时频率点的选取很关键。在转折频率附近要多取几个点因为幅值和相角变化最快。低频段和高频段可以少取几个点。5. 系统串联校正与 EF-Servo 集成从理论设计到工程落地5.1 频率法串联校正的设计步骤实验五要求用频率法设计串联校正装置。给定开环传递函数 G(s) K / [s(s1)(s2)]要求设计串联校正环节使系统满足一定的相位裕度和截止频率要求。设计步骤如下先画出未校正系统的伯德图确定当前的相位裕度和截止频率。然后根据目标相位裕度计算需要的相位超前量。再根据相位超前量确定校正装置的参数。最后验证校正后系统的性能指标是否满足要求。用 MATLAB 实现这个流程% 未校正系统 num 1; den conv([1 0], conv([1 1], [1 2])); G tf(num, den); % 绘制伯德图 bode(G); margin(G); % 设计超前校正装置 % 假设目标截止频率为 2 rad/s wc 2; % 计算未校正系统在 wc 处的幅值和相角 [mag, phase] bode(G, wc); % 计算需要的相位超前量 phi_needed 60 - (180 phase); % 根据相位超前量计算校正装置参数 alpha (1 sind(phi_needed)) / (1 - sind(phi_needed)); T 1 / (wc * sqrt(alpha)); % 校正装置传递函数 Gc tf([alpha*T 1], [T 1]); % 校正后系统 G_comp Gc * G; bode(G_comp); margin(G_comp);这段代码的关键参数是目标截止频率 wc 和期望相位裕度。wc 决定了系统的响应速度相位裕度决定了系统的超调量和稳定性。如果校正后系统的相位裕度仍然不够可以适当降低 wc或者改用滞后-超前校正。5.2 EF-Servo 集成与滤波实验六涉及 EF-Servo 集成和滤波。EF-Servo 是一种伺服控制系统集成过程中需要考虑滤波器的设计。滤波器的目的是抑制高频噪声但也会引入相位滞后影响系统的稳定性。常见做法是在反馈通道中加一阶低通滤波器传递函数为 Gf(s) 1 / (Tf·s 1)。Tf 的选择需要在噪声抑制和相位滞后之间折中。Tf 越大高频噪声抑制越好但相位滞后越严重。% 滤波器设计 Tf 0.01; Gf tf(1, [Tf 1]); % 集成到系统中 G_open G_comp * Gf; G_closed feedback(G_open, 1); % 阶跃响应 step(G_closed);如果发现系统响应出现高频振荡说明滤波器截止频率太高需要增大 Tf。如果发现系统响应变慢、超调增大说明滤波器截止频率太低需要减小 Tf。注意滤波器的引入会改变系统的开环频率特性必须重新检查相位裕度。如果相位裕度低于 45°系统可能会不稳定。6. 避坑与排查六个实验里最容易翻车的五个地方6.1 传递函数系数顺序写反现象建立的传递函数模型和预期完全不符阶跃响应曲线形状奇怪。原因num 和 den 的系数必须按 s 的降幂排列缺项要补零。很多人会把 s³ 6s² 11s 6 写成 [1 6 11 6]但把 s³ 6s² 6 写成 [1 6 6]漏掉了 s 项的零系数。解决写系数向量之前先把传递函数按降幂排列写完整缺项补零。写完后再数一下向量的长度分母向量的长度应该等于系统阶数加一。6.2 反馈连接的正负号搞反现象闭环系统的阶跃响应发散或者和预期完全相反。原因feedback() 函数默认是负反馈但如果前向通道和反馈通道的位置写反了或者系统本身是正反馈结果就会出错。解决确认 feedback(G1, G2) 中 G1 是前向通道、G2 是反馈通道。如果是正反馈查一下当前版本的 feedback() 函数是否支持正反馈参数。另外如果反馈通道的传递函数是 1直接写 feedback(G, 1) 即可。6.3 状态空间模型转换后数值误差过大现象从状态空间模型转回传递函数后系数和原始值差异很大。原因状态空间模型不唯一不同的转换算法可能得到不同的数值结果。如果系统本身是病态的或者极点非常接近数值误差会被放大。解决检查系统的极点分布如果极点非常接近考虑用 zpk 模型做中间转换。另外可以用 ss(G, minimal) 先做最小实现去掉不可控或不可观的状态。6.4 频率特性测量时频率点选取不当现象画出的奈奎斯特图或伯德图不光滑关键特征点被跳过。原因频率点选取太稀疏尤其是在转折频率附近。解决在转折频率附近加密频率点。可以用 logspace() 生成对数等间距的频率向量比如 logspace(-2, 2, 1000) 生成从 0.01 到 100 的 1000 个频率点。6.5 串联校正后相位裕度不达标现象校正后系统的相位裕度仍然低于目标值或者超调量过大。原因超前校正装置的相位超前量有限一般不超过 60°。如果需要的相位超前量超过 60°单级超前校正无法满足要求。解决改用两级超前校正或者改用滞后-超前校正。另外可以适当降低目标截止频率牺牲一点响应速度来换取更大的相位裕度。7. 从脚本到报告一套可复用的实验数据整理方法做完六个实验手头会积累大量的响应曲线和性能指标数据。如果每次都用鼠标从图上读点效率低还容易出错。我一般会写一个统一的脚本把关键指标自动算出来直接填进实验报告。以二阶系统为例下面这段代码自动计算上升时间、峰值时间、超调量和调节时间function [tr, tp, Mp, ts] step_metrics(G) % 求解阶跃响应 [y, t] step(G); % 稳态值 y_ss y(end); % 峰值时间和超调量 [y_peak, idx_peak] max(y); tp t(idx_peak); Mp (y_peak - y_ss) / y_ss * 100; % 上升时间从 10% 到 90% idx_10 find(y 0.1 * y_ss, 1); idx_90 find(y 0.9 * y_ss, 1); tr t(idx_90) - t(idx_10); % 调节时间进入并保持在 2% 误差带内 idx_settle find(abs(y - y_ss) 0.02 * y_ss, 1, last); if isempty(idx_settle) ts 0; else ts t(idx_settle); end end这个函数的输入是 LTI 模型 G输出是四个性能指标。注意调节时间的计算用的是 2% 误差带如果实验要求 5% 误差带把 0.02 改成 0.05 即可。上升时间的定义是从 10% 到 90%有些教材用 0% 到 100%需要根据实验要求调整。调用方式zita 0.5; w 5; num w^2; den [1, 2*zita*w, w^2]; G tf(num, den); [tr, tp, Mp, ts] step_metrics(G); fprintf(上升时间: %.4f s\n, tr); fprintf(峰值时间: %.4f s\n, tp); fprintf(超调量: %.2f%%\n, Mp); fprintf(调节时间: %.4f s\n, ts);把不同参数下的计算结果整理成表格对比分析就一目了然了。比如固定 ωₙ 5改变 ζ阻尼比 ζ上升时间 (s)峰值时间 (s)超调量 (%)调节时间 (s)0.50.480.7316.32.900.7070.550.974.31.851.00.78无峰值01.521.51.12无峰值02.10这张表比六张单独的响应曲线图更有说服力因为数字是精确的不依赖读图误差。我后来养成了一个习惯每次做完仿真先把性能指标算出来存成表格再画曲线。曲线用来直观展示趋势表格用来精确对比。两者配合实验报告的数据部分就非常扎实了。从那以后我每次做控制实验都强制自己先写指标计算脚本再画图。这样即使参数扫描几十组也能快速定位到最优参数。希望这套方法帮到你。本文还有配套的精品资源点击获取
阅读完成 · 觉得有帮助?
咨询建站