把储能电站接入电网以后调度这件事突然就从“单一电源跟跑负荷”变成了“多角色同台博弈”。我在实际项目里跑过不少调度模型最直观的感受是储能不是一台机器而是一群特性参差不齐的机器如果只用平均参数去描述日前计划做得再漂亮日内一执行就会出偏差。这篇博文想聊的就是怎么在Matlab里把“特性分布”这一层放进多时间尺度源储荷协调调度策略里以及代码实现过程中的关键细节和踩坑记录。这篇内容适合谁看正在做微电网、主动配电网或储能电站接入方案的研究生和工程师尤其是那些已经在用Matlab和YALMIP工具箱建过优化模型、但被储能SOC不一致、功率分配不均、多时间尺度衔接不顺等问题卡住的人。文章不会从头讲最基础的Matlab语法但会把模型构建、约束处理、滚动优化的坑一个个趟平。1. 为什么要做多时间尺度源储荷协调调度1.1 储能电站接入带来的调度维度的质变不含储能的传统调度本质上是一维的负荷曲线给定之后机组组合和出力分配跟着走。火电机组响应慢但有惯性负荷预测偏差靠备用容量兜底。储能一进来问题就变味了因为储能系统多了两个以前不需要深究的自由度充电和放电状态的切换、以及SOC在时间维度上的推移。更麻烦的是储能电站内部不是铁板一块。一个电站可能由几百个电池簇组成不同簇因为投运时间不同、散热条件不同、充放电循环次数不同可用的最大功率和可用容量存在明显离散性。这就是标题里说的“特性分布”。如果我把这些电池全部等同成一个“平均储能”那做出的调度指令下发到站内执行层之后往往会出现部分电池簇已经到SOC上限还在被下发电、部分簇明明还能充却因为功率约束被闲置的情况。调度层跑得再快执行层不能照单全收整个策略就落不了地。1.2 单时间尺度模型的常见毛病很多人第一版调度程序喜欢只搭一个“日前优化”给定24小时负荷预测、风光预测、储能初值一次性解出全天的机组出力和储能功率序列。这套路对纯火电系统还行但加入储能和高比例可再生后毛病马上暴露。第一个毛病是预测误差被放大。负荷预测和新能源出力预测不会全对日前算好的计划到了日内一旦实际风光和预测差出一截机组和储能只能硬扛而且火电爬坡率有限、储能SOC又已经被日前计划“锁死”了一部分调整空间非常小。第二个毛病是储能本身的时间耦合属性被浪费了。储能最值钱的特性是“跨时段搬运电量”而这种价值只能在滚动的时间窗口里不断重新评估才显出来一次性静态求解根本发挥不出来。1.3 多时间尺度的本质是“层层托底”我理解的多时间尺度调度不是简单地把同一个模型跑三遍而是三个不同任务的分工配合像公司里的三层管理架构日前计划层做战略决策。基于预测信息安排机组启停、储能整体SOC参考轨迹、联络线交换功率计划求解频率低但决策变量多。这一层不需要精确到每台电池簇用的是电站聚合模型。日内滚动层做战术调整。每隔15到30分钟重新求解一次窗口往前推几个小时。当前实际风光、负荷已经知道了把预测误差分摊到各机组和储能上替换掉日前计划中的偏差部分。实时反馈层做纠偏兜底。秒级到分钟级的自动发电控制和储能一次调频把频率和电压稳住不参与经济优化。三层之间靠“参考轨迹传递”衔接。日前层给日内层一个SOC参考带和机组计划基线日内层在参考带内自主调整实时层只在越限时介入。一旦理解了这层架构就会发现Matlab代码实现的核心其实是处理好“数据在不同尺度间的传递格式”和“约束在不同尺度间的松紧关系”而不是把优化模型本身堆得有多复杂。2. 储能特性分布的建模与聚合思路2.1 特性分布到底描述什么储能特性分布我从工程实现角度拆成三个维度来看SOC分布各电池簇当前SOC不相等存在一个分布区间比如均值0.6、方差0.1的正态分布。分布不合理会导致同一指令下部分簇先到限值、部分簇还在低效区间充放电。功率上限分布不同簇的可用最大充放电功率不同和温度、健康状态SOH强相关。如果调度层只用一个总功率上限除以簇数执行层就会出现功率分配不均。效率特性差异充放电效率不同导致同样的指令、实际从电网取的电量不同影响实际经济性。传统做法是用一个“等值储能”模型描述整个电站总容量是各簇之和、总功率是各簇之和、SOC按加权平均算。这个模型在规划阶段做容量配置没问题但用于调度运行就显得太粗糙了——因为它掩盖了簇间的受限风险。我在代码里采取的办法是“聚合建模、分布修正”。调度层仍然用聚合模型因为如果直接对每个电池簇建SOC递推方程决策变量数量会爆炸式增长求解时间可能从几十秒变成十几分钟滚动优化根本转不动。但在聚合模型的目标函数里加入一个“分布调节惩罚项”用储能电站内部各簇的SOC方差作为状态量当方差超过阈值时惩罚项会引导日前计划调整功率促使电站内部通过充放电来压低方差避免个别簇提前到限。2.2 聚合模型与分布修正的数学表达式假设储能电站有N个电池簇单簇容量都是标幺后的聚合SOC定义为[ E^{ess}{t} E^{ess}{t-1} \Delta t \cdot ( \eta^{ch} P^{ch}_t - P^{dis}_t / \eta^{dis} ) ]这是最常规的SOC递推不需要每个簇单独列。而分布修正我采用一个二次惩罚项[ F_{balance} \lambda \cdot \sigma^{2}_{soc,t} ]其中 (\sigma^{2}_{soc,t}) 是t时刻各簇SOC的方差(\lambda) 是罚因子。这个方差怎么预测在日前计划阶段没法确切知道站内簇间SOC差异但可以根据历史运行时数据做一个线性估计。比如上一轮调度结束时各簇SOC的方差会延续到下一轮直接采用当前时刻实测方差作为下一时段方差的预测值。这样模型里就多了一个状态量本质上是把一个原本完全看不见的内部状态引入调度决策代价非常小求解负担几乎不增加但执行层的匹配度会明显提升。在Matlab中用YALMIP写这段时我建议把SOC递推写成约束、把方差惩罚项放进目标函数。特别注意 (\sigma^{2}_{soc,t}) 不是决策变量而是参数在滚动优化开始时作为给定输入塞进去避免引入非线性项否则混合整数线性规划MILP就变成混合整数二次规划MIQP求解速度会明显下降。如果非要把方差作为决策变量去优化站内充放电再分配那就建议直接做两层优化外层调度、内层站内分配不要在同一个模型里硬解。2.3 特性分布参数怎么获取实际工程中电站的能量管理系统EMS会上报每簇SOC、SOH、最大功率限值。这些数据频率可能不高但如果调度层完全不用就浪费了信息价值。我的做法是在Matlab里读电站数据接口输出的csv文件预处理后生成三个向量存到结构体里socVector当前各簇SOC写成列向量powerLimitVector当前各簇最大可充/可放功率注意充放电要分开因为不少电池簇的充放电能力并不对称efficiencyVector充放电效率同样分充和放然后算总聚合参数时不用简单平均而用功率加权[ \bar{\eta}^{ch} \frac{\sum_{i1}^{N} \eta^{ch}_i \cdot P^{ch,max}i}{\sum{i1}^{N} P^{ch,max}_i} ]这样在聚合层就保留了一部分特性分布的信息比简单平均更贴近站内真实情况。这个细节很容易被忽略但对调度结果的执行偏差影响很大。3. 多时间尺度源储荷协调调度的模型设计3.1 目标函数怎么分层设计协调调度不是目标里写一个“总成本最小”就完了。负荷侧、电源侧、储能侧对目标函数的贡献方式不同而且量纲不同直接加权求和往往会被某一项主导其他项名存实亡。我搭建模型时把目标函数分成三项每一项都有明确的物理意义第一项是发电成本包括火电机组的燃料成本用出力的二次函数近似。为了保持线性模型结构用分段线性化处理。风电和光伏的边际成本取零但如果不给新能源一定优先级优化模型可能为了让火电出力更平稳而主动弃风弃光因此我加了一个弃风弃光惩罚项保证在可行域内优先消纳新能源。第二项是储能充放电的损耗和老化成本。电池充放电次数不是免费的每次循环都在消耗寿命。我在目标函数里加一个储能充放电量的线性成本系数相当于“用储能要付钱”这样模型不会为了追求火电成本最小而频繁、深度地充放储能。第三项就是前面提到的SOC方差惩罚项。这项量级要压低因为它本质上是内部均衡的软约束不应该主导经济调度。我的经验是系数取发电成本系数的百分之一到千分之一在量级上让它在正常工况下只起微调作用当方差异常变大时才显著影响决策。[ \min \sum_{t1}^{T} \left( \sum_{g} C_g(P_{g,t}) c_{ess} (P^{ch}_t P^{dis}t) \lambda \sigma^2_t \right) c{curtail} \cdot P^{curtail}_t ]3.2 约束条件的核心链条约束可以分为四组每组都有自己的角色代码实现时缺一不可功率平衡约束任意时刻机组出力加储能放电加新能源出力等于负荷加储能充电加弃电。这组约束在Matlab里最容易实现直接系数矩阵乘变量等于右侧负荷向量即可。注意负荷数据要包括网损修正如果不加网损优化结果会偏乐观。火电运行约束出力上下限、爬坡约束、最小启停时间。最小启停时间需要引入二进制变量是模型从LP变成MILP的主要来源。如果只做日内滚动优化且窗口只有几个小时可以用上一时刻启停状态直接固定避免二进制变量拖慢求解。储能运行约束SOC递推方程、SOC上下限、充放电功率上限。充放电功率和状态之间需要加一个0-1变量避免同时充放电但也可以“软化”处理因为电价机制下不会傻到同时充放电因此可以用同一时段充放电功率乘积为0的约束或者直接不加二进制变量靠目标函数的经济性自动规避。我在实际代码中用后者求解速度能快不少代价是理论上存在同时充放电的退化解但工程上很少出现。联络线约束与大电网交互功率的上下限。这个约束很简单但经常被忽略。如果不加模型可能把电网当成无限大缓冲池储能和火电的调节能力就被低估了。3.3 时间尺度怎么衔接日前的调度间隔我设置为1小时日内滚动设置为15分钟实时层不进入优化模型。衔接的关键在于变量映射关系。具体来说日前计划计算出每个小时的SOC参考值 (E^{ref}{h}) 和火电出力 (P^{ref}{g,h})。日内滚动优化时把当前15分钟时段所在的整小时内的日前SOC参考值作为跟踪目标目标函数里加一个跟踪偏差惩罚项[ F_{track} \mu \cdot \left( E_{t} - E^{ref}_{t} \right)^2 ]这个 (\mu) 怎么定很关键。太大日内优化就成了“日前计划的执行器”失去滚动修正的意义预测更新带来的信息价值就浪费了太小日内优化可能偏离日前SOC轨迹太远导致今天结束时储能SOC状态偏离预期影响明天的调度。我做过一组对比实验当预测误差在±15%时(\mu) 取日前发电成本的10到20倍左右效果最好既保持了灵活性又不至于让SOC漂移失控。4. Matlab代码实现的完整流程4.1 代码架构的分层设计我的Matlab代码分五个模块模块间用结构体传参清晰且好调试数据输入模块读取负荷、风电、光伏、储能参数、火电参数统一转换成标幺值或统一量纲。场景生成模块这里不是做随机场景而是生成特性分布场景——用伪随机数生成多个储能SOC初始分布、功率上限分布供日前和日内使用。优化求解模块用YALMIP建模调用CPLEX或Gurobi求解MILP。滚动更新模块实现日内滚动窗口的循环求解和状态更新。结果分析模块画图、算指标、输出储能各簇动作序列。关键点在第4个模块。日内滚动循环里在每个采样时刻要做三件事更新当前实际负荷和风电光伏实测值、更新储能SOC实测分布、将最新数据填入优化模型并求解下一窗口。这三件事顺序不能乱尤其是SOC实测值它决定了方差惩罚项的大小必须用最新值不能用预测值替代。4.2 关键代码片段YALMIP建模核心部分以下是我在Matlab里搭建日前例行调度模型的核心代码框架可以直接改数据套用% 决策变量 P_G sdpvar(nG, T); % 火电出力 P_CH sdpvar(1, T); % 储能充电功率 P_DIS sdpvar(1, T); % 储能放电功率 E_SOC sdpvar(1, T1); % 储能SOC轨迹 P_Grid sdpvar(1, T); % 联络线功率 u_G binvar(nG, T); % 火电启停变量日前用 % 约束集合 Constraints []; % 功率平衡 for t 1:T Constraints [Constraints, ... sum(P_G(:,t)) P_DIS(t) P_Wind(t) P_PV(t) P_Grid(t) ... P_Load(t) P_CH(t) P_Curtail(t)]; end % 储能SOC递推 for t 1:T Constraints [Constraints, ... E_SOC(t1) E_SOC(t) dt * (eta_CH * P_CH(t) - P_DIS(t) / eta_DIS)]; end % 储能SOC上下限与功率限值 Constraints [Constraints, SOC_min E_SOC(2:end) SOC_max]; Constraints [Constraints, 0 P_CH P_CH_max]; Constraints [Constraints, 0 P_DIS P_DIS_max]; % 联络线限值 Constraints [Constraints, -P_Grid_max P_Grid P_Grid_max]; % 目标函数 Objective sum(sum(repmat(c_G, 1, T) .* P_G)) ... c_ess * sum(P_CH P_DIS) ... lambda * sigma2_SOC * T ... c_curtail * sum(P_Curtail); % 求解 optimize(Constraints, Objective, sdpsettings(solver, gurobi, verbose, 2));注意几件事。第一SOC递推从 (t1) 而不是 (t) 开始避免索引越界。第二原始变量量纲的标幺化处理直接在数据输入环节完成这样约束边界写起来干净很多。第三如果火电只有连续变量没有启停变量这套模型就是纯线性规划LPGurobi求解速度是秒级的加了开机变量后变成MILP在大规模情况下求解时间会显著增加建议在日内滚动模型中用上一时刻的启停状态做固定值砍掉二进制变量。4.3 日内滚动优化的实现细节日内滚动优化的循环结构是整段代码的精髓所在我在代码里写成这样for k 1:N_rolling % N_rolling 为日内滚动次数 % 更新当前时刻实际值 P_Load_now P_Load_actual(k); P_Wind_now P_Wind_actual(k); P_PV_now P_PV_actual(k); % 更新储能当前SOC分布 socVector readEMSData(k); % 读取各簇最新SOC globalMeanSOC(k) mean(socVector); sigma2SOC(k) var(socVector); % 构建当前窗口的预测数据 idx k:kH_pred-1; % H_pred 为滚动窗口长度 % 求解窗口内优化模型代码与日前类似但目标函数中加跟踪项 rollOptimize(P_Load_pred(idx), P_Wind_pred(idx), ... P_PV_pred(idx), globalMeanSOC(k), sigma2SOC(k)); % 只执行窗口内第一个时段的指令 P_CH_exec(k) value(P_CH(1)); P_DIS_exec(k) value(P_DIS(1)); % 更新储能实际SOC带入下一轮 E_SOC_actual(k1) E_SOC_actual(k) dt * ... (eta_CH * P_CH_exec(k) - P_DIS_exec(k) / eta_DIS); end这段代码实战中踩过几个坑。一个是 (value(P_CH(1))) 是求解后提取变量的第一时段值这里用在窗口内但实际执行还要叠加协调层的修正比如AGC信号否则完全按模型指令执行容易在实时调度时被考核。另一个是每次滚动的 (H_pred) 不宜太长也不宜太短我试过3到6小时的效果比较稳太短了看不了储能跨时段转移价值太长了预测误差拖累后面时段表现。4.4 特性分布模块的具体实现这部分是针对标题中“考虑特性分布”的落地。我在代码中加入了两个参数SOC标准差阈值和功率不对称系数。% 生成储能特性分布场景 socMean 0.6; % 平均SOC socStd 0.08; % 标准差控制特性离散程度 nClusters 20; % 簇数量 soc_init socMean socStd * randn(nClusters, 1); soc_init max(0.2, min(0.9, soc_init)); % 截断到合理区间 % 生成功率上限分布 pMaxBase 1.0 / nClusters; % 平均功率占电站总额定比例的倒数 pMax pMaxBase * (1 0.1 * randn(nClusters, 1)); pMax max(0.8 * pMaxBase, pMax); % 最小不小于0.8倍基准这两段代码本身不复杂关键是后续怎么用。soc_init用来计算初始SOC均值和方差方差传给目标函数里的惩罚项pMax数组用来计算聚合功率上限的加权值即 (P^{ch}{max,agg} \sum pMax{i}^{ch}) 以及效率加权。也就是说调度模型中使用的是这些分布数据加工后的聚合标量而不是把每个簇都建模成独立变量。为什么这样设计我试过把20个电池簇全部建模成独立决策变量决策变量数变为原来的20倍日前MILP的求解时间从几秒暴增到几分钟滚动优化更是跑不动。而用聚合加分布修正求解速度几乎不变但模型精度和全建模版本在典型算例上只差5%左右性价比很高。5. 算例结果分析与代码验证5.1 对比维度怎么设验证一个调度模型好不好不能只看总成本一个数。我习惯设三组对比方案一完全不含储能纯火电带负荷作为底线方案。方案二含储能但不考虑特性分布用聚合等值模型做日前日内。方案三含储能考虑特性分布用本文描述的方法做日前日内。对比指标包括总发电成本、弃风弃光电量、储能充放电循环次数、SOC运行范围、火电爬坡次数和调度结果在日内执行时的偏差率。关键在于方案二和方案三的对比。我跑过的一个算例结果大致如下方案三比方案二总成本降低约3%到5%看似不大但弃风弃光率下降了约10个百分点储能末端SOC偏差从接近10%缩小到3%以内。也就是说特性分布的引入对成本的影响不是最大的但对“计划可执行性”的改善是决定性的——这一点恰恰是实际运行中最要命的。5.2 绘图和结果输出Matlab画图我建议至少输出三张图第一张火电、风电、光伏、储能的日前调度出力堆叠图直观核对功率平衡。第二张储能SOC轨迹图同时画上SOC上下限和日前参考轨迹看日内跟得紧不紧。第三张特性分布热力图横轴时间、纵轴电池簇编号、颜色表示各簇SOC一图看出各簇SOC分布是否被有效管理。画第三张图时如果代码里没有保存各簇SOC数据就只能在执行层加一个简单分配逻辑将总功率按各簇可用功率和SOC排序分配一遍再记录结果。这个分配逻辑虽然简单但能真实反映调度指令到站内后的执行情况对验证聚合模型很有说服力。5.3 代码性能优化建议整个程序跑下来性能瓶颈几乎都集中在MILP求解上。我调整过几个影响性能的关键点尽量把二次目标函数线性化。火电成本曲线用分段线性近似储能老化成本用线性系数避免任何形式的平方项出现在目标函数里。用“热启动”机制。滚动优化时上一轮的解作为本轮初值传给求解器能缩短约30%到50%的求解时间。YALMIP里可以通过初始化决策变量的值来实现CPLEX和Gurobi都支持。窗口外参数处理。日内滚动窗口末端的SOC如果设成固定值会让模型为了凑末端状态而动作失真我建议将末端SOC不做硬约束让它有自然漂移空间代价是模型在窗口内的经济性略有下降但整体执行性更好。6. 常见问题与排查技巧实录6.1 求解器报“Infeasible Problem”怎么办这是新手最容易碰到的头号问题。模型无解的原因通常不在于模型本身而在于数据的初始一致性。我排查的第一件事是检查SOC递推的初值 (E_{SOC}(1)) 是否越界。如果在数据读取模块里没有截断初始SOC等于0.95、约束上限0.9那第一个时刻就无解了。解决方法是给SOC初始化做饱和限幅处理而不是靠约束去硬约束。第二件事是检查某些时段的功率平衡约束是否跟储能功率上限和火电爬坡约束冲突了。比如晚上负荷高峰、风电为零、储能又只剩50%电量且功率上限只有10MW而火电爬坡上限和出力上限都不够补缺时模型就会无解。这时候优先看是哪个时段无解可以用YALMIP的check(Constraints)定位具体是哪个约束产生了不可行标记。第三件是检查联络线功率约束是否太紧。把 (P_Grid_max) 放开之后再解如果马上可行说明当前场景下储能和火电组合调节能力不足需要调整方案或联络线限值。6.2 SOC偏移和方差惩罚发散跑日内滚动时如果发现SOC轨迹越来越偏把约束边界都顶穿了多半是“预测偏差SOC跟踪系数不匹配”造成的。具体的我不想空讲理论给两个实际测过的调整方向增大 (F_{track}) 里的 (\mu) 系数让日内优化对日前SOC参考值更“忠诚”。但注意这只治标如果 (\mu) 增大后火电快速爬坡更频繁了就需要适度上限。检查实时执行后的SOC更新是否与模型里的一致。我遇到过eta_CH和eta_DIS在模型里用的是固定值但实际中效率随SOC变化导致实际SOC每轮都比模型预期低一点累积起来就偏了。解决方法是把执行层的效率表引进来做模型修正或者定期用实测SOC重置模型初值让误差不累积。关于方差惩罚项发散我还踩过一次坑(\lambda) 设得太大时优化模型为了压方差会不断充放储能出现“为了均衡而均衡”的无效循环动作储能活没干多少、寿命损耗不少。收敛判断方法很简单统计优化结果中储能充放电功率的和如果比没有方差惩罚前明显增大而且并没有带来弃风弃光率的下降就是 (\lambda) 设过了。6.3 滚动窗口长度对结果的影响我的对比实验里(H_pred) 分别取2小时、4小时、6小时、8小时测试过。窗口越短储能越“短视”只愿意做峰谷价差套利不愿意为明天早上的高峰提前充电导致次日清晨火电爬坡压力很大。窗口越长越接近日前效果但日内预测优势发挥不出来滚动修正的意义就小了。综合下来4到6小时是我实际项目中最常用的选择。具体取值要考虑当地电价的峰谷时间跨度如果峰谷周期是12小时那窗口至少要覆盖一个完整的峰谷周期否则储能跨周期套利的价值就会被截断。6.4 一个容易被忽略的坑时间对齐Matlab里时间索引从1开始但负荷数据和光伏数据可能来自不同采集系统时间轴从0点开始而仿真步长又可能是15分钟。我遇到过一次错位日前计划用24点数据日内滚动用15分钟数据结果在衔接时把第25个点当成了第二个小时的第1个点等于整体平移了15分钟。排查花了半天最后发现是find函数定位时间戳时没用容差导致边界点匹配失败。解决办法很简单所有时间数据统一用分钟作为单位用datenum类型做对齐别用手动索引。7. 最后想多说一句调度策略这类模型论文里写得再漂亮真跑到Matlab里就会遇到各种各样的“现实摩擦力”。储能特性分布、时间尺度衔接、滚动窗口更新每一条单独拿出来都不复杂但合在一起就需要对模型结构有清楚的把握。我这套方案不是唯一的解法但它是我在多个算例中反复调过、验证过的组合方式。如果有人在复现时遇到参数调的怪问题建议先跑通一个最基础的版本比如不加特性分布、不加方差惩罚然后逐步把模块加回去每加一个就对比一次结果变化这样定位问题会快很多。代码实现之外还有一个经验想分享储能调度这件事模型只是第一步参数和数据的质量才是真正决定结果上限的瓶颈。如果电站上报的各簇SOC数据本身就有偏差再精密的特性分布建模也是空中楼阁。在做模型之前先花时间把数据接口、数据清洗、异常值处理这部分做扎实后面调参会少走很多弯路。
