基于序贯蒙特卡洛的配电网可靠性评估：Matlab实现与避坑指南

做配电网可靠性评估，最常被人问的一句话是：这套网架一年到底会停几次电、平均停多久、折算到用户侧又是什么水平？这些问题靠拍脑袋估不准，尤其当网架里出现分段开关、联络线路、分布式电源以后，故障影响关系变得特别绕。我前段时间刚好把一套基于序贯蒙特卡洛模拟法的配电网可靠性评估流程用Matlab完整跑通了，从元件故障抽样、时序状态推演、负荷点停电记录，到指标统计与收敛判断都做了落地。这篇文章把整体思路、代码结构和调试过程中踩过的坑整理出来，给正在做可靠性方向课程设计、毕业论文或者规划方案比选的读者做个参照。

1. 配电网可靠性评估到底在算什么：从停电到指标

1.1 停电是怎么发生和传导的

配电网可靠性评估不是简单统计某台变压器坏了几次就结束了，它关心的是“故障从元件出发，经过开关动作与网络重构之后，最终影响到了哪些负荷点”。这里面的因果关系链条往往比想象中长。

举个最简单的中压辐射网例子：一条馈线带着若干个负荷点，沿线装了分段开关，馈线出口有断路器。当某一段线路发生永久性故障时，变电站出口断路器会先跳闸，导致整条馈线失压——这是“故障隔离前”的停电影响。随后保护动作完毕，分段开关通过人工或遥控方式打开，将故障段隔离出来；这时故障段上游的负荷点可以由变电站恢复供电，故障段下游的负荷点如果存在联络线，则可能倒换到相邻馈线继续供电，如果没有联络就只能等待故障修复。

这就产生了一个关键问题：同一个故障事件，不同位置的负荷点，停电时间和停电次数可能完全不一样。有的用户只停“断路器跳闸到分段开关操作完成”的几分钟，有的用户要停到故障修复完成，有的用户还要算上联络开关倒闸的时间。可靠性评估要做的，就是把每一次故障事件以及它引发的所有负荷点停录事件一一记录下来，最后汇总成系统级指标。这也是为什么可靠性评估不能只看元件故障率，而要和网络拓扑、开关配置、运行方式绑在一起分析。

1.2 核心可靠性指标体系与公式

行业里衡量配电网可靠性，绕不开下面这张表里的几个经典指标。很多时候一个项目的最终交付，就是给业主一张这样的指标汇总表。

SAIFI和SAIDI是绝大多数报告里最显眼的两个数。ASAI则更直观——99.9%的可用率听起来很美好，但折算一下就是每年约8.76小时的停电时间，对用户来说体感一点也不轻。ENS经常被用来做停电损失评估，把电量指标和电价、产损系数相乘，就能算出一个网架方案对应的停电经济损失，这也是可靠性评估经常被规划项目拉来做方案比选的原因。

理解这些指标是理解序贯蒙特卡洛的第一步，因为无论用什么算法，最终输出都在朝这套指标体系收敛。

2. 为什么是序贯蒙特卡洛：三种评估路线的横评

2.1 解析法：靠谱但活儿只能干一半

配电网可靠性评估的传统主力是解析法，尤其以故障模式与后果分析法（FMEA）为代表。思路很直接：枚举所有可能的元件故障事件，对每个事件做网络影响分析，再乘上故障率和修复时间，叠加得到系统指标。

解析法的优势是精确、可解释，对简单辐射网和规模很小的配电网非常友好。手工就能列出几条馈线的所有故障模式，再根据分段开关位置划出停电范围。但一旦网络规模变大，或者引入多电源、分布式电源、储能和复杂的运行策略，枚举组合会指数膨胀。比如一条馈线上出现两个元件同时故障，或者一个故障叠加一次计划检修，类似场景的排列组合数量直接起飞，解析法要么算不动，要么得做大量简化假设。

另外还有更现实的一点：解析法处理时变负荷非常别扭。传统FMEA通常假设负荷恒定，或者只按峰值负荷计算，这样ENS的精度会打折扣。当你需要评估“不同季节、不同时段的负荷水平对可靠性指标的影响”时，解析法的公式会越写越复杂，最后维护成本高到让人想放弃。

2.2 非序贯蒙特卡洛与序贯蒙特卡洛：差在一个“时间”

蒙特卡洛模拟法给了另一条路：不枚举，而是抽样。由随机数驱动元件状态，通过大量重复模拟来逼近真实可靠性水平。蒙特卡洛内部又分两派：非序贯蒙特卡洛和序贯蒙特卡洛。

非序贯的思路是“拍快照”：随机抽样一年中某一天甚至某个小时，判定该时刻系统状态，抽样几万次后统计平均。计算效率高，但它天然缺乏时间先后顺序，难以精确处理需要时间逻辑的细节，比如元件故障后修复持续多久、什么时候恢复供电、负荷随时间的变化曲线怎么走。很多情况下非序贯法只能给出期望值，给不出时间相关的分布信息。

序贯蒙特卡洛则是按照时间轴一步步推演。它给每个元件建立“运行-故障-修复-再运行”的状态时间序列，所有元件并行推进，从系统全年的时序响应里统计停电事件。这样负荷曲线能真实参与计算，分布式电源的出力波动也能按小时建模，检修计划、季节差异这类时间性很强的影响都可以塞进去。我所配电网项目里经常被问到“某个方案的年停电时间是多少”，序贯蒙特卡洛给出的结果可以直接落到SAIDI、SAIDI这些指标上，和实际运行统计口径更接近。

2.3 什么时候该果断放弃序贯蒙特卡洛

不过序贯蒙特卡洛不是万能药，它有两个绕不开的痛点：计算量大、收敛慢。因为要逐小时甚至逐分钟推进时序状态，一个规模稍大的配电网，仿真一万年（模拟时间）往往意味着上亿步运算，每步还要做网络影响分析，跑起来确实慢。不少人在项目里一上来就选序贯蒙特卡洛，跑到一半发现收敛判据迟迟不达标，才意识到前期估算不足。

我的判断标准是这样：如果网络比较简单、负荷模型可以取恒定值、又需要快速得到比较确定的指标，直接用解析法。如果网络中等复杂度、只关心期望指标而不强调时序特性，非序贯蒙特卡洛配合方差缩减技巧，性价比很高。只有当场景里有明显的时序逻辑——比如时变负荷、储能充放电策略、分布式电源季节性出力、复杂的转供策略——序贯蒙特卡洛才有不可替代的价值。另外，正式开跑前最好用非序贯或解析法粗算一轮，给序贯仿真提供一个初始解和参考量级，能提前发现参数错误，避免白跑几十万步才发现模型根本不对。

3. Matlab仿真主流程拆解：从拓扑参数到指标输出

3.1 元件可靠性建模与参数设定

无论是哪种蒙特卡洛，第一步都是把网络和参数变成程序能处理的数据结构。配电网可靠性仿真至少需要四类输入。

第一类是拓扑信息。常见的建模方式是把网络抽象成节点和支路，支路类型包括线路、变压器、分段开关、联络开关等。每个负荷点挂在某个节点上，并关联用户数、平均负荷和时序负荷曲线。Matlab里我习惯用结构体数组存节点，用邻接表表达连接关系，这样后续做连通性判断和影响分析都方便。

第二类是元件可靠性参数。线路和变压器需要给出故障率与平均修复时间。故障率单位一般用“次/(km·年)”，注意乘以线路长度才是单条线路的年故障率。架空线路典型故障率在0.05~0.2次/(km·年)之间，电缆会低不少，修复时间通常在2~8小时之间。这些参数可以来自电网运行统计，也可以参考经典文献的典型值。仿真结果的绝对数值完全由参数决定，所以参数来源必须记录清楚。

第三类是运行参数。比如分段开关操作时间、联络开关倒闸时间、计划检修周期、保护动作时间等。这些时间通常只有十几分钟到几小时，但对SAIDI的影响不小，不能忽略。

第四类是负荷模型。如果只有全年平均负荷数据，可以把负荷曲线设成单位比例曲线再乘以年峰值。更精细的做法是按季节、工作日/休息日、24小时分段给曲线。在Matlab里可以准备一个维数为[小时数, 负荷点数量]的矩阵，模拟推到某一小时时直接查表得到该时刻负荷。

3.2 序贯抽样的核心代码与状态推进

序贯蒙特卡洛最核心的机制是：对每个元件抽样它的无故障工作时间TTF和故障修复时间TTR，然后让所有元件的状态沿着时间轴交替变化。

抽样公式是可靠性工程里的老朋友。

% lambda: 元件的年故障率 % MTTR: 平均修复时间，单位小时 % rand: [0,1)均匀分布随机数 % 抽样无故障工作时间（单位: 年） u1 = rand; TTF = -log(u1) / lambda; % 抽样修复时间（单位: 小时） u2 = rand; TTR = MTTR * (-log(u2));

这里用到的原理是逆变换法：当随机变量服从指数分布时，通过-ln(U)/λ就能从均匀随机数得到对应的采样值。元件“运行-故障”交替过程接近指数分布假设，所以这套抽样公式在工程上是主流选择。

主循环的写法有很多种，核心思路都是“下一次事件驱动”。每轮循环找到所有元件里最先发生的状态转移，把系统时钟推进到那个时刻，然后更新发生转移的那个元件状态，再做故障影响分析。

Clock = 0; % 系统仿真时钟，单位年 SimYears = 5000; % 仿真总年数 N = length(lambda); % 元件数量 State = ones(1, N); % 1=运行, 0=故障 NextChangeTime = zeros(1, N); % 初始化: 给每个元件抽样首次状态转移时间 for k = 1:N NextChangeTime(k) = Clock + (-log(rand) / lambda(k)); end while Clock < SimYears % 找到最早发生转移的元件 [dt, idx] = min(NextChangeTime - Clock); Clock = NextChangeTime(idx); % 状态翻转: 运行->故障 或 故障->修复 State(idx) = 1 - State(idx); if State(idx) == 0 % 刚发生故障，做故障影响分析 DoOutageAnalysis(idx, Clock); % 安排下一次修复完成时间 NextChangeTime(idx) = Clock + MTTR(idx) * (-log(rand)) / 8760; else % 刚修复完成，安排下一次故障时间 NextChangeTime(idx) = Clock + (-log(rand) / lambda(idx)); end % 若有计划检修，这里还需要插入检修状态及其持续时间 end

这段代码里有个小坑：时间单位很容易混。故障率λ单位如果是“次/年”，-log(rand)/lambda得到的TTF单位就是“年”；而修复时间一般按小时给，如果MTTR的单位是小时，就得除以8760换算成年，否则系统时钟会乱掉。我见过不少人在这个单位换算上翻车，改起来倒简单，但定位起来特别费劲。

3.3 故障影响分析：谁在停电、停多久

抽样完成之后，最关键的环节是故障影响分析。这一步要做的事情是：给定某个元件在当前网络状态下故障，找出哪些负荷点失去供电，以及各自要停多久。

在辐射状配电网里，这个分析可以按“上游、故障段、下游”三个区域来划分。

首先，找到故障元件所在的馈线以及变电站母线。故障发生后，馈线出口断路器跳闸，保护动作瞬间整条馈线失压。接下来隔离故障段：分段开关打开后，故障点上游且与母线连通的部分恢复供电，恢复时间是开关操作时间。故障点下游的部分，要看有没有联络开关能够倒供。如果存在联络开关且对侧电源容量足够，带电合联络开关后，下游恢复供电，恢复时间为倒闸时间；如果下游没有联络通路，这部分负荷就只能等故障修复完成再恢复。

注意，实际项目中“容量足够”这四个字背后也要判断。联络馈线的剩余容量不够，合闸后可能过载，这时候要么削减部分负荷，要么放弃转供。严谨的仿真应该在倒供判断里加入线路容量约束，简单实现时可以查一回馈线额定容量与当前负荷之差。

故障影响分析的结果要记录到按负荷点维度的数组里：本次停电事件持续时间、停电时刻、停电期间缺供电量。由于序贯蒙特卡洛按时间轴推进，每个负荷点在仿真年份里的所有停电区间都能完整记录，天然支持后续指标统计和时间分布分析。

3.4 指标统计与收敛判据

仿真结束后，指标统计就是纯算术。把仿真期间每个负荷点的停电次数、停电时长、缺供电量累加起来，再按用户数加权，就能得到表1里的各个指标。

totalCustomers = sum(CustomerNum); % 总用户数 SAIFI = sum(InterruptionCount) / totalCustomers; SAIDI = sum(InterruptionDuration) / totalCustomers; CAIDI = SAIDI / SAIFI; ASAI = 1 - SAIDI / 8760; ENS = sum(EnergyNotSupplied); % 单位: kWh/年

除了算指标，蒙特卡洛仿真一定要做收敛判断，否则仿真结束得很随意。这里引入方差系数β，一般定义为目标指标估计值的标准差与均值之比。工程上常用的收敛标准是β小于0.02~0.05，对应的仿真年数通常在几千到几万不等，视网络规模和指标类型而定。

% 以SAIFI为例，计算逐年样本的方差系数 betaSAIFI = std(SAIFI_year) / (mean(SAIFI_year) * sqrt(nYears));

更稳妥的做法是边跑边检查，实现“自适应停”：每仿真200年检查一次β，达到阈值就停止，达不到就继续。这样比固定仿真年数更科学，因为不同网络参数下达到相同精度所需的仿真年数可能差一个数量级。如果每次跑完发现不收敛，优先考虑增加仿真年数或降低指标方差，而不是盲目修改随机数种子。

4. 实操中的坑与排查技巧

4.1 仿真年数、方差系数与随机数种子的三角关系

第一次跑序贯蒙特卡洛的人，最容易上来就设一个“看起来很大”的仿真年数，比如1000年，然后默认数据够用了。实际上，1000年对于结构比较复杂、故障率较低的网络，SAIDI的方差系数可能还在0.1以上，指标抖动非常明显。

我的经验是：先在Matlab里做一个“年数-指标”的敏感性测试。从1000年开始，每次翻倍仿真年数，观察SAIFI和SAIDI估计值的变化幅度。如果2000年的结果和4000年的结果差得不多，几千年的量级基本够用；如果差异还很大，继续翻倍。这个测试本身也花不了多少时间，但能避免两种尴尬情况：一种是不够收敛就出报告，另一种是不必要地堆仿真年数，白白让程序跑一下午。

随机数种子也需要刻意管理。MATLAB里每次启动rand的状态不同，如果不设置全局种子，同样是5000年仿真，两次结果也会有可观察到的偏差。正规做法是在主程序入口固定种子：

rng(2024); % 固定随机数种子，保证结果可复现

这在调试阶段特别重要。当开发人员修改了某段影响分析逻辑后，固定种子能保证两次运行只因为代码变化而出现差异，排除随机噪声的干扰，比对结果才有说服力。

4.2 时变负荷处理：恒定功率会算错什么

很多初版程序为了省事，把所有负荷点都设成全年恒定功率，结果ENS算出来往往和真实情况偏差很大。原因很直观：停电事件发生在凌晨低谷负荷时段，损失的电量比发生在晚高峰时段小得多。恒定功率模型相当于把每个时段都按同一个负荷水平计，必然扭曲缺供电量。

要修这个问题，我建议按三条曲线来做负荷模型：工作日负荷曲线、周末负荷曲线、节假日负荷曲线。Matlab里可以预先生成一个全年8760小时的负荷比例矩阵，每行对应一个小时的负荷倍数，然后与负荷点的年峰值功率相乘。系统时钟推进到某小时时，直接查表得到该小时每个负荷点的功率值。这个做法的代码量不大，但对ENS和SAIDI的结果影响很明显。

如果项目实在没有负荷曲线数据，至少要分“峰、平、谷”三段负荷水平，而不是一个常数走到底。我自己接手过一个评估项目，最初的恒定负荷模型ENS比修正后高出约18%，工程上这个差距足以改变方案排序了。

4.3 联络开关与孤岛逻辑的边界判断

另一类高频bug出现在联络转供逻辑里。最常见的错误是“无条件合联络”：只要故障段下游有联络开关，就认为负荷能全部转供，不管联络馈线是否过载、对侧是否还有剩余容量。改进逻辑至少应该加三个判断。

第一，联络开关对侧的馈线必须处于运行状态，不能那边也在检修。第二，转供后联络馈线的载流量不能超过上限，这个要结合当前时刻负荷算。第三，转供操作时间要合理计入，一般为几分钟到几十分钟，不是瞬时完成。盲目合联络的结果是SAIDI被低估，尤其重负荷时段故障时，实际调度员根本不敢直接合闸，程序却默认合上了。

还有边界条件：故障元件如果恰好在联络开关所在支路上，或者两个联络开关同时可用，影响分析逻辑要分清楚谁先合、合哪个，不同顺序可能影响后续容量分配。遇到这种逻辑，我建议先把情况拆成单故障事件表，手算一遍期望结果，再用代码跑同一事件对比，快速定位逻辑错误。

4.4 常见问题速查表

5. 一个算例，看指标怎么变成决策依据

5.1 典型三馈线系统的仿真结果长什么样

拿一个我实际做过的模拟项目举例。场景是某园区规划三回10kV馈线，每回馈线带8个负荷点，线路总长约18公里，全部为架空线，故障率取0.08次/(km·年)，平均修复时间4小时。分段开关配置在每回馈线中段，两回馈线之间有一处联络开关，故障后可倒供部分负荷。所有负荷点年峰值负荷之和约5200kW，全年负荷曲线按时段变化。

用上面的Matlab流程固定种子跑5000年，输出的核心指标大概是：

这类绝对值本身并没有标准答案，完全依赖参数设定和网架结构，但它能反映出方案之间的相对差异。同一模型把联络开关取消掉以后，SAIDI直接升到5.8小时以上，ENS翻了一倍多，这就是可靠性评估最有决策价值的时刻——多花多少钱加一个联络开关，能换回多少停电损失，一目了然。

5.2 怎么从负荷点指标反查网架薄弱环节

只看系统级指标不够，规划里通常还要看每个负荷点的个体指标。我每次跑完仿真，一定会输出一张负荷点指标清单，然后按SAIDI从大到小排序。排在前面的负荷点，要么处于馈线末端，要么处于转供盲区，要么紧邻高故障率线段。

这类分析对规划很有帮助。比如某负荷点SAIDI明显高于系统平均值，一看拓扑，它处于故障段下游且没有联络通路，所有上游故障都要等修复完成才能送电。对策可以是加装联络开关，也可以调整线路分段位置来缩小故障隔离范围。可靠性评估从“算个数”变成“定位短板”，靠的就是这些负荷点级别的中间结果。

我实际做项目时还会把负荷点ENS按年汇总后再除以该负荷点的年用电量，得到“损失电量占比”，这个相对量对于比较不同性质的用户（比如工业用户与商业用户）更公平，也方便进一步做停电经济损失评估。

最后说几句体会

做了这么多轮可靠性仿真，我自己最大的感受是：这套方法的难点其实不在蒙特卡洛抽样，而在故障影响分析逻辑。抽样只是把随机性送进来，判断哪些用户停、停多久、什么时候恢复，才是决定指标准不准的核心。实际项目里，我习惯把故障影响分析单独写成函数，用一组手工构造的事件逐项测试，确认无误后再接入蒙特卡洛主循环。另外再分享一个小技巧：正式跑大批次仿真之前，先固定随机数种子跑一轮，记录下指标输出，以后每次改代码都用同一批随机数做回归对比，避免算法改进半天，最后发现是随机波动在干扰判断。这套Matlab流程目前还在持续扩展，后面计划把储能充放电策略和分布式电源出力的时序模型加进去，进一步贴近真实的运行场景。