做电力系统仿真这些年我一直把序分量分析当成故障诊断的“透视镜”。电网正常运行时我们习惯用A、B、C三相来观察电压电流可一旦发生短路三相之间的对称关系被打破直接看原始波形往往是一团乱麻。这时候把三相分量拆成正序、负序、零序三组对称分量故障的特征反而清清楚楚地“浮”出水面。今天这篇就用Python完整实现这个流程从对称分量法的数学原理到四种典型短路故障的序分量计算再到用matplotlib把相量图和故障过程可视化出来。整个过程我都会给出可直接运行的代码你在自己机器上装好Python环境就能复现。这套内容适合两类人一是刚接触电力系统分析的学生如果你正在为“零序电流怎么算”“单相接地为什么是正负零序串联”这类问题挠头跟着代码走一遍会豁然开朗二是做配网运维、保护整定或新能源并网研究的工程师想把手头的算例从手工推导变成可复用工具。我不会只丢公式而是把我踩过的坑、绕过的弯都写出来尽量让每个细节都能直接“抄作业”。1. 序分量分析的本质为什么故障仿真要先拆三相1.1 对称分量法不是“高级技巧”而是故障分析的坐标系很多教材一上来就搬出变换矩阵把对称分量法搞得像线性代数练习题。我的理解更朴素对称分量法相当于给三相系统换了一套坐标系。就好比一张彩色照片可以拆成R、G、B三个通道一组不对称三相向量也可以拆成三组完全对称的向量之和正序分量A相、B相、C相依次滞后120度电机旋转方向为正向负序分量A相、C相、B相依次滞后120度对应反向旋转磁场零序分量三相大小相等、相位相同没有旋转方向。这样拆的好处是三相对称的电力元件线路、变压器、发电机在正序、负序、零序下各自表现出确定的阻抗而且三者之间互不相干。于是一个看起来不对称的故障电路被拆成三个独立、对称的序网再用故障点边界条件把三个序网重新连起来计算量立刻降了一个维度。用Python做这件事核心就是两个复矩阵正变换把三相相量变成三个序分量逆变换把序分量还原成三相相量。很多人记不住公式里的算子a其实a就是旋转120度对应复数域的e^{j120°}数值等于-0.5 j0.8660254。Python的复数可以直接表示比手工计算舒服太多。1.2 常见故障类型对序分量的约束关系搞懂序分量之后接下来要记住四种典型故障分别对序分量提出什么“约束”。这些约束全部来自故障点的边界条件比如某相电压为零、某相电流为零还有电流的流入流出关系。整理成表格如下故障类型边界条件要点序分量关系三个序网的连接方式三相短路三相电压对称跌落仍为对称系统只有正序分量负序、零序为零正序网络单独计算单相接地短路故障相电压为零非故障相电流为零I1 I2 I0正、负、零序网串联两相短路两故障相电流大小相等、方向相反非故障相电流为零I1 -I2I0 0正、负序网并联两相接地短路两故障相电压为零非故障相电流为零I1 I2 I0 0U1 U2 U0正、负、零序网并联这张表是代码实现的地图。拿到一个故障类型第一步不是写代码而是先确定序网怎么连、电流关系是什么。我在实际项目里见过不少同事跳过这一步直接拿公式硬套结果算出来的相量图自相矛盾。2. 动手搭模型一套能算短路故障的最小Python框架2.1 安装环境与基础库不会Python也能跟下来这部分的代码量并不大。你需要安装Python 3.8以上版本然后装三个库numpy做复数矩阵运算matplotlib画图pandas用来整理结果。如果你之前没装过直接执行pip install numpy matplotlib pandas如果pip下载速度慢可以换国内镜像源比如清华或阿里云的PyPI镜像。装的时候偶尔会遇到“找不到msvcp140.dll无法继续执行代码”这种报错多半是系统缺Visual C运行库去微软官网下载对应版本的vc_redist安装就行。类似这种环境坑几乎每个做数据分析的人都遇到过提前说一句能省很多事。2.2 定义系统参数和序网阻抗为了把注意力放在故障分析上我采用一个最简化的单电源系统发电机经过线路接到无穷大母线故障点在线路末端。用标幺值建模电源电动势E_a 1.0∠0°正序阻抗Z1 j0.2负序阻抗Z2 j0.2零序阻抗Z0 j0.1。之所以都用纯电抗是为了模拟短路瞬间线路的感性特征电阻忽略不计也让手工验算更方便。实际工程里Z1、Z2、Z0不是随便取的需要根据变压器接线方式、线路结构、发电机中性点接地方式来确定。比如零序阻抗会受中性点接地电阻影响单相接地时的零序电流大小几乎就是由零序网络决定。不过作为仿真演示这套参数已经足够说明问题。在Python里我习惯把系统参数集中放到一个字典里后边所有函数都从这里面读改算例的时候只动一处非常省事。代码长这样import numpy as np # 系统参数标幺值角度用弧度 a np.exp(2j * np.pi / 3) # 旋转120度算子 system { Ea: 1.0 * np.exp(0j), # A相电动势 Z1: 1j * 0.2, # 正序阻抗 Z2: 1j * 0.2, # 负序阻抗 Z0: 1j * 0.1, # 零序阻抗 Zf: 0.0, # 故障点过渡电阻先取0 }注意Ea这里用了1.0 * np.exp(0j)也就是实轴上幅度为1的向量。后面所有相量都以此为参考画图时会发现这个“参考习惯”特别重要否则角度会乱。2.3 实现对称分量变换矩阵很多教材把变换矩阵写成3x3的复数矩阵看起来有点吓人其实代码短短几行。正变换是把三相相量V_abc变成序分量V_012公式为V_012 T * V_abc其中T是T (1 / 3) * np.array([ [1, 1, 1], [1, a, a**2], [1, a**2, a] ]) T_inv np.array([ [1, 1, 1], [1, a**2, a], [1, a, a**2] ])注意正变换里有1/3逆变换里没有。很多初学Python的人在这里翻车直接用numpy矩阵乘完发现三相电流还原不回去检查半天才想起来逆矩阵才是完整的变换而正变换需要除以3。我把这两个矩阵封装成函数顺手加了一个往返校验def abc_to_012(v_abc): return T v_abc def seq_to_abc(v_012): return T_inv v_012 # 往返校验还原后和原值误差应接近0 test np.array([10j, -0.50.866j, -0.5-0.866j]) v_back seq_to_abc(abc_to_012(test)) print(np.max(np.abs(v_back - test)))校验这个步骤虽然简单却能省掉后面排查的很大一部分精力。我见过太多案例结果一乱就怀疑参数最后发现是变换矩阵写错了。3. 四种典型故障的序分量分析与代码实现3.1 单相接地短路三个序网直接串联单相接地是最常见也最容易理解的故障类型。以A相接地为例故障点边界条件是A相电压为0B相和C相电流为0。把U_A 0写成序分量形式得到U_1 U_2 U_0 0把I_B I_C 0写成序分量形式得到I_1 I_2 I_0。这正好说明正、负、零三个序网在故障点串联。串联序网的计算就简单了总阻抗等于三个阻抗相加故障点正序电流为def single_phase_ground(system): Z1, Z2, Z0 system[Z1], system[Z2], system[Z0] Ea system[Ea] I1 Ea / (Z1 Z2 Z0 3 * system[Zf]) I2, I0 I1, I1 I_012 np.array([I0, I1, I2]) # 注意这里排列是 [零序, 正序, 负序] I_abc seq_to_abc(I_012) # 故障点各序电压 U1 Ea - Z1 * I1 U2 -Z2 * I2 U0 -Z0 * I0 return I_012, I_abc, np.array([U0, U1, U2])注意我这里I_012的排列是[零序正序负序]和很多教材反着来。这是为了写矩阵变换时方便T矩阵第一行对应零序。后面所有函数我都沿用这个排列只要你在代码文件头部注释写清楚就不会乱。按前面给的参数算I1 1 / (j0.2 j0.2 j0.1) -j2.0A相短路电流是3 * I1 -j6.0。不要被这个6吓到这是标幺值不是实际电流实际电流要乘以基准值通常几千安培。3.2 两相短路负序是正序的镜像BC两相短路是另一类典型故障。边界条件A相电流为0B相和C相电流大小相等、方向相反B相和C相电压相等。写成序分量就是I_1 -I_2I_0 0U_1 U_2。正序网和负序网在故障点并联零序网不参与。计算时正序电流等于Ea除以Z1加Z2的并联严格说是两个序网串联后的总阻抗因为电流是从正序网流出去再流回负序网路径上是串联但从端点看是并联结构。不管口径怎么叫算式就是def two_phase_short(system): Z1, Z2 system[Z1], system[Z2] Ea system[Ea] I1 Ea / (Z1 Z2) I2 -I1 I0 0 I_012 np.array([I0, I1, I2]) I_abc seq_to_abc(I_012) U1 Ea - Z1 * I1 U2 U1 # 两相短路时 U2 U1 return I_012, I_abc, np.array([0, U1, U2])代入参数I1 1 / j0.4 -j2.5。B相和C相电流绝对值等于根号3倍的|I1|约4.33。这里有一个特别容易混淆的点两相短路电流并不等于三相短路电流而是比三相短路小大约13.4%当Z1Z2时。这个比例在保护整定里经常用到比如复合电压闭锁过流保护为什么能区分相间短路和三相短路原理也在这里。3.3 两相接地短路序网并联的典型两相接地短路比前面两个复杂一点因为正、负、零三个序网都要参与而且是并联关系。以BC两相接地为例边界条件是A相电流为0B相和C相电压为0。序分量关系是I_1 I_2 I_0 0且U_1 U_2 U_0。三个序网在故障点并联正序电流Ea需要除以Z1加上Z2和Z0并联的复阻抗。代码写出来也不复杂但阻抗并联的复数运算要小心def two_phase_ground(system): Z1, Z2, Z0 system[Z1], system[Z2], system[Z0] Ea system[Ea] # 负序与零序并联 Z_parallel (Z2 * Z0) / (Z2 Z0) I1 Ea / (Z1 Z_parallel) I2 -I1 * Z0 / (Z2 Z0) I0 -I1 * Z2 / (Z2 Z0) I_012 np.array([I0, I1, I2]) I_abc seq_to_abc(I_012) # 故障点序电压 U1 Ea - Z1 * I1 U2 U1 U0 U1 return I_012, I_abc, np.array([U0, U1, U2])代入例子Z2和Z0并联是(j0.2 * j0.1) / (j0.3) j0.0667I1 1 / (j0.2 j0.0667) -j3.75。I2 -(-j3.75) * j0.1 / j0.3 j1.25I0 -(-j3.75) * j0.2 / j0.3 j2.5。注意I2和I0相位都是正90度I1是负90度所以三相电流之和为0满足A相电流为0的边界条件。这种“用约束条件倒推计算结果是否正确”的思路是排查代码问题最有效的方法之一。3.4 三相短路最简单却最容易算错的特例三相短路时系统仍然对称理论上不需要序分量分析直接单相计算就行I Ea / Z1然后乘以旋转因子得到三相电流。但你如果真用对称分量法去算会发现负序和零序分量全部为0正序分量等于三相短路电流除以3这个说法绕得厉害我每次都会被自己绕进去。正确的理解是三相短路时只有正序网络所以I_012 [0, Ea/Z1, 0]还原成三相电流时A相电流等于I1B相等于a^2 * I1C相等于a * I1。这里不需要再乘3因为变换矩阵的排列我已经明确是[0, 正序, 负序]。换成人话三相短路电流的幅值就是|Ea/Z1|然后B相、C相分别旋转120度。代码最简单但也最考验你对序分量排列顺序的掌握def three_phase_short(system): Z1 system[Z1] Ea system[Ea] I1 Ea / Z1 I_012 np.array([0, I1, 0]) I_abc seq_to_abc(I_012) U1 0 # 故障点三相电压都降为0 return I_012, I_abc, np.array([0, U1, 0])3.5 故障电流与相电流的换算四种类型都算完之后我建议把所有结果放到一个表里对比。这里不是简单打印而是统一用幅值和相位来表示方便眼睛扫。我自己常用的输出函数是这样def fmt_phasor(z): mag abs(z) ang np.degrees(np.angle(z)) return f{mag:.4f}∠{ang:.2f}° def print_summary(name, I_012, I_abc, U_012): print(f\n {name} ) labels [零序, 正序, 负序] for lab, val in zip(labels, I_012): print(f序电流 {lab}: {fmt_phasor(val)}) for lab, val in zip([A相, B相, C相], I_abc): print(f相电流 {lab}: {fmt_phasor(val)})我把四种故障函数都跑一遍输出的结果再拿手算的笔算对比能对上就基本放心了。这里有个很重要的实操习惯仿真的第一步不是追求复杂模型而是用一个能手工验算的简单算例把所有代码跑通。模型太复杂、手工没法验算的话一旦出错你连错误方向都找不到。4. 让电气参数“看得见”相量图与可视化实现4.1 从一串复数到一张相量图算完故障电流光是盯着打印出来的“6.0000∠-90.00°”这种数字还是不够直观。尤其是三相电压、电流的相位关系用文字描述半天不如一张相量图直接秒懂。可视化在故障分析里的作用不只是摆着好看它更像测试工具如果A相短路时画出来的A相电流相量没有落在正确位置你一眼就能发现问题。matplotlib画相量图的核心思路是把每个复数当成二维平面上的一个点从原点画一条带箭头的线段。numpy的real和imag分别对应横纵坐标。为了看得清楚我会把坐标轴纵宽比设为等比例不然120度、240度会被拉伸得不像。下面是一个通用画相量图的函数支持传入任意多个复数相量import matplotlib.pyplot as plt def plot_phasors(phasors, labels, title, colorsNone): fig, ax plt.subplots(figsize(6, 6)) ax.axhline(0, colorgray, lw0.8) ax.axvline(0, colorgray, lw0.8) ax.set_xlim(-8, 8) ax.set_ylim(-8, 8) ax.set_aspect(equal) ax.grid(True, ls--, alpha0.6) for i, (ph, lab) in enumerate(zip(phasors, labels)): color colors[i] if colors else None ax.annotate(, xy(ph.real, ph.imag), xytext(0, 0), arrowpropsdict(arrowstyle-, colorcolor, lw1.8)) ax.text(ph.real * 1.15, ph.imag * 1.15, lab, hacenter, vacenter, fontsize10) ax.set_title(title) plt.tight_layout() plt.savefig(f{title}.png, dpi150) plt.show()4.2 单相接地故障的相量图与序分量分解拿单相接地短路为例。故障后A相电压为0B相、C相电压保持不变在我们的简化模型里不完全保持但趋势是对的。把A相电流等于-j6、B相电流约0、C相电流约0画出来你会看到只有一个很大的A相电流。单独看这个结果你可能觉得没什么意思但如果把序分量I1、I2、I0画在同一张图上会发现三个序分量幅度相同、相位关系完全一致直观体现了串联约束。我最喜欢画的是“故障电流的三相相量与序分量对比图”上半部分画三相电流下半部分画三个序电流同一比例尺一眼就能看出“不对称”是怎样由序分量叠加出来的。这个图放到汇报PPT里比满页公式说服力强得多。4.3 故障前后对比动态图静态相量图还不够过瘾我把故障瞬间的三相电流变化过程做成了动画。其实原理不复杂故障前的三相电流是三相对称故障后的稳态电流是刚刚算出来的不对称相量中间用一个过渡系数逐步从故障前插值到故障后。matplotlib.animation可以轻松实现我这里给出一个简化的动态演示版本重点在旋转相量from matplotlib.animation import FuncAnimation def animate_phasor(ph_before, ph_after): fig, ax plt.subplots(figsize(6, 6)) ax.set_xlim(-8, 8) ax.set_ylim(-8, 8) ax.set_aspect(equal) ax.grid(True, ls--, alpha0.6) line, ax.plot([], [], o-, lw2, markersize6) def update(frame): t frame / 100 # 简单线性插值也可以加旋转演化 ph (1 - t) * np.array(ph_before) t * np.array(ph_after) xs, ys [0, ph.real], [0, ph.imag] line.set_data(xs, ys) return line, ani FuncAnimation(fig, update, frames101, interval40, blitTrue) ani.save(fault_transition.gif, writerpillow, dpi120)这个动画你可以自己改着玩把三相电流都传进去就能看到从对称到不对称的突变过程。我实际测试下来用pillow保存gif最稳不需要额外装ffmpeg如果想保存成MP4就需要ffmpeg支持很多Windows用户的ffmpeg环境没配好折腾半天。所以默认我都存gif省事。5. 实操中的坑与排查清单重点5.1 短路电流结果异常偏大或偏小怎么排查仿真结果和手算对不上是最常见的情况。我总结了一套排查顺序先看复数类型有没有被强转成float再看矩阵排列顺序对不对最后看边界条件的符号有没有写反。出现过一次印象特别深的bug单相接地算出来的A相短路电流是-j2而不是-j6。排查半天发现我在构建T矩阵时把a和a^2的顺序写反了。正变换后面还跟了个逆变换两个错误互相抵消了一部分导致结果不是完全错而是错得特别隐蔽。这种问题最坑人因为你觉得“大概对了”。所以我的建议是每个故障函数写完立刻手动构造一个边界条件做断言测试比如单相接地后I_B应该接近0两相短路后I_A应该接近0用np.isclose检查错一毫秒就报错。5.2 画相量图时角度乱转相位参考与unwrap问题画相量图时角度漂移特别常见。表现是明明计算结果正确但画出来B相和C相的箭头方向全乱了。原因多半是matplotlib的横纵坐标范围不同或者用了极坐标时没有统一角度参考。我统一把参考方向设为正实轴即Ea 1∠0°所有相量都基于同一个复数平面这样图就不会乱。还有一个小坑numpy.angle返回的角度范围是(-π, π]如果相量经过复数运算后出现了相位跳变画出来的箭头会突然跳到负半平面。比如正序电流I1 -j2.5angle返回的是-90度而不是270度看起来没问题但如果算电压降时出现了符号反转角度可能从170度跳到-190度画出来就会变成反向。这时候可以先对角度做unwrap处理再画连续曲线如果是静态相量图先把每个复数的real和imag直接作为坐标也就不用关心角度范围问题了。5.3 复数运算的广播陷阱NumPy的复数矩阵运算里最容易被忽视的是形状和类型。T矩阵是3x3v_abc是3元向量直接用没问题。但如果你从某个表格里读到电压幅值和相位然后拼成字符串再解析成复数很容易触发类型转换问题。我试过直接用np.array([1, 0, 0])结果数组变成int型后面的复数乘法直接报错或者更糟不报错但精度丢失。解决方案构造复数相量时统一用.astype(complex)或者用np.exp(1j * angle)从极坐标直接生成避免中间类型转换。# 推荐从幅值相位构造复数 def polar(mag, deg): return mag * np.exp(1j * np.deg2rad(deg))5.4 零序参数测不准带来的连锁影响零序阻抗在很多实际工程里很难精确获取尤其受大地回路、线路换位、变压器绕组接线影响。你仿真时用j0.1算出来单相接地电流-j6实际工程里零序阻抗很可能在0.05到0.2之间波动短路电流跟着能差出一倍。这就要求我们在做方案设计时不能只算一种工况要做敏感性分析——把Z0从0.05逐步变到0.2观察单相接地电流的变化范围。这个敏感性分析用Python做几乎不费劲循环调用single_phase_ground函数就行。很多工程计算书里那一堆“最不利工况”的数据其实背后就是这种参数扫描。我在做配网保护定值校验时会把不同接地电阻下的短路电流曲线全部画出来保护人员看了直说“这比给一个单点数值有用多了”。5.5 别让可视化误导判断最后提一个可视化层面的坑坐标轴等比例是必须的否则120度会被压成看起来像90度或150度。此外相量图里箭头长短的感知会受到刻度范围影响如果A相电流是6而其他相电流接近0放在同一张图里小相量可能被压得看不见。这时候可以做两个子图一个保留全局视角一个对故障相做局部放大。我见过有人只看放大图以为非故障相也有大电流判断直接跑偏。6. 把仿真代码扩展成自己的工具6.1 封装成可复用的故障分析函数到这一步手写函数已经能覆盖四种典型故障。但每次都要手动改参数、手动调用再手动查表效率太低了。我建议把整个模块封装成一个类比如叫FaultAnalyzer初始化时传入系统参数然后暴露几个方法analyze(fault_type)、plot_phasors()、summary()。这样不管是教学演示还是实际工程只需要几行代码就能完成一次完整的故障评估class FaultAnalyzer: def __init__(self, system): self.system system def analyze(self, fault_type): if fault_type AG: return single_phase_ground(self.system) elif fault_type BC: return two_phase_short(self.system) elif fault_type BCG: return two_phase_ground(self.system) elif fault_type ABC: return three_phase_short(self.system) else: raise ValueError(funsupported fault type: {fault_type})当然真正的工程系统不会只有一台发电机和一条线路但底层逻辑相同。你可以把每个元件的序阻抗做成输入数据表形成关联矩阵再通过节点导纳矩阵构建更复杂的网络。这个方向一旦走通就能分析任意拓扑下的各种短路故障甚至支持故障点位置扫描。6.2 从稳态故障到暂态过程的扩展思考我前面写的都是稳态短路分析也就是假设故障已经发生了一段时间电气量进入稳定状态。实际上保护装置动作发生在故障后的几个周波内暂态过程中的非周期分量、谐波分量对保护影响很大。要把稳态分析扩展成暂态仿真一个简单方法是叠加直流衰减分量用指数衰减项模拟非周期分量。比如把短路电流从纯工频相量改成i(t) I_m * sin(ωt θ - φ) C * e^{-t/τ}其中τ由X/R决定。这个模型虽然还是近似但已经能把断路器开断电流的峰值估出来对设备选型很有帮助。Python里用numpy生成时间序列再用matplotlib画波形几行代码就能看到完整的暂态电流包络线。我在实际项目里还会把故障仿真结果接入保护定值校核流程用故障电流算保护安装处的测量阻抗判断距离保护会不会误动。说起来有点复杂但底层的相量计算和序分量分解就是这篇文章里这套东西。6.3 个人体会做故障仿真这几年最深的体会是代码本身不重要重要的是代码背后的电气直觉。序分量分析给了我从另一个角度审视不对称故障的能力而Python帮我免去了大量重复的复数运算让我能把精力放在判断“这个结果合不合理”上。你有空的话建议把代码里的参数换成自己手头真实系统的数据哪怕只是把Z0从0.1改成0.15看看电流变化了多少这种直观感受比读十篇论文都更来得实在。希望这篇分享能帮你少踩几个坑也欢迎你在这套框架上继续玩出更多花样。