多微电网、共享储能、博弈调度这几个词放在一起乍一看像三个独立方向的拼盘但真正做过配电网优化的人应该能立刻get到这是一个非常典型的“多主体利益博弈”场景。这几年分布式光伏、风电在配电网层面大规模接入负荷的峰谷差也在拉大单一微电网自己配储能容量小了不够用容量大了闲置严重投资回报率很难看。而共享储能这种模式——相当于把储能站从“私有财产”变成“公共资源”谁缺电谁租、谁有余电谁存——就成了一个很自然的解决思路。但问题也随之而来多个微电网都要用同一个储能站谁优先充、谁优先放、租金怎么定、配电网的网损由谁承担这些矛盾靠单独优化某一个主体是解决不了的必须用博弈的思想来做调度。这篇总结我把这套建模和Matlab实现思路完整梳理一遍从博弈模型设计、KKT条件转换到YALMIP求解器的代码框架再到算例结果分析和常见的坑尽量写得能直接照着复现而不是给一个花架子流程图。1. 问题背景与整体设计思路1.1 为什么单个微电网“自建储能”不划算先看一个典型场景。某园区配电网下挂了3个微电网每个微电网都装了分布式光伏容量大概在500kW到800kW之间。白天光照好的时候各微电网内部光伏出力可能超过自身负荷产生倒送功率到了傍晚负荷高峰光伏出力又快速跌到零微电网不得不从上级配电网大量购电。如果每个微电网各自装一套储能会面临两个尴尬数字按各自峰值负荷的20%配置储能3个微电网加起来要装大概900kW/1800kWh的储能初始投资接近500万但实际运行时由于各微电网的负荷曲线和光伏出力曲线不完全重合很多储能容量一天里根本用不满等效利用小时数极低。加上储能电池的循环寿命是有限的每多充放一次都是成本利用率越低折算到每度电上的成本就越高。共享储能的意义就在这里用一个集中式储能站比如600kW/1200kWh替代三个分散的小储能总容量降下来了但通过错峰调度实际效果能覆盖住三个微电网的调峰需求。储能利用率上去了单位成本自然下来。不过这个“错峰调度”不是拍脑袋排出来的而是要通过博弈让每个微电网主动调整自己的充放电计划。1.2 多主体博弈到底在博什么在这个系统里至少有三类利益主体配电网运营商DSO关注的是整个配电网的网损、电压质量、以及联络线功率峰值。它手里有电网拓扑信息和跟上级电网的购电电价可以通过调整内部电价信号、或者收取网损费用来影响微电网的行为。多个微电网MG每个微电网有自己的负荷和分布式电源运行目标是让自身一天的总运行成本最小包括向上级购电费用、向共享储能租赁的费用、以及自己分布式电源的运维成本。共享储能运营商SESO建了一个集中储能站通过向微电网出租充放电容量来赚钱。它的目标是在满足自己储能安全约束的前提下最大化租赁收益。这三个主体的目标函数是互相冲突的微电网希望租金越低越好储能运营商希望租金越高越好微电网希望多从配电网买便宜电配电网希望联络线功率不要太冲。这种多主体冲突没法用一个统一的全局优化目标来表达所以需要用博弈论来建模。标题里的“博弈优化调度”最常用、也最适合配电网场景的是主从博弈Stackelberg game配电网运营商或共享储能运营商作为领导者先发布价格信号多个微电网作为跟随者在这个价格下各自独立优化自己的用电计划。领导者根据跟随者的响应再调整价格迭代到双方都不愿意单方面改变策略就得到了Stackelberg均衡。1.3 为什么选择这种“先定价格、再定功率”的主从结构我在最初设计这个问题时也纠结过是直接用纳什均衡多个微电网同时博弈、地位对等还是用主从博弈。实际做下来在配电网这种分层结构里主从博弈更符合物理本质。原因是微电网之间通常不存在直接的电力交易通道它们都要通过公共配电网连接。真正掌握网络信息和定价权的是上级配电网或储能运营商它天然处于“领导者”的位置。如果硬把微电网设计成对等的纳什博弈网络约束比如节点电压、馈线容量就不好处理因为网络约束是全局性的不属于任何一个微电网。而主从博弈天然分层领导者可以把网络约束放在自己的上层问题里下层微电网只需要响应价格信号结构非常清晰。从Matlab实现的角度主从博弈也比纳什均衡好求得多。纳什均衡往往需要求解一个变分不等式通用解法不强而主从博弈可以通过KKT条件把下层问题转成上层的约束最终变成一个单层优化问题直接用现成的混合整数线性规划求解器就能解。2. 博弈优化模型的数学建模2.1 上下层决策变量的设计先明确各主体手里有哪些“牌”可以打。领导者共享储能运营商或配电网的决策变量售电/购电价格储能向微电网放电时的售电价格、从微电网充电时的购电价格或者可分配的储能容量每时段给每个微电网分配的充放电功率上限。跟随者每个微电网MG_i的决策变量自身各种分布式电源的出力与配电网的交换功率买电/卖电向共享储能申请的充电功率和放电功率。这里有一个常见的设计选择如果领导者的决策变量是价格那么下层微电网是一个价格接受者求得的解是均衡价格下的最优用电计划如果领导者直接决策容量分配那更像是一个资源分配问题博弈色彩会弱一些。实际文献里常见的是“储能运营商决策租赁价格微电网决策租赁电量”的组合所以我也按这个思路来建模。2.2 微电网下层优化模型每个微电网MG_i的目标函数是自身总成本最小化min ∑(购电费用 储能租赁费用 分布式电源运行成本 - 向电网卖电收益)约束条件包括功率平衡约束光伏出力 储能放电 上级购电 负荷 储能充电 上级卖电联络线功率上下限微电网与配电网的交换功率不能超过线路容量储能充放电功率约束充放电不能同时进行且有最大功率限制储能SOC约束每个时段SOC在安全范围内且一天结束时SOC要恢复初值防止把储能当“无底洞”用。这个下层问题本身是一个线性规划LP因为目标函数和约束都是线性的。要做成博弈的话还需要往目标函数里加入对价格信号的响应项也就是“被领导者定价引导”的部分。2.3 共享储能上层优化模型储能运营商的目标是收益最大化max ∑(售电收入 充电服务费 - 购电成本 - 储能折旧成本)它的约束包括储能本身的物理约束充放电功率上限、SOC动态方程、SOC上下限、日始日终SOC平衡。同时它的价格决策会影响下层微电网的用电行为所以它不是简单的单层优化而是“带约束的下层问题”的领导者问题数学上是一个MPECMathematical Program with Equilibrium Constraints带均衡约束的数学规划问题。这里要特别说一下MPEC是不满足常规凸优化条件的因为下层问题的最优性条件引入了互补松弛约束直接扔给Gurobi是解不了的。标准的处理办法是用KKT条件把下层问题“提升”到上层然后对互补松弛约束做线性化。2.4 KKT条件替换下层问题对于一个线性规划的下层问题它的KKT条件包含三类稳定性条件stationarity拉格朗日函数对决策变量求导等于零原始可行性原问题约束成立对偶可行性 互补松弛对偶变量非负且约束和对偶变量的乘积等于零。把这三类条件全部写成上层问题的约束原来的双层问题就变成了单层问题。步骤很清楚写出下层问题的拉格朗日函数对所有决策变量求偏导得到稳定性条件保留原约束作为原始可行性约束把互补松弛条件用大M法或强对偶定理线性化。互补松弛约束的形式是“λ ⊥ (Ax - b)”也就是λ≥0、Ax-b≥0、λ·(Ax-b)0这三者同时成立。最常用的线性化方法是引入二进制变量和大M参数λ ≤ M·zAx - b ≤ M·(1 - z)其中z是0/1变量M是一个足够大的正数。这个M的取值有讲究取得太小可能切掉可行解取得太大可能导致数值病态一般取该约束物理量最大可能值的10到100倍。2.5 从MPEC到MILP最终交给求解器的形态下层问题如果只是LP经过KKT替换后上层问题变成带互补线性化约束的MILP。整条技术路线总结成一句话上层领导者决策价格下层微电网响应功率用KKT把下层卷进上层用大M把互补松弛摊平最后得到一个可以用商业求解器解的混合整数线性规划。不过这里还有一个备选方案如果下层问题数量多、规模大整体转成MILP之后整数变量会爆掉每个互补约束需要1个0/1变量下层约束一多就成百上千个整数变量求解速度会很慢。这种情况下用分布式的迭代求解更实用上层先给定价格下层并行优化各自问题把结果返回上层上层根据响应更新价格迭代到收敛。这种方法虽然不能保证收敛到全局最优但工程上完全够用而且每一步只需要解一个LP速度优势巨大。后面Matlab实现部分我会讲这两种方式的代码结构。3. Matlab代码整体结构与核心函数设计3.1 工程化文件组织我建议把代码按功能拆成几个文件不要全都堆在一个脚本里。一方面是好调试另一方面是后面换算例场景时不用大改。main_stackelberg.m % 主程序负责数据加载、循环调用、结果汇总 init_system_data.m % 初始化微电网、配电网、储能参数 build_upper_model.m % 构建储能运营商/配电网的上层优化模型 build_lower_model.m % 构建单个微电网的下层优化模型 solve_mpec.m % 用KKT转换求解单层MILP方式一 solve_iterative.m % 迭代式主从博弈求解方式二 plot_results.m % 绘图功率曲线、SOC曲线、收敛曲线、成本对比这种结构下换场景只需要改init_system_data.m里的数据逻辑代码不用动。算例规模大的时候也方便把求解部分改成parfor并行。3.2 数据准备别让数据成为复现的拦路虎算例数据是整个复现过程中最容易翻车的地方。很多论文的代码下载下来跑不通80%的原因是数据没给全比如负荷曲线和光伏出力曲线不匹配、储能参数缺失、电价曲线单位不对。我用的这组典型参数如下参数数值微电网数量3个调度周期24h单位时段1h光伏装机容量各MG800kW / 600kW / 500kW微电网峰值负荷1200kW / 900kW / 700kW共享储能额定容量1200kWh储能额定功率600kW储能SOC下限/上限0.1 / 0.9储能充放电效率0.95 / 0.95分时电价峰/平/谷1.2 / 0.8 / 0.4 元/kWh微电网购电价上级电网分时电价 * 1.05微电网售电价上级电网分时电价 * 0.9负荷曲线我一般用实际采集的典型日数据如果没有可以按“早晚高峰、中午低谷”的形态用正弦函数叠加随机噪声生成。光伏出力用“晴天日”曲线从早上6点开始上升中午12点到1点达到峰值下午5点降到零。注意要让光伏峰值出现在中午低谷电价时段这样才有“午间存电、晚峰放电”的经济动机。3.3 主程序循环逻辑迭代式主从博弈的流程是这样的初始化价格储能租赁电价取分时电价的某个倍数作为初值每个微电网以当前价格独立求解自己的LP得到充电/放电功率需求储能运营商汇总所有微电网的申报功率加上网络约束求解自己的优化问题更新价格检查前后两轮微电网申报功率的偏差是否小于阈值比如1e-4如果满足就停止否则回到步骤2。这个流程非常像拍卖微电网报需求储能运营商调整报价来来回回几次就收敛了。关键点在于步骤3里储能运营商的模型必须包含“微电网的响应函数”也就是前面说的KKT条件。用迭代的方式做就不需要把KKT显式写出来而是靠外部循环来逼近均衡点代码简单很多代价是收敛速度会受初值影响。3.4 YALMIP建模核心代码片段下面给一段建立下层微电网优化模型的YALMIP代码变量声明和约束写法都以可直接运行为目标function [model] build_lower_model(mg_data, price, storage) % mg_data: 微电网负荷、光伏、联络线限值 % price: 储能租赁电价(领导者给定) % storage: 储能参数 % 决策变量 Pbuy sdpvar(24, 1); % 向上级电网购电 Psell sdpvar(24, 1); % 向上级电网售电 Pch sdpvar(24, 1); % 储能充电功率 Pdis sdpvar(24, 1); % 储能放电功率 SOC sdpvar(24, 1); % 储能荷电状态 % 目标函数 Obj sum(Pbuy .* mg_data.grid_buy_price) ... - sum(Psell .* mg_data.grid_sell_price) ... sum((Pch - Pdis) .* price) ... % 租赁费用 sum(Pdis .* storage.degradation_price); % 储能磨损折算 % 约束条件 Constraints []; % 功率平衡光伏购电放电 负荷售电充电 Constraints [Constraints, ... mg_data.pv Pbuy Pdis mg_data.load Psell Pch]; % 联络线约束 Constraints [Constraints, 0 Pbuy mg_data.Plink_max]; Constraints [Constraints, 0 Psell mg_data.Plink_max]; % 储能充放电约束同一时段不能同时充放 Constraints [Constraints, 0 Pch storage.Pmax * (1 - is_charging)]; Constraints [Constraints, 0 Pdis storage.Pmax * is_charging]; % 储能SOC动态 Constraints [Constraints, SOC(1) storage.SOC_init Pch(1)*storage.eta_ch - Pdis(1)/storage.eta_dis]; for t 2:24 Constraints [Constraints, SOC(t) SOC(t-1) Pch(t)*storage.eta_ch - Pdis(t)/storage.eta_dis]; end Constraints [Constraints, storage.SOC_min SOC storage.SOC_max]; Constraints [Constraints, SOC(24) storage.SOC_init]; % 日始日终一致 % 求解 options sdpsettings(solver, gurobi, verbose, 0); optimize(Constraints, Obj, options); model.Pch value(Pch); model.Pdis value(Pdis); model.SOC value(SOC); model.Obj value(Obj);这里面有个容易写错的点充放电同时性的约束。初学者经常写成Pch * Pdis 0这是一个非线性约束求解器会报错或者求解极慢。正确做法是引入二进制变量is_charging把约束线性化。上面的代码里is_charging是一个24维的二进制变量当时段t处于充电状态时等于1放电状态时等于0。3.5 上层储能运营商的MPEC目标上层模型在YALMIP里本质上和下层很像只不过多了一大堆KKT条件。我的建议是先不要手写所有的KKT太容易出错。可以直接用YALMIP的kkt函数来生成下层问题的KKT条件然后append到上层问题的约束里。% 下层微电网问题用优化模型对象表示 lower_opt optimizer(Constraints, Obj, options, {u}, {Pch, Pdis, SOC}); % 生成KKT条件 [KKT_system, details] kkt(lower_opt, {Pch, Pdis}); % 上层问题约束拼接 Constraints_upper [Constraints_upper, KKT_system];虽然YALMIP的kkt函数在旧版本里有一些坑比如有些非线性算子识别不了但在我们这种纯线性模型里是没问题的。生成完KKT之后还需要手动处理互补松弛约束通常是用binvar声明二进制变量辅助大M法这部分建议封装成单独的函数方便复用。3.6 求解器的选择与参数配置这个问题用到的求解器主要是Gurobi或CPLEX因为转成MILP之后整数变量规模不小。如果没有商业求解器也可以先用intlinprog顶着但性能会差不少。我的配置建议options sdpsettings(solver, gurobi, ... gurobi.MIPGap, 0.001, ... gurobi.TimeLimit, 600, ... verbose, 0); % 如果用CPLEX % options sdpsettings(solver, cplex, ... % cplex.mip.tolerances.mipgap, 0.001, ... % verbose, 0);Gurobi求解MILP的默认MIPGap是1e-4对于这种调度问题其实可以放宽到1e-3求解速度快很多而目标函数值差异几乎可以忽略。TimeLimit设置很重要——尤其是在调试阶段如果不加时间限制遇到一个病态模型可能几分钟都算不完严重影响调试效率。我习惯先把TimeLimit设成120秒能算出来再逐步加大。4. 两种求解路线的对比与代码实现说明4.1 单层MILP路线一次求解理论最优把KKT条件代入上层整个问题变成一个单层的MILP好处很明显一次求解就能得到博弈均衡的全局最优解不需要关心初值、收敛速度这些麻烦事。缺点是整数变量多规模大了求解时间呈指数增长。我做过一个3微电网、24时段、每时段300个约束左右的模型转成MILP后大概有800多个二进制变量Gurobi求解时间在3到10分钟之间。这个速度拿来发论文是没问题的但要反复调参时就很痛苦。4.2 迭代路线先求可行再求最优迭代路线代码量少很多逻辑直观每一步只解线性规划LP速度快到飞起。3个微电网并行求解每一步耗时基本在0.5秒以内迭代30次左右就能收敛总共不到20秒。缺点是需要手动设置价格更新规则和收敛判据而且价格更新策略直接影响收敛质量。常用的更新方式是次梯度法price_new price_old alpha * (申报功率总需求 - 储能可提供功率)步长alpha的取值是门手艺活取大了来回震荡取小了收敛太慢。我的经验是先取一个较大的步长快速逼近等目标函数值变化开始变慢时再缩小步长。4.3 两种路线的代码对比例子这是迭代法的核心循环for iter 1:max_iter % 步骤1: 给定价格求解所有微电网下层问题 for k 1:num_mg mg_model(k) build_lower_model(mg_data(k), price, storage); end % 步骤2: 汇总申报功率判断收敛 demand_charge sum([mg_model.Pch], 2); demand_discharge sum([mg_model.Pdis], 2); gap max(abs([demand_charge - demand_charge_prev; ... demand_discharge - demand_discharge_prev])); if gap tol break; end % 步骤3: 更新价格次梯度法 price price alpha * (demand_charge - storage.Pmax * ones(24,1)); price max(price, price_min); price min(price, price_max); demand_charge_prev demand_charge; demand_discharge_prev demand_discharge; end其中price_min和price_max是价格的上下限根据分时电价来确定不能放开了让价格乱跑否则很容易发散。我在实际项目中更推荐一个混合做法先用迭代法快速得到一个近优解作为MILP的初值热启动再切到单层MILP做精细求解。这样既保证了求解质量也避免了MILP冷启动时前几百秒的探索过程。Gurobi的Start参数可以直接传入这个热启动解。5. 算例设计、结果分析与出图技巧5.1 典型算例场景设计为了把博弈效果展示清楚我设置了3个对比场景场景A无储能。每个微电网只从上级电网购电光伏多余电量直接上网卖电。场景B各微电网自建储能。储能容量按各自峰值负荷的20%配置不参与共享。场景C共享储能。3个微电网共同租赁一个集中储能站。从运行结果看场景C的总运行成本最低。关键在于储能利用率这个指标自建储能场景里每个微电网储能的日循环次数只有0.6到0.9次而共享储能站全天循环次数可以达到1.6次以上。电池的度电成本跟循环次数直接挂钩循环次数上去了单位储能成本就降下来了。5.2 结果图表要怎么画才说明问题Matlab绘图的几个关键图件储能SOC曲线3个微电网共享一个储能站时多时段内不同微电网在交替充放电SOC曲线呈现明显的“锯齿形”说明储能被高频利用。而自建储能场景下每条SOC曲线都比较平放电深度也不深。各微电网购电功率曲线共享储能场景下中午光伏大发时微电网会大量充电晚高峰时放电从配电网购电的峰值会被削掉一截。这个削减幅度是论文里最有说服力的指标。博弈收敛曲线画出每一轮迭代后微电网申报总功率的偏差值横轴是迭代次数纵轴是gap可以看到偏差指数级下降。这张图能证明你用的是博弈收敛而不是随便迭代了几次就停了。成本对比柱状图横轴是3个场景纵轴是总运行成本。重点标注出共享储能相比自建储能的成本下降百分比。5.3 从结果中能读出什么不只是“省了多少钱”算例做完之后还需要留意一个更深层的现象博弈均衡点下微电网的成本未必比“只优化自己而不考虑别人”时更低。这跟囚徒困境有点类似——每个微电网个体理性地追求成本最小但所有微电网同时这么做反而可能造成储能资源挤兑大家都不好过。共享储能的博弈调度本质上是用一个价格信号来协调多个理性主体的行为让个体最优和全局最优尽量接近。我在调试时多次观察到如果把租赁价格定得过高微电网宁可多从配电网买电也不用储能储能收益率反而下降价格过低又会出现多个微电网同时抢着租用储能形成新的峰谷矛盾。找到那个微电网愿意用、储能运营商又有利可图的“中间价格区间”是调参时最有意思的部分也是这类问题的核心难点。6. 常见问题与排查技巧实录6.1 算不出来、求解器报错怎么办我自己踩坑最多的是YALMIP的kkt函数报错。最常见的情况是最新版YALMIP对某些变量类型的约束生成KKT条件时会缺少某几行对偶变量导致上层模型约束不完整求解器返回infeasible。排查方法是把KKT条件逐一打印出来跟手写的理论推导对照一遍看看哪条稳定性条件缺失了。另外大M的取值也是个高频坑。我见过有人设M1e6结果Gurobi直接报数值问题。大M的正确取值原则是“刚刚大于该约束在物理上可能的最大取值即可”。比如联络线功率上限是500kW那么M取5000就足够了完全没必要取1e6。M太大会把MILP的LP松弛边界推离真实可行域求解器内部的branch-and-cut算法效率暴跌。6.2 迭代不收敛或结果震荡用次梯度法更新价格最常见的现象是价格在上限和下限之间来回震荡。解决手段有三个方向减小步长alpha这最简单有效对价格做“惯性”处理即新价格 0.7 * 新价格 0.3 * 旧价格能显著抑制震荡不用次梯度改用“最佳响应”更新每次把所有微电网的申报需求汇总后直接求解一个不考虑互补约束的储能优化问题用它的对偶变量作为新价格信号。第三个方法收敛性更好但实现复杂度稍高。如果只是做课程作业或验证方法可行性前两个就够用了。6.3 MATLAB版本和工具箱的坑kkt函数依赖YALMIP的优化模型接口有些旧版本2020年以前对sdpvar和optimizer的交互支持不佳建议直接把YALMIP升级到最新版。另外YALMIP本身是第三方工具箱不是Matlab官方自带的安装时要把路径添加到MATLAB的搜索路径里很多报错“Undefined function or variable sdpvar”就是路径没设置好。求解器方面Gurobi和CPLEX都需要单独安装并配置license。安装完一定记得在MATLAB里执行gurobi_setup或者cplex_setup然后输入yalmiptest验证YALMIP能正常调用求解器。这一步不做好后面所有代码都会卡在“No solver found”这一行上。6.4 常见问题速查表现象可能原因解决办法sdpvar未定义YALMIP路径未添加检查setpath重启MATLAB求解器提示无输出许可证未激活运行yalmiptest诊断模型infeasible大M太小或约束冲突核查KKT条件、放宽SOC约束求解时间过长整数变量过多增加MIPGap到1e-3设TimeLimit迭代结果震荡步长设置不合理缩小alpha加惯性系数SOC曲线不闭合日末SOC约束缺失添加SOC(24)SOC_init功率平衡不满足光伏/负荷曲线时段错位检查单位确认都是24维向量6.5 调试技巧先跑通最小用例再上规模这个过程我踩了很多次坑。一开始直接在24时段、3个微电网的规模上调试出了问题根本不知道是模型逻辑错还是数据错排查效率极低。后来学聪明了先用2个微电网、6个时段的最小用例跑通确认博弈收敛、成本计算结果符合直觉再逐步扩展到24时段。具体做法是把光伏和负荷曲线截取前6个点储能SOC约束也按6个时段设置。跑通之后检查一个关键指标——储能只有在这6个时段里有充放行为且充放电量完全对得上SOC变化。这个验证一旦通过基本可以确定模型逻辑没问题后面的扩展只是数据规模变化。还要强调一个容易被忽略的细节每个微电网的日末SOC约束最好都设为初值。如果不加这个约束优化结果会出现储能一天结束时电量被放空的“偷能量”现象——它在利用日与日之间的边界漏洞虽然单日成本好看但第二天的运行没法继续。我们做调度时边界条件必须闭环。7. 一些个人实操心得这套“多微电网共享储能博弈调度”的Matlab实现这两年从理论到落地我已经跑了不下十来组算例最大的两个感悟是博弈建模的成败不在于算法多复杂而在于下层模型是否真正反映了每个微电网“个体理性”的行为逻辑以及价格信号的设计能不能同时引导个体和整体的利益趋同。只改KKT推导不动价格更新规则模型很容易变成一个花瓶。再说一个实用建议——如果只是做单时段的博弈分析不需要24时段这么重可以用3个典型时段峰、平、谷代表时段来搭建模型计算量小出图也更直观。等把博弈过程彻底跑通了再换24时段精算。最后分享一个小技巧在跑大规模算例前先用profile查看一下代码的耗时分布。很多时候90%的时间都花在求解器内部而不是建模过程本身。如果是迭代法可以试试把各微电网下层的LP放到parfor循环里并行求解在3微电网时效率提升有限但如果算例扩展到10个以上微电网并行求解几乎是必须的。共享储能博弈的魅力就在于微电网越多、负荷互补性越强这套方案的收益就越明显——不过要是你只想跑个demo3个微电网、单台储能站就已经能讲清楚整个故事了。