储能参与现货-调频市场的双层交易决策与Matlab实现

储能参与现货电能量-调频辅助服务市场的双层交易决策研究(Matlab代码实现)

最近在复盘储能参与电力市场的项目，发现不少同学一上来就盯住“现货价差套利”不放，却忽略了调频辅助服务这条同样重要的收益路径。现在很多省份已经允许储能同时参与现货电能量市场和调频辅助服务市场，两个市场的出清逻辑不同、结算规则不同，储能要在两个市场之间分配容量、申报价格，本质上是一个典型的双层决策问题：上层是储能运营商的收益最大化，下层是电力市场的出清优化。用Matlab做双层建模与求解，是我个人比较推荐的一条技术路线，Yalmip + Gurobi组合基本能覆盖大部分研究场景。这篇文章我会把整个思路、模型设计、代码实现框架和踩坑记录都整理出来，适合正在做储能市场决策、调频辅助服务、双层优化相关课题的研究生和工程师参考，拿到思路就能自己动手改。

1. 储能参与现货电能量-调频辅助服务市场的思维起点

1.1 为什么储能要把两个市场放在一起决策

储能参与单一现货市场时，收益模型其实比较简单：低电价时段充电，高电价时段放电，赚的是峰谷价差。但实际项目里，储能设备通常具有快速响应的优势，调频辅助服务市场对这类资源的奖励往往比单纯峰谷套利更可观。问题也随之而来：储能容量有限，如果一部分容量用于调频，另一部分用于电能量套利，两者之间就产生了竞争关系。

更麻烦的是，两个市场不是独立运行的。储能申报电能量市场的充放电曲线，会影响现货出清价格；储能申报的调频容量和调频价格，又会影响调频辅助服务市场的出清结果。反过来，两个市场的出清价格会决定储能的最终收益。这个闭环用传统单层优化描述不了，必须用双层模型：储能作为领导者（Leader）先做容量分配与报价决策，市场出清作为跟随者（Follower）在给定报价下优化出清结果。

我见过不少论文把这两个市场拆分优化，先算现货最优充放电，再算调频最优申报，最后简单相加。这种做法在理论上不严谨，因为调频备用容量的机会成本来自电能量市场，而电能量市场的充放电计划又受到调频约束的影响，两个问题天然耦合。把一个耦合系统拆成两个独立的子问题，计算出来的“最优策略”基本不是全局最优。

1.2 现货电能量市场和调频辅助服务市场的运行逻辑

现货电能量市场目前普遍以日前市场加实时市场为主，储能主要参与日前市场，按节点边际电价（LMP）结算。市场出清的核心任务是安全约束经济调度，也就是在满足负荷、网络安全约束的前提下，让发电成本最小化。储能作为可充电可放电的灵活资源，在出清模型中通过充放电变量参与功率平衡。

调频辅助服务市场的机制各省级有差异，一般有两种：一种按调频容量付费，一种按调频里程付费，更多的地区是“容量+里程”组合。储能申报调频容量后，在调度周期内需要实时跟随AGC指令，所以储能的可放电功率会被调频容量占用，这部分功率不能再参与电能量市场的放电报价。

在双层模型里，两个市场的耦合约束主要体现为：

• 储能总功率约束：充放电功率与调频备用容量不能超过额定功率；
• 储能能量（SOC）约束：调频动作会消耗能量，所以调频容量会影响SOC变化；
• 竞价策略耦合：储能申报电能量价格和调频价格会影响市场出清价格，进而影响收益。

1.3 双层交易决策到底在优化什么

用一句话概括：上层决策储能“怎么用”——每时段申报多少充放电功率、申报什么价格、预留多少调频容量；下层市场出清决定储能“被怎么分配”——在给定申报值的情况下，储能实际上中标多少、获得什么结算价格。

上层目标函数是储能的总收益最大化，包括电能量收益和调频收益，同时扣除电池循环损耗成本或者折旧成本。下层目标函数是市场出清目标，通常是最小化系统总购电成本，或者最大化社会福利。上下层之间的关联通过储能申报的报价曲线和市场中其他机组报价共同参与出清来实现。

这类模型最常见的形式是MPEC（带均衡约束的数学规划），而下层线性规划问题可以通过KKT条件或者强对偶理论等价转换到上层约束中，再用大M法或SOS1处理互补松弛条件，最终转化为一个混合整数线性规划（MILP），Matlab配合Yalmip+Gurobi就能求解。这也是我接下来要重点展开的实现路径。

2. 下层市场出清模型与上层决策目标的建模细节

2.1 下层市场出清的数学表达

做储能参与市场的研究，我建议先把下层市场出清模型写熟练。以一个简化但保留核心机制的现货出清模型为例，目标函数是最小化系统总成本：

目标函数： 系统总成本 = 常规机组发电成本 + 储能放电成本（或报价成本） + 调频容量购买成本 + 调频里程成本

常规机组发电成本可以用二次函数或者分段线性函数近似。为了把下层模型保持为线性规划，通常用分段线性逼近，或者直接用线性报价曲线。储能放电成本可以设为一个很小的正数，代表电池放电的边际损耗。

系统约束至少包括：

• 功率平衡约束：每个时段的负荷等于常规机组出力加储能放电功率减去储能充电功率；
• 机组出力上下限约束；
• 储能充放电功率上限约束，其中放电上限要扣除调频容量；
• 调频容量平衡约束：系统调频需求等于各机组和储能中标调频容量之和；
• 调频容量上限约束；
• 线路潮流约束（如果采用直流潮流模型）。

下层模型用线性规划表示时，变量包括各机组出力、储能充放电功率、储能调频容量、调频里程中标量等。所有约束都是线性的，KKT条件可以严格写出。

2.2 上层储能运营商的收益函数

上层的核心是收益表达式，我实际建模时用的收益公式如下：

总收益 = 电能量市场收益 + 调频容量收益 + 调频里程收益 - 电池充放电损耗成本

其中：

• 电能量市场收益 = 各时段放电功率 × 现货节点电价 - 各时段充电功率 × 现货节点电价；
• 调频容量收益 = 各时段调频中标容量 × 调频容量出清价格；
• 调频里程收益 = 各时段调频里程中标量 × 调频里程出清价格。

注意，现货电价、调频容量价格、调频里程价格在单层优化里都是外生参数，但在双层模型里它们是下层市场出清的对偶变量，是内生变量。这就是双层模型“难”的本质：上层决策变量出现在下层约束中，下层对偶变量又出现在上层目标函数里。

2.3 两个市场之间怎么通过容量约束耦合

容量耦合是最容易写错的地方，我重点说一下。储能系统额定功率为P_max，假设储能参与电能量放电和调频容量申报，那么在任意时段t，必须满足：

P_dch(t) + R(t) <= P_max

其中P_dch(t)是放电功率，R(t)是调频容量。充电功率单独由P_ch(t)表示：

P_ch(t) <= P_max

同时，充电和放电不能同时进行： P_ch(t) + P_dch(t) <= P_max

我最初建模时图省事，只写了P_dch(t) + R(t) <= P_max和P_ch(t) <= P_max，结果出现充电和放电同时为正的荒谬结果。后来加了防同时充放电的约束才正常。在MILP中，这个约束需要引入二进制变量，但在很多学术简化模型里也可以用P_ch(t) + P_dch(t) <= P_max来近似，虽然严格说这个约束在物理意义上不是完全等价于“不能同时充放电”，但在目标函数激励下，同时充放电没有经济意义，所以可以接受。

SOC耦合约束也很关键：

SOC(t+1) = SOC(t) + (P_ch(t) * eta_ch - P_dch(t) / eta_dch - beta * R(t)) * dt

这里的beta表示调频容量对SOC的等效消耗系数。严格说，调频动作是双向的，上下调节会部分抵消，但考虑到调节过程存在能量损耗，简化模型中用一个线性系数折算到SOC变化量上，是论文里很常见的近似。

2.4 为什么必须用双层而不是加权单层

我经常被问到一个问题：既然下层出清是线性规划，能不能直接把两个市场的目标函数加权成一个单层优化？理论上，如果没有上下层之间的策略性交互，单层确实可以做。但储能的市场力体现在它可以调整申报策略来影响出清结果，这时候单层优化就无法表达“储能预测市场出清结果并响应”的逻辑。

举个例子：储能知道夜间申报低价充电会导致自己小额度中标，白天申报略低于边际机组价格可以保证放电优先中标，这种策略性报价行为天生就是Stackelberg博弈，储能是领导者，市场是跟随者。单层优化等价于把市场价格当常数，根本反映不出储能报价对市场出清价格的反馈影响。这就是双层模型的价值所在。

3. 双层模型的KKT转换与求解思路

3.1 KKT条件的工程化处理

下层是线性规划，满足Slater条件，所以最优解等价于其KKT条件。把下层KKT条件加入上层模型，得到MPEC，然后用大M法将互补松弛条件线性化，就可以交给Gurobi求解MILP。

具体来说，下层每个不等式约束都会对应一个对偶变量和一个互补条件。互补条件的标准形式是：

0 <= 对偶变量 ⊥ (约束右端项 - 约束左侧表达式) >= 0

这个“⊥”符号表示两者至少有一个为零。要转化为线性条件，引入二进制变量和大M常数：

对偶变量 <= M * z (约束右端项 - 约束左侧表达式) <= M * (1 - z)

M的取值很关键。取太小会把最优解排除在外，取太大会导致数值不稳定。我的经验是先求解一次不含互补条件的松弛问题，观察对偶变量和约束松弛量的取值范围，再据此设置M，一般取10^3到10^4之间比较稳妥。

3.2 Yalmip建模与Gurobi求解的代码框架

我在Matlab里用Yalmip建模，根据是否引入了二进制变量分别调用gurobi或求解器处理MILP。下面的代码框架是我调试过的版本，逻辑完整，可以直接在此基础上改算例。

%% 参数定义 % 系统参数 T = 24; % 调度时段数 load_demand = [40;35;...]; % 各时段负荷，24x1 gen_max = [80;50;...]; % 常规机组容量上限 gen_cost = [25;30;...]; % 常规机组边际成本 % 储能参数 Pmax = 20; % 额定功率 MW Emax = 100; % 额定容量 MWh eta_ch = 0.95; % 充电效率 eta_dch = 0.95; % 放电效率 SOC_init = 0.2 * Emax; % 初始SOC SOC_min = 0.1 * Emax; SOC_max = 0.9 * Emax; %% 变量定义 % 上层决策变量 P_ch = sdpvar(T,1); % 充电功率 P_dch = sdpvar(T,1); % 放电功率 R = sdpvar(T,1); % 调频容量申报 % 下层变量 P_gen = sdpvar(T,size(gen_cost,1)); % 机组出力 R_freq = sdpvar(T,1); % 系统调频中标变量（简化为单资源） % 对偶变量由Yalmip自动生成，通常在KKT转换时手动写出

严格的双层MPEC求解比较繁琐，我建议分两步：第一步写清楚下层的拉格朗日函数；第二步用KKT替换下层约束。Yalmip里有一个功能是可以用kkt(Constraints,Objective)直接得到KKT系统，但实际用下来发现变量多的时候容易出错，我反而更喜欢手动写出KKT条件，逻辑更清楚，也方便调试。

手动写KKT的一个代码片段如下：

%% 下层KKT条件 % 对下层变量P_gen求导 Stationarity = [...]; % 导数表达式 % 互补条件线性化 M = 10000; % 大M常数 z = binvar(size(ConstraintMatrix,1),1); % 二进制变量 Complementarity = []; for i = 1:num_constraints Complementarity = [Complementarity, dual_var(i) >= 0]; Complementarity = [Complementarity, constraint_residual(i) >= 0]; Complementarity = [Complementarity, dual_var(i) <= M * z(i)]; Complementarity = [Complementarity, constraint_residual(i) <= M * (1 - z(i))]; end

3.3 强对偶法的替代方案

KKT+大M法是经典做法，但有一个明显缺点：引入了大量二进制变量，求解速度会随着节点数和时段数增加而迅速恶化。我在做24时段、96时段的算例时，明显感到求解时间从几秒上升到几十秒甚至几分钟。

一个替代方案是强对偶法。利用强对偶定理，可以用下层原始目标值等于对偶目标值来替代互补松弛条件。强对偶等式不引入二进制变量，模型规模小很多，得到的是一个非线性优化问题，需要精心处理双线性项（对偶变量乘上层变量）。如果上层价格变量只出现在对偶目标中，可以采用线性化的技巧处理。

实际工程中，我还用过一种迭代逼近方法：先固定上层申报，求解下层出清；根据下层出清价格，回代求解上层优化；反复迭代。这个方法不保证收敛到全局最优，但算得快，适合大规模系统初步筛选方案。做学术研究建议还是用KKT+强对偶的严格方法，工程预研可以考虑迭代法。

4. 算例设计与Matlab实现过程

4.1 算例系统的构建

我用的算例是一个简化火电-储能测试系统，包含2台常规火电机组、1台储能、1个聚合负荷节点。这种简化足够验证储能参与两个市场后的行为变化，又不至于让下层出清的KKT条件过于复杂。

系统参数设置如下：

负荷曲线设置为典型的两峰一谷形状，这样储能可以在低谷充电、两个高峰放电，同时在高峰时段预留部分容量参与调频，能明显看出容量分配的变化。

4.2 仿真结果的关键观察点

调好代码后，我最先关注三个输出：

• 各时段储能充放电曲线；
• 各时段调频容量申报量；
• 两个市场分别贡献的收益。

实际上，储能参与调频后，电能量市场的充放电量会明显减少。以我上面参数的系统为例，如果不参与调频，储能会在两个高峰时段各放电约15MWh；参与调频后，因为每个时段都要预留5MW的调频容量，放电能力被压缩，高峰放电量降到了10MWh左右。但总收益反而是增加的，因为调频容量收益补偿了电能量套利减少的损失，总收益提升约15%到20%。这个结论符合直觉：储能响应速度快的天生优势，在调频市场上比单靠峰谷价差更有价值。

我还对比了不同电价场景。如果峰谷价差超过0.6元/kWh，电能量套利收益占比会上升；如果价差小于0.4元/kWh，储能更倾向多分一些容量给调频市场。这说明双层模型的自适应能力是有价值的，固定比例分配容量的经验方法是不可靠的。

4.3 代码运行效率优化

算例规模不大时双层模型很快，但想扩展到IEEE 30节点或者96时段，求解效率就成了硬伤。我常用的优化手段：

一是减少二进制变量。互补条件里的二进制变量数量 = 下层不等式约束数。如果某些不等式约束在最优解中不可能取等号，可以在建模时就剔除，避免白白增加变量。

二是对大M做灵敏度分析。先跑一个初步结果，观察各互补条件实际的对偶变量最大量级，再调整大M值。大M从10^6降到10^3后，Gurobi的求解速度能提升30%以上。

三是对称性消除。如果下层出清模型中存在多个同质机组，可以用聚合机组替代，减少约束数量。

我实际测试过，不加这些优化，24时段模型求解耗时约45秒；做了变量削减和M调整后，降到8秒左右。对需要批量跑蒙特卡洛或灵敏度分析的研究任务，这个提升非常关键。

5. 常见问题与排查技巧实录

5.1 求解器报“无法求解非线性模型”

这个是我最常遇到的问题。用Yalmip建双层模型时，如果直接把下层对偶变量和上层变量相乘，模型就变成了非线性，Gurobi无法直接处理。报错通常提示“Nonlinear constraints detected”。

解决办法是先看目标函数里是否有双线性项。上层收益中的价格（对偶变量）乘以储能功率（上层变量）就是双线性项。如果价格是市场出清价格，这个乘法必然存在。这时要么用强对偶方法消去双线性项，要么在KKT框架下把价格变量对应的互补约束线性化后，再用目标函数中的“收益”替换为线性表达式。

一个实用技巧：如果坚持要用Gurobi求解MILP，那么在目标函数里不要出现对偶变量与上层变量的乘积，考虑用下层最优目标的对偶等价表达式替代收益，这样可以把双线性项消掉。

5.2 模型可行但结果不符合物理直觉

比如储能莫名其妙地在高电价时段充电，或者在低电价时段放电。我排查的顺序是：

1. 检查SOC约束是否写反了充放电方向；
2. 检查P_dch和P_ch是否用了同一个功率上限，导致同时充放电；
3. 检查调频容量R是否错误地减了充电功率上限；
4. 检查对偶变量符号是否与约束方向匹配。

其中对偶变量符号问题最隐蔽。下层约束是“<=”和“>=”两种方向，对应不同符号的对偶变量，一旦写反，求出“价格”可能为负，储能就会反向套利。我的建议是先在无储能报价参与的情况下，用Matpower或者通用优化求解器跑一遍下层出清，检查节点电价是否为正且随负荷升高而升高，确认下层模型正确后再加入上层模型。

5.3 Matlab代码报错排查

代码层最常见的错误集中在Yalmip变量维度不匹配、sdpvar拼接成大矩阵时索引越界、以及kkt函数对包含二进制变量的约束不适用。

我在写代码时习惯每个约束组合都用size()打印维度，确认左端和右端维度一致再append到总约束集。另外，遇到“Index exceeds array bounds”时，第一时间不是看循环变量而是看Yalmip的变量矩阵维度，因为很多表达式展开后维度超出了预期。

关于Yalmip本身的安装，不少同学卡在setup没反应。这个组件下载依赖网络，可以手动下载对应版本放入工具箱路径，再用addpath(genpath(...))加载。Matlab版本太旧也容易出兼容问题，建议至少在R2020a以上跑Yalmip+Gurobi，我自己在R2022b和R2023a环境里都测试过，没有遇到兼容性问题。

5.4 大M取值不当导致的数值问题

大M的取值不能拍脑袋。M太小，模型切掉了真实的最优解，结果总是与原问题不一致；M太大，Gurobi求解时会出现numerical trouble，甚至报“infeasible or unbounded”。

我的排查方法：先写出下层约束的残差范围，比如功率平衡约束的残差可能是机组出力累计值，通常不超过500。对应的大M取5000就足够；调频容量约束残差不超过50，大M取500就行。不要全局统一一个M，最好每个互补条件单独设置一个差别不大的M，能显著改善数值稳定性。

另外，使用SOS1约束替代大M也是常见做法，Gurobi对SOS1支持良好，且不需要手动调M，只是建模时要Yalmip中手动指示变量组成SOS1集合。

6. 项目扩展与经验总结

储能参与两个市场的双层交易决策，我做完一遍的最大感受是：模型结构不难，难在把市场机制转换成数学表达式时不要失真。国内各省电力市场的规则还在快速迭代，比如有的省份现货市场采用节点边际电价出清，调频市场又有性能指标和里程补偿系数，这些都直接影响储能参与策略。做研究时可以先抓主要矛盾，用简化模型验证方法，再逐步叠加地区规则细节。

如果后续想扩展，我有几个建议方向：

• 引入不确定性的两阶段或分布鲁棒优化，考虑新能源出力波动对现货价格的影响；
• 多储能主体博弈，上层从单领导者变成多领导者，模型从MPEC变成EPEC，求解复杂度明显上升；
• 在约束中加入调频里程的动态响应过程，把现在简化的线性折算改成更精细的动态模型；
• 结合储能电池寿命模型，把循环老化成本内生到上层收益函数中，这样储能频繁充放电的行为会被自动抑制，结果也更有工程参考价值。

我个人的习惯是把Matlab代码模块化：数据参数文件、下层模型文件、上层模型文件、求解主文件分开写，算例调整只需要改参数文件。遇到模型改版，比如从确定性改成鲁棒优化，只需要在约束层加场景变量，不需要大规模改动主程序。做储能市场决策方向的课题，建议一开始就按这个结构来，否则反复调试会非常消耗时间。希望这篇拆解能帮你绕开我踩过的坑，直接跑通自己的储能双层交易决策模型。