做配电网规划项目的朋友应该都有同感:分布式电源选址与定容这个问题,看着就是把DG挂到IEEE33节点这样的测试网络上算一算,实际上背后牵扯多目标优化、潮流计算、智能算法调度和Matlab代码实现一整条技术链。我最初接触这个项目时以为只是跑一个粒子群优化算法(PSO)而已,真正动手才发现,最花时间的不是算法本身,而是怎么把离散的节点选择、连续的容量配置和配电网潮流模型合理地耦合在一起。这篇文章我会完整拆解我做这个项目的思路、模型、代码结构和踩过的坑,适合正在做分布式电源规划方向研究的学生,也适合刚入手智能算法在电力系统中应用、想找一个能快速上手且可复现的Matlab算例的工程师。
1. 这个项目解决什么问题:从拍脑袋选址到科学定容
1.1 选址定容的本质是一个组合优化问题
传统的配电网规划里,分布式电源的位置和容量往往依赖经验判断:看看哪条线路末端电压低、哪段线路重载,就在附近装一台光伏或者燃气轮机。这种做法在DG渗透率低的时候问题不大,但DG装多了以后,系统的潮流方向不再是从变电站单向流向负荷。如果DG位置选得不好、容量配得过大,会带来两个典型问题:一是局部电压越上限,尤其是轻负荷时段,DG出力送不出去,节点电压被抬得很高;二是线路出现反向潮流,配电保护的整定和线损都会变得很麻烦。
把这个问题抽象出来,就是一个选址与定容的组合优化问题:在配电网的若干个候选节点中,选出一组节点作为DG接入点,同时为每个接入点确定最优的安装容量。目标可以是网损最小、电压质量最好、DG投资运行费用最低,或者是这几个指标的综合最优。IEEE33节点系统作为经典的辐射状配电网测试算例,承担了这个问题的标准验证平台角色:规模适中、参数公开、潮流计算简单,非常适合测试不同的优化算法和策略。
1.2 为什么不能靠穷举法,必须用智能优化算法
有朋友可能会问:33个节点而已,穷举所有DG接入位置和容量组合不就行了?这里需要算一笔账。假设系统里规划接入3台DG,从33个节点中选3个作为接入点(实际中还要排除变电站根节点),组合数是C(32,3)=4960种;如果每个DG的容量在0到500kW之间,按10kW步长离散化,就有51个档位,三台DG的容量组合是51的三次方,约13万种。乘起来就是6亿多个方案。更关键的是,每评估一个方案都要做一次潮流计算,这个计算量在单机环境下完全不可接受。
所以这个问题的标准解法,是让优化算法在搜索空间中用启发式规则快速逼近较优解。粒子群优化算法(PSO)就是其中非常适合的一种:结构简单、参数少、收敛快,而且连续变量优化的底层机制天然适配容量配置这类实数编码问题。把粒子位置中的一部分映射到候选节点编号,另一部分映射到DG容量,然后用潮流计算的结果去评估粒子的适应度值,PSO就能自动找到网损低、电压质量好的DG接入组合。
2. IEEE33节点系统:经典算例为什么值得反复使用
2.1 系统结构盘点
IEEE33节点配电网系统是1989年提出的经典测试系统,也是目前分布式电源规划、配电网重构、无功优化等方向用得最多的算例之一。它的结构是典型的辐射状配电网:首端通过一个变电站节点(通常编号为0号或1号,各版本略有差异)向下供电,经过32条支路连接33个节点,系统基准电压12.66kV,基准功率通常取10MVA。系统总有功负荷3.715MW,总无功负荷2.3Mvar,峰值时段线路末端电压偏低,这正好给了分布式电源“大展拳脚”的空间。
不同文献里IEEE33系统的节点编号习惯不太一样,有的从0开始编号,有的从1开始,这个问题在实现代码的时候非常容易踩坑,后面我会专门讲。系统的基本参数我整理在下面:
| 参数 | 数值 |
|---|---|
| 节点总数 | 33 |
| 支路总数 | 32 |
| 基准电压 | 12.66 kV |
| 总有功负荷 | 3.715 MW |
| 总无功负荷 | 2.3 Mvar |
| 网络拓扑 | 单电源辐射状 |
| 最大允许电压偏差 | ±5%(即0.95~1.05 p.u.) |
2.2 为什么用它而不是更大的实际配电网系统
实际工程中的配电网动辄上百节点,拓扑复杂,数据也可能涉及保密不能公开。IEEE33的好处在于:第一,数据完全公开,网上很容易找到标准参数表,复现成本极低;第二,规模合适,前推回代潮流计算一次只需要几毫秒,即便粒子种群设成100、迭代200次,总计算时间也就几十秒到几分钟,很适合算法迭代调试;第三,它本身的电压问题明显——重载情况下末端节点电压会跌到0.9左右,接入合适的DG后电压改善效果非常直观,这就让算法的优化效果有了很强的对比度。
我个人的建议是,如果你做DG选址定容方向的研究,先用IEEE33把算法调通、把结果跑稳,再迁移到IEEE69节点、IEEE118节点或者自己手头的实际馈线上去。算法逻辑是通用的,换了系统只需要换节点参数表和支路参数表,再把候选节点约束改一下就行。
3. 多目标优化模型构建:网损、电压偏差与DG成本的三方权衡
3.1 三个目标函数分别衡量什么
标题里写的是多目标优化,这个“多目标”在DG选址定容问题里通常不是指用MOPSO或者NSGA-II去求解帕累托前沿,而是把多个性能指标通过加权方式合成一个综合适应度。至于选哪些指标,工程上最常用的就是下面三个。
第一个是系统有功网损。DG接入后如果位置合适,相当于在负荷中心附近注入了功率,减少了功率在长线路上的流动,网损会下降;但如果位置不合适,反而会增大网损。计算方法是跑完潮流后累加所有支路的功率损耗,优化目标就是让这个总和最小。
第二个是节点电压偏差。配电网规划十分关注电能质量,通常把各节点电压与额定电压的偏差绝对值累计起来作为指标。DG接入能起到电压支撑作用,特别是接在末端节点附近时,对抬升线路末端电压效果显著。但如果DG容量过大,末端电压又会越上限,因此这个指标能倒逼算法找到容量合适的方案。
第三个是DG的年投资运行成本。光伏、风电、燃气轮机的单位装机成本不同,运行维护费用也不同。优化时需要对每种DG类型的单位容量投资成本和年运行费用做估算,折算成年值后与网损费用一起构成经济性指标。多目标加权中给这个指标一个合理的权重,能防止算法为了降网损和调电压而无限加大DG容量。
3.2 指标归一化与权重设置
三个目标函数的量纲完全不同。网损的单位是kW,电压偏差的单位是p.u.,成本是万元,如果直接相加,数值大的那一项会完全主导适应度,其他指标就无法参与优化。所以必须先做归一化处理,把每个指标除以它对应的基准值。常见做法是取未接入DG时原始系统的网损、电压偏差和估算成本作为基准值,优化后各指标除以这些基准,得到一个无量纲的比值。
综合适应度函数的形式可以写成:
f = w1 × (Ploss / Ploss_base) + w2 × (ΔU / ΔU_base) + w3 × (C / C_base)
三者的权重w1、w2、w3之间满足和等于1。默认情况下可以取均匀权重0.33/0.33/0.34,如果项目更关注经济性,可以适当调高w3;如果项目重点解决末端低电压问题,则把w2调大。需要特别提示的一点是:权重本质上体现了规划决策者的偏好,没有绝对正确的取值。我在测试中发现,如果w3设得过高(超过0.7),算法会倾向于完全不装DG或者只装很小的DG,因为不投资是最省钱的;而w2过低时,又可能出现电压越限的方案被当成最优解。这个平衡要靠多次试验来找到合适的区间。
3.3 约束条件:容易被忽略的四个边界
优化模型除了目标函数,还必须带约束条件。第一个是潮流等式约束,即优化得到的接入方案必须能通过潮流计算得到收敛的系统状态,这一步在代码里就是每次评估粒子时调用潮流函数。第二个是节点电压约束,正常情况下所有节点电压应在0.95~1.05 p.u.之间。第三个是DG总渗透率约束,所有DG总容量不能超过系统总有功负荷的一定比例,一般取20%-40%,防止反向潮流过重。第四个是单节点DG容量上限,比如1000kW,避免把过多电源堆在同一个节点上。
这些约束的处理方式我在第5章详细展开,这里先说结论:在Matlab代码实现中,最稳妥的方法是“潮流不收敛或越限就直接给一个很大的惩罚适应度”,配合“对粒子边界进行物理约束”,两者结合可以兼顾搜索效果和稳定性。
4. 粒子群算法落地到选址定容的三个核心设计
4.1 粒子编码:连续位置与离散节点的双向映射
PSO本身是面向连续优化的,但DG接入节点是离散整数变量。这是整个代码实现中最核心的一个设计点。我采用的方案是把一个粒子设计成两段式编码,假设规划方案最多接入NDG台DG,那么粒子的位置向量长度为2×NDG,前NDG维表示各台DG接入的节点编号,后NDG维表示对应节点的DG容量。
问题是PSO迭代过程中粒子的位置分量是连续实数,无法直接作为节点编号使用。解决办法是维护一个候选节点列表,比如[3, 5, 7, 9, 12, 15, 18, 22, 25, 28, 30, 32],粒子位置的前半段分量是候选列表的下标索引,取整到[1, 12]的范围,再通过索引查表得到真正的节点编号。后半段容量分量直接用连续实数表示,范围限定在[Pmin, Pmax]之间。
这种做法的好处是PSO的速度更新公式完全不需要改动,只在解码环节做一次round取整和索引映射即可,逻辑清晰且不易出错。
4.2 解码与潮流计算的衔接
粒子解码之后,下一步就是把DG接入方案写入潮流计算的输入数据中。IEEE33系统的原始数据用两个矩阵表示:节点数据矩阵(节点编号、负荷有功、负荷无功)和支路数据矩阵(首端节点、末端节点、支路电阻、支路电抗)。接入DG后,对应节点的注入有功功率需要减去DG发出的有功功率,即净负荷有功 = 原始负荷有功 − DG有功输出。
这里有一个容易被忽略的细节:DG配电网潮流中,一般把DG简化为PQ节点,即给定有功出力和功率因数(通常取1.0,即纯有功注入)。如果要考虑DG的无功调节能力,可以设为PQ节点并给定无功出力为负值或零。对于光伏和风电这类通过逆变器并网的DG,功率因数设为0.98或1.0都是合理的。在我这个项目中,为了简化并突出优化算法的核心逻辑,按功率因数1.0处理。
4.3 速度更新公式与边界处理
PSO的核心更新公式是经典的速度-位置模型:
v(i+1) = w × v(i) + c1 × r1 × (pbest(i) − x(i)) + c2 × r2 × (gbest(i) − x(i))
x(i+1) = x(i) + v(i+1)
这里的w是惯性权重,c1和c2是学习因子,r1和r2是[0,1]之间的随机数,pbest是个体历史最优位置,gbest是全局最优位置。惯性权重w的取值直接影响全局搜索和局部开发能力的平衡:w较大时粒子飞行速度快,全局探索能力强;w较小时粒子趋于收敛,局部精细搜索能力强。常用做法是线性递减,从0.9逐步降到0.4,让迭代前期尽量探索,后期压榨精度。
边界处理上我试过三种方法:直接截断、随机重置、反弹修正。直接截断最简单,越界就拉回边界值;随机重置适合位置分量,越界后在可行域内重新随机生成;反弹修正在一些测试中收敛更平滑,但代码量稍大。对容量分量建议用直接截断,对节点索引分量建议在四舍五入后再检查一遍候选列表范围,越界就用随机重置,这样粒子多样性损失较小。
5. Matlab代码实现与代码结构解析
5.1 整体流程
这个项目的Matlab实现我建议按模块化方式组织,主程序、潮流计算函数、PSO核心函数、结果可视化函数分开写,方便调试和后续移植。整体流程可以概括为下面几步:
- 加载IEEE33节点系统的节点参数、支路参数;
- 设置PSO参数(种群规模、最大迭代次数、惯性权重范围、学习因子);
- 设置DG参数(最大接入台数、候选节点集合、容量上下限、单位投资成本);
- 初始化粒子群,对每个粒子的位置向量做随机初始化;
- 进入主循环:对每个粒子解码、修改节点注入功率、调用潮流计算函数、计算多目标适应度、检查约束并施加惩罚;
- 更新个体最优pbest和全局最优gbest,更新粒子速度和位置;
- 判断是否达到最大迭代次数,输出最优方案及对应各目标值;
- 绘制收敛曲线和接入DG前后的节点电压分布对比图。
主循环的骨架结构大致是:
for iter = 1:maxIter for i = 1:Npop % 解码粒子 [dgNode, dgPower] = decodeParticle(position(i, :), candNode, NDG); % 修改节点注入功率 loadData = baseLoadData; for k = 1:NDG loadData(dgNode(k), 2) = loadData(dgNode(k), 2) - dgPower(k); end % 潮流计算 [V, Ploss] = powerFlow(branchData, loadData); % 计算适应度 fitness(i) = calcFitness(Ploss, V, dgPower, w1, w2, w3); % 约束检查并施加惩罚 fitness(i) = fitness(i) + penalty(V, dgPower); end % 更新pbest、gbest,更新速度与位置 end5.2 前推回代潮流计算的工程细节
IEEE33节点是辐射状配电网,潮流计算用前推回代法比牛顿-拉夫逊法更合适。前推回代法的原理可以概括为两步:第一步从末端节点向首端节点回推,根据节点功率和末端电压初值逐段计算各支路电流或功率流;第二步从首端节点向末端节点前推,由首端电压和各支路电流逐段计算各节点电压。如此反复迭代,直到两次迭代的电压差小于收敛阈值。
这个方法的优势在于不需要求雅可比矩阵,不需要做矩阵分解,对初值不敏感,而且对于辐射状网络天然收敛。代码实现中需要特别注意支路数据的存储顺序问题——前推回代要求能够从末端逐支路向上回溯,因此支路表最好按照“首端节点靠近变电站、末端节点远离变电站”的方向排列,或者维护一个父节点索引数组。
下面是一个简洁版的前推回代核心代码思路:
function [V, Ploss] = powerFlow(branchData, loadData) V = ones(33, 1); % 电压初值 tolerance = 1e-6; maxIter = 50; for iter = 1:maxIter V_old = V; % 回推:从线路末端向首端计算支路功率 S_branch = zeros(32, 1); for k = 32:-1:1 endNode = branchData(k, 2); S_branch(k) = complex(loadData(endNode, 2), loadData(endNode, 3)); % 叠加下游支路功率... end % 前推:从首端向末端计算节点电压 for k = 1:32 startNode = branchData(k, 1); endNode = branchData(k, 2); V(endNode) = V(startNode) - (S_branch(k) / V(startNode)) * ... complex(branchData(k, 3), branchData(k, 4)); end % 计算网损 Ploss = real(sum(S_branch .* conj(S_branch)) .* branchData(:, 3)); if max(abs(V - V_old)) < tolerance break; end end end这里为了展示核心逻辑做了简化,实际使用中回推阶段特别要注意支路功率是从末端往首端逐级累加的,前推阶段电压是从首端往末端逐步更新的,两个方向不能搞反。
5.3 PSO核心更新模块
PSO模块本身代码量不大,最重要的是确保速度和位置数组的维度与粒子编码长度一致。假设NDG取3,粒子位置向量长度为6。速度初始化为0到1之间的小随机数,速度上限Vmax设为位置范围宽度的20%左右。
速度更新和位置更新的核心代码如下:
w = wMax - (wMax - wMin) * (iter / maxIter); % 惯性权重线性递减 v(i, :) = w * v(i, :) ... + c1 * rand(1, dim) .* (pbest(i, :) - position(i, :)) ... + c2 * rand(1, dim) .* (gbest - position(i, :)); % 速度限幅 v(i, :) = max(min(v(i, :), Vmax), -Vmax); % 位置更新 position(i, :) = position(i, :) + v(i, :); % 边界处理 for d = 1:dim if d <= NDG % 节点索引维度:向上四舍五入并限制在候选索引范围 position(i, d) = round(position(i, d)); if position(i, d) < 1, position(i, d) = 1; end if position(i, d) > numel(candNode), position(i, d) = numel(candNode); end else % 容量维度:截断到[Pmin, Pmax] position(i, d) = max(min(position(i, d), Pmax), Pmin); end endpbest和gbest的更新逻辑需要注意一点:如果新一代粒子的适应度值更小,则更新个体最优;全局最优gbest取所有个体最优中适应度最小的那个。在实际测试中,因为PSO是随机算法,单次运行结果不稳定,建议循环运行多次取最优或者对多次运行结果做统计分析。
5.4 结果输出与可视化
结果输出部分我认为至少应该包含三张图:一是适应度收敛曲线,横轴是迭代次数,纵轴是每次迭代的全局最优适应度值,用来判断算法是否收敛;二是IEEE33节点接入DG前后的节点电压分布对比图,这是最直观体现优化效果的图;三是DG接入位置和容量方案表。
收敛曲线的画法很简单,每次迭代结束后记录一次gbest的适应度值,最后用plot画出来即可。电压分布对比图需要把无DG时各节点电压和优化方案下各节点电压放到同一个坐标系里画,纵轴范围建议设为0.90到1.06,这样能清楚看到DG对末端电压的抬升作用和最大值有没有越限。
我自己的经验是,如果收敛曲线在前20代就快速下降然后进入平缓,说明算法收敛速度不错,但如果曲线一直上下波动不下降,先检查罚函数和边界处理是否正确,其次再考虑调整PSO参数。
6. 实测结果与参数调优:从收敛曲线看算法在干什么
6.1 一组典型优化结果
我用Matlab在IEEE33节点系统上跑了一组测试,条件设置为:种群规模50,最大迭代次数100,规划接入3台DG,候选节点为[8, 12, 15, 18, 22, 25, 28, 32],单个DG容量上限500kW,总渗透率上限不超过系统总负荷的30%。权重设置为网损0.4、电压偏差0.4、成本0.2。
得到的典型优化结果大致如下:
| 指标 | 无DG原始系统 | PSO优化接入DG后 |
|---|---|---|
| 系统有功网损 | 约202 kW | 约95 kW |
| 最低节点电压 | 约0.90 p.u. | 约0.96 p.u. |
| DG接入方案 | 无 | 节点18接入约480kW,节点32接入约420kW,节点25接入约200kW |
| 总DG渗透率 | 0 | 约29.6% |
可以看到,优化后网损下降了一半多,最低电压从接近越限的0.90抬升到0.96以上。这个结果再次验证了DG选址和定容的重要性:同样是装了1100kW的DG,如果位置选得不对,网损可能不仅不降反而上升,电压也可能越限。
6.2 为什么参数微调会导致结果翻车
PSO对参数比较敏感,这个项目里最影响结果稳定性的三个参数是惯性权重、学习因子和种群规模。惯性权重如果用固定值0.4,算法容易过早收敛到局部最优,结果每次跑都不一样;如果用线性递减0.9到0.4,收敛稳定性和解的质量都会好不少。c1和c2取2.0是经典配置,但我在测试中发现c1=c2=1.5配合线性递减w,后期的局部搜索能力更细腻,不过全局探索弱一点。种群规模从30提高到80,解的质量有明显提升,但计算时间线性增长,对于IEEE33这种小系统,50到60的种群规模是比较经济的选择。
还有一个影响结果的因素是惩罚系数的设置。惩罚太轻,越限方案可能在适应度上反而占优;惩罚太重,又会严重压缩粒子在边界附近的搜索空间。我的做法是对电压越限量做平方惩罚,即超出0.95~1.05范围越多,惩罚值按平方量级增长,这样既能让算法主动避开越限区域,又不至于一棒子打死边界附近的可行解。
7. 排坑实录:复现和修改这个项目时踩过的六个坑
7.1 节点编号偏移导致的潮流结果对不上
我在第一次实现潮流计算时发现,算出来的初始网损和文献对不上。排查了很久发现是节点数据里0号节点对应Matlab下标1,1号节点对应下标2,而我直接在nodeData矩阵里用原始节点编号做索引,导致负荷数据错位了一个节点。这个问题在用IEEE33系统时最容易出现,因为不同来源的原始数据节点编号起始不同。建议做法是读入数据之后统一生成一个从原始节点编号到Matlab数组下标的映射表,后续所有操作都走映射表,不要直接拿节点编号当数组索引。
7.2 DG容量过大导致潮流不收敛
优化算法在初始阶段会随机生成大量的粒子,其中难免出现DG总容量偏大的情况。当某个节点注入的DG有功功率超过负荷很多时,潮流计算可能出现不收敛,表现为前推回代迭代次数达到上限但电压差仍大于收敛阈值。如果此时程序直接报错中断,整个优化就没法继续了。处理办法是在潮流函数中增加收敛标志位,返回成功或失败状态,然后在适应度计算时对失败粒子直接返回一个大惩罚值。
7.3 权重分配不当把算法带偏
有一组测试我把成本权重w3设成了0.5,结果算法给出的“最优方案”竟然是什么DG都不装或者只装很小容量的DG。原因很简单:不装DG,投资成本就是0,即使网损高一点、电压差一点,综合适应度可能仍然最小。这说明多目标加权优化中,权重不仅反映偏好,还直接影响解的走向。实际工程里如果必须考虑经济性,建议给成本项加上“投资回收期”或“年收益”之类的表达,把投资和收益挂钩,而不是只算支出。
7.4 单次运行结果波动大
PSO是随机搜索算法,初始粒子群和随机数的变化都会影响最终结果。如果只用单次运行的最优解作为结论,很可能会被一次偶然的糟糕收敛误导。建议对同样的参数组合至少运行10次,记录最优值、平均值和标准差。一组参数如果10次运行的最优适应度标准差很大,说明算法稳定性差,优先排查边界处理和惯性权重设置,而不是纠结单次结果好坏。
7.5 两代粒子都选了同一个节点
当NDG大于1时,解码后可能出现两台DG映射到同一个候选节点的情况,这在实际规划中没有意义。最简单的处理办法是在解码环节后做一个查重,把重复的节点用附近未被选择的候选节点替换。更稳妥的做法是在初始化粒子时保证每个粒子的前NDG维分量互不相同,但这样会增加代码复杂度。我在项目中用的是后置查重,因为它的逻辑独立且不容易破坏PSO的速度更新。
7.6 传统PSO后期收敛变慢
标准PSO在迭代后期容易出现种群多样性下降、粒子聚集在局部最优附近的问题。如果发现收敛曲线在后期基本走平但解的网损不算很理想,可以尝试三个改进:一是在惯性权重递减到0.4之后保持一段时间,让粒子在局部精细搜索;二是对gbest做小范围随机扰动,相当于给全局最优加一点变异;三是引入压缩因子c,将速度更新公式中的w、c1、c2统一缩放。最简单的验证方式是把迭代次数从100增加到200,观察收敛曲线是否还有明显下降空间。
8. 最后一点实操心得
项目做到后面我最大的感受是:粒子群优化算法本身不是这个项目的难点,难点在于把配电网潮流、多目标评价、约束处理和算法搜索空间这四个模块组合成一个稳定可靠的闭环。我建议第一次复现的朋友不要急着追求复杂版本,先把单目标(比如只优化网损)跑通,确认潮流计算无误、PSO能够收敛,再加入电压指标和成本指标做多目标加权,最后再考虑惩罚函数和边界优化。每一步都验证过,出问题的时候才能快速定位是算法的问题还是模型的问题。另外一个小技巧:在调试阶段把种群规模和迭代次数设小一些,比如20个粒子迭代30次,这样跑一组测试只要几秒钟,调参效率会高很多。等参数和逻辑都稳定了,再放大到正式规模。这个习惯帮我节省了大量时间,也避免了在错误代码上反复浪费时间。