简介：面向配电网中分布式发电机（DG）接入场景的Matlab脚本，聚焦前推回代潮流计算与网损分析，适合电力系统专业学生、科研人员及配电网规划人员用于掌握DG并网对电压分布和网络损耗的影响规律。压缩包内仅1个.m文件，压缩后约3KB，代码量不大但完整呈现潮流计算核心流程，便于阅读和二次开发。目前已有181人学习下载，说明该主题在配电网研究与应用中具有一定关注度。通过运行脚本并调整DG接入位置与容量，可观察不同配置下的节点电压和网损变化，理解分布式电源就地消纳与馈线电流变化对网损的影响机制，为DG并网选址定容和配电网降损优化提供数据参考。前推回代法从电源节点逐级推算至负荷端，再反向迭代修正，兼具计算效率与内存占用优势，尤其适用于辐射状配电网结构；尽管包体袖珍，但基础潮流方程到迭代求解的关键环节均有体现，适合作为课程设计、毕业设计或工程问题排查的入门示例。
1. 配电网网损算不准，分布式发电机的并网评估就没法谈
前阵子处理一个台区光伏接入评估，业主拿着两份报告来问：同一台变压器，A 单位算出来装 800kW 分布式发电机后线损能降 12%，B 单位算出来反而升 3%。两边都说自己做了潮流计算，差别只在于 B 单位把 33 节点的配电网直接塞进了输电网用的牛顿法程序，结果不收敛。这就是配电网网损计算的典型困境：方法选错，结论完全反过来。DGdS.zip 这类工程代码包，解决的正是“分布式发电机建模、潮流计算求解、配电网网损统计”这条链路怎么在同一套模型里跑通。适合做配网规划、光伏并网评估、台区线损分析的工程师。下面会从方法选型、最小可复现代码、网损对比、收敛优化和常见坑依次讲透。
2. 配电网潮流计算选型：前推回代法为什么比牛顿法更配分布式发电机
2.1 配电网的“病态”结构让牛顿-拉夫逊法经常翻车
输电网和配电网的潮流计算，表面上是同一个问题，实际条件差得很远。输电网强调环网、长距离、高电压，线路电抗 X 远大于电阻 R，节点间电气距离大，牛顿-拉夫逊法在这种条件下收敛性好。配电网是典型的辐射状结构，中低压线路 R/X 比值经常到 2 到 3 甚至更高，再加上单相负荷、单相光伏、多分段多联络，雅可比矩阵的条件数很差。用输电网那一套直接来算，最典型的现象就是迭代发散，或者电压算出来在 0.3 p.u. 到 1.8 p.u. 之间来回跳。
牛顿法并不是不能在配电网里用，只是需要做很多修正：改坐标、加阻尼因子、做潮流解初值估计。但 DGdS.zip 这种专门处理配电网的代码包，主流的核心求解器往往不会首选牛顿法，而是选前推回代法。前推回代法的思路特别贴合辐射状配网：从末端开始往回代计算支路电流，再从根节点向前推算节点电压，两个方向交替迭代。它不需要组装大的雅可比矩阵，每轮迭代只做一次前推和一次回代，内存占用小，对 R/X 高的问题不敏感，速度也够快。
2.2 前推回代法的核心假设与分布式发电机的节点类型
前推回代法能成立，依赖于一个前提：网络拓扑是树状或者弱环网。也就是说每条支路只有一个父节点、一个或多个子节点，潮流方向从根节点流向末端。配电网正常运行时基本满足这个条件，和分布式发电机接入不冲突。分布式发电机在潮流计算里一般按节点类型分两类看待：恒功率因数控制的逆变器，通常建模成 PQ 节点，给定有功出力和无功出力；电压支撑型的机组，比如同步发电机或部分储能变流器，建模成 PV 节点，给定有功出力和节点电压幅值。
PV 节点是配电网潮流计算里的麻烦来源。输电网里 PV 节点无功是没有上限的，可以随便发，配电网里一个 500kW 的分布式发电机，逆变器无功能力可能就限制在 ±200kvar。若把它的电压目标定得太高，迭代时无功会一路越限，结果网损算出来完全失真。我在用 DGdS.zip 跑案例前，第一件事就是把所有分布式发电机按“PQ 节点 + 无功限幅”处理，只有在接入点电压支撑要求严格时，才开 PV 模式并做无功越限回退。
另外要注意，前推回代法在遇到分布式发电机接入时会改变局部功率流向。比如某个节点负载轻、光伏容量大，支路电流可能反灌，这个时候回代阶段单纯按净负荷累加电流就会出错。工程上处理方式是：回代阶段仍然把“负荷减去分布式电源出力”作为净注入功率，用复数运算自然体现反向电流，而不是提前判断方向。这样不仅代码简单，反向潮流造成的网损变化也能被正确算出来。
2.3 最小可跑代码：一版前推回代的迭代骨架
下面这段代码是我复现 DGdS.zip 核心流程时常驻在工程笔记里的骨架，去掉了读文件部分，保留潮流主循环。它要求节点编号满足父节点编号小于子节点编号，DGdS.zip 自带的 33 节点算例都是这个顺序。
import numpy as np def backward_forward_sweep(bus, branch, dg, max_iter=50, tol=1e-6): """ 前推回代法主循环。 bus: [[bus_id, p_load, q_load], ...] 单位 p.u. branch: [[from, to, r, x], ...] 单位 p.u. dg: {bus_id: (p_dg, q_dg)} 单位 p.u. """ n = len(bus) idx = {int(b[0]): i for i, b in enumerate(bus)} V = np.ones(n, dtype=complex) # 平启动，所有节点电压 1.0∠0 S_load = np.array([b[1] + 1j * b[2] for b in bus]) S_dg = np.zeros(n, dtype=complex) for k, (p, q) in dg.items(): S_dg[idx[k]] = p + 1j * q Z = {} parent = {} children = {i: [] for i in range(n)} for f, t, r, x in branch: fi, ti = idx[f], idx[t] Z[(fi, ti)] = complex(r, x) parent[ti] = fi children[fi].append(ti) # 按拓扑生成回代顺序：子节点先算 order = [] stack = [0] while stack: u = stack.pop() order.append(u) stack.extend(children[u]) backward_order = order[::-1] for it in range(max_iter): I = np.zeros(n, dtype=complex) # 回代：从末端往根节点方向算支路电流 for node in backward_order: if node == 0: continue S_net = S_load[node] - S_dg[node] + V[node] * np.conj(I[node]) I_node = np.conj(S_net / V[node]) I[node] = I_node I[parent[node]] += I_node V_new = V.copy() # 前推：从根节点往末端更新电压 for node in range(1, n): p = parent[node] V_new[node] = V_new[p] - Z[(p, node)] * I[node] diff = np.max(np.abs(V_new - V)) V = V_new if diff < tol: return V, I, it + 1 return V, I, max_iter逻辑说明：回代阶段从叶子节点算起，每个节点先叠加子支路电流，再算出流向父节点的电流；前推阶段用父节点电压减支路压降得到子节点电压；两次迭代电压差小于容差则收敛。分布式发电机以负的净功率形式出现在回代里，因此当光伏出力大于负荷时，I 的相位会自然翻转，反向潮流就体现在这里。
参数说明：max_iter 设 50，实际 33 节点算例一般 5 到 8 次迭代就收敛；tol 设 1e-6，对应电压精度约 0.0001 p.u.，算配电网网损足够。如果接入的分布式发电机容量接近所在馈线短路容量，需要把 tol 收紧到 1e-8，防止网损在小数点后第三位还在抖。
3. 在 IEEE 33 节点上跑通 DGdS.zip：分布式发电机接入前后的网损对比
3.1 先准备算例：为什么选 IEEE 33 作为基准
做配电网潮流计算验证，IEEE 33 节点系统几乎是默认的起点。它是一条 12.66kV 的辐射状馈线，总负荷约 3.7MW 加 2.3Mvar，有 32 条支路，拓扑结构简单到适合手算校核，又包含了足够多的分段负荷，能让分布式发电机接入后的网损变化体现得清楚。相比 IEEE 123 节点系统的三相不平衡和复杂负荷模型，33 节点更适合先验证代码正确性。
用 DGdS.zip 做这类算例，常见做法是把母线表和支路表整理成两个 CSV，或者直接复用代码包里现成的 case33 数据结构。母线表至少要有节点编号、有功负荷、无功负荷；支路表要有首端节点、末端节点、电阻、电抗。基准容量一般取 10MVA，基准电压取 12.66kV，这样所有数据都落在 0.0001 到 0.01 这个量级，迭代过程比较稳定。
以下是构造七个节点的示意数据，用来快速验证代码流程。IEEE 33 节点就是把 bus 和 branch 换成完整表，逻辑完全一致：
import numpy as np # 7节点链式馈线：节点0为变电站出口，负荷单位kW/kvar，换算成p.u.后喂给潮流函数 bus = [ [0, 0, 0], [1, 100, 60], [2, 120, 80], [3, 90, 40], [4, 60, 30], [5, 40, 20], [6, 30, 15], ] # 支路单位欧姆 branch = [ [0, 1, 0.12, 0.05], [1, 2, 0.10, 0.04], [2, 3, 0.08, 0.03], [3, 4, 0.10, 0.04], [4, 5, 0.12, 0.05], [5, 6, 0.11, 0.04], ] # 分布式发电机：节点4接入一台100kW光伏，无功0 dg = {4: (100 / 10000.0, 0.0)} V, I, iters = backward_forward_sweep(bus, branch, dg) print("迭代次数:", iters)逻辑说明：这段代码把实际工程数据换算成 p.u. 后交给第 2 章的潮流函数，节点 4 接入的 100kW 光伏以 0.01 p.u. 参与回代计算。运行后如果迭代次数在 10 次以内，说明拓扑和单位没问题。
3.2 分布式发电机接入参数怎么设：有功容量、功率因数、接入节点
接入参数直接决定配电网网损结论。第一是有功容量，不能只看光伏板标称容量，要看逆变器有功限值；第二是功率因数，恒功率因数模式下需要给定无功；第三是接入节点，负荷重还是轻，位置在馈线前端还是末端，实际效果差异很大。
功率因数这里我一般设 0.95 到 1.0 之间，并保留无功上下限。DGdS.zip 这类工具通常会提供“恒功率因数”和“恒电压”两种模式。做网损对比评估时，我习惯全部用恒功率因数模式，因为并网逆变器最常见的控制策略就是这样，结果也更保守。如果某个节点接入的是同步发电机，再单独开恒电压模式。
# 把光伏出力按功率因数换算成无功 p_dg_kw = 500.0 pf = 0.95 q_dg_kvar = p_dg_kw * np.tan(np.arccos(pf)) # 或者直接用给定的无功范围 dg_config = { 18: (500 / 10000.0, 0.0), # 节点18：500kW光伏，功率因数1.0 33: (300 / 10000.0, q_dg_kvar / 10000.0), # 节点33：300kW光伏，功率因数0.95 }参数说明：q_dg_kvar 的计算用了功率三角形，如果功率因数设为 0.95，无功约为有功的 0.33 倍。实际项目中逆变器的无功能力一般在 ±0.4 倍容量以内，超过这个范围的设定值会被限幅，需要在下发到 DGdS.zip 前先做一次截断。
3.3 跑潮流并统计配电网网损：核心代码与输出解释
潮流跑完拿到各支路电流和节点电压，网损计算就变得直接。三相配电网的支路有功损耗公式是 I²R，用 p.u. 值计算后乘以基准容量就能换算回 kW。下面的代码把结果统一整理成表格，并计算网损率：
def calc_loss(branch, I, base_mva=10.0): """ 计算全网有功损耗，单位kW。 branch: [[from, to, r, x], ...] I: 各支路电流复数数组，下标对应branch顺序 """ total_loss = 0.0 for k, (f, t, r, x) in enumerate(branch): loss_pu = np.abs(I[k]) ** 2 * r total_loss += loss_pu * base_mva * 1000.0 return total_loss # 接入分布式发电机前 V0, I0, _ = backward_forward_sweep(bus, branch, {}) loss0 = calc_loss(branch, I0) # 接入分布式发电机后 V1, I1, _ = backward_forward_sweep(bus, branch, dg_config) loss1 = calc_loss(branch, I1) print(f"基态网损: {loss0:.2f} kW") print(f"接入后网损: {loss1:.2f} kW") print(f"变化量: {loss1 - loss0:.2f} kW")逻辑说明：calc_loss 遍历每一条支路，用该支路电流幅值的平方乘以电阻得到损耗。注意这里的 I 并不是“节点注入电流向量”，而是每条支路首端到末端的电流，和前推回代里 I 的更新逻辑一致。如果直接用节点电流模值套 I²R，结果必然偏大。
参数说明：base_mva 必须与前期换算 p.u. 时的基准容量一致。常见翻车点是在 10MVA 基准下算出结果，输出时当成 MVA 直接写 kW，导致损耗数值相差 1000 倍。IEEE 33 节点基态网损用 DGdS.zip 算出来一般在 200kW 上下，如果超出这个范围很多，优先检查电阻单位是不是用了毫欧。
4. 网损升高还是降低？最优因子法与收敛性排查
4.1 分布式发电机接入网损的两种相反效应
分布式发电机接入后网损并不是一定下降。道理不复杂：负荷节点附近接入电源，馈线首段电流减小，网损降低；但光伏大发时多余功率从末端倒送，靠近变电站的支路电流在低谷时段反而大于无光伏时的电流，网损被抬高。两种效应同时存在，谁占主导取决于光伏容量和负荷曲线的匹配度。
这就是为什么接入评估不能只看容量和年发电量，必须跑典型日曲线。工程上常用的做法是取春、夏、秋、冬四个典型日，每个典型日按 24 小时逐时潮流计算，再把网损累加。DGdS.zip 这类工具很多时候是单时段的快照分析，做规划够用，做运行评估就要自己在外面包一层循环。算完 96 次潮流后，网损率通常能区分出分布式发电机是“助减损”还是“助增损”。
4.2 最优因子法：当潮流收敛严重振荡时的“后悔药”
前推回代法在大部分配电网算例里几轮就收敛，但如果分布式发电机无功设定得离谱，或者接入点电压支撑特别弱，迭代曲线会出现等幅振荡：前一次电压偏高，下一次偏低，总有某个节点跳不出收敛判据。这种时候我一般不会换求解器，而是给迭代加一个阻尼因子，也就是把电压修正量乘上一个小于 1 的系数。所谓潮流计算最优因子法，就是把这个系数从拍脑袋改成自动寻优。
最优因子的思路很简单：每次前推回代得到新的电压 V_new 后，不直接用它替换 V，而是在 V 和 V_new 的连线上找一点，使该点的最大节点电压误差最小。一维搜索用黄金分割就能解决，下面的代码就是阻尼牛顿法的实现骨架，换到前推回代上也一样适用：
def golden_search(f, a, b, tol=1e-4): """黄金分割法找一维最小值。f 是步长到误差范数的映射。""" phi = (np.sqrt(5) - 1) / 2 c = b - phi * (b - a) d = a + phi * (b - a) while abs(c - d) > tol: if f(c) < f(d): b = d else: a = c c = b - phi * (b - a) d = a + phi * (b - a) return (a + b) / 2 def damped_forward_backward(bus, branch, dg, max_iter=50): """ 带最优阻尼因子的前推回代。 每次迭代后，用黄金分割求最优步长 alpha，再更新电压。 """ alpha = 1.0 # ... 初始化 V 与拓扑，同第2章 ... for it in range(max_iter): # 常规回代 + 前推，得到 V_new V_new, I, _ = backward_forward_sweep_once(bus, branch, dg, V) # 把 V_new 与 V 的差作为搜索方向 dx dx = V_new - V # 定义步长 alpha 到误差范数的函数 def err_func(alpha): V_trial = V + alpha * dx return np.max(np.abs(V_trial - forward_sweep_result(bus, branch, V_trial))) # 在 [0.2, 1.5] 之间寻找最优因子，范围下界不为0是为了避免原地踏步 alpha = golden_search(err_func, 0.2, 1.5) V = V + alpha * dx if np.max(np.abs(dx)) < 1e-6: break return V, I, it + 1逻辑说明：damped_forward_backward 把一次前推回代的结果作为搜索方向，再用黄金分割找到使误差最小的步长 alpha。alpha 大于 1 表示当前方向值得加速，小于 1 表示需要压步子防振荡。对配电网高 R/X 比场景，alpha 落在 0.5 到 0.9 之间是常态，接近 1.5 的情况多出现在重负荷率很低的馈线。
参数说明：golden_search 容差 1e-4 足够，因为 alpha 本身对网损结果不敏感；搜索范围 [0.2, 1.5] 是经验值。如果反复出现 alpha 贴着 0.2 的情况，说明模型本身有问题，优先检查无功限幅和支路参数，而不是继续调因子。
4.3 最优因子的经验参数与调参建议
| 场景特征 | 建议初始因子 | 迭代上限 | 说明 |
|---|---|---|---|
| 轻负荷、光伏容量占比低 | 1.0 | 20 | 常规前推回代即可 |
| 光伏容量接近馈线负荷 | 0.8 | 30 | 收敛变慢，适度阻尼 |
| 末端弱电压、长线路 | 0.6 | 50 | 阻尼太强会迭代过慢 |
| 出现等幅振荡 | 0.3 起自动寻优 | 50 | 先排除无功越限 |
阈值设定没有统一公式，核心经验是：先不加阻尼跑一次，如果发散，把初始因子设 0.8 再跑；还发散就检查模型参数，不要用阻尼因子硬扛。
5. 配电网网损计算避坑：五个让结果直接报废的常见问题
5.1 潮流迭代发散，电压一路跌到 0.6
现象：前推回代跑了十几轮，末端电压从 1.0 一路跌到 0.6，且没有回升趋势。原因多见于分布式发电机的无功功率设成了正值且远超逆变器能力，叠加到重负荷末端后电压持续下滑。我遇到最离谱的一次，是配置文件里把无功单位 kvar 写成了 Mvar，等于给节点硬塞了 10 倍无功。解决：先检查单位换算，再给分布式发电机配无功限幅表，PV 节点模式下无功越限后立刻回退为 PQ 节点重新迭代。
5.2 网损是负数？基准值单位不统一是头号原因
现象：接入分布式发电机后网损变化量为负，但算出来的绝对损耗也是负的。原因：母线有功负荷用 MW，分布式发电机出力用 kW，在 p.u. 换算时一个除了 100，另一个除了 1000，净功率出现负数。解决：所有有功无功统一换算到 p.u. 后再进入潮流函数，不要在函数内部混合单位。我的习惯是在读 CSV 时直接用一列“单位倍率”做转换，例如 kW 除以 base_kva 再除以 1，MW 除以 base_mva 乘以 1，逻辑一目了然。
5.3 装了分布式发电机网损反而升了
现象：同一算例，光伏接入后网损率从 2.8% 升到 4.1%。原因：分布式发电机装在了馈线末端，且光伏容量远超末端负荷，午后反送电流使首段支路重载。解决：用网损灵敏度判断接入位置，不要只看消纳能力。末端装小容量光伏降损效果明显，末端装大容量光伏很容易帮倒忙。
5.4 单相光伏当三相平衡算
现象：三相潮流里，基态网损正常，接入单相光伏后结果波动剧烈，且某一相电压偏高。原因：单相分布式发电机接入 A 相，却按三相平衡容量给每一相都注入功率，等于实际容量放大了 3 倍。解决：用三相潮流模型，不满足条件时把单相容量除以 3 近似为三相平衡电源，同时要在报告中注明这是近似。真要评估单相光伏密集接入，建议直接用支持三相不平衡前推回代的工具，不要做近似。
5.5 理论网损与关口表对不上
现象：计算出的网损是 180kW，关口电能表统计出来是 230kW。原因：模型里没有计及变压器铁损、线路对地导纳和老化电阻，也没有考虑实际电压不是额定电压。解决：把变压器铁损作为固定损耗叠加，线路电阻按实际温度修正，再对 24 个时段的网损做积分，而不是取一个时段乘 24。对比关口表时，至少留 10% 到 15% 的偏差余量。
6. 进阶：多点分布式发电机接入下的网损灵敏度排序
单点接入的网损结论可以靠一次潮流对比得到，但实际项目很少只装一台。多点接入时，每个候选节点的网损灵敏度差异很大。我现在的做法是：先对每个候选节点分别做一次分布式发电机出力变化测试，比如统一加 500kW，记录全网有功损耗变化量，再按灵敏度从低到高排序，优先选择网损下降幅度大的节点。公式上就是近似偏导：ΔP_loss / ΔP_DG。
下面这段差分法代码是我在每个项目里固定要跑一遍的：
def loss_sensitivity(bus, branch, dg_config, candidate_nodes, step_kw=500): """ 计算候选节点的网损灵敏度。 返回值：{节点: 网损变化量 kW} """ base_mva = 10.0 # 基态网损 _, I0, _ = backward_forward_sweep(bus, branch, {}) loss0 = calc_loss(branch, I0, base_mva) result = {} step_pu = step_kw / (base_mva * 1000.0) for node in candidate_nodes: dg_test = {node: (step_pu, 0.0)} _, It, _ = backward_forward_sweep(bus, branch, dg_test) losst = calc_loss(branch, It, base_mva) result[node] = loss0 - losst # 正值表示降损 return result # 候选节点：支路末段、中段、重负荷点各取一个 candidate_nodes = [18, 25, 33] sens = loss_sensitivity(bus, branch, dg_config, candidate_nodes) for n, s in sorted(sens.items(), key=lambda x: x[1], reverse=True): print(f"节点 {n}: 接入500kW后网损变化 {s:.2f} kW")参数说明：step_kw 取 500 是折中方案。取得太小，网损差量可能淹没在潮流收敛误差里；取得太大，灵敏度数值受非线性影响偏离真实值。差分结果只用于排序，不用于定量预测。最终确定方案后，再按典型日曲线做逐时潮流复核。
这种做法帮我在至少三个项目里避免了“选了最短电气距离节点结果却增损”的尴尬。验证灵敏度结论也简单，找两个灵敏度差异最大的节点，分别跑 24 时段的潮流积分对比，结果和排序基本一致。希望这个思路能帮你在配电网网损分析上少走弯路。
本文还有配套的精品资源，点击获取
