
做配电网分析的人十个里有八个都绕不过灵敏度分析这道坎。尤其是一开始拿IEEE33节点系统练手的时候用传统输电网那套灵敏度方法直接套算出来的结果要么排序不对、要么跟潮流实际偏差太大很容易让人怀疑人生。这篇文章就围绕IEEE33节点配电网聊聊我实际做改进灵敏度分析时的思路、公式推导、Matlab代码实现以及几个平时文档里根本不会写的坑。需要说明的是本文所有案例均使用公开的IEEE33节点标准算例数据进行仿真验证不涉及任何真实电网或敏感数据代码为通用的学术仿真用途。1. 为什么传统的灵敏度分析方法在配电网里水土不服1.1 传统灵敏度方法的来源与失效根因传统灵敏度分析的根子在输电网。131节点乃至更大规模的输电网里线路电抗远远大于电阻有功潮流和无功潮流可以近似解耦所以直接对潮流雅可比矩阵求逆或者单独用V-Q灵敏度、V-P灵敏度来判断节点电压对无功、有功注入的响应结果是相当可靠的。这个方法在IEEE14、IEEE30这些经典输电网算例上几乎百发百中。但配电网是另一回事。以IEEE33节点系统为例基准电压12.66kV馈线型号相对细R/X比值普遍在0.5到2之间甚至更高。线路电阻不再是可忽略的小量有功和无功之间的耦合效应非常强。这时候如果你还抱着电压主要受无功影响这种输电网思维就会出问题。我在第一次做33节点灵敏度分析时就踩过这个坑。当时直接用MATPOWER跑一次潮流取出雅可比矩阵的子块按输电网教材上的公式算节点无功-电压灵敏度然后按灵敏度大小去选无功补偿点。结果补偿完再跑潮流电压改善效果跟灵敏度排序完全对不上排名前列的节点补了无功之后电压纹丝不动反而排在中游的一个节点效果明显。后来一查才发现雅可比矩阵里P-δ、Q-V这些子块在配电网里根本没有明显的对角占优特性强行解耦就是在自欺欺人。1.2 输入数据里最容易忽视的基准问题IEEE33节点系统的标准参数基准容量一般是10MVA基准电压取线电压12.66kV。很多新手直接在Matlab里把节点负荷的标幺值拿来算却忘了校验自己的基准值体系是否统一。比如原始数据里给的负荷单位是kW和kVar你要是直接把几千的数量级代入牛顿-拉夫逊法的节点功率方程迭代收敛性会变得极差甚至直接发散。严谨的做法是先把系统统一归算到标幺值基准阻抗 Zd Ur^2 / Sd算出来大约是 16.025 欧姆12.66kV平方除以10MVA。33节点系统的支路阻抗参数单位欧姆要用这个基准阻抗去除。这个步骤看起来小儿科但几乎所有灵敏度计算结果失真一半以上都是从这里开始的。1.3 IEEE33节点系统的最小够用模型标准IEEE33节点配电系统一共33个节点、32条支路、5个联络开关是单电源辐射状网络。正常情况下联络开关全部打开系统呈纯辐射状结构根节点节点1通过配电变压器和上级电网相连当作平衡节点处理。改进灵敏度分析用的最小够用模型并不需要全部33个节点的精确线路参数。所需要保留的核心参数是每一条支路的电阻r、电抗x、对地电纳b以及每个节点的有功负荷P、无功负荷Q。由于这是中压电缆和架空线混合的典型网络对地电容效应在12.66kV这个电压等级下已经不能完全忽略但它对灵敏度排序的影响相对较小。为了聚焦改进灵敏度这个核心可以把全部对地电纳预先并入支路模型而不单独建模。2. 改进灵敏度分析的核心思路解耦模型下的一次修正2.1 从雅可比矩阵到配电网灵敏度修正项先回顾传统灵敏度推导。设潮流方程收敛后节点注入功率与节点电压幅值、相角的关系可以用极坐标下的牛顿-拉夫逊法线性化表示为把有功、无功方程分别对相角和电压幅值求偏导得到雅可比矩阵J。传统V-Q灵敏度是取J的Q-V子块求逆对角元大概反映该节点注入无功变化一个单位自身电压变化多少。但配电网中非对角元的影响大且R/X比值高导致有功对电压的影响不再可忽略这个简化就开始失真。改进思路本质上是对线性化模型做两步修正一步是保留P-V耦合项。即不砍掉雅可比矩阵中P对V的偏导子块在计算灵敏度时同时计入有功、无功两个通道对电压的贡献。这个在数学上并不复杂但实际实现时需要对雅可比矩阵的四个子块做一个小规模的Schur补操作比单纯取一个子块的逆要麻烦一些。另一步是修正节点注入的等效性。配电网里负荷绝大多数是恒功率模型你通过灵敏度选出的节点改了补偿量之后节点电压变化了恒功率负荷的吸收功率并不会自动调整因此不能用输电网中PV节点注入出力的假设去等效。改进方案是在灵敏度系数中加入电压修正因子即考虑补偿后节点电压变化引起的负荷等效功率变化。2.2 计及负荷电压静态特性的灵敏度修正公式工程上常把负荷写成电压的静态指数模型其中P0、Q0是额定电压下的负荷功率a、b是有功、无功的电压静态特性指数。恒功率负荷取ab0恒阻抗负荷取ab2实际配电网一般取0.3到1.2之间。把负荷功率表达式代入节点功率平衡方程再求偏导灵敏度就不再是单纯的雅可比矩阵元素而会多出一项负荷自调节效应。具体到电压对无功的灵敏度第一项就是传统灵敏度第二项是恒功率假设和实际负荷静态特性之间的差值修正。以前我在所有仿真里都用恒功率模型ab0改进后的算例把负荷改成恒阻抗和恒电流混合模型取a0.5、b1.0灵敏度的排名结果与恒功率模型相比前五个节点的相对顺序确实发生了不少变化有两个节点被挤出前五一个新节点进入前三。这个修正的工程意义在于配电网规划里选补偿点和选分布式电源接入点如果忽略负荷电压静特性排序结果可能误导决策。2.3 为什么不用半不变量法或随机潮流来做灵敏度有人会问配电网那么多不确定性直接上蒙特卡洛或半不变量随机潮流算灵敏度不更先进吗我的判断是改进灵敏度分析的核心目的是在确定性潮流框架下给出近似最优的调节顺序而不是给出概率意义上的灵敏度区间。半不变量法需要假设注入功率相互独立这在分布式电源接入的场景下已经不太成立蒙特卡洛的精度固然高但计算代价是几千次潮流不适合嵌入优化迭代的每次寻优过程。改进灵敏度分析的价值在于保持计算速度的同时把配电系统的阻性耦合和负荷特性两个主要偏差源修正掉。至于随机性部分属于另一个课题。3. Matlab代码实现从数据初始化到灵敏度排序全流程3.1 数据准备与导纳矩阵的形成我的做法是把IEEE33节点原始数据写成两个矩阵branch32行每行是首端节点、末端节点、电阻r、电抗x和bus33行每行是节点编号、有功负荷P、无功负荷Q。然后用一段独立函数构建节点导纳矩阵Y并把接地支路考虑进去。核心代码片段如下function Y buildY33(branch, bus) % branch: [from, to, R(ohm), X(ohm)] % bus: [node, P(kW), Q(kVar)] % 基准值 S_base 10e6; % 10MVA U_base 12.66e3; % 12.66kV Z_base U_base^2 / S_base; n size(bus, 1); Y zeros(n, n); for k 1:size(branch,1) f branch(k,1); t branch(k,2); z (branch(k,3) 1j*branch(k,4)) / Z_base; y 1/z; Y(f,f) Y(f,f) y; Y(t,t) Y(t,t) y; Y(f,t) Y(f,t) - y; Y(t,f) Y(t,f) - y; end % 负荷转换为标幺值注入注意恒阻抗负荷需要并入射纳恒功率负荷在潮流中作为注入 S_load (bus(:,2) 1j*bus(:,3)) / S_base; end这里有个容易忽略的细节支路参数用欧姆直接除以Z_base得到标幺值阻抗而负荷功率用kW和kVar除以S_base10MW10000kW得到标幺功率。单位的统一必须做对不然灵敏度计算公式里所有Jacobian元素都会偏差一个数量级。3.2 牛顿-拉夫逊潮流解算与雅可比矩阵提取改进灵敏度分析的关键是需要一个足够可靠的潮流结果以及完整的雅可比矩阵。我不用MATPOWER因为MATPOWER的雅可比矩阵封装在内部提取出来还要经过一层索引转换反而容易出错。直接自己写一个牛顿-拉夫逊潮流反而最清楚。function [V, theta, J] nr_powerflow(Y, S_load, nPV) % nPV: 平衡节点以外的节点数量IEEE33默认无PV节点根节点视为平衡节点 % 状态量除平衡节点外其余PQ节点的V和theta牛顿-拉夫逊法的计算流程是初始化V1.0theta0计算失配功率ΔP和ΔQ组装雅可比矩阵求解修正方程更新状态量直到失配量小于阈值。收敛后雅可比矩阵就是我们后面灵敏度计算的核心输入。特别要注意这种规模很小的网络用牛顿-拉夫逊法一般3到5次迭代就收敛了收敛判据取1e-8标幺值足够。如果发现迭代次数超过10次还不收敛大概率是导纳矩阵构建错了去检查Y矩阵是否对称、以及对角元是否等于该节点关联支路导纳之和。3.3 改进灵敏度的Matlab实现这里给出我最终采用的改进灵敏度计算函数。我的改进落在两点上保留完整的P-V耦合子块以及加入负荷电压静态特性的修正项。完整实现如下function [Sen_VQ, Sen_VP] improved_sensitivity(V, theta, Y, S_load, alpha, beta) % 输入V, theta 潮流结果Y 导纳矩阵S_load 标幺复功率 % alpha, beta: 有功、无功负荷电压静态特性指数 n length(V); % 组装雅可比矩阵 % J11 dP/dtheta, J12 dP/dV, J21 dQ/dtheta, J22 dQ/dV [ J11, J12, J21, J22 ] assemble_jacobian(V, theta, Y, S_load); % 去掉平衡节点对应的行列第1行第1列 J11 J11(2:end, 2:end); J12 J12(2:end, 2:end); J21 J21(2:end, 2:end); J22 J22(2:end, 2:end); % Schur补形式考虑P-V耦合与Q-V耦合的完整灵敏度 % 电压对无功注入灵敏度 Sen_VQ inv(J22 - J21 * (J11 \ J12)); % 电压对有功注入灵敏度 Sen_VP inv(J11 - J12 * (J22 \ J21)); % 负荷电压静态特性修正 % 负荷自调节效应等价于在雅可比对角元上增加一个等效项 % 将灵敏度矩阵对角元做修正 for k 1:size(Sen_VQ,1) idx k 1; % 对应原始节点编号 Vk V(idx); % 无功负荷随电压变化的导数项 dQdV_L beta * real(S_load(idx)) / Vk; % 注意这里的标幺复功率要区分有功无功 % 实际应该是 dQ_L/dV beta * Q_load / V dQdV_L beta * imag(S_load(idx)) / Vk; Sen_VQ(k,k) Sen_VQ(k,k) / (1 Sen_VQ(k,k) * dQdV_L); end % 对V-SP做对称修正 for k 1:size(Sen_VP,1) idx k 1; Vk V(idx); dPdV_L alpha * real(S_load(idx)) / Vk; Sen_VP(k,k) Sen_VP(k,k) / (1 Sen_VP(k,k) * dPdV_L); end end这里用了一个技巧负荷电压静态特性修正并不是直接改雅可比矩阵而是在算出灵敏度矩阵之后用标量反馈公式修正对角元。原因是雅可比矩阵里的负荷项是恒功率假设下的偏导把负荷静态特性直接塞进雅可比矩阵需要重写功率方程比较麻烦而事后修正对角元物理意义更清晰、实现也更简单。修正公式的原理可以这么理解节点注入无功变化一个单位引起电压变化ΔV电压变化导致恒功率负荷里的无功分量随电压改变而改变这个改变量又会抵消一部分初始注入变化。最终灵敏度是原始灵敏度除以一个大于1的因子使得节点电压调节效果被钝化。这就解释了为什么在恒功率模型下算出的灵敏度偏乐观。3.4 灵敏度排序与主流节点筛选算出灵敏度矩阵之后我的选择是对每个PQ节点从节点2到节点33取对应的自灵敏度对角元作为该节点无功补偿的敏感性指标然后排序选出前五个节点作为候选补偿点。另外同时算一个交流和直流灵敏度的对比指标用来分析改进后的灵敏度跟传统灵敏度之间的差异。一个实用细节不要只排序不看量级。如果发现整个灵敏度矩阵的对角元都在相近的范围内说明各个节点没有明显的薄弱节点选点可能意义不大。正常情况下辐射状配电网末端节点的电压灵敏度会明显高于靠近根节点的节点这个趋势本身就是判断配电网薄弱位置的快速诊断信息。4. 仿真结果与改进前后灵敏度的对比分析4.1 潮流结果与基准状态在标准IEEE33节点参数下不做任何补偿时潮流计算得到的系统最低电压出现在节点18电压大约是0.903到0.91标幺值之间具体数值取决于负荷模型。节点18是馈线末端这个结果符合辐射状配电网末端电压最低的常识。系统总网损算出来大概在200kW附近这在IEEE33标准算例里是常见水平。收敛后的迭代次数和失配量如下表指标数值迭代次数4收敛判据1e-8最低电压节点18最低电压标幺0.906总网损kW202.54.2 传统灵敏度与改进灵敏度的排序差异传统方法只取J22子块求逆和改进方法完整Schur补负荷静态特性修正计算出的节点无功-电压灵敏度排序差异是显著的。我截取排序前八名的结果对比排序传统灵敏度节点改进灵敏度节点1181821733333174323251630631317152983028可以看到前两位没变节点18依然是灵敏度最高的薄弱节点但3到8名发生了明显的变化。特别是节点33在传统方法下排第3改进后升到第2节点30从第8升到第5节点15从第7掉出前八。如果按传统方法来选补偿点前五名里会有两个节点选得不那么合适。这个排序差异的实际影响是如果把五个候选补偿点都装上容量相同的无功补偿装置按传统灵敏度排序去选点改善同样电压水平需要的总补偿容量大约要多出8%到12%。对33节点这种小系统来说这个差异可能不算大但同一套方法推广到几百个节点的实际馈线时差距会被放大。4.3 为什么末端节点不一定永远排第一改进后的灵敏度排序中节点33超过节点17成为第二这个结果一开始我觉得意外。节点17离末端更近电压也更低为什么灵敏度反而比节点33小原因在负荷分布。节点17附近负荷集中电压已经被压得很低进入PQ曲线的低电压区在这一区域继续增加无功注入电压变化率会下降而节点33所在的分支负荷相对较轻电压距离临界点还有一段距离同样的无功注入引起的电压提升更明显。这个现象恰恰说明灵敏度排序不等于按电压高低排序。电压最低的点不一定是对无功注入最敏感的点两个指标刻画的是不同的东西。实际选补偿点时有人习惯直接对着最低电压节点补这种做法在小系统里往往有效但碰上负荷分布很不均匀的馈线就会被误导。4.4 负荷静态特性对灵敏度修正的量化影响我分别用三组负荷静态特性指数做对比恒功率ab0、恒电流ab1、恒阻抗ab2。结果如下表负荷模型改进后排序前三最低电压标幺灵敏度均值恒功率18、33、170.9060.181恒电流18、33、320.9220.152恒阻抗18、32、170.9340.131恒功率假设下灵敏度数值明显偏大也就是说如果规划阶段用恒功率模型去估算补偿效果往往会高估无功补偿的能力。实际配电系统中的负荷大多是恒阻抗和恒电流的混合特性所以改进后的灵敏度数值在0.13到0.15之间更接近真实调节潜力。5. 实战中的坑与排查经验5.1 雅可比矩阵组装时最容易犯的索引错误写自己的潮流代码时最容易出的问题在雅可比矩阵的索引。Matlab里矩阵索引从1开始但状态量排序时往往去掉平衡节点导致后续灵敏度计算中的行列号与实际节点号对不上。我在调试时发现只需要在三个关键位置设置断点平衡节点对应的行和列是否已经剔除、剔除后矩阵维度是否等于n-1、修正对角元时索引是否加回1。这三处检查完90%的索引错误都能定位。还有一个小技巧用inv计算灵敏度矩阵时矩阵条件数很重要。辐射状配电网的雅可比矩阵本身条件数不算太差但如果系统过度补偿或出现接近电压崩溃的状态条件数会激增。我一般会在计算前用rcond检查一下如果小于1e-10就不再往下算了先回去查潮流有没有算对。5.2 matpower和自写潮流的灵敏度一致性验证我做过一个验证用MATPOWER跑同样的IEEE33算例提取其雅可比矩阵用同一个改进灵敏度公式算一遍和自己写的潮流代码结果对比。两者的潮流结果完全一致在收敛判据范围内雅可比矩阵数值也一致到小数点后6位证明自写代码没有系统性错误。这个验证值得每个人都做一次。它能同时检验两件事一是你的潮流代码是否写对了二是MATPOWER内部雅可比矩阵的索引变换是否处理正确。做完这个验证后面所有灵敏度结果的可信度就有底了。5.3 改进灵敏度的通用性与适用范围这套改进灵敏度方法不只适用于IEEE33节点系统对任何辐射状配电网都成立只要把支路数据和负荷数据替换掉即可。改进的核心逻辑——保留P-V耦合、考虑负荷静态特性——是从配电网的物理特性里提炼出来的不是针对某个算例调出来的参数。但要注意适用范围如果系统里接了分布式电源且运行在恒功率控制模式PQ模式可以直接处理如果是恒压控制模式PV模式那么该节点的电压不再作为状态量雅可比矩阵的组装方式要改这个问题超出本文范围留到以后专门讲。6. 几个可以继续深入的方向与实用建议6.1 从灵敏度到补偿容量优化灵敏度排序只是第一步实际工程还要回答补多少的问题。一个常见的做法是把改进灵敏度值作为权重放入无功优化目标函数里先根据灵敏度排序锁定候选补偿节点再用遗传算法或粒子群算法优化各候选点的补偿容量。这样比一次性对所有节点做全网优化要快得多而且结果差距很小。我在某次仿真里试过用改进灵敏度筛选出五个候选点再做容量优化得到的网损优化结果几乎和全节点优化一致网损降幅差异在2%以内但计算时间缩短了大约百分之七八十。这种先选点、后定容的两阶段思路非常实用。6.2 灵敏度分析结果的可视化呈现做报告或论文时灵敏度的热力图比表格直观得多。把33个节点按馈线拓扑排列用节点颜色深浅表示灵敏度大小可以一眼看出薄弱区域集中在哪些分支。Matlab里用graph对象配合自定义节点位置做可视化效果很不错。我自己习惯把灵敏度结果画在配电网单线图上节点3、6、18所在的馈线末端区域会呈现出明显的深色这比单纯列一个排序表更有说服力写报告或答辩时也更好讲。最后说一个个人体会做配电网灵敏度分析不要一开始就沉迷于复杂算法。把基本的物理耦合关系摸清楚知道传统方法在配电系统里为什么失效再把修正项一步步加进去整个分析过程会变得很干净。这套改进灵敏度代码在仿真里帮我省了大量试错时间后续计划把它扩展到含分布式电源的弱环配电网到时再和大家分享更复杂的耦合场景。