还在学校里啃电力系统分析的时候最头疼的事情之一就是“暂态稳定”这一块。课本上把等面积法则、摇摆曲线、极限切除时间解释得清清楚楚可一旦让你回答“三相短路到底比两相接地严重多少”“故障发生在线路首端和末端差别有多大”这类问题光靠手算根本算不过来。到后来做仿真才知道电力系统暂态稳定性仿真分析——用 Matlab 编程配合 Simulink 模型去模拟各类故障对系统的冲击——是把这个抽象问题落到实处的唯一靠谱路径。这篇内容就是我个人在这个方向上的完整复盘从建模思路、故障模拟、代码实现到结果解读把能踩的坑和能抄的作业都整理出来希望对正在做相关课程设计、毕业设计或者入门研究的你有点实际帮助。1. 为什么要专门研究“暂态稳定”这个方向1.1 三类稳定性问题先把概念拎清楚在展开仿真之前必须先分清电力系统稳定性研究里的几个层次。静态稳定关注的是系统在某个稳态运行点附近能否抵抗微小扰动而不失去同步动态稳定一般指的是在控制器励磁调节器、调速器等作用下系统在受到扰动后恢复到稳定运行状态的能力暂态稳定则是针对大扰动——也就是各类短路故障、断线、大容量机组跳闸——发生后的几秒甚至几十秒内发电机转子之间能否保持同步的问题。这个区分极其重要因为很多初学者一上来就在 Simulink 里搭了个模型加了故障模块结果出来的波形乱七八糟就开始怀疑模型错了。其实大概率是把研究范畴搞混了你做暂态稳定分析关心的核心物理量是发电机的功角功率角随时间的变化轨迹而不是母线电压的谐波畸变率也不是潮流分布的变化。暂态稳定研究的本质是“转子动力学”问题。发电机转子是一个具有惯性的物理实体当系统发生短路时电磁功率输出瞬间骤降而原动机输入的机械功率来不及变化于是转子上出现不平衡转矩转子开始加速。如果故障切除得足够快加速过程中积累的动能可以在故障切除后被减速区域消化系统回到稳定如果故障切除得太慢转子加速过度系统就失稳了。这个“快慢”的临界值就是极限切除时间CCT, Critical Clearing Time。1.2 故障如何一步步把系统推向失稳边缘从物理过程来看一个典型的三相短路故障对电力系统的影响可以拆解成三个阶段。故障发生瞬间系统网络结构突变发电机端电压急剧下降发电机的电磁功率从故障前的约 1.0 标幺值掉到接近零机械功率不变转子开始加速功角开始增大。故障切除后系统网络结构恢复或切换到备用线路发电机电磁功率回升但只要功角已经摆得过大减速面积不足转子仍可能越过不稳定平衡点进入失步状态。不同故障类型对系统的冲击程度差异很大。三相短路是金属性对称短路扰动最强功率传输能力几乎完全被破坏两相短路和两相接地短路属于不对称故障故障期间部分功率还能传输冲击相对小一些单相接地短路虽然发生概率最高但扰动强度最低。实际系统中单相接地占全部短路故障的百分之八十以上但其对暂态稳定的威胁反而小于概率低得多的三相短路。这就形成了一个很有意思的研究点概率和风险如何权衡也是故障影响研究里值得深挖的内容。这些问题用手算只能处理单个简单场景要做“各类故障影响研究”这种系统性对比就必须依靠仿真工具了。2. 仿真方案怎么定Matlab 编程与 Simulink 模型各有分工2.1 为什么选 Matlab Simulink 这套组合电力系统的暂态仿真工具有不少专业的机电暂态程序有很多商业软件电磁暂态方向也有专门的工具。但我的选择是 Matlab 加 Simulink原因有三个。第一Matlab 的矩阵计算能力让“批量仿真”这件事变得非常舒服。故障类型有四五种故障位置可以沿线路上每隔一段设置一个点故障持续时间从 0.1 秒到 0.5 秒扫一遍再加上不同的运行方式组合出来的仿真场景轻松上百个。这种批量处理的工作流用脚本驱动 Simulink 模型自动运行再自动收集结果整套流程在同一个环境下闭环极其顺手。第二Simulink 里的 SimPowerSystems现在叫 Simscape Electrical模块库提供了同步电机、变压器、线路、负荷、故障发生器等现成模块。你不需要自己推导和编写全部的微分方程而是像搭积木一样把系统拓扑搭出来这对于验证理论模型和快速获得直观结果非常重要。尤其是可视化方面功角曲线、转速偏差、电磁功率、机端电压等波形可以直接在 Scope 里看到调整参数后对比非常直观。第三Matlab 编程本身可以做纯数值仿真。用经典二阶模型摇摆方程写一个 m 文件几十行代码就能对单机无穷大系统做完整个暂态稳定计算这个过程对理解暂态稳定的物理本质、对编写自己的仿真工具很有帮助。两条路线互为验证得到的结论更有说服力。我个人的建议是两种方式都要做先写代码理解算法再搭模型验证工程细节。2.2 两套方案各自适合的场景需要明确的是Matlab 编程和 Simulink 模型不是替代关系它们的侧重点完全不同。纯 Matlab 编程m 文件适合做核心算法验证和批量参数扫描。它的优势在于你完全掌控每一步计算过程可以方便地修改微分方程模型、加入不同的控制环节励磁系统模型、PSS 等。坏处是如果系统复杂到多机多节点编程模型会迅速变得繁琐光是对应节点导纳矩阵的修改就够你调试一阵子。Simulink 模型适合做系统级仿真和“看得见”的波形分析。模块化搭出来的系统直观故障模块设置方便可以模拟各种不同类型的故障包括各种短路和断线。它的内部算法比如变步长求解器已经比较成熟数值稳定性好。但正因为模块封装程度高一旦仿真发散或者结果不合理排查问题的难度也比较大需要一定的数值计算功底才能定位问题根源。我的建议是两套方案同时推进互相验证。编程模型算出一个极限切除时间Simulink 模型再去验证一下两者的结果应该非常接近误差在 5% 以内是正常的。如果差得太多要么模型参数不一致要么哪个环节的设置出了问题。2.3 系统模型单机无穷大系统是研究起点不管是编程还是 Simulink第一套模型都建议用单机无穷大系统Single Machine Infinite Bus, SMIB。无穷大系统的意思是母线电压幅值和频率恒定相当于一个容量无限大的系统。实际中不可能有真正的无穷大系统但当发电机容量远小于受端系统容量时比如一台 300MW 机组接入一个几千 MW 的区域电网近似为单机无穷大系统是足够准确的。SMIB 模型的优势在于所有核心概念都有解析解做对照。等面积法则可以手算极限切除时间的近似值经典二阶模型的摇摆方程可以数值求解所以仿真结果的正确性容易验证。把 SMIB 搞透了再把单机扩展成双机系统、三机九节点系统思路完全一致。仿真参数的设定要注意标幺值系统。发电机的暂态电抗 xd 通常取 0.2 到 0.3 标幺惯性时间常数 H 取 2 到 8 秒汽轮发电机组典型值约 3 到 6 秒水轮机组约 5 到 8 秒变压器漏抗约 0.1 标幺线路电抗视长度而定220kV 线路典型值在 0.4 欧姆每公里左右折算到标幺值后一般在 0.05 到 0.2 之间。这些参数直接决定了仿真结果是否合理也需要和你研究的具体对象匹配。3. Simulink 建模实操从零搭出带故障注入的暂态仿真系统3.1 模型搭建全过程以单机无穷大系统为例打开 Simulink 新建一个模型先放一个 Synchronous Machine 模块在 Simscape Electrical 的 Machines 库下面参数设置里选择 pu 单位模式设定额定功率比如 100MVA、额定电压比如 13.8kV、惯性时间常数 H 和阻尼系数。阻尼系数 D 是个容易忽略的参数它代表转子运动中的机械阻尼和电气阻尼的综合效果取 0 到 10 之间的值太小了曲线振荡衰减慢太大了会掩盖系统的真实失稳特性。然后是汽轮机调速器或者简化的原动机输入。在暂态稳定仿真里机械功率 Pm 通常可以假设为恒定值因为调速器的动作时间常数秒级远比暂态过程亚秒级慢所以直接用常量模块给一个恒定机械功率即可。如果要做得更精细可以加一个 Hydraulic Turbine and Governor 或 Steam Turbine and Governor 模块但第一版不建议加先把主脉络跑通。变压器用 Three-Phase TransformerTwo Windings连接发电机出口的 13.8kV 母线升压到 220kV 或 500kV。线路用 PI Section Line 模块一条或两条并联。无穷大母线用 Three-Phase Source 模块代替设好电压幅值和频率220kV50Hz内阻抗设一个很小的值比如 0.001 欧姆来模拟理想电压源。负荷可以用静态阻抗负荷或恒定功率负荷暂态分析里一般简化为恒阻抗。这里有一个新手特别容易忽视的点所有模块的基准频率必须设为一致50Hz 或 60Hz各模块的基准容量也要协调。Powergui 模块一定要放在模型里它相当于整个仿真环境的管理员。在 Powergui 里把仿真类型设为 Continuous求解器建议用 ode23tb 或 ode15s这类刚性求解器适合电力电子和电机模型因为时间常数差异很大的情况。3.2 各类故障的模拟方法与参数选择Simulink 里模拟故障靠的是 Three-Phase Fault 模块。这个模块可以设置各相的状态A、B、C 相的开断/接地组合用来模拟各种类型。关键设置包括Transitional 时间故障开始时间和结束时间、Transitional 电阻通常设为 0.001 欧姆代表金属性短路、接地电阻。具体来说模拟三相短路把 A、B、C 三相全部选中并且选中 Ground 选项表示三相同时接地短路。两相短路比如只选中 A 相和 B 相不选 Ground代表两相之间的相间短路如果同时选 Ground就是两相接地短路。单相接地短路只选 A 相加 Ground对应最常见的单相接地故障。断线故障的模拟要稍微绕一点用 Three-Phase Breaker 模块或者手动修改线路连接来实现某相断开的效果。故障位置对仿真结果的影响很大这一点在建模时必须认真考虑。故障点与发电机之间的电气距离越近故障对发电机的影响越大。所以建议在模型中设置多个故障位置参考点比如发电机机端母线最严重、升压变压器高压侧出口、线路首端离发电机 10% 处、线路中点50% 处、线路末端90% 处。在每个位置各自放一个故障模块或者通过编程动态修改故障模块的位置参数。故障持续时间的设置是研究极限切除时间的关键操作。对于固定故障场景采用二分法或者扫描法先设一个较大的故障持续时间比如 0.3 秒观察系统是否失稳再设一个较小的值比如 0.1 秒观察是否稳定然后不断缩小范围逐步逼近临界值。这个过程用手动操作太痛苦后面的章节会讲怎么用 Matlab 脚本自动完成。3.3 从功角曲线中读出故障特性的门道仿真跑完后的第一件事是把发电机的功角曲线Rotor angle delta 随时间的变化调出来看。正确建立模型后你会看到几条典型的响应形态。第一条形态是“先上扬后回落再衰减振荡”。功角在故障发生瞬间开始上升故障切除后爬升斜率变缓到达最大值后回落然后围绕新的稳态功角做衰减振荡振幅越来越小。这是稳定状态的曲线说明故障在极限切除时间内被清除了系统能够回到稳态。第二条形态是“单调攀升直至失步”。功角持续增加越过 180 度后继续上升或者直接翻转转速偏差持续扩大。这种形态说明故障切除时间超过了极限切除时间发电机与无穷大系统失去同步。在观察功角曲线时有一个细节需要注意功角是“相对角度”在单机无穷大系统里可以直接取发电机转子角度和无穷大母线角度之差。在多机系统里通常取各发电机转子角度相对于系统惯性中心角度的偏差。如果直接看绝对角度可能因为参考点选取不当而得出一堆无法解释的曲线。故障类型对功角曲线形状的影响也很有意思。三相短路引起的功角初始上扬斜率最大曲线最陡单相接地引起的上扬相对平缓。这意味着同样的故障持续时间下三相短路更容易导致失稳——这也是为什么继电保护对三相短路要求切除时间更短的原因之一。4. Matlab 编程数值仿真的实现思路4.1 经典二阶模型与程序框架除了 Simulink 的模块化仿真用 m 文件写一个纯数值的暂态稳定计算程序能让你把整个物理过程掌握得更彻底。单机无穷大系统的经典二阶模型由两个一阶微分方程组成d(delta)/dt omega_s * (S_m - S_m0)d(S_m)/dt (P_m - P_e) / (2H)其中 delta 是功角S_m 是转子转速标幺值偏差P_m 是机械功率P_e 是电磁功率H 是惯性时间常数。电磁功率 P_e 在三种情况下分别取不同的表达式故障前P_e EV / (X1 X2) * sin(delta)故障期间P_e EV / (X1) * sin(delta)其中 X1 是故障时网络等值电抗比正常运行电抗大得多故障切除后P_e EV / (X3) * sin(delta)其中 X3 是故障切除后的网络等值电抗。这个模型之所以叫“经典二阶模型”是因为它把发电机等效为暂态电抗之后的恒定电压源忽略了励磁调节和凸极效应。精度虽然不高但对理解暂态稳定的核心机制绰绰有余。程序框架用 ode45 就能解决。核心逻辑是先正常仿真一小段时间让系统到达稳态比如 0.5 秒然后在故障时刻切换电磁功率表达式故障持续期间用故障时的网络参数计算电磁功率故障切除后切换回切除后的网络参数。积分到仿真结束时间比如 5 到 10 秒最后画出功角曲线。4.2 故障矩阵批量仿真的工作流用编程做批量仿真最大的优势是循环。我把所有需要研究的工况组织成一个“故障矩阵”每一行代表一个仿真场景列包含故障类型、故障位置、故障起始时间、故障持续时间。然后用一层循环逐行调用核心仿真函数把每一行的功角最大值或者稳定状态标记保存下来。一段建议的脚本构架% 定义故障矩阵 % 第1列: 故障类型 (1三相短路, 2两相短路, 3单相接地) % 第2列: 故障位置 (pu, 0机端, 0.5线路中点) % 第3列: 故障持续时间 (s) caseList [ 1, 0.0, 0.15; 1, 0.5, 0.15; 2, 0.0, 0.15; 2, 0.5, 0.15; 3, 0.0, 0.15; 3, 0.5, 0.15; ]; % 遍历所有工况 for k 1:size(caseList, 1) faultType caseList(k, 1); faultLoc caseList(k, 2); tClear caseList(k, 3); [delta, t] smibTransientSim(faultType, faultLoc, tClear); status(k) checkStability(delta, t); end这样跑完后你可以很方便地把每个工况的稳定性状态和故障类型、位置、持续时间关联起来做进一步的定量分析。4.3 关键步骤的代码逻辑拆解这里给出一个经过调试验证的核心函数的关键实现片段帮大家理解每一步的计算逻辑。function [delta, t] smibTransientSim(faultType, faultLoc, tClear) % 系统参数标幺值 H 4.0; % 惯性时间常数, 秒 V 1.0; % 无穷大母线电压 Ep 1.1; % 发电机暂态电动势 Xd 0.3; % 发电机暂态电抗 Xt 0.1; % 变压器电抗 Xl 0.4; % 线路总电抗 Pm 0.8; % 机械功率 % 三种状态下的等值电抗 X1 Xd Xt Xl; % 正常运行 X2 Xd Xt faultLoc * Xl; % 故障期间 X3 Xd Xt Xl; % 故障切除后本模型不考虑线路切换 % 初始稳态功角 Pmax1 Ep * V / X1; delta0 asin(Pm / Pmax1); % 时间参数 t0 0.5; % 故障发生时刻 tEnd 10; % 仿真结束时间 % 用 ode45 分段积分 % 第一段: 故障前, 0 ~ t0 % 第二段: 故障期间, t0 ~ tClear % 第三段: 故障切除后, tClear ~ tEnd % 第一段: 稳态, 功角不变 % 第二段 [tFault, yFault] ode45((t, y) swingEq(t, y, Pm, Ep, V, X2, H), ... [t0, tClear], [delta0; 0]); % 第三段 [tPost, yPost] ode45((t, y) swingEq(t, y, Pm, Ep, V, X3, H), ... [tClear, tEnd], yFault(end, :)); % 拼接结果 t [tFault; tPost]; delta [yFault(:, 1); yPost(:, 1)]; end function dydt swingEq(~, y, Pm, Ep, V, Xeq, H) delta y(1); omega y(2); Pe Ep * V / Xeq * sin(delta); ddelta omega; domega (Pm - Pe) / (2 * H); dydt [ddelta; domega]; end有几个细节值得特别说明。第一初值功角 delta0 asin(Pm / Pmax) 是稳态功角必须通过求解潮流或者功率平衡方程得到不能随便拍脑袋给一个。第二故障期间的等值电抗 X2 用了 faultLoc 表示故障位置这个处理是把故障点与发电机之间的电抗和故障点之后的电抗并联后的等效结果每类故障的具体表达式都不一样需要单独推导这是编程求解中最容易出错的位置务必小心核对。第三积分段必须分开用故障切除时刻的末值作为下一阶段的初值这正是数值计算里“事件驱动”的基本思路。5. 结果分析仿真数据里最有价值的几个指标5.1 极限切除时间CCT的求取方法极限切除时间是暂态稳定研究里最重要的单一指标。它的物理意义是从故障发生到故障切除的最大允许时间超过这个时间系统就会失稳。用仿真的方式求 CCT比较实用的做法是“二分法扫描”。举个例子三相短路故障下先试 0.2 秒稳定再试 0.3 秒失稳那就把区间缩小到 0.2 到 0.3 秒试 0.25 秒稳定再试 0.28 秒失稳继续二分下去直到区间宽度小于你想要的精度比如 0.001 秒。每一轮试一次仿真几十次迭代就能得到满意的精度。这个流程完全可以自动化因为稳定性判据可以用程序判断功角最大值是否超过 180 度或者功角是否持续单调增加超过某个阈值。对于单机无穷大系统还有一个验证方法用等面积法则从解析角度估算 CCT。原理是故障期间的加速面积等于故障切除后的最大减速面积时对应的切除时间就是极限切除时间。实测下来数值求解的结果和等面积法则的估算误差通常在几个百分点以内这个一致性本身就是对你仿真模型正确性的有力验证。5.2 不同故障类型对稳定性的影响排序跑完故障矩阵后你会得到一个非常清晰的结论同样的故障位置和持续时间下故障严重程度的排序是确定的。三相短路最严重两相接地次之两相短路再次单相接地最轻。我跑过的一组典型数据某模拟系统机端故障线路全长 100km大致是这样的三相短路 CCT 约 0.17 秒两相接地约 0.22 秒两相短路约 0.26 秒单相接地约 0.32 秒。这个排序可以用故障期间电磁功率的传输能力来解释三相短路几乎完全切断了送端到受端的功率传输路径两相接地和两相短路还有健全相可以部分传输功率单相接地有两相健全相传输能力损失最小。故障位置对 CCT 的影响同样显著。同一故障类型下故障点越靠近发电机CCT 越小。还是上面那组数据三相短路在线路末端90% 处的 CCT 大约是机端故障的 1.5 倍左右。这个结论对继电保护整定有直接参考价值距离保护在近区故障时要求更快动作因为留给切除的时间窗口更短。另外值得关注的是“故障持续时间对系统恢复形态的影响”。我试过把故障持续时间设在 CCT 的 90% 附近功角曲线的最大摆动角会明显变大恢复时间明显拉长。这说明即使系统最终稳定接近极限切除时间切除也会让系统承受很大的暂态冲击。实际运行中继电保护和断路器配合的切除时间往往需要留有一定的裕度比如至少低于 CCT 的 70%就是为了避免这种“勉强稳定”的状态。5.3 从振荡曲线反推系统薄弱点把多个场景的功角响应曲线放在同一张图上对比是寻找系统薄弱环节的有效手段。具体来说我做过的操作是固定故障类型为三相短路固定切除时间为 0.15 秒然后把故障位置从 0%机端一直扫描到 100%线路末端每 10% 取一个点画出功角最大值的曲线。这条曲线会呈现明显的单调下降趋势曲线中部如果有异常的“凸起”或者“凹陷”往往对应网络拓扑中的特殊点比如变压器出口、线路换位点或者无功补偿装置的接入点。还有一个十分重要的观察维度是“多摆稳定性”。有些工况在故障切除后的第一摆看起来是稳定的功角没有越过 180 度但在第二摆甚至第三摆时失稳了。这种情况在阻尼较弱的系统里比较常见。所以稳定性判据不能只看最大功角是否超过某个阈值还要结合转速偏差曲线如果转速偏差的振荡幅值在数秒内无法衰减即使功角没有越过失稳界限实际运行中也可能因为持续振荡而被迫停机。这类分析对于理解和改善系统阻尼很有价值。如果你的研究中包含 PSS电力系统稳定器或者励磁系统模型可以用同样的仿真框架对比加装前后的暂态响应曲线这又是一个比较完整的扩展方向。6. 常见问题与排查技巧实录6.1 仿真不收敛与步长设置的坑Simulink 仿真发散是最常见的问题。我在刚接触这套模型时把故障模块的 Transitional 电阻设成 0仿真直接报错。后来才知道纯理想短路会导致数值奇异必须给故障模块设置一个很小的非零电阻比如 0.001 欧姆或者在故障点并联一个小电阻接地。还有一个典型情况是故障切除瞬间仿真步长突然变小、计算时间骤增甚至提示“Stepsize becomes too small”。这通常意味着模型中存在时间常数差异过大的环节——发电机转子的时间常数是秒级而故障开关动作的时间是微秒级两者叠加导致数值刚性。解决的方案是用刚性求解器ode15s、ode23tb替代默认的 ode45同时把相对容差Relative tolerance从默认的 1e-3 改到 1e-4 或更小。对于纯 Matlab 编程实现步长问题同样存在。ode45 是变步长求解器但如果你的故障切换逻辑是在积分函数内部用 if 语句判断当前时刻来切换模型数值稳定性会比较差。更稳妥的做法是像 4.3 节展示的那样把整个积分过程拆成三段故障前、故障中、故障切除后每段独立积分用上一个阶段的末值初始化下一个阶段。6.2 初始化与参数匹配的雷区仿真开始前必须完成潮流或稳态初始化。Simulink 的 Powergui 模块提供了 Load Flow 工具可以用来初始化同步电机。刚接触时我试图直接跳过初始化运行仿真结果波形一开始就剧烈振荡完全得不到有意义的功角曲线。正确流程是先设置好发电机内电势和机械功率用 Load Flow 工具算出合理的初始状态再开始暂态仿真。参数不匹配的问题就更隐蔽了。同步电机模块里的额定容量、额定电压变压器模块里的变比、容量线路模块的长度和单位参数这些必须和 Powergui 里的基准值保持一致。我遇到过一次变压器容量设成了 50MVA而发电机是 100MVA结果功角曲线出现了一个不明所以的偏置排查了半天最后发现是容量不匹配导致的标幺值基准混乱。对于 m 文件编程最常见的错误是“万用标幺值”误导。在写代码前必须统一整个系统的基准功率和基准电压所有参数折到同一套标幺体系里再参与计算。如果一半用有名值、一半用标幺值混着算结果必然是一笔坏账。我建议从开始就只用标幺值并且在代码注释里写明每套参数的基准值。6.3 排查技巧速查表现象可能原因排查方法仿真开头波形剧烈振荡未进行稳态初始化使用 Powergui 的 Load Flow 工具初始化波形发散或出现 NaN故障电阻为 0数值奇异设置 0.001 欧姆级故障电阻功角曲线出现异常偏置各模块基准容量/电压不一致核对全部模块额定参数与基准值求解器报步长过小刚性过大或容差太紧改用 ode15s/ode23tb调整容差批量仿真中途异常中断某个工况参数超出有效范围在循环体内加 try-catch 记录异常场景故障切除后功角持续振荡不减阻尼系数 D 设置过小增大阻尼系数 D 或加入 PSS 模型编程模型与 Simulink 结果偏差大等值电抗推导不一致重查故障期间网络等值电抗表达式最后再分享一个我踩过几次坑之后的习惯每跑一组成果把工况且关键参数记录在表格里仿真结束后及时归档。每次研究的工况设定都做到可溯源——哪些改变了哪些没变为什么要这么改——这一步虽然不起眼但会让后续对比分析省下大量返工时间。这也是我好几次解决“结果说不出理由”类问题的关键手段。这个课题做到后面你会发现最有意思的事情不只是把所有故障类型排序一遍——更值得深挖的是故障位置、系统电抗、惯性常数与极限切除时间之间的定量关系。掌握了两套工具之后完全可以沿这个方向扩展到多机系统、考虑励磁模型甚至新能源接入的暂态稳定分析思路是相通的。希望这篇内容能帮你减少一些起步阶段的试错成本让你把精力集中在真正核心的研究问题上。