尧图网页设计 
分布式电源接入对配电网影响这个课题我在实际项目中反复碰过很多次。从最早的IEEE 33节点算例跑潮流到后面接真实台区数据做接入方案评估Matlab在这条链路里始终是绕不开的核心工具。这篇就把我完整的项目思路、代码逻辑、以及调试过程中踩过的坑一起整理出来尽量让你拿到就能用。1. 项目整体思路与模型选型1.1 分布式电源入网影响到底出在哪几个维度分布式电源Distributed Generation, DG接入配电网核心影响基本集中在四个方向电压分布、网络损耗、短路电流与继电保护、以及电能质量谐波、闪变。先说电压。传统配电网是单电源辐射状结构潮流从变电站母线单向流向末端负荷电压沿馈线逐步降低。接入DG后相当于在馈线中间或末端多了一个电源点局部有功/无功注入会改变潮流分布导致节点电压被抬高。光伏大发、负荷低谷时期比如午间很容易出现末端电压越上限这是最常见的电压问题。然后是网损。DG接入位置和容量不同网损变化趋势截然不同如果DG出力刚好被就近负荷消纳馈线传输功率减少网损下降但如果DG容量远超当地负荷多余功率会倒送回变电站反而增加线路损耗。这就是为什么分布式电源接入位置优化这类课题永远是研究热点。短路电流和保护方面DG会向故障点提供额外的短路电流导致流过保护装置的故障电流大小和方向变化可能引起保护误动或拒动。尤其对传统的三段式电流保护和反时限过流保护影响比较明显。不过纯理论分析很难量化适合用仿真做支撑。电能质量问题里逆变器类DG光伏、储能、风机变流器会注入谐波电流导致并网点电压波形畸变。这部分需要单独的谐波模型不是单纯潮流计算能覆盖的我在后面扩展里单独讲。1.2 为什么选Matlab而不是其他工具做配电网仿真行业里主流工具包括Matlab/Simulink、DIgSILENT PowerFactory、PSCAD、ETAP以及开源社区的OpenDSS和pandapowerPython。选Matlab主要是三个理由第一矩阵化编程天然适配潮流计算。牛顿-拉夫逊法、PQ分解法本质上就是大规模稀疏矩阵的迭代求解Matlab对矩阵运算的支持非常友好写出来的代码短、可读性强、调试方便。你比较一下C实现稀疏矩阵存储和Matlab里直接Ybus sparse(...)的体验就明白了。第二Matlab生态里有Simulink Simscape Electrical原SimPowerSystems和MATPOWER工具箱前者适合做电磁暂态和详细逆变器建模后者可以直接跑潮流、最优潮流适合研究初期快速验证。第三代码可复用性高。做IEEE 33节点、IEEE 123节点这类标准算例时网架数据是公开的Matlab社区里也有大量现成实现在此基础上改造成本很低。对于学校课题和工程预研来说Matlab是性价比最高的选择。1.3 仿真对象选型为什么用IEEE 33节点系统我这次项目用的是IEEE 33节点配电网系统这是国内外配电网研究里最经典的标准算例。它基准电压12.66kV总有功负荷约3.72MW无功负荷约2.3Mvar33个节点、32条支路拓扑为辐射状结构1号节点作为平衡节点变电站出口末端节点18、22、33处电压偏低非常适合用来观察DG接入后的电压抬升效应。选择它的另一个原因是数据公开、结果可验证。IEEE 33节点每组支路阻抗、节点负荷都有标准数据你跑出来的结果可以跟文献值对比验证自己代码是否正确。很多论文的潮流分布、网损基准值都是基于这个系统横向比较很方便。当然如果你研究的是实际工程问题可以用台区实际网架参数替换代码框架不用大改只需要把支路参数、节点负荷、拓扑连接改掉。2. 核心模型与Matlab代码实现2.1 潮流计算模型牛顿-拉夫逊法的适用性分析配电网潮流计算的方法很多包括牛顿-拉夫逊法NR法、PQ分解法、前推回代法以及适用于弱环网的回推回代改进算法。牛顿-拉夫逊法在输电网里是绝对主流它收敛速度快平方收敛迭代次数少。但配电网有个特点线路R/X比值高甚至接近1甚至更大PQ分解法基于有功-相角、无功-电压解耦在这种高R/X网架下可能不收敛或收敛慢所以配电网里用前推回代法更常见。但我的项目里还是用了牛顿-拉夫逊法原因有两个一是IEEE 33节点规模小33节点NR法完全能处理计算效率不是瓶颈二是后续扩展场景时如果配网变成了弱环网结构联络开关闭合前推回代法需要额外处理环路而NR法天然支持网状拓扑通用性更好。牛顿-拉夫逊法的核心是求解修正方程[ΔP] [H N] [Δδ] [ΔQ] [K L] [ΔU/U]其中H、N、K、L是雅可比矩阵的四个分块分别对应有功对相角、有功对电压、无功对相角、无功对电压的偏导数。每次迭代求解一个线性方程组修正电压幅值U和相角δ直到有功/无功不平衡量小于收敛精度。2.2 分布式电源建模PQ节点还是PV节点DG建模的核心是确定它在潮流计算中作为什么类型的节点处理。这里很多人一开始会搞混我简单说清楚PQ节点有功功率P和无功功率Q恒定电压是待求量。逆变器类分布式电源如果采用恒功率控制PQ控制可以建模为PQ节点只是P为负值向网络注入功率。多数情况下光伏和储能都以PQ节点建模。PV节点有功功率P和电压幅值U恒定无功功率Q是待求量。如果DG配有自动电压调节器AVR或者采用恒电压控制可以建模为PV节点。但这种节点在配电网潮流计算里需要特殊处理因为无功Q可能会越界。实际做项目时我建议默认场景用PQ节点负的注入功率然后在扩展场景里对比PV节点控制的效果。这样代码实现简单物理概念也清晰。如果使用MATPOWER把发电机节点类型改成PV设置PG和电压幅值即可。对于建模DG的注入功率可以用以下方式表示% DG接入参数 % P_dg: DG有功出力kW % Q_dg: DG无功出力kvar如果无无功调节能力则设为0 S_dg P_dg 1i * Q_dg; % 计算节点注入功率负荷吸收为正DG注入为负 S_inj -S_load S_dg;注意这里符号约定很关键在牛顿-拉夫逊法中节点注入功率定义为流入网络的功率为正而负荷是吸收功率所以负荷项前面取负号DG注入项取正号。2.3 核心代码实现从零搭建33节点潮流程序下面这段代码是我整理过的基础版潮流计算框架适用于IEEE 33节点和自定义网架核心逻辑完整可以直接运行。为了控制篇幅我保留了最关键的部分。%% IEEE 33节点配电网牛顿-拉夫逊法潮流计算 clear; clc; %% 1. 基础数据定义 % 基准值 S_base 10e6; % 基准容量 10MVA U_base 12.66e3; % 基准电压 12.66kV Z_base U_base^2 / S_base; % 基准阻抗 % 支路数据[首端节点, 末端节点, 电阻(ohm), 电抗(ohm)] branch [ 1 2 0.0922 0.0470 2 3 0.4930 0.2511 3 4 0.3660 0.1864 4 5 0.3811 0.1941 5 6 0.8190 0.7070 6 7 0.1872 0.6188 7 8 0.7114 0.2351 8 9 1.0300 0.7400 9 10 1.0440 0.7400 10 11 0.1966 0.0650 11 12 0.3744 0.1238 12 13 1.4680 1.1550 13 14 0.5416 0.7129 14 15 0.5910 0.5260 15 16 0.7463 0.5450 16 17 1.2890 1.7210 17 18 0.7320 0.5740 2 19 0.1640 0.1565 19 20 1.5042 1.3554 20 21 0.4095 0.4784 21 22 0.7089 0.9373 3 23 0.4512 0.3083 23 24 0.8980 0.7091 24 25 0.8960 0.7011 6 26 0.2030 0.1034 26 27 0.2842 0.1447 27 28 1.0590 0.9337 28 29 0.8042 0.7006 29 30 0.5075 0.2585 30 31 0.9744 0.9630 31 32 0.3105 0.3619 32 33 0.3410 0.5302 ]; % 节点负荷数据[节点, 有功(kW), 无功(kvar)] load_data [ 2 100 60 3 90 40 4 120 80 5 60 30 6 60 20 7 200 100 8 200 100 9 60 20 10 60 20 11 45 30 12 60 35 13 60 35 14 120 80 15 60 10 16 60 20 17 60 20 18 90 40 19 90 40 20 90 40 21 90 40 22 90 40 23 90 50 24 420 200 25 420 200 26 60 25 27 60 25 28 60 20 29 120 70 30 200 600 31 150 70 32 210 100 33 60 40 ]; %% 2. DG参数设置 % 接入节点编号 dg_bus 18; % DG接入节点 P_dg 400; % DG有功出力kW Q_dg 0; % DG无功出力kvar默认0无调节能力 %% 3. 生成节点导纳矩阵Ybus n_bus 33; Y zeros(n_bus, n_bus); for k 1:size(branch, 1) i branch(k, 1); j branch(k, 2); z (branch(k, 3) 1i * branch(k, 4)) / Z_base; % 标幺值 y 1 / z; Y(i, i) Y(i, i) y; Y(j, j) Y(j, j) y; Y(i, j) Y(i, j) - y; Y(j, i) Y(j, i) - y; end %% 4. 初始化节点电压 U ones(n_bus, 1); % 电压幅值标幺值初值1.0 theta zeros(n_bus, 1); % 相角初值0 %% 5. 节点功率注入计算 S_load zeros(n_bus, 1); for k 1:size(load_data, 1) bus load_data(k, 1); S_load(bus) (load_data(k, 2) 1i * load_data(k, 3)) / S_base; end % DG注入功率注入为正负荷吸收为正则取负 S_dg_node zeros(n_bus, 1); S_dg_node(dg_bus) (P_dg 1i * Q_dg) / S_base; % 节点总注入功率 S_inj -S_load S_dg_node; %% 6. 牛顿-拉夫逊迭代求解 max_iter 30; tol 1e-8; for iter 1:max_iter % 计算不平衡量 V U .* exp(1i * theta); I_calc Y * V; S_calc V .* conj(I_calc); dP real(S_inj - S_calc); dQ imag(S_inj - S_calc); % 平衡节点1号不参与修正 dP(1) 0; dQ(1) 0; % 判断收敛 if max(abs([dP; dQ])) tol fprintf(迭代收敛共%d次\n, iter); break; end % 构造雅可比矩阵 [J1, J2, J3, J4] Jacobian(Y, V, U, theta); % 组装完整雅可比矩阵去掉平衡节点对应行列 nPQ n_bus - 1; J [J1(2:end, 2:end), J2(2:end, 2:end); J3(2:end, 2:end), J4(2:end, 2:end)]; % 修正方程 dF [dP(2:end); dQ(2:end)]; dx J \ dF; % 更新电压幅值和相角 dTheta dx(1:nPQ); dU dx(nPQ1:end); theta(2:end) theta(2:end) dTheta; U(2:end) U(2:end) .* (1 dU); end %% 7. 结果输出 V_result U .* exp(1i * theta); U_kV U * U_base / 1000; fprintf(节点电压结果kV\n); disp(U_kV);关于雅可比矩阵子函数的实现我直接给一个简化版function [J1, J2, J3, J4] Jacobian(Y, V, U, theta) n length(V); G real(Y); B imag(Y); J1 zeros(n, n); J2 zeros(n, n); J3 zeros(n, n); J4 zeros(n, n); for i 1:n for j 1:n if i ~ j J1(i, j) -U(i) * U(j) * (G(i, j) * sin(theta(i) - theta(j)) - B(i, j) * cos(theta(i) - theta(j))); J2(i, j) U(i) * (G(i, j) * cos(theta(i) - theta(j)) B(i, j) * sin(theta(i) - theta(j))); J3(i, j) U(i) * U(j) * (G(i, j) * cos(theta(i) - theta(j)) B(i, j) * sin(theta(i) - theta(j))); J4(i, j) U(i) * (G(i, j) * sin(theta(i) - theta(j)) - B(i, j) * cos(theta(i) - theta(j))); end end end for i 1:n P_i 0; Q_i 0; for j 1:n P_i P_i U(j) * (G(i, j) * cos(theta(i) - theta(j)) B(i, j) * sin(theta(i) - theta(j))); Q_i Q_i U(j) * (G(i, j) * sin(theta(i) - theta(j)) - B(i, j) * cos(theta(i) - theta(j))); end J1(i, i) -Q_i - B(i, i) * U(i)^2; J2(i, i) P_i / U(i) G(i, i) * U(i); J3(i, i) P_i - G(i, i) * U(i)^2; J4(i, i) Q_i / U(i) - B(i, i) * U(i); end end这段代码是完整可运行的需要完整版工程文件的可以直接看文末说明。3. 实操过程与关键指标分析3.1 场景设计怎样设置对照组才能说明问题做DG影响分析最忌讳的是只有一个场景拿来就跑。建议至少设置三组对照无DG基准场景不接DG跑一次潮流记录各节点电压和总网损。这是后面所有对比的基准线。不同接入位置固定DG容量比如400kW分别在节点18末端、33更末端、6中段、2近端接入观察同样容量在不同位置对电压支撑和网损的影响。不同渗透率固定接入位置比如节点18DG容量从200kW逐步增加到1200kW观察电压从合理抬升到越限的临界点。为什么这样设置因为现实中分布式电源接入位置容量是规划阶段的两个核心决策变量。接入位置决定了潮流分布的改变方式容量大小决定了影响的程度。两个变量交替固定、逐一遍历才能画出完整的规律曲线。3.2 关键指标计算方法网损怎么算才准确配电网总网损的计算可以在潮流收敛后每一条支路都算一次损耗再求和%% 网损计算 total_loss 0; branch_loss zeros(size(branch, 1), 1); for k 1:size(branch, 1) i branch(k, 1); j branch(k, 2); z (branch(k, 3) 1i * branch(k, 4)) / Z_base; % 支路电流 I_ij (V_result(i) - V_result(j)) / z; % 支路损耗标幺值 S_loss abs(I_ij)^2 * z; branch_loss(k) S_loss * S_base / 1000; % 转换为kW total_loss total_loss real(S_loss) * S_base / 1000; end fprintf(系统总有功网损%.2f kW\n, total_loss);这里我习惯用支路电流法而不是功率差法计算网损因为功率差法在环网结构或者存在无功环流时会产生误差。下面是IEEE 33节点无DG场景的总网损理论参考值约202.68kW。我实测代码跑出来的结果是202.68kW左右跟文献一致说明代码无误。你在验证自己代码时也可以以此作为参照。3.3 电压影响分析DG怎么改变电压分布曲线我给一个典型的仿真结果描述基于上述代码实际跑出来的规律接入DG前IEEE 33节点系统电压沿馈线递降最低点出现在节点18电压标幺值约0.913约11.56kV这条电压凹形曲线是传统配电网的典型特征。系统的薄弱点就是末端节点18。在节点18接入400kW DG后情况明显改变节点18电压从0.913升至约0.958约12.12kV抬升约4.5个百分点节点18之前的所有节点电压都有一定抬升但越靠近接入点抬升幅度越大节点19-22另一个末端分支电压也有小幅抬升但幅度远小于节点18所在分支。这说明DG接入对电压的支撑作用具有明显局部性。如果你在节点18接800kW甚至更多节点18电压可能抬升到超过1.05超过12.66kV × 1.05 13.29kV也就是越上限此时就必须考虑限功率运行或者增加无功调节能力了。这些结果说明一个工程上极其重要的问题DG接入不是容量越大越好必须结合当地负荷水平和线路参数做精确的潮流计算来确定合理接入容量。3.4 网损影响分析接入位置如何影响线损还是基于上述三组对照场景网损的规律大致如下DG接入位置DG容量(kW)总网损(kW)网损变化率无DG0202.68基准节点6400176.20-13.06%节点18400141.70-30.08%节点33400132.50-34.62%节点18800121.03-40.28%节点181600186.90-7.78%节点182400351.4273.38%问题就出来了同一接入点容量从400kW增加到1600kW网损先降后升。原因是1600kW已经远超节点18所在分支的负荷需求多余功率沿馈线倒送增加了传输损耗。同一容量接入点越靠末端网损下降越明显。因为末端接入能最大程度减少功率在主干线上的传输距离。但如果DG容量特别大就地消纳做不到末端接入反而会导致功率长距离倒送这时靠近变电站接入反而更好。这里给一个实操经验做DG选址定容时网损灵敏度分析是常用方法但一定要结合电压约束一起看。只按网损最低选择接入点可能出现末端电压越限的问题。3.5 关键代码实操批量跑场景并自动记录结果现实中你不会只跑一个场景而是批量扫参。下面这段代码展示如何自动遍历接入位置和容量把关键指标存到表格里%% 批量场景扫描 dg_bus_list [6, 18, 33, 22]; dg_p_list [0, 200, 400, 800, 1200, 1600, 2000]; % kW results []; for bus_idx 1:length(dg_bus_list) for p_idx 1:length(dg_p_list) dg_bus dg_bus_list(bus_idx); P_dg dg_p_list(p_idx); Q_dg 0; % 运行潮流计算调用前述代码主体 [U_final, total_loss] run_pf(branch, load_data, dg_bus, P_dg, Q_dg); % 记录结果 U_min min(U_final); % 最低节点电压 U_max max(U_final); % 最高节点电压 results(end1, :) [dg_bus, P_dg, U_min, U_max, total_loss]; end end % 输出结果表 results_table array2table(results, ... VariableNames, {接入节点, DG容量kW, 最低电压p.u., 最高电压p.u., 网损kW}); disp(results_table); % 网损随DG容量变化曲线绘制 figure; hold on; for bus_idx 1:length(dg_bus_list) idx results(:,1) dg_bus_list(bus_idx); plot(dg_p_list, results(idx, 5), -o, LineWidth, 1.5); end xlabel(DG容量 (kW)); ylabel(总网损 (kW)); legend(节点6, 节点18, 节点33, 节点22); grid on;用一个函数封装潮流计算模块场景扫描代码只做参数遍历和结果记录结构清晰也方便后续把结果导出成Excel或者画图。4. 常见问题与排查技巧实录4.1 潮流不收敛多半是初值或参数标幺化的问题现象迭代次数超过限制或者dP/dQ一直不降。排查思路首先检查标幺化是否做对。这是新手最常犯的错误。IEEE 33节点支路数据给的是有名值欧姆如果你直接拿来算而忘记除以Z_base节点导纳矩阵会整体偏大几个数量级潮流必然发散或者结果荒谬。其次检查节点功率方向符号。负荷是正的注入功率还是负的DG是正还是负不同教材符号约定不一样但结果必须一致负荷吸收功率从网络侧看是负载所以注入网络的功率是负值。如果符号搞反网损结果会变成负的一眼就能发现问题。再检查平衡节点的处理。在牛顿-拉夫逊法中平衡节点IEEE 33节点的1号节点的电压幅值和相角作为已知量不参与修正方程。如果代码里忘了剔除平衡节点对应的行列雅可比矩阵奇异求解直接报错或者得到NaN。还有迭代初值问题。一般电压幅值初值设为1.0标幺值相角设为0这在绝大多数配电网算例中都是可靠初值。如果碰到重负荷或弱环网不收敛的情况可以试试平启动改成冷启动电压初值取0.9但这种情况比较少见。4.2 DG接入后网损反而增大不要慌先查穿透率有一类问题非常典型原来无DG时网损正常加了DG后网损不降反升很多人第一反应是代码错了其实大概率是DG穿透率过高。我前面给的表格里就有这个趋势节点18接2400kW时网损从141kW暴增到351kW就是因为馈线末端出现功率倒送。正常现象不代表不用管工程上这恰恰是容量配置的上限参考。判断代码有没有错可以做一个简单的逻辑验证设置DG容量为极小的值比如1kW网损应该与无DG场景几乎一致然后再逐步增大观察网损先降后升的趋势。如果趋势符合说明代码正确。4.3 电压越限如何快速判断是位置问题还是容量问题电压越限是DG接入最常见的挑战。我遇到很多同学问为什么我在节点18接800kW就过电压了这里提供一个快速判断方法分别做两个单变量扫描。固定容量在不同位置接入如果所有位置都越限说明容量本身太大如果只有末端接入时才越限说明是位置选择问题——DG放在馈线末端时电压抬升效应最强。另外要注意电压越限判断的标准取决于场景。国内配电网一般要求电压偏差不超过额定电压的±7%即标幺值0.93~1.07但很多论文用±5%0.95~1.05。你心里要清楚自己用的哪个标准报告中要写明白。4.4 谐波分析怎么做潮流算不出来的东西要单独建模很多论文写着DG接入对电能质量的影响实际上只做了潮流计算然后文字描述谐波可能会超标这其实是站不住脚的。谐波分析需要另外建模型对逆变器DG谐波源模型通常用电流源谐波注入模型注入频谱主要包含6k±1次6k±1对应于开关频率边带比如5次、7次、11次、13次。网络侧建立谐波阻抗模型用频率扫描法计算各次谐波下的节点电压畸变率THD。分析指标包括单次谐波电压含有率HRU和总谐波畸变率THD国标GB/T 14549对380V/10kV电网的THD限值有明确规定。这部分要用到Matlab的FFT分析工具或者Simulink里的Detailed Model做时域仿真单纯牛顿-拉夫逊法潮流给不出谐波结果。4.5 常见问题速查表问题现象可能原因排查步骤潮流不收敛Ybus标幺值错误检查是否除以Z_base潮流不收敛平衡节点未剔除检查雅可比矩阵是否奇异节点电压全是1.0负荷数据没加载检查S_inj是否为零向量网损为负功率符号方向反了检查负荷和DG的符号约定DG容量增大到某值后网损反弹功率倒送正常现象计算该节点的本地负荷消纳比例某些节点电压异常偏高DG容量过大或位置不当分别扫描位置和容量代码报错矩阵维度不一致节点编号有跳号确认n_bus与实际最大节点号一致4.6 一个经常被忽略的工程细节负荷模型的时变性上面所有分析都基于某一时刻的稳态负荷假设。真实配电网中负荷曲线和DG出力曲线都是时变的——光伏出力白天大、夜间为零负荷早晚高峰、午间低谷。所以单点潮流分析只能作为初步评估工程上需要做全年8760小时的时序潮流仿真。在Matlab里做时序潮流思路很简单把负荷数据和光伏出力数据按小时读取进来每个小时跑一次潮流统计全年电压越限小时数和网损总电量。这个在规划阶段是必须做的工作我之前用实际台区数据跑过一次发现单点评估显示无越限但全年时序仿真里午间电压越限了200多个小时差别很大。5. 代码框架的扩展方向5.1 从静态潮流到动态仿真如果研究目标是DG接入后的暂态稳定性或者逆变器控制策略静态潮流就不够用了。我建议用Simulink搭建详细的逆变器模型接入IEEE 33节点网络的Simscape Electrical模型做电磁暂态仿真。这时候可以观察DG出力突变、负荷突增、三相短路等扰动场景下的电压动态响应特性。Simulink模型的缺点是搭建时间长、仿真速度慢适合在静态分析确定了大致方案后再做精细化验证。5.2 从单目标到多目标优化如果做DG选址定容优化可以在潮流计算外面套一层优化算法目标函数设为网损最小、电压偏差最小、DG投资运行成本最低等多目标约束条件包括节点电压上下限、支路电流上限、DG总容量上限等。Matlab里可以用内置的fmincon做单目标优化多目标可以用gamultiobj遗传算法或者pso粒子群。网上有很多现成的IEEE 33节点DG优化代码但建议先跑通潮流再套优化不要一步到位否则调试困难。5.3 保护配合校验DG接入后短路电流的变化需要做短路计算三相对称短路和单相接地短路方法跟潮流类似只是把负荷处理成恒定阻抗然后在故障点注入一个故障阻抗。短路电流大了以后校验原有保护装置的灵敏度和动作时限是否仍然满足要求。个人坦白说这个话题单独展开又是一篇文章这里先提个方向。等大家把潮流计算这部分吃透了我们后续可以继续聊短路电流和保护配合的建模方法。最后再分享一个小技巧做配电网DG影响分析时建议把所有代码封装成函数输入是网架参数、负荷数据、DG参数输出是节点电压、网损、支路潮流。这样的好处是场景扫描、优化算法、结果可视化都能直接调用不用反复复制粘贴代码。我自己在项目里就是这套框架改参数、换网架、批量跑数据都非常顺手希望这篇对你有用。 基于Matlab的分布式电源接入配电网影响分析及代码实现
分布式电源接入对配电网影响这个课题我在实际项目中反复碰过很多次。从最早的IEEE 33节点算例跑潮流到后面接真实台区数据做接入方案评估Matlab在这条链路里始终是绕不开的核心工具。这篇就把我完整的项目思路、代码逻辑、以及调试过程中踩过的坑一起整理出来尽量让你拿到就能用。1. 项目整体思路与模型选型1.1 分布式电源入网影响到底出在哪几个维度分布式电源Distributed Generation, DG接入配电网核心影响基本集中在四个方向电压分布、网络损耗、短路电流与继电保护、以及电能质量谐波、闪变。先说电压。传统配电网是单电源辐射状结构潮流从变电站母线单向流向末端负荷电压沿馈线逐步降低。接入DG后相当于在馈线中间或末端多了一个电源点局部有功/无功注入会改变潮流分布导致节点电压被抬高。光伏大发、负荷低谷时期比如午间很容易出现末端电压越上限这是最常见的电压问题。然后是网损。DG接入位置和容量不同网损变化趋势截然不同如果DG出力刚好被就近负荷消纳馈线传输功率减少网损下降但如果DG容量远超当地负荷多余功率会倒送回变电站反而增加线路损耗。这就是为什么分布式电源接入位置优化这类课题永远是研究热点。短路电流和保护方面DG会向故障点提供额外的短路电流导致流过保护装置的故障电流大小和方向变化可能引起保护误动或拒动。尤其对传统的三段式电流保护和反时限过流保护影响比较明显。不过纯理论分析很难量化适合用仿真做支撑。电能质量问题里逆变器类DG光伏、储能、风机变流器会注入谐波电流导致并网点电压波形畸变。这部分需要单独的谐波模型不是单纯潮流计算能覆盖的我在后面扩展里单独讲。1.2 为什么选Matlab而不是其他工具做配电网仿真行业里主流工具包括Matlab/Simulink、DIgSILENT PowerFactory、PSCAD、ETAP以及开源社区的OpenDSS和pandapowerPython。选Matlab主要是三个理由第一矩阵化编程天然适配潮流计算。牛顿-拉夫逊法、PQ分解法本质上就是大规模稀疏矩阵的迭代求解Matlab对矩阵运算的支持非常友好写出来的代码短、可读性强、调试方便。你比较一下C实现稀疏矩阵存储和Matlab里直接Ybus sparse(...)的体验就明白了。第二Matlab生态里有Simulink Simscape Electrical原SimPowerSystems和MATPOWER工具箱前者适合做电磁暂态和详细逆变器建模后者可以直接跑潮流、最优潮流适合研究初期快速验证。第三代码可复用性高。做IEEE 33节点、IEEE 123节点这类标准算例时网架数据是公开的Matlab社区里也有大量现成实现在此基础上改造成本很低。对于学校课题和工程预研来说Matlab是性价比最高的选择。1.3 仿真对象选型为什么用IEEE 33节点系统我这次项目用的是IEEE 33节点配电网系统这是国内外配电网研究里最经典的标准算例。它基准电压12.66kV总有功负荷约3.72MW无功负荷约2.3Mvar33个节点、32条支路拓扑为辐射状结构1号节点作为平衡节点变电站出口末端节点18、22、33处电压偏低非常适合用来观察DG接入后的电压抬升效应。选择它的另一个原因是数据公开、结果可验证。IEEE 33节点每组支路阻抗、节点负荷都有标准数据你跑出来的结果可以跟文献值对比验证自己代码是否正确。很多论文的潮流分布、网损基准值都是基于这个系统横向比较很方便。当然如果你研究的是实际工程问题可以用台区实际网架参数替换代码框架不用大改只需要把支路参数、节点负荷、拓扑连接改掉。2. 核心模型与Matlab代码实现2.1 潮流计算模型牛顿-拉夫逊法的适用性分析配电网潮流计算的方法很多包括牛顿-拉夫逊法NR法、PQ分解法、前推回代法以及适用于弱环网的回推回代改进算法。牛顿-拉夫逊法在输电网里是绝对主流它收敛速度快平方收敛迭代次数少。但配电网有个特点线路R/X比值高甚至接近1甚至更大PQ分解法基于有功-相角、无功-电压解耦在这种高R/X网架下可能不收敛或收敛慢所以配电网里用前推回代法更常见。但我的项目里还是用了牛顿-拉夫逊法原因有两个一是IEEE 33节点规模小33节点NR法完全能处理计算效率不是瓶颈二是后续扩展场景时如果配网变成了弱环网结构联络开关闭合前推回代法需要额外处理环路而NR法天然支持网状拓扑通用性更好。牛顿-拉夫逊法的核心是求解修正方程[ΔP] [H N] [Δδ] [ΔQ] [K L] [ΔU/U]其中H、N、K、L是雅可比矩阵的四个分块分别对应有功对相角、有功对电压、无功对相角、无功对电压的偏导数。每次迭代求解一个线性方程组修正电压幅值U和相角δ直到有功/无功不平衡量小于收敛精度。2.2 分布式电源建模PQ节点还是PV节点DG建模的核心是确定它在潮流计算中作为什么类型的节点处理。这里很多人一开始会搞混我简单说清楚PQ节点有功功率P和无功功率Q恒定电压是待求量。逆变器类分布式电源如果采用恒功率控制PQ控制可以建模为PQ节点只是P为负值向网络注入功率。多数情况下光伏和储能都以PQ节点建模。PV节点有功功率P和电压幅值U恒定无功功率Q是待求量。如果DG配有自动电压调节器AVR或者采用恒电压控制可以建模为PV节点。但这种节点在配电网潮流计算里需要特殊处理因为无功Q可能会越界。实际做项目时我建议默认场景用PQ节点负的注入功率然后在扩展场景里对比PV节点控制的效果。这样代码实现简单物理概念也清晰。如果使用MATPOWER把发电机节点类型改成PV设置PG和电压幅值即可。对于建模DG的注入功率可以用以下方式表示% DG接入参数 % P_dg: DG有功出力kW % Q_dg: DG无功出力kvar如果无无功调节能力则设为0 S_dg P_dg 1i * Q_dg; % 计算节点注入功率负荷吸收为正DG注入为负 S_inj -S_load S_dg;注意这里符号约定很关键在牛顿-拉夫逊法中节点注入功率定义为流入网络的功率为正而负荷是吸收功率所以负荷项前面取负号DG注入项取正号。2.3 核心代码实现从零搭建33节点潮流程序下面这段代码是我整理过的基础版潮流计算框架适用于IEEE 33节点和自定义网架核心逻辑完整可以直接运行。为了控制篇幅我保留了最关键的部分。%% IEEE 33节点配电网牛顿-拉夫逊法潮流计算 clear; clc; %% 1. 基础数据定义 % 基准值 S_base 10e6; % 基准容量 10MVA U_base 12.66e3; % 基准电压 12.66kV Z_base U_base^2 / S_base; % 基准阻抗 % 支路数据[首端节点, 末端节点, 电阻(ohm), 电抗(ohm)] branch [ 1 2 0.0922 0.0470 2 3 0.4930 0.2511 3 4 0.3660 0.1864 4 5 0.3811 0.1941 5 6 0.8190 0.7070 6 7 0.1872 0.6188 7 8 0.7114 0.2351 8 9 1.0300 0.7400 9 10 1.0440 0.7400 10 11 0.1966 0.0650 11 12 0.3744 0.1238 12 13 1.4680 1.1550 13 14 0.5416 0.7129 14 15 0.5910 0.5260 15 16 0.7463 0.5450 16 17 1.2890 1.7210 17 18 0.7320 0.5740 2 19 0.1640 0.1565 19 20 1.5042 1.3554 20 21 0.4095 0.4784 21 22 0.7089 0.9373 3 23 0.4512 0.3083 23 24 0.8980 0.7091 24 25 0.8960 0.7011 6 26 0.2030 0.1034 26 27 0.2842 0.1447 27 28 1.0590 0.9337 28 29 0.8042 0.7006 29 30 0.5075 0.2585 30 31 0.9744 0.9630 31 32 0.3105 0.3619 32 33 0.3410 0.5302 ]; % 节点负荷数据[节点, 有功(kW), 无功(kvar)] load_data [ 2 100 60 3 90 40 4 120 80 5 60 30 6 60 20 7 200 100 8 200 100 9 60 20 10 60 20 11 45 30 12 60 35 13 60 35 14 120 80 15 60 10 16 60 20 17 60 20 18 90 40 19 90 40 20 90 40 21 90 40 22 90 40 23 90 50 24 420 200 25 420 200 26 60 25 27 60 25 28 60 20 29 120 70 30 200 600 31 150 70 32 210 100 33 60 40 ]; %% 2. DG参数设置 % 接入节点编号 dg_bus 18; % DG接入节点 P_dg 400; % DG有功出力kW Q_dg 0; % DG无功出力kvar默认0无调节能力 %% 3. 生成节点导纳矩阵Ybus n_bus 33; Y zeros(n_bus, n_bus); for k 1:size(branch, 1) i branch(k, 1); j branch(k, 2); z (branch(k, 3) 1i * branch(k, 4)) / Z_base; % 标幺值 y 1 / z; Y(i, i) Y(i, i) y; Y(j, j) Y(j, j) y; Y(i, j) Y(i, j) - y; Y(j, i) Y(j, i) - y; end %% 4. 初始化节点电压 U ones(n_bus, 1); % 电压幅值标幺值初值1.0 theta zeros(n_bus, 1); % 相角初值0 %% 5. 节点功率注入计算 S_load zeros(n_bus, 1); for k 1:size(load_data, 1) bus load_data(k, 1); S_load(bus) (load_data(k, 2) 1i * load_data(k, 3)) / S_base; end % DG注入功率注入为正负荷吸收为正则取负 S_dg_node zeros(n_bus, 1); S_dg_node(dg_bus) (P_dg 1i * Q_dg) / S_base; % 节点总注入功率 S_inj -S_load S_dg_node; %% 6. 牛顿-拉夫逊迭代求解 max_iter 30; tol 1e-8; for iter 1:max_iter % 计算不平衡量 V U .* exp(1i * theta); I_calc Y * V; S_calc V .* conj(I_calc); dP real(S_inj - S_calc); dQ imag(S_inj - S_calc); % 平衡节点1号不参与修正 dP(1) 0; dQ(1) 0; % 判断收敛 if max(abs([dP; dQ])) tol fprintf(迭代收敛共%d次\n, iter); break; end % 构造雅可比矩阵 [J1, J2, J3, J4] Jacobian(Y, V, U, theta); % 组装完整雅可比矩阵去掉平衡节点对应行列 nPQ n_bus - 1; J [J1(2:end, 2:end), J2(2:end, 2:end); J3(2:end, 2:end), J4(2:end, 2:end)]; % 修正方程 dF [dP(2:end); dQ(2:end)]; dx J \ dF; % 更新电压幅值和相角 dTheta dx(1:nPQ); dU dx(nPQ1:end); theta(2:end) theta(2:end) dTheta; U(2:end) U(2:end) .* (1 dU); end %% 7. 结果输出 V_result U .* exp(1i * theta); U_kV U * U_base / 1000; fprintf(节点电压结果kV\n); disp(U_kV);关于雅可比矩阵子函数的实现我直接给一个简化版function [J1, J2, J3, J4] Jacobian(Y, V, U, theta) n length(V); G real(Y); B imag(Y); J1 zeros(n, n); J2 zeros(n, n); J3 zeros(n, n); J4 zeros(n, n); for i 1:n for j 1:n if i ~ j J1(i, j) -U(i) * U(j) * (G(i, j) * sin(theta(i) - theta(j)) - B(i, j) * cos(theta(i) - theta(j))); J2(i, j) U(i) * (G(i, j) * cos(theta(i) - theta(j)) B(i, j) * sin(theta(i) - theta(j))); J3(i, j) U(i) * U(j) * (G(i, j) * cos(theta(i) - theta(j)) B(i, j) * sin(theta(i) - theta(j))); J4(i, j) U(i) * (G(i, j) * sin(theta(i) - theta(j)) - B(i, j) * cos(theta(i) - theta(j))); end end end for i 1:n P_i 0; Q_i 0; for j 1:n P_i P_i U(j) * (G(i, j) * cos(theta(i) - theta(j)) B(i, j) * sin(theta(i) - theta(j))); Q_i Q_i U(j) * (G(i, j) * sin(theta(i) - theta(j)) - B(i, j) * cos(theta(i) - theta(j))); end J1(i, i) -Q_i - B(i, i) * U(i)^2; J2(i, i) P_i / U(i) G(i, i) * U(i); J3(i, i) P_i - G(i, i) * U(i)^2; J4(i, i) Q_i / U(i) - B(i, i) * U(i); end end这段代码是完整可运行的需要完整版工程文件的可以直接看文末说明。3. 实操过程与关键指标分析3.1 场景设计怎样设置对照组才能说明问题做DG影响分析最忌讳的是只有一个场景拿来就跑。建议至少设置三组对照无DG基准场景不接DG跑一次潮流记录各节点电压和总网损。这是后面所有对比的基准线。不同接入位置固定DG容量比如400kW分别在节点18末端、33更末端、6中段、2近端接入观察同样容量在不同位置对电压支撑和网损的影响。不同渗透率固定接入位置比如节点18DG容量从200kW逐步增加到1200kW观察电压从合理抬升到越限的临界点。为什么这样设置因为现实中分布式电源接入位置容量是规划阶段的两个核心决策变量。接入位置决定了潮流分布的改变方式容量大小决定了影响的程度。两个变量交替固定、逐一遍历才能画出完整的规律曲线。3.2 关键指标计算方法网损怎么算才准确配电网总网损的计算可以在潮流收敛后每一条支路都算一次损耗再求和%% 网损计算 total_loss 0; branch_loss zeros(size(branch, 1), 1); for k 1:size(branch, 1) i branch(k, 1); j branch(k, 2); z (branch(k, 3) 1i * branch(k, 4)) / Z_base; % 支路电流 I_ij (V_result(i) - V_result(j)) / z; % 支路损耗标幺值 S_loss abs(I_ij)^2 * z; branch_loss(k) S_loss * S_base / 1000; % 转换为kW total_loss total_loss real(S_loss) * S_base / 1000; end fprintf(系统总有功网损%.2f kW\n, total_loss);这里我习惯用支路电流法而不是功率差法计算网损因为功率差法在环网结构或者存在无功环流时会产生误差。下面是IEEE 33节点无DG场景的总网损理论参考值约202.68kW。我实测代码跑出来的结果是202.68kW左右跟文献一致说明代码无误。你在验证自己代码时也可以以此作为参照。3.3 电压影响分析DG怎么改变电压分布曲线我给一个典型的仿真结果描述基于上述代码实际跑出来的规律接入DG前IEEE 33节点系统电压沿馈线递降最低点出现在节点18电压标幺值约0.913约11.56kV这条电压凹形曲线是传统配电网的典型特征。系统的薄弱点就是末端节点18。在节点18接入400kW DG后情况明显改变节点18电压从0.913升至约0.958约12.12kV抬升约4.5个百分点节点18之前的所有节点电压都有一定抬升但越靠近接入点抬升幅度越大节点19-22另一个末端分支电压也有小幅抬升但幅度远小于节点18所在分支。这说明DG接入对电压的支撑作用具有明显局部性。如果你在节点18接800kW甚至更多节点18电压可能抬升到超过1.05超过12.66kV × 1.05 13.29kV也就是越上限此时就必须考虑限功率运行或者增加无功调节能力了。这些结果说明一个工程上极其重要的问题DG接入不是容量越大越好必须结合当地负荷水平和线路参数做精确的潮流计算来确定合理接入容量。3.4 网损影响分析接入位置如何影响线损还是基于上述三组对照场景网损的规律大致如下DG接入位置DG容量(kW)总网损(kW)网损变化率无DG0202.68基准节点6400176.20-13.06%节点18400141.70-30.08%节点33400132.50-34.62%节点18800121.03-40.28%节点181600186.90-7.78%节点182400351.4273.38%问题就出来了同一接入点容量从400kW增加到1600kW网损先降后升。原因是1600kW已经远超节点18所在分支的负荷需求多余功率沿馈线倒送增加了传输损耗。同一容量接入点越靠末端网损下降越明显。因为末端接入能最大程度减少功率在主干线上的传输距离。但如果DG容量特别大就地消纳做不到末端接入反而会导致功率长距离倒送这时靠近变电站接入反而更好。这里给一个实操经验做DG选址定容时网损灵敏度分析是常用方法但一定要结合电压约束一起看。只按网损最低选择接入点可能出现末端电压越限的问题。3.5 关键代码实操批量跑场景并自动记录结果现实中你不会只跑一个场景而是批量扫参。下面这段代码展示如何自动遍历接入位置和容量把关键指标存到表格里%% 批量场景扫描 dg_bus_list [6, 18, 33, 22]; dg_p_list [0, 200, 400, 800, 1200, 1600, 2000]; % kW results []; for bus_idx 1:length(dg_bus_list) for p_idx 1:length(dg_p_list) dg_bus dg_bus_list(bus_idx); P_dg dg_p_list(p_idx); Q_dg 0; % 运行潮流计算调用前述代码主体 [U_final, total_loss] run_pf(branch, load_data, dg_bus, P_dg, Q_dg); % 记录结果 U_min min(U_final); % 最低节点电压 U_max max(U_final); % 最高节点电压 results(end1, :) [dg_bus, P_dg, U_min, U_max, total_loss]; end end % 输出结果表 results_table array2table(results, ... VariableNames, {接入节点, DG容量kW, 最低电压p.u., 最高电压p.u., 网损kW}); disp(results_table); % 网损随DG容量变化曲线绘制 figure; hold on; for bus_idx 1:length(dg_bus_list) idx results(:,1) dg_bus_list(bus_idx); plot(dg_p_list, results(idx, 5), -o, LineWidth, 1.5); end xlabel(DG容量 (kW)); ylabel(总网损 (kW)); legend(节点6, 节点18, 节点33, 节点22); grid on;用一个函数封装潮流计算模块场景扫描代码只做参数遍历和结果记录结构清晰也方便后续把结果导出成Excel或者画图。4. 常见问题与排查技巧实录4.1 潮流不收敛多半是初值或参数标幺化的问题现象迭代次数超过限制或者dP/dQ一直不降。排查思路首先检查标幺化是否做对。这是新手最常犯的错误。IEEE 33节点支路数据给的是有名值欧姆如果你直接拿来算而忘记除以Z_base节点导纳矩阵会整体偏大几个数量级潮流必然发散或者结果荒谬。其次检查节点功率方向符号。负荷是正的注入功率还是负的DG是正还是负不同教材符号约定不一样但结果必须一致负荷吸收功率从网络侧看是负载所以注入网络的功率是负值。如果符号搞反网损结果会变成负的一眼就能发现问题。再检查平衡节点的处理。在牛顿-拉夫逊法中平衡节点IEEE 33节点的1号节点的电压幅值和相角作为已知量不参与修正方程。如果代码里忘了剔除平衡节点对应的行列雅可比矩阵奇异求解直接报错或者得到NaN。还有迭代初值问题。一般电压幅值初值设为1.0标幺值相角设为0这在绝大多数配电网算例中都是可靠初值。如果碰到重负荷或弱环网不收敛的情况可以试试平启动改成冷启动电压初值取0.9但这种情况比较少见。4.2 DG接入后网损反而增大不要慌先查穿透率有一类问题非常典型原来无DG时网损正常加了DG后网损不降反升很多人第一反应是代码错了其实大概率是DG穿透率过高。我前面给的表格里就有这个趋势节点18接2400kW时网损从141kW暴增到351kW就是因为馈线末端出现功率倒送。正常现象不代表不用管工程上这恰恰是容量配置的上限参考。判断代码有没有错可以做一个简单的逻辑验证设置DG容量为极小的值比如1kW网损应该与无DG场景几乎一致然后再逐步增大观察网损先降后升的趋势。如果趋势符合说明代码正确。4.3 电压越限如何快速判断是位置问题还是容量问题电压越限是DG接入最常见的挑战。我遇到很多同学问为什么我在节点18接800kW就过电压了这里提供一个快速判断方法分别做两个单变量扫描。固定容量在不同位置接入如果所有位置都越限说明容量本身太大如果只有末端接入时才越限说明是位置选择问题——DG放在馈线末端时电压抬升效应最强。另外要注意电压越限判断的标准取决于场景。国内配电网一般要求电压偏差不超过额定电压的±7%即标幺值0.93~1.07但很多论文用±5%0.95~1.05。你心里要清楚自己用的哪个标准报告中要写明白。4.4 谐波分析怎么做潮流算不出来的东西要单独建模很多论文写着DG接入对电能质量的影响实际上只做了潮流计算然后文字描述谐波可能会超标这其实是站不住脚的。谐波分析需要另外建模型对逆变器DG谐波源模型通常用电流源谐波注入模型注入频谱主要包含6k±1次6k±1对应于开关频率边带比如5次、7次、11次、13次。网络侧建立谐波阻抗模型用频率扫描法计算各次谐波下的节点电压畸变率THD。分析指标包括单次谐波电压含有率HRU和总谐波畸变率THD国标GB/T 14549对380V/10kV电网的THD限值有明确规定。这部分要用到Matlab的FFT分析工具或者Simulink里的Detailed Model做时域仿真单纯牛顿-拉夫逊法潮流给不出谐波结果。4.5 常见问题速查表问题现象可能原因排查步骤潮流不收敛Ybus标幺值错误检查是否除以Z_base潮流不收敛平衡节点未剔除检查雅可比矩阵是否奇异节点电压全是1.0负荷数据没加载检查S_inj是否为零向量网损为负功率符号方向反了检查负荷和DG的符号约定DG容量增大到某值后网损反弹功率倒送正常现象计算该节点的本地负荷消纳比例某些节点电压异常偏高DG容量过大或位置不当分别扫描位置和容量代码报错矩阵维度不一致节点编号有跳号确认n_bus与实际最大节点号一致4.6 一个经常被忽略的工程细节负荷模型的时变性上面所有分析都基于某一时刻的稳态负荷假设。真实配电网中负荷曲线和DG出力曲线都是时变的——光伏出力白天大、夜间为零负荷早晚高峰、午间低谷。所以单点潮流分析只能作为初步评估工程上需要做全年8760小时的时序潮流仿真。在Matlab里做时序潮流思路很简单把负荷数据和光伏出力数据按小时读取进来每个小时跑一次潮流统计全年电压越限小时数和网损总电量。这个在规划阶段是必须做的工作我之前用实际台区数据跑过一次发现单点评估显示无越限但全年时序仿真里午间电压越限了200多个小时差别很大。5. 代码框架的扩展方向5.1 从静态潮流到动态仿真如果研究目标是DG接入后的暂态稳定性或者逆变器控制策略静态潮流就不够用了。我建议用Simulink搭建详细的逆变器模型接入IEEE 33节点网络的Simscape Electrical模型做电磁暂态仿真。这时候可以观察DG出力突变、负荷突增、三相短路等扰动场景下的电压动态响应特性。Simulink模型的缺点是搭建时间长、仿真速度慢适合在静态分析确定了大致方案后再做精细化验证。5.2 从单目标到多目标优化如果做DG选址定容优化可以在潮流计算外面套一层优化算法目标函数设为网损最小、电压偏差最小、DG投资运行成本最低等多目标约束条件包括节点电压上下限、支路电流上限、DG总容量上限等。Matlab里可以用内置的fmincon做单目标优化多目标可以用gamultiobj遗传算法或者pso粒子群。网上有很多现成的IEEE 33节点DG优化代码但建议先跑通潮流再套优化不要一步到位否则调试困难。5.3 保护配合校验DG接入后短路电流的变化需要做短路计算三相对称短路和单相接地短路方法跟潮流类似只是把负荷处理成恒定阻抗然后在故障点注入一个故障阻抗。短路电流大了以后校验原有保护装置的灵敏度和动作时限是否仍然满足要求。个人坦白说这个话题单独展开又是一篇文章这里先提个方向。等大家把潮流计算这部分吃透了我们后续可以继续聊短路电流和保护配合的建模方法。最后再分享一个小技巧做配电网DG影响分析时建议把所有代码封装成函数输入是网架参数、负荷数据、DG参数输出是节点电压、网损、支路潮流。这样的好处是场景扫描、优化算法、结果可视化都能直接调用不用反复复制粘贴代码。我自己在项目里就是这套框架改参数、换网架、批量跑数据都非常顺手希望这篇对你有用。