做配电系统规划这个方向也有好几年了最常被问的一句话是到底怎么在投资预算和供电可靠性之间找一个既能向董事会交代、又能让运维少挨骂的方案只算经济账往往抠门到最后故障频发投诉工单一堆只追可靠性冗余投资又会撑爆预算。所以这两年我越来越倾向于把经济性和可靠性当成两个平行的优化目标放进同一个框架里而不是先只做单目标再回头校核。这篇文章就来拆一套基于经济与可靠性双目标的混合配电系统规划及可靠性评估 Python 实现从双目标建模、DG 与储能建模、可靠性仿真、NSGA-II 求解到常见的坑位全部按我实际跑代码的思路来讲。适合正在做配电网规划、分布式电源接入研究或者可靠性评估的电力专业学生也适合想快速把原型跑通的工程师。1. 项目整体拆解双目标规划与可靠性评估在解决什么问题1.1 为什么把“经济”和“可靠性”放进同一个优化框架很多人习惯的做法是先以投资和运行成本最小为目标算一版方案然后把可靠性指标作为校核比如要求 SAIDI 不超过某个值不满足就在个别节点加装设备加到满足为止。这种思路的问题在于可靠性并不是均匀分布的你补上校验值不等于所有薄弱点都被修复而且单目标优化天然倾向于去“省成本”很可能得到一个刚好擦线、抗风险能力很弱的方案。反过来如果只追求可靠性DG、储能、开关和线路冗余会越加越多边际成本一路飙升最后投资效益极差。经济性和可靠性本质上是一对互相拉扯的指标所以把它做成双目标优化更符合工程直觉最后输出一条 Pareto 前沿前沿上的每个点都代表“在某种投资水平下能做到的最好可靠性”决策者再根据自己的风险偏好挑选方案。可以用生活类比理解 Pareto 前沿买手机时预算和配置互相打架预算有限又想拍照好不可能同时满足两个极端但每个价位都有一个“该价位性价比最高”的机型把这些机型连起来就是 Pareto 前沿。双目标规划输出的就是配电系统里的这一条“机型曲线”。1.2 混合配电系统的“混合”体现在哪三个层面标题里的“混合配电系统”不是说配电网多了几台风机而是系统组成方式有了本质变化。我自己习惯把它拆成三个层面来建模电源层面的混合常规网络电源、光伏、风电、储能、微型燃气轮机可能同时存在。光伏和风电是不可控、间歇性电源储能是可控的功率型/能量型元件微型燃气轮机虽然运行成本高但可调度性强。规划时不能简单用一个“等值电源”替代必须按各自出力模型单独处理。网架形态层面的混合传统配电网是单电源辐射状网络而含 DG 的配电网在故障后可以转入孤岛运行甚至刻意划分微网。可靠性和潮流计算都得考虑“故障后联络开关合闸”、“DG 孤岛支撑”这类动态过程。规划变量层面的混合变量里有连续量光伏容量、储能能量、离散量节点编号、台数、甚至整数量变压器分接头档位、电容器组数。不加处理直接丢给求解器通常是一个混合整数非线性规划MINLP收敛性非常差。所以我的做法是把 MINLP 拆成两层外层用进化算法搜索选址定容方案内层用确定性潮流和可靠性评估计算每个方案的目标值。这样外层是“黑箱多目标”内层是“工程计算引擎”结构清楚也方便逐步加复杂约束。1.3 规划与评估闭环Python 在这里到底扮演什么角色规划不是算一遍账就结束的事。外层的每一个候选方案都要经过潮流计算和可靠性评估然后一次次迭代形成“方案生成→评估→进化→再评估”的闭环。这个闭环用 Python 写特别合适——numpy 处理矩阵运算pandas 管理节点和支路参数表networkx 做网络拓扑连通性分析scipy 可以做简单优化和数值积分matplotlib 画 Pareto 前沿多目标优化直接用 pymoo 库避免自己从零实现 NSGA-II。我见过不少课题组还在用 Matlab 写老代码也不是不行但 Python 在数据接入、后续扩展上优势明显负荷数据可能是 Excel、CSV甚至直接从 SCADA 导出的表pandas 一行就能读进来如果要接入项目数据库或者做 Web 展示Python 生态也顺得多。如果你还没配好环境建议先把 numpy、pandas、networkx、scipy、matplotlib、pymoo 这六个库装上版本冲突基本不会遇到直接 pip install 就行。2. 核心模型与关键参数哪些细节决定了计算结果2.1 经济性目标怎么折算成年值经济目标我习惯取“年综合费用”单位是万元/年。它通常由四块组成投资等年值、运行维护费、网损费用、以及与主网交互的电费。投资不是一次性算完而是用资本回收因子 CRF 折算成每年成本否则建设期投资和十年运行成本放在一起会严重失真。CRF r * (1 r)^n / ((1 r)^n - 1)r 是折现率一般取 6%8%n 是设备寿命光伏、风机按 2025 年储能按 10 年变压器按 30 年。这样一台 1000 kW、造价 4000 元/kW 的光伏在 r6%、n20 年时CRF 约 0.087年投资成本就是 4000×1000×0.087 ≈ 34.8 万元。如果直接把一次性投资 400 万丢进目标函数算法会被建设成本主导完全不考虑后期运行收益这是新手最常犯的错误。运行维护费一般按初始投资的固定比例估算光伏 1%2%风机 2%3%储能 2% 左右。网损费用等于年网损电量乘以综合电价。和主网交互的电费则要考虑上网与购电价差分布式电源多发时可以减少购电。表里给一组我在课题里常用的参数参考值参数典型取值说明折现率 r6%8%反映资金时间价值企业内部收益率不同可调光伏单位造价30005000 元/kW含组件、逆变器、安装近年已明显下降风电单位造价60008000 元/kW受机位和地形影响明显储能系统造价8001200 元/kWh磷酸铁锂电芯加 PCS 加安装综合购电价0.50.8 元/kWh按全年峰谷加权平均光伏年利用小时8001400 h地区差异大南方低、西部高需要强调的是参数不追求绝对精确关键是口径一致。规划模型里所有设备成本都用等年值可靠性评估里的停电损失也折成年值才能放到同一个目标函数中比较。2.2 可靠性目标为什么选 EENS而不是只盯 SAIFI/SAIDI可靠性指标的基准在 1.1 里提到要有最优权衡空间工程中最常用的是 SAIFI、SAIDI、CAIDI 和设备可用率。SAIFI 衡量平均每年停电多少次SAIDI 衡量平均每年累计停多少小时CAIDI 是 SAIDI 除以 SAIFI表示每次停电平均持续多久。这三个指标反映的是“频率”和“时长”但都没有回答最关键的问题到底因为这几次停电损失了多少电量。在双目标规划里第二个目标我倾向于用 EENS期望缺供电量单位万 kWh/年。它等于每个负荷点停电概率与该点平均负荷的乘积再求和物理意义直观也能直接与电价挂钩换算成停电损失费用。用 EENS 做第二目标还有一个好处它能反映 DG 和储能在停电期间的“兜底”作用。SAIDI 只管用户停电累计时长不管负荷是 1 kW 还是 1000 kW而 EENS 天然加权了负荷重要性重负荷用户在可靠性模型里自动获得更高权重。可靠性评估本身有两种路线。解析法基于故障模式后果分析或者最小割集思路是把每个元件故障对负荷点的影响算出来再按元件可用率和修复时间加权求和。优点是计算快、结果稳定缺点是一旦要考虑时序性负荷、光伏出力随机变化、储能充放电策略模型复杂度会迅速爆炸。蒙特卡洛模拟法正好相反它把整个年份的元件状态、负荷曲线、DG 出力一条一条抽出来仿真复杂度与元件数近似线性建模灵活。混合配电系统里不确定性来源多我最终选的是序贯蒙特卡洛为主、典型日聚类加速后面第 3 章会细讲。2.3 DG 出力与负荷的不确定性建模光伏出力受辐照度影响标准做法是把某时段辐照度用 Beta 分布描述再乘上光伏方阵面积和逆变器效率得到出力风机出力则先对风速用 Weibull 分布建模再通过风机功率曲线映射到电功率。看起来数学模型都不复杂但实际规划中最大的坑是“时序相关性”光伏晚上出力为 0风机后半夜容易满发负荷晚高峰和光伏出力高峰错位。如果图省事只做非序贯抽样不考虑某时段内“光伏出力低而负荷高”这种相关性新能源的电量替代效应会被系统性高估可靠性结果也会偏乐观。我通常的做法是对历史数据做 K-means 聚类抽出每个季节的典型日辐照、典型日风速和典型日负荷曲线组成 12 个“典型日场景”作为全年代表。规划迭代时先用这 12 个场景快速算目标优化出最终 Pareto 解集后再用完整 8760 小时数据做一次精确校验。这样既保证了模拟法对时序性的要求又把计算量压缩了一个数量级。储能模型相对标准化一些规划变量是额定功率与容量运行模型里要满足 SOC 上下限、充放电功率上限以及充放电效率。储能的成本要单独计寿命通常按 10 年或 6000 次循环折算千万别和光伏风机一样用 25 年否则成本会严重虚低。3. Python 实现流程从潮流计算到多目标优化的代码骨架3.1 整体程序结构我的代码目录一般长这样data 文件夹放网架参数 CSV 和负荷曲线src 文件夹里分几个模块core 模块负责读数据和定义网络类powerflow 模块放前推回代潮流reliability 模块放蒙特卡洛可靠性评估opt 模块放 pymoo 的 Problem 子类和求解主程序。data/ ieee33_bus.csv ieee33_branch.csv load_profile.csv src/ core.py powerflow.py reliability.py opt.py主程序的大循环是读入 IEEE 33 节点或自定义配电网网架参数建立 networkx 图对象。定义候选 DG 安装节点和容量上限。用 pymoo 初始化种群比如 40 个个体。对每个个体调用潮流函数算网损、电压再调用可靠性函数算 EENS。两个目标值返回给优化器进行非支配排序和进化。迭代若干代后输出 Pareto 前沿并对最终解做 8760 小时全时段精确验证。这个结构里最关键的设计原则是目标函数计算函数一定要是纯函数输入是容量向量输出是两个目标值和可能的约束违反量内部不要碰全局随机种子和外部状态。这样后面加并行计算、加存档都能平稳过渡。3.2 用前推回代法做配电网潮流配电网潮流我强烈建议用前推回代法而不是常规牛拉法。配电网是辐射状结构R/X 比值高牛拉法初始化和迭代反而容易出问题前推回代按“从负荷端往根节点回推支路电流、从根节点往末端前推电压”的方式交替迭代简单又稳。以 IEEE 33 节点为例负荷数据和支路阻抗表都是现成的格式一般是bus 表节点编号、有功负荷、无功负荷、电压初始值。branch 表首端节点、末端节点、电阻、电抗、是否联络开关。核心迭代公式不复杂。设节点 i 的注入功率为 S_i P_i jQ_i负荷为正、DG 注入为负支路 b 连接父节点 p 和子节点 c。回代时从末端向根节点逐条计算支路电流 I_b conj(S_c / V_c)前推时从根节点向末端更新电压 V_c V_p - Z_b × I_b。两个方向交替直到前后两次迭代的电压差最大值小于阈值比如 1e-5 pu。DG 接入时只需要在对应节点注入功率处减掉出力。光伏、风机在规划模型里通常按 PQ 节点处理有功按出力曲线取无功要么按固定功率因数要么给电压支撑策略。我给一个简化的函数骨架注意这里假设 branch 列表已经按根到叶的层次顺序排好import numpy as np def backward_forward_sweep(n_bus, branch, S_load): V np.ones(n_bus, dtypecomplex) # 初始平启动 max_iter, tol 50, 1e-5 for _ in range(max_iter): V_old V.copy() # 回代从末端向根节点求支路电流 I_branch np.zeros(len(branch), dtypecomplex) for k in range(len(branch) - 1, -1, -1): child branch[k][1] I_branch[k] np.conj(S_load[child] / V[child]) # 前推从根节点向末端更新电压 for k in range(len(branch)): parent branch[k][0] child branch[k][1] V[child] V[parent] - (branch[k][2] 1j * branch[k][3]) * I_branch[k] if np.max(np.abs(V - V_old)) tol: break return V这个函数对单电源辐射网可直接用DG 出力已经折算进 S_load 里。需要处理联络开关时可以先做拓扑生成树把联络开关打开再做前推回代算完再评估闭环后的负载均衡。3.3 把可靠性评估封装成目标函数可靠性评估函数是整套代码里最费时间的一块。我的实现思路是对每个规划方案先对全年小时序列做蒙特卡洛抽样每一小时根据元件可用率抽样出“哪些支路故障、哪些变压器故障”然后调用一个恢复分析函数判断故障后哪些负荷能通过联络开关转移、哪些能靠 DG 孤岛恢复剩下的就是切负荷量。把所有小时的切负荷电量累加得到年 EENS。关键点是恢复分析不能做得太理想化也不能太粗糙。太理想化会得出“任何故障都能秒速恢复”的错误结论太粗糙则比如把所有下游负荷全部切掉DG 的可靠性价值就完全体现不出来。我一般按三段式判断第一主网侧能否通过闭合联络开关带起故障区域的负荷带得动就全部转移第二如果带不动看区域内 DG 总出力能否支撑孤岛能支撑就按可发功率与负荷的比值分摊恢复第三仍无法恢复的负荷按切负荷计入 EENS时长等于故障修复时间。核心骨架给在下面def reliability_evaluation(dg_capacity, net, load_curve, rng, repair_hour4): ens 0.0 for h in range(len(load_curve)): failed_branches sample_faults(net, rng) # 抽样故障支路 if not failed_branches: continue demand load_curve[h] * net.load_base dg_out dg_output_curve(h, dg_capacity) # 依据时序出力 shed restoration_analysis(net, failed_branches, demand, dg_out) ens shed * repair_hour return ens / len(load_curve) * 8760实际工程里 repair_hour 只是简化写法更精确的做法是按元件平均修复时间 MTTR 抽样每次修复时长。但这个函数骨架已经足够说明评估闭环的逻辑。可靠性评估对同一个方案跑多次结果会有波动所以优化迭代里一定要固定随机种子或者用公共随机数否则两个基因几乎一样的个体可能因为随机波动被排到不同层级NSGA-II 会非常不稳定。3.4 用 NSGA-II 跑双目标优化多目标优化我直接用了 pymoo 库它内置 NSGA-II比自己动手写快速非支配排序和拥挤距离省很多事。安装就是pip install pymoo。规划变量定义成候选节点处的安装容量长度为候选节点数量每个变量范围为 0 到单节点容量上限0 表示不装。这样变量是连续实数算法友好。自定义 Problem 类很直观import numpy as np from pymoo.core.problem import Problem class DGPlanningProblem(Problem): def __init__(self, net_data, n_candidate_nodes, max_cap1000): self.net net_data self.max_cap max_cap super().__init__(n_varn_candidate_nodes, n_obj2, xl0.0, xu1.0) # 归一化到0~1 def _evaluate(self, X, out, *args, **kwargs): f1, f2 [], [] for x in X: cap x * self.max_cap # 还原真实容量 f1.append(annual_cost(cap, self.net)) f2.append(reliability_evaluation(cap, self.net)) out[F] np.column_stack([f1, f2])求解主线同样简洁种群 40、进化 100 代作为起点。实际跑完后画 Pareto 前沿横轴年综合费用、纵轴年 EENS你会看到一个典型的 L 形曲线右下角对应“花大钱买高可靠”左上角对应“省钱但每年都要承受较多停电电量损失”。中间膝盖点附近通常就是工程上最推荐的方案。需要注意pymoo 的版本更迭偶尔改 API跑之前看一下本地版本的 docstring不要照搬老博客代码。4. 踩坑实录常见问题、排查思路与提速技巧4.1 蒙特卡洛仿真收敛慢、目标值抖动这个坑几乎每个人都会踩。直接把 8760 小时塞进规划优化40 个种群跑 100 代就是 4000 次方案评估每次评估都要遍历数千小时、做网络连通性分析和恢复判断算一次可能要几十秒整个优化跑几天都未必收敛。更麻烦的是不同方案的 EENS 都带随机噪声噪声一大NSGA-II 的非支配排序就会把本应淘汰的个体误判成 Pareto 优前沿会抖得不像话。我的解法是把评估拆成两档优化迭代阶段用 12 个典型日场景替代全年 8760 小时每个场景代表一类天气和负荷组合得到最终 Pareto 前沿后再对前沿上的方案跑完整全年仿真精确验证。两档之间的计算时间通常能差出 20 倍以上。为了进一步降低方差所有个体共用一套随机数也就是固定 rng 种子按顺序使用这叫公共随机数法试验设计里非常常见实现成本几乎为零。4.2 电压越限却找不到切负荷点我调试时遇到过一种诡异现象潮流算出来好多节点电压低于 0.9 pu但可靠性评估函数报告“没有切负荷”。原因是我在恢复分析里默认故障后的网络还是连通网而实际上某些支路故障后下游区域与主网的连接完全断开了没有给它安排 DG 或联络开关恢复前推回代在非连通网里会算出毫无意义的电压值切负荷判断自然失联。之后我在可靠性函数开头就先用 networkx 求连通分量先判断故障后的图里哪些负荷节点与主网节点仍连通哪些不连通对不连通的孤岛再判断里面 DG 总出力能否覆盖岛内负荷能覆盖才允许继续带负荷。这一改电压越限和切负荷的逻辑才对得上。代码里nx.connected_components或nx.is_connected都能直接干这个事但注意要先把光伏、风机在故障状态下的出力也按该时刻的实际出力算不能拿额定容量拍脑袋。4.3 规划结果里 DG 容量过度集中跑出来的 Pareto 前沿有时候会特别“极端”某个节点容量顶到上限周围节点零安装。原因是目标函数里网损费用占了主导算法倾向于把 DG 堆在电气距离最有利的一个点希望用最小投资降低全网络损。这在单纯网损优化里是“正确”的但在实际工程里完全不可行——单点接入容量过大会造成局部电压抬升、短路电流超标、逆变器容量不够等一系列问题。解决思路是在模型里加上单节点容量上限和全网渗透率上限前者按馈线载流量来定后者按变压器容量来定如果想让结果更均匀还可以在目标函数里加一项节点电压偏差的惩罚。加上这些约束后 Pareto 前沿会明显更“温和”不会再出现单点堆砌的疯方案。4.4 Pareto 前沿两端稀疏、中段跳跃优化跑完后前沿不是一条光滑曲线中间还有大段空缺这种情况我见得不少。原因主要有三个一是目标函数噪声导致的排序混乱二是种群多样性不够早熟收敛到局部区域三是约束处理不当让可行域被割裂。处理方法我按优先级排序先固定随机数把噪声压下来再把种群规模提高一倍进化代数从 100 代加到 200 代同时把模拟二进制交叉的分布指数调得稍大一些让子代更容易跳出局部。如果还不行就在支配关系上做文章。pymoo 里可以用 epsilon 支配设置一个适当的 epsilon 值允许目标值差异小于 epsilon 的个体互不支配这样会在前沿上保留更多中间个体曲线的连续性会好很多。这个技巧对工程问题特别管用因为在物理世界里成本和可靠性本来就做不到绝对意义上的“连续”。4.5 提速技巧能向量化就别写 for 循环可靠性评估里最耗时的部分是对全年小时做恢复分析。如果每个小时都走一遍 python 循环加 networkx 图遍历速度会非常难看。我在工程上把全过程优化过的经验是先用np.random.Generator一次性生成全年元件状态矩阵再用 numpy 布尔索引找出所有故障小时只对故障小时调用恢复分析函数无故障小时直接跳过。这样 8760 小时里可能只有几十个小时需要真算速度至少提升一个量级。如果还需要更快恢复分析函数本身就是热点可以用 numba 的njit装饰器把它编译成机器码或者用多进程并行评估种群。pymoo 通过设置 evaluate 时的并行 runner把种群内个体分给多个进程同时算核多的时候 4 倍提速是能轻松拿到的。备一个坑要小心并行评估时随机种子管理必须严格否则每个子进程会拿到不同随机流结果比串行时更抖。最后说一点我自己的实际体会。这套规划与可靠性评估闭环最早我用 Matlab 写后来为了让学生复现、方便接各种格式的电网数据才逐步迁移到 Python。迁移过程中最头疼的不是优化器也不是潮流算法而是把可靠性评估从“能跑”改到“算得快、结果稳”前前后后花了快两周时间。如果你现在也在做类似的规划课题我建议先别急着上 NSGA-II先把网架数据清洗干净、把前推回代潮流跑稳、把恢复分析的逻辑图手画明白再去碰多目标优化。模型堆得再高评估逻辑错了Pareto 前沿再漂亮也是空中楼阁。上面这些代码骨架和参数参考可以直接照抄进你的工程替换成你自己的网架数据就能用。