基于粒子群算法的主从博弈综合能源调度Matlab实现
最近我一直在调一套园区综合能源交易与调度的仿真程序核心是“基于粒子群优化算法的主从博弈优化模型”。这个方向在综合能源系统、微电网和电力市场研究里很常见上层售能主体负责制定价格策略下层用户根据价格调整用能行为双方各自追求自己利益最大化。粒子群优化算法在这里不是单纯去做一个调度求最小值而是被用来搜索上层策略下层还要嵌套一层最优化响应——整套逻辑比平时做常规机组组合或者经济调度要绕一些。这篇文章我想把这套模型的建模思路、Matlab代码框架和调试中容易踩的坑梳理一遍。适合正在做综合能源系统、微网交易、需求响应方向并且需要落地一套可运行仿真代码的读者。无论你是刚接触主从博弈还是已经用粒子群做过常规优化我觉得这套“外层粒子群搜索结构 内层独立优化”的套路都有一定参考价值。1. 这套模型到底在解决什么问题1.1 为什么不能直接把所有主体当成一个整体来优化我先说个最容易产生误区的地方很多人拿到这类题目第一反应是把所有参与方的收益加在一起或者做一个带权重的综合目标函数然后调一次粒子群就完事。这种做法在纯数学上可能得到一个让整体福利最优的解但在真实能源交易中无法落地。原因是不同主体的利益目标常常冲突。售电方希望售电单价高一点用户希望购电成本低一点储能运营商希望在峰谷价差大的时段多充放。如果强行合成一个目标等于让所有主体放弃自己的独立诉求去服从一个“计划型”的整体最优这不符合实际交易里面的决策顺序。主从博弈解决的核心问题就是处理这种“有先后行动顺序、且目标不一致”的优化场景。1.2 三方三层到底指什么标题里的“三方三层”在不同论文里解释略有差异。我这里用一个比较通用的角度来拆解。所谓三方指的是三个不同利益主体参与方在交易链中的身份主要目标上游供能方/批发代理提供基础能源给中层一个批发购电成本曲线保障自身售电收益希望走量且减少峰时压力中层综合能源服务商从上游购能制定面向用户的零售价格同时管理储能/光伏等设备零售收入减购能成本最大化通过调节储能赚取价差终端用户负荷集群根据服务商给出的零售电价调整可转移负荷和用能计划购电成本最低同时尽量不牺牲用电舒适度所谓三层则对应决策结构和信息流动的层级第1层是价格决策层中层服务商制定面向用户的销售电价第2层是运行调度层服务商优化内部储能、新能源和购能计划第3层是用户响应层用户根据价格重新安排自己的负荷曲线。从博弈关系来看中层服务商是价格制定者即领导者用户是价格接受者即跟随者。上游供能方对整个模型的影响主要是通过批发价参数和电能交互约束传导进来的。所以在具体写程序时我会把上游处理成外部边界条件把中层的价格策略作为外层粒子群优化变量把第三层的负荷响应作为内层优化问题求解。1.3 “先动者定价后动者响应”的互动关系为了理解这个结构可以类比一个社区团购团长和团员之间的互动。团长提前公布鸡蛋价格团员根据价格决定今天买两斤还是一斤团长根据所有团员报上来的需求量再决定是否调整下一轮报价。这里的核心是团长做决策时必须考虑团员会如何反应不能要求团员“为了整体利益”多买鸡蛋。主从博弈里的数学表达是先把上层的价格策略固定住下层求解出在这个价格下的最优负荷曲线然后把这条负荷曲线反馈给上层计算收益。我第一次实现这套模型时最花时间的不是粒子群本身的代码而是把三层之间的变量关系和求解顺序理清楚。只要顺序错了粒子群再强也容易搜出一个莫名其妙的结果。2. 为什么选用粒子群优化算法而不是其他方案2.1 主从博弈模型给传统算法出了什么难题主从博弈的关键特征是上层目标函数里的某个变量比如用户用电量并不是直接给定的常数而是下层优化问题的解。这带来一个很麻烦的问题上层收益对价格策略没有显式的解析梯度。如果我想用梯度下降法去优化价格向量就需要对每个价格变量做一次摄动每摄动一次又得重新求解一遍下层优化计算量非常大。而且下层优化问题含有储能充放电功率、负荷转移比例等变量使得整个映射关系是非线性甚至不连续的梯度算法在这种“黑箱式”映射上很容易卡在局部极值或者压根收敛不到一个有经济意义的解。相比之下粒子群这类启发式算法不依赖目标函数的梯度信息。它通过一群候选粒子在价格变量空间里同时搜索不断把表现好的粒子的历史最优和全局最优信息传递给后续粒子。这种种群机制天然适合处理“上层策略评价一下层优化响应”的黑箱问题。2.2 粒子群优化算法的三个具体优势第一个优势是容易实现。粒子群的速度和位置更新公式很短对Matlab这类支持向量化运算的环境尤其友好。写一个通用的粒子群主循环半小时内就能搞定。第二个优势是并行度好。粒子群每一代里的每个粒子在计算适应度时只依赖粒子自身位置相互之间没有耦合。这意味着可以放心地做矩阵化批量求解或者用并行计算工具箱加速。相比之下很多序列化算法在计算下层响应时要一个个跑循环速度差异会很明显。第三个优势是超参数相对少。粒子群需要调的参数主要是惯性权重、个体学习因子和社会学习因子、种群规模和迭代次数。只要把边界条件、约束处理好使用一组经典参数往往就能得到一个比较合理的搜索结果。即便后续换成多目标粒子群算法也只是在外部加一个档案集和Pareto支配关系主体逻辑不用推翻。2.3 和遗传算法、差分进化算法的对比我也用遗传算法和差分进化算法做过对比。遗传算法因为有选择、交叉、变异三个操作代码量会多出一截而且离散编码和连续价格变量之间的转换稍显繁琐。差分进化算法在连续优化上表现不错但是它的变异策略引入的缩放因子和交叉概率需要根据不同变量维度调参不一定比粒子群省事。需要说明的是粒子群也不是万能药。它同样存在早熟收敛、后期收敛速度变慢等毛病。所以在这套模型里我并没有把内层用户侧响应也交给粒子群去求解而是采用确定性的优化器来做。这样既发挥粒子群在离散型策略空间上的搜索优势又保证了每个粒子适应度评价是稳定可复现的。这一点放到后面代码部分细讲。3. 数学模型与关键设计细节3.1 决策变量与优化目标的对应关系在搭建Matlab代码之前首先要把变量归属搞清楚。不同主体掌握的决策变量完全不同不能混在一个变量向量里。层级决策者典型控制变量目标函数形式中层价格策略综合能源服务商24小时零售电价向量售电收入减购电成本与设备维护成本最大中层运行计划综合能源服务商储能充放电功率、向上游购电功率配合价格策略实现收益最大化下层用户响应终端用户可转移负荷投入量、可削减负荷比例用电效用最高或购电成本最小在主从博弈的求解框架中中层既决定向用户卖电的价格又决定自身储能、购电等运行方案。为了让模型不至于过于复杂通常把价格变量作为主从博弈的主变量把运行调度问题作为给定价格后的次变量处理。3.2 上层目标为什么要“嵌套”下层响应如果用一句话描述这个模型的核心公式可以写成上层收益 函数(上层价格策略, 下层最优负荷响应)其中下层负荷响应并不是随便给出的而是在每一个给定的价格策略下通过求解下层最小化购电成本问题得到的。也就是说上层粒子每评价一次就必须调用一次下层优化。这个嵌套关系决定了整个程序的结构不会是单层粒子群而是“PSO外循环内层优化器”的组合结构。很多初学者容易犯的一个错误是把价格策略和用户负荷都放进同一个粒子群变量里去优化。这样做的结果一般是粒子群发现如果把负荷变量调得很低用户侧成本会显得很低但售电侧收益也会被压低。由于两者目标冲突粒子群很难找到让各方都能接受的解。而且这种做法会丢失“用户会理性响应价格”这个现实约束仿真结果经济上不可解释。3.3 用户侧约束与储能约束怎么处理下层用户响应模型我一般会设置三类约束。第一类是功率平衡约束用户某个时段的用电量等于基础刚性负荷、可转移负荷投入量、光伏出力和电网购入电量之间要满足平衡关系。第二类是储能运行约束包括充放电功率的上限、荷电状态SOC的动态更新约束、以及SOC上下限保护。第三类是用户舒适度边界比如可转移负荷总量一天必须满足特定次数不能为了省电把所有负荷全都挪走否则用户效用过低会导致模型失真。储能约束里我特别强调SOC更新表达式的写法SOC(t1) SOC(t) 充电功率*充电效率 - 放电功率/放电效率这段式子要按时间步长迭代展开最后会形成一组包含耦合变量的线性约束。在Matlab里建议用稀疏矩阵去表达这组约束否则24个时段可能还好如果扩展到8760小时用稠密矩阵会严重拖慢内层求解速度。3.4 主从博弈求解的收敛含义在主从博弈里“收敛”意味着什么需要先说清楚。理论上如果模型存在Stackelberg均衡那么在均衡点处中层服务商单方面改变零售价格策略无法提高自身收益用户单方面改变负荷响应计划也无法降低自身购电成本。由于整个模型是非线性双层规划很难保证一定有唯一的全局均衡所以工程上通常采取的策略是用粒子群在多个随机初始种群下多次求解观察最优价格向量是否稳定在一个小范围内。在实际Matlab调试中我会额外设置一个“稳定性检验”环节取粒子群输出的最优价格向量把它固定住然后从不同初始点重新求解下层优化看得到的负荷曲线是不是同一个结果。如果下层优化器因非凸问题得到不同结果那整个上层收益计算的可信度就会大打折扣需要先处理下层优化而不是继续调上层粒子群。4. Matlab代码实现与实操过程4.1 程序模块到底应该怎么划分这套模型的Matlab代码我建议按五个模块来组织不要全部堆在一个脚本里参数初始化模块设置各主体的成本参数、储能参数、负荷参数、粒子群参数粒子群主循环模块负责粒子速度位置更新、越界处理、适应度排序适应度计算模块把每一个粒子的价格向量传给下层求解函数回收用户负荷曲线后计算上层收益下层用户响应优化模块求解给定价格下的用户最优负荷计划结果分析与绘图模块输出最优价格曲线、负荷曲线、储能充放电曲线和收益构成。这里最重要的划分边界是粒子群主循环只负责和上层价格策略打交道它绝不能直接去调整用户侧负荷变量。下层用户响应优化模块则只接收当前价格向量返回用户侧最优结果。两个模块通过输入输出参数松耦合调试时可以单独测试下层优化是否正确。4.2 粒子群参数设置的经验取值范围粒子群参数不是越大越好也不是越小越好。对于典型的24小时价格优化问题我常用的参数组合如下参数推荐范围备注种群规模30到50粒子太少容易早熟太多增长计算负担最大迭代次数100到300每次迭代要跑Np次内层优化不宜盲目调大惯性权重w0.9线性递减到0.4前期全局搜索强后期局部精细搜索学习因子c1, c21.49到2.0c1太大粒子只顾自己c2太大过早向全局集中速度限制变量范围的10%到20%防止粒子在解空间里“飞”出边界过多在实际调试里我习惯把惯性权重写成线性递减的形式每代更新一次而不是全程固定。刚才说的c1、c2我一般都用1.49附近的值这样的参数组合虽然多了一套衰减逻辑但整体收敛速度比固定权重更平稳。如果发现结果反复震荡我首先会检查速度限制是不是设得太松了尤其是电价变量这种有明确经济含义的变量。4.3 粒子群主循环的核心代码框架下面给出一个可以直接扩展的粒子群主循环框架。这里假设外部已经通过函数get_price_lb_ub获得每个时段电价的上下界内层函数solve_user_response会返回用户侧最优负荷以及对应的上层收益。rng(2024); % 固定随机种子保证结果可复现 % 基础参数 Np 40; % 粒子数量 MaxIt 100; % 迭代次数 dim 24; % 价格变量维度这里默认24小时 c1 1.49; c2 1.49; w_max 0.9; w_min 0.4; % 获取电价上下界 lb_price get_price_lb(); ub_price get_price_ub(); range_price ub_price - lb_price; % 初始化粒子位置和速度 X repmat(lb_price, Np, 1) rand(Np, dim) .* repmat(range_price, Np, 1); V zeros(Np, dim); Vmax 0.15 * range_price; pbest_X X; % 个体历史最优位置 pbest_F inf(Np, 1); % 个体最优适应度 gbest_X zeros(1, dim); gbest_F inf; % 主循环 for it 1:MaxIt w w_max - (w_max - w_min) * it / MaxIt; for i 1:Np price X(i, :); [user_load, ~] solve_user_response(price); profit calc_upper_profit(price, user_load); fit -profit; % 因为粒子群默认求最小收益取负 if fit pbest_F(i) pbest_F(i) fit; pbest_X(i, :) X(i, :); end if fit gbest_F gbest_F fit; gbest_X X(i, :); end end % 更新速度与位置 r1 rand(Np, dim); r2 rand(Np, dim); V w * V c1 * r1 .* (pbest_X - X) c2 * r2 .* (gbest_X - X); % 限制速度 V max(V, -Vmax); V min(V, Vmax); X X V; % 边界约束处理直接吸收到可行边界 X max(X, repmat(lb_price, Np, 1)); X min(X, repmat(ub_price, Np, 1)); % 输出迭代过程 fprintf(迭代次数 %03d, 当前最优收益 %.4f\n, it, -gbest_F); end fprintf(最优价格策略\n); disp(gbest_X);这段代码结构上把粒子群和上层定价问题耦合在一起但没有把下层的约束细节塞进主循环方便后续维护。外层粒子群的适应度函数里把profit取了负号因为标准的粒子群更新逻辑都是朝目标函数最小的方向走收益最大就等价于负收益最小。这里有个容易忽略的细节在我实际调试的模型里上层收益的数量级可能在几万到几十万之间而价格变量的数量级只有零点几到一点几。粒子群速度更新公式是用当前粒子位置和最优位置的差值来驱动如果速度限制设成与价格变量同尺度粒子一步就可能跨过整个决策区间。所以Vmax不能太随意我一般取价格变量范围的0.1倍。4.4 下层用户响应优化怎么求下层用户响应在实际中可能是简单的价格弹性负荷也可能是带储能的复杂模型。我在这里给出一个带储能和可削减负荷的简化版本供参考。下层问题的目标是最小化用户购电成本同时要维持合理的负荷曲线。变量至少包括各时段购电功率、各时段储能充电功率、放电功率和可削减负荷量。如果完全不考虑整数变量用linprog就能求解。比如把储能充放电、可削减负荷等都近似成连续变量只是加一个同一时刻充放电不能同时进行的逻辑约束。在Matlab的linprog框架里这一步表达成两个变量之和不超过一个极大值就可以了。但有些模型里可转移负荷会涉及“开/关”状态这种就需要用intlinprog。如果模型里必须用整数变量而你又用linprog优化器会直接报错这也是新手常遇到问题的地方。我个人的建议是第一版先不要急着把可转移负荷建成整数变量。可以把可转移负荷建模成“分时段的连续替代量”也就是允许它部分转移并设置一个总用电量约束和按小时可以削减的倍率上下限。这样下层问题保持线性规划算法稳定且跑得快。后面如果需要更细致的行为经济学建模可以再替换成intlinprog版本但那时候要对求解耗时做好心理准备。4.5 内层优化函数的输入输出设计内层优化函数我会写成类似下面的形式function [user_load, total_cost, exitflag] solve_user_response(price, params) % price: 1x24 当前零售电价 % user_load: 1x24 用户从电网购入的功率曲线 % total_cost: 用户购电总成本 % exitflag: 内层优化器返回的状态 n 24; x_init zeros(4 * n, 1); % 这里假设变量顺序为 [购电功率, 储能充电, 储能放电, 可削减负荷] Aieq []; bieq []; Aeq []; beq []; lb zeros(4 * n, 1); ub []; [fval, ~, exitflag] ... linprog(cost_coeff, Aieq, bieq, Aeq, beq, lb, ub, options); if exitflag 0 warning(内层优化未收敛exitflag%d, exitflag); end end这个函数每次被外层粒子群调用时都会以当前价格向量为基础重新构建目标系数和约束。这里最影响性能的是重复构建大型稀疏矩阵的操作建议把与价格无关的约束矩阵提前算好只把目标系数和常数项在每次调用时更新。另外如果循环里频繁调用linprog建议预先设置好options特别是关闭迭代显示否则命令窗口会被刷屏刷到没法看。4.6 结果输出与合理性验证程序跑完以后不要只盯着最终收益数字。我会把这几张图画出来看最优零售电价曲线和用户原始负荷曲线的叠加图如果价格高峰时段对应的是负荷低谷那收益模型大概率有问题储能充放电功率曲线正常情况下储能应该表现为低价时段充电、高价时段放电迭代收敛曲线粒子群全局最优收益应当随迭代次数上升并逐渐平稳而不是断崖式跳来跳去。如果发现储能行为不对优先检查储能参数里的效率设置。一个常见的现象是充放电效率过低导致储能没有经济性可言所以优化结果里储能几乎不动作。这时候模型不是错了而是经济参数本身让储能没有获利空间。这种情况要区分清楚是参数导致的不运行还是算法没搜到储能的可行解。5. 常见问题与排查技巧实录5.1 一眼定位问题在先检查方向这部分是我反复调试后总结的基本覆盖了初学者会遇到的大部分情况表现可能原因处理方向内层优化报无解或不可行约束过紧或上下限取值不合理先放宽储能SOC终值约束查看内层单独求解是否可行粒子群收益曲线前期乱跳后期立刻平稳粒子过早集中到局部最优增大种群数量或调大惯性权重上限最优价格永远卡在边界上电价上下限边界不合理或者收益随价格单调变化扩大电价边界检查目标函数符号是否写反收益迭代曲线很平滑但结果不具备解释性粒子群是在“背数字”而不是在找规律检查上层收益是否对用户响应存在过度依赖内层是否稳定多次运行结果差异大没有固定随机种子或下层优化存在多解固定rng并让内层优化多次随机重启程序非常慢内层优化器被频繁调用且每次重建约束矩阵提前抽取常数矩阵关掉linprog迭代显示必要时并行计算适应度5.2 惩罚函数里最容易踩的坑很多主从博弈模型在实现时会在上层目标里加入惩罚项用来处理某些非严格约束。粒子群在搜索时会尝试大量极端价格组合比如把某个时段的电价抬到很高这会导致用户侧购电量为0或者非常小在没有惩罚的情况下上层收益反而会虚高。我遇到过的典型情况是价格上限没有约束好粒子群发现把白天所有时段电价都推高、夜间降低总收益比一个平滑报价策略更高。从纯数学上看没有错但经济上不合理因为运行商不可能给出明显背离供需走势的价格。解决办法不是单纯靠惩罚系数死压而是先给价格变量设置一个较紧的上下界让粒子群只能在可能被市场接受的报价范围内搜索。如果确实需要加入惩罚项惩罚系数不能拍脑袋定。我习惯先跑一次无惩罚的模型记录上层收益的数量级再把惩罚项系数设置为收益数量级的0.1到1倍之间。否则会出现两个问题惩罚系数太小约束根本不起作用惩罚系数太大粒子群的搜索方向完全被约束项牵着走失去优化目的。5.3 为什么模型总卡在局部最优主从博弈的目标函数并不是一个光滑凸函数粒子群即使设置得当也只能做到“在多次随机重启下得到稳定最优”无法在数学上保证全局最优。所以每完成一轮优化我会用不同的随机种子重新跑几遍。如果几次跑下来的最优价格曲线大致相近只有个别时段有毛刺可以取多次结果里收益最高的一组再对最优解附近的邻域做一次精细搜索。具体做法是把最优价格向量叠加一个很小的随机扰动生成若干新粒子继续迭代二三十代往往能进一步微调结果。这种做法也被很多论文称为局部精搜或者重启策略。相比之下单纯把迭代次数调大并不一定有效因为粒子后期如果没有多样性会在同一个区域反复绕圈。5.4 固定随机种子到底是好事还是坏事在调试期间固定随机种子很有必要。如果不固定粒子群的初期种群每次都不同你很难判断算法改动到底是改善了性能还是只是随机运气好。我会在调试阶段用固定的随机种子跑通流程确认每一处逻辑调整都符合预期比如收益提升、收敛曲线变化等一系列结果都对得上。但最终提交结论或者写论文的时候就不能只依赖固定的那一次随机结果了。更稳妥的做法是设置一组不同的随机种子比如从2021到2030这10个值分别跑完整模型然后取最优值、平均值和方差。只报单次运行结果审稿人或者导师很容易怀疑程序的稳定性。5.5 内层优化无解时先不要怀疑粒子群这部分我想单独拎出来强调。每次外层粒子群生成一个新的价格向量下层都会重新求解一次用户响应。如果下层约束里SOC初值和末值设定得过于严格某些价格组合下就会造成下层优化无解最终导致粒子群报错或收益被错误地赋成一个很大的惩罚值。遇到这种情况首先应该将某个极端价格向量单独提出来直接调用下层求解函数看linprog返回的exitflag是多少。如果单独调用还是无解说明下层问题本身约束矛盾。例如储能初始SOC为0.5、要求结束时SOC也是0.5但充电功率上限过低且负荷曲线又要求在末端强制放电就可能推不出可行解。先把这些约束松一松确认内层在各种边界价格下都能返回合理最优解后再回头调粒子群。6. 进阶扩展方向与复盘建议6.1 从单目标粒子群到多目标粒子群算法这套主从博弈框架目前是在单目标意义下做的也就是中层服务商收益最大化。如果后面要同时考虑用户用能满意度、碳排放总量或者系统峰谷差单目标粒子群就不够了需要切换成多目标粒子群优化算法。多目标粒子群的底层速度更新逻辑和单目标一致差别主要在于两个地方一是适应度变成了多个目标构成的向量粒子之间的优劣要通过Pareto支配关系来判断二是维护一个外部档案集用来保存当前已经找到的非支配解。如果只在目标函数里把多个指标加权成一个总指标然后用单目标粒子群求解这种做法虽然简单但权重系数很难选。不同权重会得到完全不同的报价策略而且权重的微小变化可能导致最优解跳变。改用多目标粒子群以后得到的是一整条Pareto前沿收益和碳排或者用户满意度的权衡关系可以直接在图上看到不用反复调权重。6.2 从确定性场景到随机鲁棒场景现实中光伏出力和用户基础负荷都有很强的不确定性。当前很多主从博弈模型已经往前跨了一步把外层的目标函数改写成多个典型场景下的期望收益。算法结构上变化不大外层粒子群每评价一个粒子就对多个随机场景分别调用下层优化最后求期望收益或者最差场景收益。场景数增加以后计算量成倍增长。这时候可以利用粒子群天然并行的特点用parfor把不同粒子的适应度计算分配到多个worker上能明显缩短运行时间。我在扩展阶段就把光伏出力的随机场景做成了24×N个矩阵每次粒子群迭代时同一个粒子对应多个场景曲线内层求解函数里用循环或批量linprog来处理。只要矩阵预先分配好而不是在循环里不停动态增加行数效率是可以接受的。6.3 用户行为模型可以做得更贴近现实如果论文后面想强化创新性可以考虑把用户的响应机制从“完全理性经济主体”改成更贴近实际的行为模型比如参考消费者心理学中的价格阈值概念。用户不会因为电价比上一时刻贵了0.01元就立刻转移负荷而是会有一个触发阈值认为价格变化幅度达到一定程度才会响应。这种模型会让下层优化变成更复杂的非线性问题可能要用到fmincon或者外部近似求解。但要注意一旦下层模型非线性程度上升粒子群外层嵌套的求解稳定性也会受到冲击调试前最好先把下层优化单独测稳。6.4 复盘做这一整套模型我最大的体会如果让我重新做一遍这套代码我会先把主体之间的逻辑用一张简单的数据结构图画清楚再动手写粒子群。很多人拿到题目后第一件事就是找粒子群代码模板结果接进程序以后才发现主从关系和决策变量边界没定义清楚最后只能反复返工。其实粒子群实现本身五分钟能写完真正决定模型能不能出结果的是内层优化和上层收益计算之间的接口设计。另外参数敏感性分析一定要做。电价上下界、储能初始SOC、内层优化器终值容差这些参数对最终结果的影响远大于粒子群自身参数。每次拿到一组新结果我会随手把价格上下界放宽或收窄20%观察最优解是否变化剧烈。如果变化太剧烈说明模型边界条件设置有问题不能靠调粒子群掩盖。如果这篇内容能帮你少踩几个坑那我花在这些调试经验上的时间也算没白费。希望你的粒子群算法和主从博弈模型也能跑出一个既收敛合理、经济上又说得通的均衡结果。