
简介面向电力系统与配电网研究者的MATLAB实现资源聚焦最小路法和非序贯蒙特卡洛算法在配电网可靠性评估中的应用。资源包内含IEEE RBTS可靠性测试系统的原始参数文档PDF与Excel格式、IEEE33节点系统参数以及RBTS BUS6的MATLAB数据文件同时提供两套完整程序基于最小路法的可靠性评估程序以及采用节点影响分析法判断受影响负荷的蒙特卡洛主程序。压缩包体积仅1.34MB主要为PDF、Excel和MATLAB脚本便于离线查阅与直接运行。已有3816人学习下载适合电气工程专业学生、电网科研人员及从事可靠性分析的工程师。文件组织清晰PDF参数文档方便核对原始数据Excel表格可直接修改导入MATLAB脚本注释完整既可用于计算可靠性指标也能作为理解两类算法原理、改造适应自有网络拓扑的参考模板。整体精简而完整适合作为课程设计或科研起步的实用参考。1. 配电网可靠性评估的两种算法最小路法和蒙特卡洛法是不是只能二选一我去年调一个乡镇10kV配网的可靠性数据SAIFI怎么算都跟真实统计对不上。后来发现问题不在数据脏而在方法选错了一条带联络开关的网格状线路我硬用最小路法去解析结果永远解释不了负荷转移那几小时。配电网可靠性评估的核心问题其实很朴素——一年里一个用户平均停几次电、每次停多久。最小路法把拓扑拆成每条负荷到电源的必经之路用枚举方式算出期望值蒙特卡洛法则模拟几十年甚至上百年的随机故障再统计停电事件。两者不是二选一的竞争关系而是不同网架、不同数据条件下的互补工具。这篇文章写给要做配网规划、供电可靠性管理或新能源接入方案比选的工程师读完你能判断该用哪个并且能直接跑通一个最小可用的评估脚本。2. 最小路法到底怎么算从找最小路到输出SAIFI的全过程2.1 最小路法的建模思路先画出每个负荷点的“命脉”最小路并不是图论里的最短路径而是电力实际流通的那条必经之路从变电站母线出发经过各段线路、开关、配变直到某个负荷点。对辐射状配电网来说这条路径通常唯一所以最小路法才特别适合。方法的核心是给每个负荷点分类元件在最小路上的元件故障负荷点只能等它修好停电时长等于修复时间不在最小路上的元件故障负荷点不一定完全停电可能只经历一次开关操作期间的短时停电甚至因为熔断器或分支保护而完全不受影响。这个分法实际就是故障模式影响分析FMEA的变体。它的好处是模型透明每个元件故障影响谁、影响多久都能摊开来查。坏处也很明显遇到多电源、有联络开关或者能转供的闭环网架最小路的“唯一性”被打破必须额外定义转供时间否则算出来的指标会偏大。所以在跑最小路法之前我会先花半小时把网架图和保护配置画清楚尤其是隔离开关、熔断器、联络开关的位置这些直接决定非最小路元件的故障影响范围。下面的表格是影响时间分类我每次建模都会对照它过一遍元件位置对负荷点的影响影响时间取值负荷点最小路上的元件停电至修复完成修复时间 r分支上的元件且分支有熔断器/保护不影响其他分支负荷0分支上的元件没有可靠保护上游开关动作后恢复开关操作时间 s馈线入口与上游主网元件整条馈线失压随场景取 r 或 s这个表不是死的具体取决于保护配置。写代码前把每个元件的影响负荷点和影响时间列出来是最小路法最重要的准备工作比写代码本身更容易出错。2.2 用Python实现FMEA式的最小路法一个可以直接抄的脚本最小路法的常规落地做法是枚举所有元件故障逐个判断它对每个负荷点的影响。下面这段代码我用一个带两个负荷点的简化辐射网演示数据写在字典里便于你改成自己的网架。import networkx as nx # 拓扑S为电源LP1、LP2为负荷点 G nx.Graph() G.add_nodes_from([S, A, B, C, LP1, LP2]) edges [ (S, A, 0.20, 5.0), # 进线故障率0.2次/年修复5小时 (A, B, 0.10, 4.0), # 主干线段 (B, LP1, 0.10, 3.0), # 分支到LP1 (A, C, 0.30, 5.0), # 分支到C (C, LP2, 0.15, 3.0), # 分支到LP2 ] for u, v, lam, r in edges: G.add_edge(u, v, lamlam, rr) # 每条边故障时对每个负荷点的影响时间小时 # 影响时间为0表示不影响该负荷点 impact_time { (S, A): {LP1: 5.0, LP2: 5.0}, (A, B): {LP1: 4.0, LP2: 4.0}, (B, LP1): {LP1: 3.0, LP2: 0.0}, (A, C): {LP1: 0.0, LP2: 5.0}, (C, LP2): {LP1: 0.0, LP2: 3.0}, } loads {LP1: 120, LP2: 80} # 每个负荷点的用户数 bd {lp: 0.0 for lp in loads} # 年停电小时数期望 U ld {lp: 0.0 for lp in loads} # 年停电次数期望 lambda for (u, v) in impact_time: lam dict(G[u][v])[lam] for lp, hours in impact_time[(u, v)].items(): if hours 0: ld[lp] lam bd[lp] lam * hours # 系统指标SAIFI / SAIDI / CAIDI total_customers sum(loads.values()) SAIFI sum(loads[lp] * ld[lp] for lp in loads) / total_customers SAIDI sum(loads[lp] * bd[lp] for lp in loads) / total_customers CAIDI SAIDI / SAIFI if SAIFI 0 else 0 for lp in loads: print(f{lp}: lambda{ld[lp]:.3f} 次/年, U{bd[lp]:.3f} 小时/年) print(fSAIFI{SAIFI:.3f} 次/(用户·年), SAIDI{SAIDI:.3f} 小时/(用户·年), CAIDI{CAIDI:.3f} 小时/次)这段代码的关键不是networkx画图而是impact_time这个字典。它就是你把2.1节的保护配置表翻译成程序的过程。lam是元件故障率单位是“次/年”hours是影响时间单位是“小时”。把每条边的lam * hours累加到负荷点上就得到该负荷点的年停电小时期望值 U把所有lam累加得到年停电次数期望 lambda。最后用负荷点的用户数做加权就是标准的SAIFI和SAIDI。跑这组参数LP1的停电频率约0.40次/年、年停电时长约1.7小时LP2约0.75次/年、3.35小时/年系统SAIFI约0.54次/(用户·年)SAIDI约2.36小时/(用户·年)。读者可以拿这个结果做基准验证自己写的版本有没有对齐。实际工程里网架会比这个demo大得多。我一般不会手写impact_time而是从电网GIS导出馈线拓扑再用一段脚本自动生成“元件—影响负荷点”矩阵。但自动化的前提还是要有一套明确的开关逻辑否则矩阵生成出来也是错的。2.3 最小路法的适用边界辐射网好用遇到转供要小心最小路法最大的优势是快、透明、结果确定。同样的网络跑一百次结果都一样不需要担心随机误差。所以它特别适合规划阶段的方案比选十个接线方案每个都算一遍SAIFI/SAIDI几分钟就能排个高低。但边界条件也很硬。第一多电源或环网运行时负荷点有两条及以上供电路径最小路不再唯一传统最小路法会低估可靠性。第二网架中存在联络开关时非故障段可能通过转供恢复供电这段恢复时间既不是修复时间也不是0必须在故障影响分析里显式建模。第三分布式电源接入后可能形成孤岛运行这种状态用最小路法表达非常别扭。所以我在实际项目里的判断是纯辐射状、没有分布式电源、转供关系简单的馈线直接上最小路法出现环网转供、微网孤岛、电气距离复杂的情况就换蒙特卡洛法或者在最小路法里手动加入“转供恢复时间”这个修正项。别把最小路法当作万能工具它在简单场景下是利器在复杂场景下就是坑。3. 蒙特卡洛法模拟配电网序贯抽样的Python实现与收敛判据3.1 序贯与非序贯蒙特卡洛怎么选蒙特卡洛法的本质是用大量随机样本来逼近真实的故障-修复过程。配电网里最常用的两种是序贯和非序贯。非序贯方法对每个元件按“当前时刻处于故障状态的概率”独立抽样生成一个系统状态再统计该状态下的停电情况。它不考虑这些事件在时间上的先后计算快但很难处理时变负荷、天气变化和修复顺序等带时间约束的因素。序贯方法则沿着时间轴推进先抽样每个元件的正常运行时间再抽样故障持续时间形成一条连续的状态演化序列。它保留了故障发生的先后顺序也天然兼容“多个元件同时故障”“故障后转供”“负荷随时间波动”这些场景。代价是计算量大得多尤其当元件数量多、故障率低的时候可能跑几千“系统年”才能让指标稳定下来。我的选择原则是只算年期望值、网架规模大、需要快速筛选场景用非序贯要做精细化评估、研究某个具体故障序列或验证转供策略用序贯。实际工作中可靠性数据往往带季节性故障率和修复时间在不同天气下差别很大所以我个人更常写序贯版本宁可多等几分钟。下面用一个带两负荷点的辐射网演示序贯蒙特卡洛的完整流程这段代码可以直接改成你自己的网络。3.2 一个可运行的序贯蒙特卡洛脚本从元件状态序列到系统指标这个脚本的思路是对每条边用指数分布抽样正常运行持续时间和故障修复时间生成若干年的故障事件列表把所有边的事件按时间排序然后逐个事件把故障边从网络里移除用图的连通性判断每个负荷点是否失电。为了方便展示我假设故障事件在时间上不重叠这对故障率低的小型网络是合理的简化。import numpy as np import networkx as nx # 与2.2节相同的网络和参数 edges [ ((S, A), 0.20, 5.0), ((A, B), 0.10, 4.0), ((B, LP1), 0.10, 3.0), ((A, C), 0.30, 5.0), ((C, LP2), 0.15, 3.0), ] G nx.Graph() G.add_nodes_from([S, A, B, C, LP1, LP2]) for (u, v), lam, r in edges: G.add_edge(u, v, lamlam, rr) loads {LP1: 120, LP2: 80} source S years 20 rng np.random.default_rng(2024) # 固定种子方便复现正式评估不要固定 def simulate_component(lam, r, horizon_hours, rng): 生成单个元件的故障时段单位统一为小时 t 0.0 events [] while t horizon_hours: t rng.exponential(8760 / lam) # 正常工作时间 if t horizon_hours: break down rng.exponential(r) # 故障修复时间 events.append((t, min(down, horizon_hours - t))) t down return events horizon years * 8760 all_events [] for (u, v), lam, r in edges: for start, down in simulate_component(lam, r, horizon, rng): all_events.append((start, start down, (u, v), down)) all_events.sort(keylambda x: x[0]) # 统计负荷点停电次数和年停电小时 outage_count {lp: 0 for lp in loads} outage_hours {lp: 0.0 for lp in loads} for start, end, edge, down in all_events: u, v edge G.remove_edge(u, v) for lp in loads: # 故障期间负荷点是否与电源不再连通 if not nx.has_path(G, source, lp): outage_count[lp] 1 outage_hours[lp] down G.add_edge(u, v) # 输出结果与最小路法结果对拍 for lp in loads: avg_freq outage_count[lp] / years avg_hours outage_hours[lp] / years print(f{lp}: 年均停电 {avg_freq:.3f} 次, 年停电时长 {avg_hours:.3f} 小时) safi sum(loads[lp] * outage_count[lp] for lp in loads) / (sum(loads.values()) * years) saidi sum(loads[lp] * outage_hours[lp] for lp in loads) / (sum(loads.values()) * years) print(f模拟SAIFI{safi:.3f}, 模拟SAIDI{saidi:.3f})这段代码里最需要注意的单位统一lam是次/年r是小时而抽样时需要把故障率换算成小时尺度所以正常工作时间用8760 / lam。rng.exponential抽的都是指数分布变量这是标准做法默认元件寿命和修复时间服从指数分布。运行后你会看到模拟值接近最小路法的结果但不会完全一致。比如我曾用同一组参数跑到SAIFI0.55左右、SAIDI2.2~2.5小时这个波动就是抽样误差。如果你把years改成200误差会缩小但计算时间也会线性增加。对辐散状的网络序贯蒙特卡洛真正的费用在多次连通性判断上节点多了以后不建议每次都全图搜索。3.3 收敛判据与仿真年限抽多少次才能信很多新手拿到蒙特卡洛结果就问“这个数能不能用”。我一般不会只看一次运行结果而是看两个东西年度指标的变异系数以及置信区间。变异系数 CV 是年指标标准差除以均值CV 越小说明结果越稳定。常见做法是把 CV 控制在0.05以内否则就继续增加仿真年限。一段简单的判据代码可以这样写import numpy as np # yearly_values 是一个数组每个元素是某一年的SAIDI或SAIFI def check_convergence(yearly_values, target_cv0.05): mean np.mean(yearly_values) std np.std(yearly_values, ddof1) cv std / mean if mean 0 else float(inf) return cv target_cv, cv这里的yearly_values怎么来在序贯模拟中可以把整个长模拟切分成一年一段每一段单独统计一次停电次数和停电时长在非序贯模拟中每次抽样生成一个系统状态相当于一年。把这些年值积累下来CV就自然算出来了。要特别提醒别为了CV好看而用固定随机种子把结果卡死。固定种子只能用来复现脚本不能证明模型正确。正确的做法是跑不同年数、不同随机种子看结果是否在合理区间内波动。如果SAIDI每次变化超过20%说明仿真年限不够而不是方法有问题。对故障率偏低的配电网模拟几千年系统年仍然不收敛的情况我都遇到过这种时候不要盲目加年限先检查是不是有某个小概率高影响事件没建模再考虑用重要抽样这类方差缩减技术。4. 两法对拍用最小路法校准蒙特卡洛的五个关键参数4.1 用最小路法做基准给蒙特卡洛结果套置信区间最小路法算出的是一个确定性期望值蒙特卡洛算的是一堆随机样本的均值。两者既然面对同一张网络和同一组参数结果就不该差太远。所以我的习惯是每次跑完蒙特卡洛都拿最小路法结果当靶子看模拟结果是否落在置信区间内。置信区间可以简单用经验公式算import numpy as np def conf_interval(values, alpha0.05): mean np.mean(values) std np.std(values, ddof1) n len(values) z 1.96 if alpha 0.05 else 2.58 return mean, mean - z*std/np.sqrt(n), mean z*std/np.sqrt(n)如果你把模拟年份切分成多个观测样本这里的values就是每年的SAIFI或SAIDI。以2.2节最小路法得到的SAIFI0.54、SAIDI2.36小时为基准如果蒙特卡洛的95%置信区间是0.48~0.60和2.1~2.6小时那这组模拟是可信的如果区间完全偏离靶心就要回头查代码或参数。对拍还有一个用处当最小路法因为网架复杂而无解时它可以作为蒙特卡洛结果的一个近似参考点。反过来蒙特卡洛也能暴露最小路法没有建模的那些随机波动。所以两种方法在我的项目里总是成对出现不是互斥。4.2 五个关键参数取值量级与对指标的影响不管用哪种方法下面这五个参数决定成败。我把常用量级和影响方向整理成表供参考参数常见取值量级主要影响线路故障率 lambda架空线0.05~0.5次/km·年电缆0.01~0.1次/km·年近似正比影响SAIFI和SAIDI故障修复时间 r架空线3~6h电缆4~8h正比影响SAIDI和CAIDI隔离开关/联络开关操作时间 s人工0.5~2h自动开关2~10min影响转供场景的SAIDI负荷点用户数与负荷分布按配变容量和低压台区折算改变系统指标的加权权重仿真年限/抽样数序贯至少1000~2000系统年非序贯更高决定蒙特卡洛抽样误差故障率数据我一般优先用当地供电公司发布的可靠性统计没有就按设备型号和运行年限经验估算。修复时间也要看抢修资源偏远地区普遍偏高。开关操作时间则要区分“人工到现场操作”和“遥控自动操作”两者能差一个数量级最容易被拍脑袋带偏。4.3 故障率和转供时间的影响差异为什么不能只看平均值先说故障率敏感性如果某条馈线的故障率从0.1次/年变成0.2次/年而其他参数不变SAIFI近似翻倍SAIDI也接近翻倍。这是因为两个指标都跟故障次数线性相关。但如果你在北方的冬天把架空线路故障率调高到0.4修复时间从4h变成8hSAIDI的变化会比SAIFI更夸张CAIDI也会明显升高。这个差异说明SAIFI主要反映设备故障频率SAIDI还受抢修资源和开关策略影响做敏感性分析时不能只看一个数。转供参数的影响方向则相反它几乎不影响SAIFI因为停电次数已经发生但能显著减少SAIDI。一个原本要停电6小时的非故障段如果通过联络开关在1小时内转供SAIDI就降了5小时。所以在有联络开关的网架里对拍时如果发现蒙特卡洛SAIDI比最小路法小很多先别急着说代码错了很可能是因为最小路法把转供时间设成了0而蒙特卡洛的开关模型把它放进了流程里。我建议每个项目至少做三组敏感性一组把故障率整体上下调20%一组把修复时间上下调50%一组把转供时间在“人工2h/自动10min/不转供”三种模式中切换。这样能快速找到指标波动的最大来源也能判断手上数据够不够支撑评估精度。5. 避坑指南可靠性评估中常见的6个翻车现场5.1 现象最小路法算出的SAIFI比统计值小一截常见做法是只把负荷点最小路上的元件或者只把主干道的故障率累加忽略了分支线故障和馈线出口断路器对短时停电的贡献。我接过一个项目进线故障率明明有0.2次/年但团队只统计了馈线主干段结果SAIFI少了将近两成。原因不是公式错而是故障元件的影响范围没列全。解决办法是回到FMEA表把每个元件影响哪些负荷点、影响多久重新过一遍尤其是变电站出线以前的元件和分支线路必须显式写进impact_time。5.2 现象蒙特卡洛每次跑出来的SAIFI差很多第一次跑0.50第二次跑0.62第三次又回到0.53。如果你只跑了500个系统年这种现象在故障率低的配电网里非常正常因为年故障事件本来就稀疏偶发一次大故障就能把均值拉高。解决思路是先算变异系数低于0.05再定数据再就是增加仿真年限或者改用分层抽样让偶发大故障更平均地出现在样本里。不要试图用固定随机种子消除波动那是自欺欺人。5.3 现象把年故障率当成单次概率结果直接指数爆炸在非序贯蒙特卡洛里有人会写if np.random.uniform() lam: 故障。这行代码看似合理实际上完全错误因为lam是“每件每年故障多少回”不是“某次抽样的概率”。用0.1次/年去做单次概率相当于默认每年发生十几次故障得到的SAIFI会比真实值大一个数量级。正确做法是用泊松分布抽样一年内的故障次数np.random.poisson(lam)或者干脆按时间轴做指数分布抽样。我见过太多翻车都出在这一行值得反复提醒。5.4 现象变压器低压侧故障被计入系统指标配电网可靠性评估的统计边界一般是变电站10kV出线开关以下到配变低压侧出口。如果把低压线路、用户内部故障也算进来SAIDI会高得离谱。原因是这些部分不属于中压配网运维单位的管理范围数据口径不一致。解决办法是在建模一开始就建立边界规则只统计10kV馈线、分段开关、联络开关、配变高压侧低压侧一律不参与系统指标汇总。如果用户要求评估到户就单独做低压可靠性模块不要混在中压模型里。5.5 现象模拟了几百年收敛判据还是不达标这种情况往往不是样本量不够而是模型里有小概率高影响事件比如一条主干线一年发生一次台风倒杆一停就是十几个小时。这种事件对SAIDI方差贡献极大普通随机抽样需要极多年份才能稳定它。解决方法是把天气状态单拎出来正常天气用一个故障率恶劣天气用更高的故障率和更长的修复时间再做状态抽样。也可以引入重要抽样但工程上先把天气状态建好基本就能解决问题。5.6 现象联络开关的转供时间没进模型有联络开关的网架非故障段往往可以在隔离后转供复电停电时间不是修复时间而是“隔离合联络开关”的操作时间。我在一个新区网架上吃过亏最小路法模型里没写转供结果SAIDI算出来2.8小时实际运行只有1.9小时。反过来如果把自动化开关的转供时间设成0又会把SAIDI算得过于乐观。正确做法是给每条联络开关设置一个转供时间在故障影响分析里对每个负荷点取min(修复时间, 隔离时间转供时间)。这个逻辑要写清楚别一个开关时间套全部线路。6. 把评估做成参数驱动的工具批量场景跑法与置信区间报告前面几张图的代码都是一次性脚本但实际做规划时往往要比较几十组方案。我现在的做法是把所有参数放到一张CSV里一行动一个场景脚本循环读取后输出指标表。核心思路是拓扑结构写在图构建函数里所有故障率、修复时间、操作时间、负荷用户数都从参数表读入这样换网架只需要改一行拓扑换数据只需要改CSV。给一个最小的参数驱动骨架import pandas as pd import numpy as np # params.csv包含列: scene_name, lam, r, switch_time, customers df pd.read_csv(params.csv) results [] for _, row in df.iterrows(): lam row[lam] r row[r] # 调用你自己的最小路法或蒙特卡洛函数 saifi run_min_path(...) # 具体函数见第2、3章 saidi run_monte_carlo(...) results.append({scene: row[scene_name], SAIFI: saifi, SAIDI: saidi}) result_df pd.DataFrame(results) result_df.to_csv(reliability_report.csv, indexFalse)批量场景跑完还要给每个SAIDI补一个置信区间用4.1节的conf_interval就行。报告里我会把每次蒙特卡洛的均值、置信区间、样本年限写在一起不让别人只看到一个干巴巴的点估计。那种“最小路法算出来0.54蒙特卡洛算出来0.55”的结论只有带上置信区间才有意义。我现在的习惯是拿到任务先不写代码先画包含联络开关、分段开关和保护配置的网架图。画完再决定用最小路法还是蒙特卡洛法然后才动手写脚本。这个顺序让我少翻了好几次车。做配电网可靠性评估最重要的不是算法多高深而是边界、参数和转供逻辑都清清楚楚。希望这个习惯也能帮到你。本文还有配套的精品资源点击获取