去年年底整理实验室研究材料的时候导师在小组会上问了我一个问题咱们这个示范区的配电网可靠性到底能打几分我当时的反应是愣了几秒——手头有拓扑图、有设备台账、有负荷数据但真要我说清楚一年下来平均每户会停几次电、每次停多久、供电可用率是多少我居然给不出一个数。也就是从那时候起我开始认认真真做了一个配电网可靠性评估的Matlab程序一路踩坑也一路把以前课本上模糊的概念全部落地成了能算数的代码。这篇文章就把我这段实现之路完整记录下来包括指标怎么定义、方法怎么选、代码结构怎么搭、算例结果对不对以及调试过程中踩过的那些坑。如果你是电力系统方向的研究生、配网规划工程师或者正在被各种可靠性课程设计折磨的本科生这篇应该能直接当成一份可参考的路线图。1. 动手前先把指标吃透可靠性评估到底在算什么1.1 五个绕不开的指标从定义到物理意义很多人一上来就急着写代码结果算完之后根本不知道结果合不合理。我建议先把指标这关过了因为后面程序里每一个变量、每一次累加本质上都是在为这几个指标服务。配电网可靠性评估用得最多的是这么几个指标SAIFI系统平均停电频率指标、SAIDI系统平均停电持续时间指标、CAIDI用户平均停电持续时间指标、ASAI平均供电可用率指标以及EENS期望缺供电量。它们的定义和公式我整理成了表格方便后面写代码时对照。指标物理含义计算公式典型量级SAIFI每个用户在一年内平均停电次数用户停电次数总和 / 总用户数0.5~2次/户·年SAIDI每个用户在一年内平均停电小时数用户停电持续时间总和 / 总用户数1~5小时/户·年CAIDI每次停电平均持续多久SAIDI / SAIFI1~3小时/次ASAI一年内供电可靠时间占比1 - SAIDI / 876099.9%以上EENS因故障而少送的电量各场景停电功率×停电时间之和视系统规模而定写代码时最容易混淆的是SAIFI和SAIDI。SAIFI只看次数所以某个负荷点只要经历了一次停电哪怕只停1分钟也要计入频率而SAIDI看的是时长停1分钟和停3小时在SAIDI里的权重完全不同。这两个指标一个管频次、一个管时长对应到配网运维上就是两条不同的改进思路一个靠减少故障次数比如线路绝缘化改造、防雷改造另一个靠缩短修复时间比如加装故障指示器、完善转供能力。EENS比较特殊它需要知道负荷点的平均功率。比如一个节点有200户居民户均峰值功率按2kW估算那么该节点负荷就是400kW。如果这个节点在一年里累积停电3小时那EENS就是1.2MWh。这个指标在做经济性评估时特别有用因为可以给它乘上电价算出少送了多少电费。1.2 三种时间参数背后的运行逻辑程序里反复用到的时间参数有三个故障修复时间、故障隔离时间、负荷转供时间。很多初学的人把它们混为一谈这是导致计算结果错误的最常见原因。我先讲一下配电网故障处置的真实流程。假设一条10kV馈线的某一段线路被外破打断了变电站出线断路器会先跳闸整条馈线瞬间失电。这时候调度员并不知道具体故障点在哪需要根据保护动作信息、故障指示器信号去定位故障区段。找到故障点之后如果能遥控操作就把故障点两侧的隔离开关拉开这叫故障隔离。隔离完成之后如果故障点靠电源侧的区域没有损坏就可以合上出口断路器先恢复送电如果故障点后端还有网架结构、存在联络开关可以通过合联络开关把下游负荷转到另一条馈线上供电这叫转供。最后维修人员到现场把故障元件修好或换掉整个系统恢复到初始运行方式。对应到程序里这几种时间要分开处理故障修复时间指故障元件从损坏到修好重新投运的时间通常是几小时到十几小时线段越长、交通越不方便修复时间越长。故障隔离时间指从故障发生到隔离操作完成的时间如果是自动化开关可能是秒级到分钟级如果是人工操作可能是0.5~2小时。负荷转供时间指合上联络开关、把负荷倒到另一条馈线的时间自动化的按分钟记手动的按0.5~1小时记。我在程序里设置了一个底层逻辑所有非故障区域的停电时间都不应该超过隔离时间和转供时间的组合而故障元件所在区域下游、且无法转供的那部分负荷停电时间才是修复时间。这个逻辑直接决定了SAIDI算得对不对。后面在第3节我会详细展示这段代码怎么实现。2. 方案选型为什么我用故障枚举而不是蒙特卡洛2.1 解析法与模拟法的取舍配电网可靠性计算方法大体分两类解析法和蒙特卡洛模拟法。解析法的典型代表是故障模式后果分析法英文缩写FMEA基本思路就是枚举每一个元件故障分析后果然后汇总指标。蒙特卡洛法则完全不同它是通过生成随机数来抽样模拟元件故障和修复过程跑成千上万次仿真统计出可靠性指标的期望值。这两种方法怎么选我的经验是对大多数配电网规划和研究场景FMEA足够了而且优势很明显。第一FMEA结果稳定同一套参数跑多少遍结果都一样不会出现蒙特卡洛那种这次跑和上次跑差几个百分点的情况第二FMEA计算速度快一个几十条支路的配电网枚举一遍所有元件故障在Matlab里几秒钟就能出结果第三FMEA天然包含灵敏度信息——你能明确看到哪条线路故障对SAIFI贡献最大这对改造决策非常有价值。蒙特卡洛的优势在于处理复杂系统比如考虑时序变化、分布式电源出力波动、元件多重故障等场景。但它的收敛速度比较慢误差按抽样次数的平方根下降想提高一位精度就要多100倍计算量。对于一篇课程设计或者一个中等规模的配电网规划项目我建议优先上FMEA把模型本身吃透之后再考虑蒙特卡洛扩展。2.2 网络拓扑怎么存邻接矩阵、支路表和图对象程序的数据结构设计决定了后面代码的复杂度和可维护性。我第一次写的时候用邻接矩阵后来发现很痛苦——邻接矩阵适合描述两节点之间直接相连的关系但要表达支路上串联了开关这条支路的类型是架空线还是电缆这种工程属性它就力不从心了。我最终采用的是支路表结构这是工程上最常用、也最直观的存储方式。每一行代表一条支路包含支路编号、首端节点、末端节点、元件类型、故障率、修复时间、是否含开关等字段。这种结构的好处是你可以把任意复杂的工程属性塞进对应的列里而且排查数据错误时非常直观。节点之间怎么连通的我建议同时用一个Matlab的graph对象来维护纯拓扑关系。它的好处是调用shortestpath、conncomp这些图论函数特别方便判断某个节点是否和电源节点连通只需要一行代码。我在程序里把graph对象当成快速查询工具把支路表当成权威数据源两者配合兼顾了灵活性和准确性。2.3 分类讨论故障影响范围FMEA的核心FMEA的核心是故障后果分析。对于一个辐射型配电网某个元件发生故障后它的影响范围可以分为三个层次第一类故障元件所在馈线的变电站出口断路器跳闸导致该馈线所有负荷瞬间失电。这不是最终结果后续还要看隔离和转供的执行情况。第二类通过拉开故障段两端的隔离开关故障段靠电源侧的区域可以快速恢复这个区域的停电时间只有隔离开关操作时间。第三类故障段下游的区域如果存在联络开关且能闭合则停电时间是转供操作时间如果不存在联络开关则只能等故障修复停电时间就是修复时间。我在程序里把每个负荷点的停电时间处理成了一个三元组停电频率停电时长缺供电量。例如某负荷点在故障发生后被隔离操作恢复了那么频率加1时长加隔离时间如果该负荷点位于故障下游且无转供路径频率加1时长加修复时间。这里有一个细节容易错即使某个负荷点在故障后很快恢复了供电只要它经历了停电SAIFI里就要计入一次。很多新手在写代码时只统计最终失电时间大于0的负荷点导致SAIFI被算小这是要特别警惕的。3. Matlab代码怎么写从数据结构到核心函数3.1 数据结构设计元件表、负荷表、系统参数我习惯先把所有数据定义集中在一个脚本或函数里方便修改和调试。元件表用Matlab矩阵来表示每一行是一个元件列的含义通过注释写清楚。下面是我在测试算例里用的一组结构示例% 元件表列定义 % 编号 | 类型(1线路/2变压器) | 首端节点 | 末端节点 | 故障率(次/年) | 修复时间(h) | 是否含开关(0/1) elements [ 1, 1, 1, 2, 0.100, 5.0, 0; 2, 1, 2, 3, 0.085, 5.0, 0; 3, 1, 3, 4, 0.120, 5.0, 1; 4, 1, 2, 5, 0.070, 5.0, 0; 5, 1, 5, 6, 0.095, 5.0, 1; 6, 2, 6, 7, 0.015, 8.0, 0; ];负荷表单独存储每一行是一个负荷节点包含节点编号、用户数量、平均负荷功率比如loads [4, 150, 300; 7, 220, 440]表示节点4有150户用户、平均负荷300kW。我还单独定义了几个全局时间参数isoTime故障隔离时间小时、tieTime联络开关转供时间小时、annualHours8760小时。这些参数独立于元件表是为了方便做敏感性分析——比如我想看看如果自动化改造把隔离时间从1小时压缩到0.2小时SAIDI能改善多少只需要改一个变量程序重跑一遍就行。3.2 写一个能跑通的故障枚举主循环主循环的逻辑其实不复杂就是遍历元件表逐个故障逐个分析。核心流程是对每一个元件把它从拓扑图中去掉模拟故障状态然后判断每个负荷点是否失电再根据隔离、转供逻辑计算失电时长最后累加到总指标里。下面是主循环的框架代码% 初始化累计变量 totalInterruptions 0; % 用户停电次数总和 totalOutageHours 0; % 用户停电小时总和 totalMissingEnergy 0; % 缺供电量总和 totalUsers sum(loads(:,2)); for k 1:size(elements, 1) % 获取当前故障元件的两端节点 n1 elements(k, 3); n2 elements(k, 4); % 从拓扑图中临时删除故障支路 G_temp rmedge(G, n1, n2); % 对每个负荷点判断失电情况 for j 1:size(loads, 1) loadNode loads(j, 1); users loads(j, 2); power loads(j, 3); % 判断电源节点与负荷节点是否连通 if ~isempty(shortestpath(G_temp, sourceNode, loadNode)) continue; % 仍连通不受影响 end % 到这里说明负荷点失电了先累计频率 totalInterruptions totalInterruptions users; % 计算失电时长需要调用自定义函数见3.3节 outageHours calcOutageDuration(elements, G_temp, k, loadNode); totalOutageHours totalOutageHours users * outageHours; totalMissingEnergy totalMissingEnergy power * outageHours; end end % 最后计算系统级指标 SAIFI totalInterruptions / totalUsers; SAIDI totalOutageHours / totalUsers; CAIDI SAIDI / max(SAIFI, 1e-9); ASAI 1 - SAIDI / 8760; EENS totalMissingEnergy;这里有一个重要的编程坑rmedge之前必须先把原图备份因为每个故障场景都要基于原始完整拓扑去删边。我第一版代码就是在循环里直接改了G本身第二圈循环就出错了。建议在主循环开始前把原始图存一份比如G_base G;每次循环用G_temp rmedge(G_base, n1, n2)。3.3 故障影响范围判断的逻辑实现这一节是最核心、也最容易出错的。上面的主循环里我只是判断了负荷点是否失电但实际的停电路径还得分三种情况。calcOutageDuration函数要做的事情就是判断失电负荷点处于哪个区域进而决定停电时长的取值。这个函数可以拆成三个步骤第一步判断故障元件是否在电源点到负荷点的唯一供电路径上。因为辐射型配电网里电源到每个负荷点的路径是唯一的故障元件只要不在这条路径上删掉故障支路根本不会导致负荷失电。这一步其实已经被主循环里的shortestpath查询覆盖了失电的负荷点故障元件必然位于其上游路径上。第二步判断负荷点在故障点的哪一侧。先搞清楚故障点两侧有没有可操作的隔离开关。方法是从故障支路的上游端点出发向电源方向搜索找到最近的含开关支路从故障支路的下游端点出发向负荷方向搜索找到最近的含开关支路。如果负荷点在上游隔离开关和电源之间那它属于第一类隔离后可恢复如果负荷点在故障支路下游那它属于第二类或第三类。第三步对故障支路下游的负荷点判断是否具备转供条件。方法是从下游隔离开关往下搜索看有没有联络开关与其他电源点相连。如果存在联络开关且开关闭合后能构成新的供电路径那它停电时长取tieTime如果不存在任何联络路径停电时长取故障元件的修复时间repairTime。这个搜索逻辑用图论函数实现非常简洁。我写了一个简化版的核心判断代码function outageHours calcOutageDuration(elements, G_temp, faultIdx, loadNode) global isoTime tieTime; % 故障元件的基本信息 n1 elements(faultIdx, 3); n2 elements(faultIdx, 4); repairTime elements(faultIdx, 6); % 如果负荷点可以从电源连通不需要计算 % 这个判断在主循环里已经做过这里省略 % 情况1负荷点在故障点上游电源侧 % 判断方法删除故障支路后如果电源节点连通负荷点则为上游区域 % 这里由于已经知道失电说明不连通所以需进一步判断上游隔离开关位置 % 如果在电源到故障点之间存在隔离开关则负荷点可在隔离后恢复 if hasUpstreamSwitch(elements, n1) outageHours isoTime; return; end % 情况2负荷点在故障点下游尝试寻找联络开关转供路径 if hasTieSwitchPath(G_temp, n2, loadNode) outageHours tieTime; return; end % 情况3无法转供只能等修复 outageHours repairTime; end实际工程里判断hasUpstreamSwitch和hasTieSwitchPath要写不少辅助函数比如要能够识别某个节点到电源路径上是否经过开关类型元件某个故障区域能否通过闭合联络开关连通到备用电源。这里我建议给自己留一个调试后门——把每个故障元件的元件编号、受影响负荷点、各类时间打印到日志文件里在开发阶段逐场景核对远比调试了半天找不到错强得多。4. 一个完整的小算例3馈线系统手把手验证4.1 测试系统参数与期望结果算法写得再好不落到具体算例上验证一遍始终不放心。我自己搭了一个三馈线的测试系统三条馈线都从同一个110kV变电站的10kV母线引出每条馈线带2个负荷节点馈线1和馈线2的末端通过联络开关相连馈线2和馈线3之间也有一个常开的联络点。我把元件参数列在下面各位可以照着录入自己的程序元件类型首端节点末端节点故障率(次/年)修复时间(h)是否含隔离开关E1线路120.105否E2线路230.085是E3线路240.125否E4线路150.095否E5线路560.115是E6线路570.075否E7联络线380.013是常开E8联络线790.013是常开负荷数据节点4有200户用户、平均负荷400kW节点3有150户、300kW节点6有180户、360kW节点7有120户、240kW。节点编号中8和9是联络线另一侧的备用电源节点分别由其他馈线供电。隔离时间取0.5小时转供时间取0.5小时。4.2 单次故障的推演过程以元件E5节点5到节点6的线路故障为例手推一遍程序应该给出什么结果。故障发生后电源节点1到节点6的路径断裂节点6失电。节点6位于E5下游检查转供路径E5下游的节点6通过E7联络线可以连到节点8节点8是备用电源所以节点6的停电时间应该是转供时间0.5小时而不是修复时间5小时。节点7仍然由原路径供电因为E5在节点7的上游支路分支点之前实际经过判断后节点7并没有失去电源连通性所以不受影响。节点3、节点4也都不受影响。再看E2节点2到节点3的线路故障。节点3失电且节点3在故障下游。节点3可以通过E6联络线连到节点8所以停电时间是0.5小时。但如果我把E7的联络线从程序里删掉节点3就会变成无转供路径停电时间立刻变成5小时修复时间。这就是为什么网架结构对可靠性影响那么大——同样一条线故障有联络线时SAIDI只增加0.5小时没联络线时SAIDI就要增加5小时差别是10倍。对E1节点1到节点2的线路故障情况比较特殊。节点1是电源E1故障相当于整条馈线1的电源进线断了节点2、3、4都失电。节点2在故障点下游但节点2到节点5之间有没有其他连接因为我这里没有设置馈线1和馈线2之间的联络所以节点2、3、4全都无法转供停电时间均等于E1的修复时间5小时。这说明辐射型配电网越是靠近电源侧的线路故障时影响范围越大这就是所谓瓶颈元件——也是规划改造时最值得优先加装自动化开关或联络通道的位置。4.3 汇总指标与合理性检查按上面的参数和逻辑程序跑出来的结果大致是SAIFI在1.0次/户·年左右SAIDI在1.5~2.0小时/户·年之间ASAI在99.98%左右。这个量级和国内许多城市配电网的实际统计值是大体吻合的。拿到结果之后我自己会习惯做一个合理性检查把某一个比较关键元件的故障率调成0看看SAIFI和SAIDI是不是相应地减少。比如把E1故障率清零SAIFI应当明显下降因为E1是馈线主进线故障影响户数多如果程序结果显示SAIFI没变化那一定是拓扑判断逻辑出了bug或者这个元件根本没被纳入枚举范围。这种变量扰动检查法是验证程序正确性的有效手段比单纯把结果和文献数字对比更能说明问题。5. 我踩过的坑常见问题与排查技巧5.1 结果异常先找这几个原因调试过程中我总结出了一个异常结果排查表每次算出来的指标不合理就按表里的顺序逐项检查现象可能原因排查方法SAIFI比预期大很多负荷点短暂停电也被重复计次检查是否所有失电负荷点都计了一次频率是否把隔离后恢复的负荷点误加了修复时间SAIDI比预期小很多把修复时间错误地替换成隔离时间重点检查无联络开关的下游负荷点取的是不是修复时间EENS为0负荷节点的平均功率数据没有正确读取打印每个负荷节点的功率值看看是否被初始化为0计算时间过长每轮循环都在重复构建最短路径考虑预计算拓扑关系或用节点邻接表替代graph对象最短路径查询结果对某个元件故障率不敏感该元件可能根本没有被枚举进故障集合检查主循环范围是否漏掉了变压器和联络线等非线路元件还有一个非常隐蔽的问题开关操作时间与修复时间使用了不同的单位。我曾经在一个工程里把隔离时间定义成分钟修复时间定义成小时结果SAIDI莫名其妙小了60倍查了半天才发现是单位不一致。建议从一开始就统一用小时作为时间单位并在代码注释里大写加粗提醒。5.2 提高运算速度的几个方法FMEA本身计算量不大但当系统规模变大——比如几百条支路、上千个节点时如果代码写得粗糙跑起来也会很卡。我有三个提速经验。第一不要在主循环里用shortestpath去判断每个负荷点和电源的连通性。更高效的做法是删除故障支路后用conncomp函数一次找出所有连通分量然后查电源节点所在分量其他分量里的节点就是失电节点。这样一次图遍历能代替对所有负荷点的逐一遍历。第二预计算每个负荷点的供电路径上都有哪些元件。因为辐射型网络路径唯一我可以一次性用DFS或shortestpath求出所有负荷点到电源的路径存成稀疏矩阵。之后判断某元件故障是否影响某负荷点就只需要查这个矩阵的一行判断故障元件在不在路径集合里复杂度从O(N)降到O(1)。第三多使用向量化操作。比如累计停电时长时不要在一个循环里一条条累加而是先记录每个故障场景下各负荷点的停电时长矩阵最后用sum(A(:))一次性汇总。Matlab的循环开销比较大能用矩阵运算解决的就别用for。5.3 避免看着对实际错的校验手段程序写完之后最怕的不是报错而是结果看起来合理、实际内部逻辑有错。我有几种校验手段推荐给大家。第一种是极小系统手工校验法构造一个只有一条线路、一个负荷点的最简单系统手算出SAIFI和SAIDI再让程序跑一遍。比如线路年故障率0.2次/年修复时间5小时负荷点100户。手算SAIFI就是0.2次/户·年SAIDI就是0.2×51.0小时/户·年。如果程序连这个都算不对那后面复杂系统肯定有问题。第二种是双实现交叉验证法我自己写程序时先用一个最简单的状态搜索方式实现一版笨方法比如完全枚举所有负荷点的失电状态再用优化后的高效算法实现第二版。两版结果做对比如果不同说明优化过程引入了逻辑错误。这个方法虽然笨但特别稳。第三种是参数扰动敏感性检查把某个元件的修复时间从5小时改到10小时SAIFI应该不变SAIDI应该变大把故障率从0.1改到0.01SAIFI和SAIDI应该近似缩小10倍。如果这些简单的比例关系不成立那程序内部八成存在错误逻辑。这种检查不会花很多时间但能帮你提前暴露很多问题。说到底配电网可靠性评估程序真正难的地方不在Matlab语法而在于把电力系统运行逻辑准确地翻译成代码。我刚写的时候也总是想一步到位结果越急越容易写出一坨逻辑混乱的代码。后来学乖了先把一条馈线、一个简单开关模型写清楚再逐步扩展成完整系统每一步都用上面的校验方法确认没问题后再往前走。这个思路不管你是做课程设计还是做实际工程都值得参考。