简介:这是一份广州大学自动控制原理课程中“二阶系统阶跃响应及性能分析”实验的完整报告,面向自动化、电气等专业本科生及需要掌握时域分析的初学者。报告围绕二阶系统闭环传递函数、阻尼比ζ与自然频率ωn对动态性能的影响,以及开环增益和时间常数对稳定性的作用展开,详细记录了利用MATLAB绘制单位阶跃响应曲线、求解系统参数K与a、对比无速度反馈与带速度反馈系统瞬态性能的全过程。资源为单个PDF文件,共483KB,包含实验目的、实验内容、程序代码、图表结果及详细分析,适合对照课程要求进行预习、复习或撰写实验报告参考。已有3032人学习下载。通过阅读可清晰理解超调量、调节时间、上升时间等指标的计算方法,并掌握借助step、hold on、legend等MATLAB指令完成系统仿真的实用技能,对巩固自控原理核心概念和提升控制系统的分析与设计能力具有实际帮助。
1. 为什么同一个二阶系统,曲线形态会被两个参数“掰”成几种极端
拿到这份广州大学自动控制原理实验报告时,我先翻到结论页,发现学生对“阻尼比增大后超调量变大”的判断写反了。这类写反在二阶系统实验里很常见,根子在于没把 ζ 和 ωn 对应到特征方程的系数上。这份报告的主题是典型二阶系统单位阶跃响应与性能指标分析,全部用 MATLAB 的 step() 函数完成:先固定 ωn=6 扫 ζ,再固定 ζ=0.7 扫 ωn,然后由 σp=10%、ts(5%)=2s 反推开环参数 K、a,最后对比无速度反馈与带速度反馈系统的动态性能。适合正在做自控实验、课程设计,或者想搞明白速度反馈为什么能“加阻尼”的读者。
2. 复现典型二阶系统响应:ζ 与 ωn 对曲线的影响如何用 MATLAB 量化
2.1 标准闭环传递函数与批量绘图代码
原报告的程序写得比较原始:阻尼比用变量名 a 存,每条曲线复制一段代码,再用 hold on 叠图。如果只是复现一次,这个做法没问题;但要做参数扫描时,写成循环效率更高,也更不容易抄错。先写出二阶系统标准闭环传递函数:
Φ(s) = ωn² / (s² + 2ζωn·s + ωn²)
分母中 s 项系数是 2ζωn,常数项是 ωn²。原报告实验一第一组固定 ωn=6,把 ζ 依次取 0.1、0.4、0.7、1.0、1.3,用下面的循环代码一次画出:
wn = 6; % 自然频率 omega_n zeta_list = [0.1 0.4 0.7 1.0 1.3]; t = 0:0.01:15; % 仿真时间网格 0~15s,步长 10ms figure(1); hold on; for z = zeta_list num = wn^2; % 分子 = omega_n^2 den = [1 2*z*wn wn^2]; % s^2 + 2*zeta*wn*s + wn^2 sys = tf(num, den); step(sys, t); end grid on; xlabel('t / s'); ylabel('c(t)'); title('单位阶跃响应:\omega_n=6,\zeta 变化'); legend(arrayfun(@(z) sprintf('\\zeta=%.1f', z), zeta_list, 'UniformOutput', false)); hold off;tf(num, den) 按降幂顺序构造闭环传递函数。step(sys, t) 用给定时间网格做单位阶跃仿真,所有曲线共用同一个 t,方便直接互读峰值时间和调节时间。legend 用 arrayfun 和 sprintf 批量生成“ζ=0.1”这类曲线名,之后想多加一组阻尼比,只需要往 zeta_list 里添一个数。原报告手写五条 legend 字符串的方式也能跑,但改参数时要动好几处。
注意:原报告把阻尼比写成变量 a,分母写成 [1 12a 36],这是 ωn=6 时 2ζωn=12ζ 的特例。换成其它 ωn 后这个写法会出错,始终建议按 [1 2zwn wn^2] 显式构造。
2.2 ζ 变化时读出的规律:原报告有一处结论是反的
这组曲线的理论值可以在画图前先算出来,用如下三个公式:
σp = e^(-ζπ/√(1-ζ²)) tp = π/(ωn·√(1-ζ²)) ts(5%) ≈ 3/(ζωn)
把 ωn=6 代进去,得到:
| ζ | σp(理论) | tp(理论) | ts(5%)(理论) | 闭环极点 |
|---|---|---|---|---|
| 0.1 | 72.9% | 0.53s | 5.0s | -0.6 ± 5.97j |
| 0.4 | 25.4% | 0.57s | 1.25s | -2.4 ± 5.50j |
| 0.7 | 4.6% | 0.73s | 0.71s | -4.2 ± 4.20j |
ts≈3/(ζωn) 只在 0<ζ<0.9 的欠阻尼区近似成立。ζ=1 是临界阻尼,ζ=1.3 是过阻尼,这两种情况响应没有超调,曲线平滑逼近稳态值,但调节时间反而比 ζ=0.7 更长,因为闭环极点变成两个实极点,其中一个离虚轴更近,拖慢了收敛。
对照曲线可以看出正确规律:随着 ζ 增大,超调量从 72.9% 单调降到 0,峰值时间从 0.53s 缓慢增加到 0.73s,上升时间也变长。原报告结论里“超调量变大”是反的。唯一看起来“对”的是调节时间在 ζ 从 0.1 到 0.7 这一段确实从 5s 缩短到 0.71s,但这是因为起点选在了强欠阻尼区,如果继续增大 ζ,调节时间不会再明显缩短。把 ζ 调到 0.7 附近,超调量只有 4.6%,调节时间也最短,这就是教材常说“ζ=0.7 附近综合性能最好”的来由。
2.3 ζ=0.7 固定、ωn 从 2 到 10:时间轴被均匀压缩
第二组实验固定 ζ=0.7,ωn 取 2、4、6、8、10。代码结构与上面几乎一致,只要把 wn 变成循环变量:
z = 0.7; wn_list = [2 4 6 8 10]; t = 0:0.01:15; figure(2); hold on; for w = wn_list sys = tf(w^2, [1 2*z*w w^2]); step(sys, t); end grid on; xlabel('t / s'); ylabel('c(t)'); title('单位阶跃响应:\zeta=0.7,\omega_n 变化'); legend(arrayfun(@(w) sprintf('\\omega_n=%d', w), wn_list, 'UniformOutput', false)); hold off;这组曲线有一个容易被忽略的性质:ζ 相同,超调量就完全相同,都是 4.6%。ωn 只决定时间尺度。ωn=2 时峰值出现在约 2.2s,调节时间约 2.1s;ωn=10 时峰值出现在约 0.44s,调节时间约 0.43s,曲线形状只是整体“压扁”。原报告里写“自然频率变大,阻尼比变大”并不成立,ζ 和 ωn 是两个独立变量,增大 ωn 不会改变阻尼比,只会让系统对所有频率分量的响应都更快,代价是初段上升速率更猛,对噪声和延迟更敏感。这个区别在下一章的指标反推里会体现得更清楚。
3. 用性能指标反推参数:σp=10%、ts(5%)=2s 时如何定 K 和 a
3.1 指标到 ζ、ωn 的公式推导
实验二的原题是单位反馈系统,开环传递函数为 G(s)=K/[s(s+a)],闭环后:
Φ(s)=K/(s²+as+K)
与标准式对比,得到两个映射关系:K=ωn²,a=2ζωn。也就是说,题目里的 K 和 a 不是随便取的,它们直接决定了 ζ 和 ωn。反推参数的正确顺序是:先用性能指标解出 ζ、ωn,再映射回 K、a。
性能指标给了两个约束。超调量约束:
σp = e^(-ζπ/√(1-ζ²)) = 0.10
两边取对数后整理:
ζ = sqrt( (ln 0.1)² / (π² + (ln 0.1)²) ) = sqrt(5.30/15.17) ≈ 0.591
调节时间约束用 5% 误差带公式:
ts(5%) = 3/(ζωn) = 2s
所以 ζωn=1.5,ωn=1.5/0.591≈2.54。接着就能算出 K=ωn²≈6.45,a=2ζωn≈3.0。
注意:ts≈3/(ζωn) 对应 5% 误差带;如果题目要求 ts(2%),公式要换成 4/(ζωn)。两种误差带的常数从 3 变成 4,结果相差 33%,对参数设计影响很大。
3.2 用 step 绘制响应曲线,并标注 σp、ts
得到 K 和 a 后,直接构造闭环传递函数并画阶跃响应。下面的代码不用 stepinfo,而是手动从响应数据计算超调量和调节时间,逻辑透明一些:
K = 6.45; % omega_n^2,由指标反推得到 a = 3.00; % 2*zeta*omega_n,由指标反推得到 sys = tf(K, [1 a K]); t = 0:0.01:6; [y, t] = step(sys, t); plot(t, y, 'LineWidth', 1.2); grid on; hold on; [ymax, idx] = max(y); yss = y(end); sigma = (ymax - yss) / yss * 100; plot(t(idx), ymax, 'ro'); text(t(idx)+0.05, ymax, sprintf('\\sigma_p=%.1f%%', sigma)); outside = abs(y - yss) >= 0.05 * yss; idx_last_out = find(outside, 1, 'last'); idx_settle = idx_last_out + 1; ts_val = t(idx_settle); plot([ts_val ts_val], [0 y(idx_settle)], 'r--'); text(ts_val+0.05, y(idx_settle), sprintf('t_s=%.2fs', ts_val)); hold off;step 返回响应数组 y 和时间数组 t;max(y) 找到峰值,yss 取末值作为稳态值,超调量按定义算。outside 是一个逻辑数组,标记所有离开 5% 误差带的点;从后往前找最后一个越界点,它的下一个采样点就是系统最后进入误差带的时刻,也就是调节时间。这里从后往前扫是有意的:有些欠阻尼曲线在包络收敛前可能短暂穿进误差带又穿出,直接从前往后找会误判。
理论上 σp 正好是 10%,ts 是 2.0s,曲线实测值一般会在 9.8%~10.3%、1.9s~2.1s 之间波动,这是离散仿真步长和包络近似造成的。原报告实验二用的参数是 K=63.2、a=3.5。把这两个数代回公式验算:ωn=√63.2≈7.95,ζ=3.5/(2×7.95)≈0.22,对应的 σp≈49%、ts(5%)≈1.7s,和题目要求的 10%、2s 完全对不上。这是复现这份报告时最值得注意的地方:画图之前必须先用公式反解参数,否则曲线和性能指标标注会对不上。
| 数据来源 | ζ | ωn | σp(理论) | ts(5%)(理论) |
|---|---|---|---|---|
| 报告中 K=63.2, a=3.5 | 0.220 | 7.95 | 49.3% | 1.71s |
| 指标反推 K=6.45, a=3.0 | 0.591 | 2.54 | 10.0% | 2.00s |
4. 速度反馈改变阻尼:同一对象从 60% 超调压到 16%
4.1 无速度反馈系统 (a) 的阶跃响应
实验三给出两个结构:系统 (a) 是无速度反馈的单位反馈系统,前向通道是 10/[s(s+1)];系统 (b) 在它基础上增加了速度反馈。先看 (a) 的闭环传递函数:
Φa(s) = 10/(s²+s+10)
和标准式对比:ωn=√10≈3.16,s 项系数 1=2ζωn,所以 ζ=1/(2√10)≈0.158。这是一个典型的强欠阻尼系统,超调量和调节时间都很大。绘制并标注:
sys_a = tf(10, [1 1 10]); t = 0:0.01:15; step(sys_a, t); grid on; hold on; text(1.2, 1.30, '\\sigma_p=60.4%', 'FontSize', 11); text(5.5, 1.10, 't_s\\approx7s', 'FontSize', 11); hold off;理论上 σp=60.4%,峰值时间 tp=π/(3.16×√(1-0.025))≈1.0s,调节时间按 5% 误差带约为 6.0s。原报告在曲线上手动标了 ts=7s,这是因为图形标注时通常以曲线最后一次穿越某个阈值的位置为准,比包络公式偏大,另外 5% 带的包络判定和实际数据采样也有差异。工程上遇到这种标注不一致时,应该以代码计算的数值为准,图上标值只做可视化参考。
4.2 带速度反馈系统 (b):求解 K1 并对比指标
系统 (b) 的速度反馈等效于在闭环特征方程里增加一个虚拟阻尼项。假设速度反馈系数为 K1,前向通道仍为 10/[s(s+1)],反馈支路对输出做微分后反馈到输入端,得到闭环:
Φb(s) = 10/(s² + (1+K1)s + 10)
对比标准式,ωn 依然是 √10,但 2ζωn=1+K1。题目要求 ζ=0.5,所以:
K1 = 2×0.5×√10 - 1 = 2.16
原报告里 K1 取的是 2.2,程序里分母写成 [1 3.2 10],对应的 ζ=3.2/(2×3.16)≈0.506,非常接近 0.5,属于允许的分析误差。
把两个系统画在同一张图上对比:
K1 = 2.16; sys_b = tf(10, [1 (1+K1) 10]); step(sys_a, 0:0.01:5); hold on; step(sys_b, 0:0.01:5); grid on; legend('无速度反馈(a)', '带速度反馈(b)');代码里的分母第二项 (1+K1) 就是速度反馈带来的阻尼贡献。因为前向增益没变,稳态值都是 1,开环增益对应的稳态误差也不变,功劳全部体现在动态指标上。理论对比结果如下:
| 指标 | 系统(a) 无速度反馈 | 系统(b) 带速度反馈 |
|---|---|---|
| 阻尼比 ζ | 0.158 | 0.506 |
| 自然频率 ωn | 3.16 | 3.16 |
| 超调量 σp | 60.4% | 16.3% |
| 峰值时间 tp | 1.0s | 1.15s |
| 调节时间 ts(5%) | ≈6.0s | ≈1.9s |
可以看到,速度反馈没有动 ωn,只是把 ζ 从 0.158 提到 0.506,超调量就下降了 44 个百分点,调节时间缩短到三分之一。这正是工程上通过测速发电机或编码器微分反馈改善动态性能的原理:系统结构上多了一个与输出变化率成正比的负反馈,等效于向特征方程注入了更多阻尼,而不用降低开环增益去换稳定。相比单纯调增益,这种做法的代价是增加一个反馈环节,但换来的性能和稳态精度都要好。
顺带提一个扩展点:把阻尼比从 0.5 提到 0.7,速度反馈系数应变成 K1=2×0.7×√10-1≈3.42,超调量会降到 4.6%,但峰值时间和上升时间会再往后推。也就是说,速度反馈系数不是越大越好,K1 调得太大会让系统过阻尼、反应迟钝。这个“调 K1 观察超调与速度变化”的过程,可以直接用前面的循环脚本改成参数扫描,是检验自己是否理解阻尼作用的简单方式。
5. 用 stepinfo 自动核对指标,并避开复现中的四个坑
5.1 stepinfo 一行读出全部指标
前面手动找峰和误差带的方法适合教学演示,但做工程验证时要更快更准确。MATLAB 的 stepinfo 专门干这件事:
sys_design = tf(6.45, [1 3.0 6.45]); S = stepinfo(sys_design, 'SettlingTimeThreshold', 0.05); fprintf('Overshoot = %.1f%%\n', S.Overshoot); fprintf('SettlingTime = %.2fs\n', S.SettlingTime); fprintf('PeakTime = %.2fs\n', S.PeakTime);注意 SettlingTimeThreshold 指定误差带宽度,传 0.05 表示 5% 带。如果不传,stepinfo 默认按 2% 误差带算调节时间,即 4/(ζωn),和题目要的 ts(5%)=3/(ζωn) 直接对比会差一截。返回的结构体字段主要有:
| 字段 | 含义 |
|---|---|
| RiseTime | 从 10% 到 90% 稳态值的上升时间 |
| PeakTime | 峰值时间 tp |
| Overshoot | 超调量百分比 |
| SettlingTime | 进入并保持误差带的时间 |
| Peak | 响应峰值 |
把 stepinfo 的结果和理论公式放在一张表里逐项比对,能快速定位“是公式用错还是参数映射错”。例如第三章正确参数算出来 Overshoot≈9.97%,SettlingTime≈1.98s,与理论值 10%/2.0s 吻合。
5.2 复现原报告时的四个坑
系数位对应错。常见把分母写成 [1 2wn zetawn wn^2],少乘一个 2。判断方法是:ζ=0.7、ωn=6 时,s 项系数必须是 8.4,如果写错,曲线就会呈现虚假的阻尼比。代码里用显式表达式 2zw 而不是手算常数,可以避免这类笔误。
误差带混淆。5% 带是 3/(ζωn),2% 带是 4/(ζωn)。拿系统 (a) 说,ζ=0.158、ωn=3.16 时两者分别是 6.0s 和 8.0s。报告里标注 7s 说明是按曲线实测读的,没有严格指定误差带。
仿真时间长度不够。强欠阻尼系统调节时间接近 5s,若 t 只取 0:0.01:3,曲线还没进误差带就结束了,会误判系统“反应慢”。先按 3/(ζωn) 估算调节时间的 2 倍作为 t 的上限,再画图。
roots(c) 只证明极点位置。原报告每段程序都算了 p=roots(c),但没输出。roots(c) 得到所有闭环极点,判断“稳定”只需全部负实部;但工程系统还要看阻尼比是否落在 0.4~0.8。ζ=0.1 的极点也在左半平面,系统稳定却超调 73%,不能用“稳定”二字交差。
把这套验证逻辑收成一个函数:输入 num、den、误差带,输出阶跃响应图和 stepinfo 指标表。之后处理三阶系统、带零点系统、延迟系统时,先化成标准二阶或其他可参照形式,再套同一套验证流程,比每次重新猜参数高效得多。
本文还有配套的精品资源,点击获取