虚拟电厂调度中的鲁棒优化：考虑光伏与负荷不确定性的工程实践

1. 问题建模为什么虚拟电厂调度必须考虑不确定性最近在做虚拟电厂调度项目把光伏出力和负荷波动的鲁棒优化跑通之后最大的感受是确定性模型在工程里基本只能当教学示例用。真实场景下光伏出力跟着云层跑负荷跟着气温和人流跑如果调度方案不考虑这些波动日内实时调整时很容易出现功率越限或者储能过充过放的问题。先说清楚我们处理的对象。虚拟电厂调度本质上是一个多设备协调问题需要把分布式光伏、储能系统、可调负荷整合起来在满足电网交互功率约束的前提下尽可能降低运行成本。目标函数通常是购电费用加上储能老化成本再加上弃光惩罚。约束条件包括功率平衡、光伏出力上下限、储能SOC递推关系、储能充放电功率限制、与电网交互功率限制。核心难点在于光伏出力预测值和负荷预测值都有误差而且这两个误差还是非独立的——比如多云天光伏骤降时空调负荷往往在高位。如果优化模型把预测值当成确定量处理实际运行时经常要推翻重来所以必须引入鲁棒优化框架。我们用的是盒式不确定集加预算约束的方案。不确定参数是光伏出力和负荷功率它们的标称值来自预测系统输出偏差范围根据历史预测误差的统计分位数设定。预算约束用来控制最坏情况下同时偏离标称值的参数个数避免所有不确定参数同时取到极端值导致方案过度保守。需要说明的是鲁棒优化解决的是最坏情况下的可行性问题而不是期望值优化问题。它牺牲一部分经济性换取调度方案在预测误差范围内一定可行。这个权衡在实际项目中非常关键——电网考核的是功率越限次数和储能SOC越限不是平均成本有多低。2. 不确定集设计从统计特性到鲁棒对等模型2.1 不确定参数的统计建模光伏出力和负荷预测误差通常服从有偏分布不能简单套用正态分布假设。我这边处理方式是取项目所在地过去90天的预测数据和实测数据计算每个调度时段的分位数差。具体做法每个时段分别统计预测误差的5%和95%分位数作为该时段不确定区间的上下界。这个方案比固定比例比如取预测值的±15%更贴合实际因为不同时段光伏出力的波动特性差异很大——中午时段的偏差通常比早晚大。预算约束的取值也有讲究。理论上不确定参数个数是N个预算Γ的取值从0到N0对应确定性模型N对应完全最坏情况。实际项目中Γ取总不确定参数个数的30%~50%比较合理。我们调试下来这个区间内方案既不保守到成本失控也保证了实盘运行中不会频繁越限。2.2 对等模型推导的关键步骤带预算约束的盒式不确定集对应的鲁棒对等模型核心推导思路是引入对偶变量把内层max问题转化为约束。这里直接给出我代码里用的关键公式。对每个含不确定参数的约束原始约束形式为∑(a_i * x_i) ≥ b ∑(u_i * z_i)其中u_i是不确定参数。鲁棒对等转化后等价于∑(a_i * x_i) ≥ b ∑(u_i^min * y_i u_i^max * (z_i - y_i)) Γ * w ∑(r_i)其中w和r_i是对偶变量相关的辅助项。实际编程时不需要手动推导这么多用Gurobi的广义不确定性建模可以自动处理但理解原理对调试很有帮助——你才能判断求解器报的infeasible到底是模型逻辑错还是不确定集参数冲突。2.3 关键代码不确定集定义与约束构建我用的求解器是Gurobi 10.0Python接口。定义不确定参数和添加鲁棒约束的代码大致如下import gurobipy as gp from gurobipy import GRB # 创建模型 m gp.Model(VPP_robust_dispatch) # 不确定参数定义 # 光伏出力预测标称值 pv_forecast [12.5, 14.2, 13.8, 12.1] # 单位MW4个时段示意 # 负荷预测标称值 load_forecast [18.3, 19.5, 20.1, 18.9] # 不确定偏差范围 pv_dev_low [-1.8, -2.1, -2.5, -1.5] pv_dev_up [1.5, 1.9, 2.3, 1.4] load_dev_low [-1.2, -1.5, -1.8, -1.0] load_dev_up [1.3, 1.6, 1.9, 1.1] # 添加不确定参数 pv_unc m.addMVar(4, lb-GRB.INFINITY, ubGRB.INFINITY, namepv_uncertainty) load_unc m.addMVar(4, lb-GRB.INFINITY, ubGRB.INFINITY, nameload_uncertainty) # 设置不确定参数的偏差范围 for t in range(4): m.addConstr(pv_unc[t] pv_dev_low[t]) m.addConstr(pv_unc[t] pv_dev_up[t]) m.addConstr(load_unc[t] load_dev_low[t]) m.addConstr(load_unc[t] load_dev_up[t]) # 预算约束 # 电池储能、可调负荷的决策变量 # p_bat_ch, p_bat_dis, p_grid_buy, p_grid_sell, p_load_shift # 功率平衡约束含不确定参数 for t in range(4): m.addConstr( p_bat_dis[t] pv_forecast[t] pv_unc[t] p_grid_buy[t] load_forecast[t] load_unc[t] p_bat_ch[t] p_load_shift[t], namefpower_balance_{t} )这里有个很重要的工程细节Gurobi处理带不确定参数的线性约束时默认会为每个不确定参数生成最坏情况保护。如果你想手动控制保护强度就要用setPWLObj或直接展开对等形式。实际项目中我更倾向于展开写对等模型一方面求解速度更快另一方面你可以单独输出每个约束的鲁棒保护成本这对向领导汇报方案保守性很有用。2.4 为什么不用随机优化而选鲁棒优化很多同行会问场景法随机优化不是更精细吗我这边的判断标准是调度周期和场景数量决定选择。如果日内滚动调度每15分钟跑一次随机优化需要几百个场景做蒙特卡洛采样单次求解时间轻松超过3分钟很难满足实时性要求。鲁棒优化单次求解就是求解一个确定性的线性规划或混合整数规划即使加上整数变量比如储能充放电状态Gurobi在两三秒内基本都能解出来。对于日内滚动调度的应用场景求解速度是第一优先级鲁棒优化的计算优势非常明显。而且从工程交付角度看鲁棒优化的结果更容易解释——调度员能明确知道这个方案在光伏和负荷偏差不超过±15%或其他设定范围时一定可行。这种确定性保证在项目验收和并网考核中有实际价值。3. 完整代码实现从数据输入到结果输出3.1 储能模型与线性化处理储能是虚拟电厂调度的核心执行单元建模时必须处理好两个关键问题SOC递推关系的时序耦合以及充放电状态互斥。我们用的是以下建模方式SOC递推SOC(t1) SOC(t) η_ch * p_ch(t) * Δt / E_rated - p_dis(t) * Δt / (η_dis * E_rated)其中η_ch取0.95η_dis取0.92Δt在日内调度中取15分钟即0.25小时。充放电互斥通过引入二元变量实现这一块是混合整数规划求解速度的关键。我用的是大M法但M值不要取太大取储能额定功率的1.2倍就足够过大的M值会导致数值病态问题。# 储能系统参数 E_rated 10.0 # 额定容量MWh SOC_init 0.5 # 初始SOC SOC_min 0.1 SOC_max 0.9 eta_ch 0.95 eta_dis 0.92 P_bat_max 5.0 # 最大充放电功率MW # 决策变量 p_ch m.addVars(n_periods, lb0, ubP_bat_max, namep_ch) p_dis m.addVars(n_periods, lb0, ubP_bat_max, namep_dis) u_ch m.addVars(n_periods, vtypeGRB.BINARY, nameu_ch) u_dis m.addVars(n_periods, vtypeGRB.BINARY, nameu_dis) soc m.addVars(n_periods 1, lbSOC_min, ubSOC_max, namesoc) # 充放电互斥 for t in range(n_periods): m.addConstr(p_ch[t] P_bat_max * u_ch[t]) m.addConstr(p_dis[t] P_bat_max * u_dis[t]) m.addConstr(u_ch[t] u_dis[t] 1) # SOC递推注意时序 m.addConstr(soc[0] SOC_init) for t in range(n_periods): m.addConstr( soc[t1] soc[t] eta_ch * p_ch[t] * dt / E_rated - p_dis[t] * dt / (eta_dis * E_rated) )这里踩过一个坑SOC上下限约束和充放电功率约束是耦合的如果SOC达到上限但p_ch仍有机会取正值求解器会为了满足SOC约束而牺牲经济性导致充电功率被压缩。解决办法是在目标函数中给储能老化成本加一个合理的权重让求解器自动权衡。3.2 目标函数与多目标权衡我们的目标函数包括四项向电网购电成本分时电价峰段1.2元/kWh平段0.8元/kWh谷段0.4元/kWh向电网售电收益上网电价固定0.45元/kWh储能老化成本按充放电循环折算每MWh吞吐量对应折旧成本弃光惩罚强制弃光成本设得较高为了避免模型为了省钱而主动弃光代码实现时用系数加权的方式写成一个线性表达式权重系数基于项目财务测算给定。这段代码的核心技巧在于把分时电价通过时段判断加进去# 分时电价数据这里简化为峰平谷三个值 price_buy [0.8, 1.2, 1.2, 0.4, 0.4, 1.2, 0.8, 0.4] # 示例实际按96点曲线 obj_expr 0 for t in range(n_periods): # 购电成本 obj_expr price_buy[t] * p_grid_buy[t] * dt # 售电收益 obj_expr - price_sell * p_grid_sell[t] * dt # 储能老化成本线性近似 obj_expr c_aging * (p_ch[t] p_dis[t]) * dt # 弃光惩罚 obj_expr c_curtail * pv_curtail[t] * dt m.setObjective(obj_expr, GRB.MINIMIZE)3.3 与电网交互功率约束与电网的交互功率约束是实际项目中容易出问题的环节。电网调度要求虚拟电厂作为整体对外呈现可控的功率曲线所以交互功率不能随意波动。我们有两种约束模式一种是设定最大交互功率上限另一种是设定爬坡速率限制。实际项目里用的最多的是爬坡约束防止功率大幅跳变对电网造成冲击# 交互功率平衡含鲁棒不确定参数 # p_import[t] - p_export[t] load_t[t] p_ch[t] - p_dis[t] - pv_t[t] # 注意这里p_import和p_export互斥需要引入二元变量或分段约束 # 爬坡约束 ramp_limit 3.0 # MW/15min for t in range(1, n_periods): m.addConstr( p_import[t] - p_import[t-1] ramp_limit * dt ) m.addConstr( p_import[t-1] - p_import[t] ramp_limit * dt ) # 同样处理p_export的爬坡这里有一个容易忽略的问题如果你定义p_import和p_export是两个独立变量必须在约束里加入互斥条件否则求解器可能同时让两者取正值造成从电网买电的同时又向电网卖电的荒谬结果。加二元变量可以解决但会增加求解时间。如果允许这类情况偶尔出现也可以用罚函数法处理——但我在实际项目中不建议这么做因为电网调度那边看到同时买卖的数据会直接质疑调度策略的合理性后续并网验收很被动。3.4 完整求解流程与结果输出整体求解流程如下读取预测数据、初始化模型、添加不确定参数、构建决策变量与约束、设置目标函数、求解、输出结果。# 模型求解与结果处理 m.optimize() if m.status GRB.OPTIMAL: # 提取结果 soc_result [soc[t].X for t in range(n_periods 1)] p_ch_result [p_ch[t].X for t in range(n_periods)] p_dis_result [p_dis[t].X for t in range(n_periods)] p_grid_result [p_grid_buy[t].X - p_grid_sell[t].X for t in range(n_periods)] # 计算总成本 total_cost m.ObjVal print(f最优调度成本: {total_cost:.2f} 元) # 输出96点调度曲线供仿真验证 import pandas as pd schedule_df pd.DataFrame({ 时段: range(1, n_periods 1), 充电功率MW: p_ch_result, 放电功率MW: p_dis_result, 电网交互MW: p_grid_result, SOC: soc_result[:-1] }) schedule_df.to_csv(vpp_schedule.csv, indexFalse, encodingutf-8-sig) else: print(f求解失败状态码: {m.status})关于结果输出多说一句一定要把SOC曲线导出来人工检查一遍。我们遇到过几次SOC曲线前后跳跃的情况比如明明15分钟步长SOC在相邻时段之间跳变了5%以上这种结果看起来不违反约束但物理上完全不合理——大概率是储能功率单位换算出了问题。4. 实验设计与结果分析鲁棒优化的代价与收益4.1 不同预算水平下的成本对比为了向项目方展示鲁棒优化的保守性代价我通常会跑一组参数扫描Γ取不同值时对比总成本和方案可行性。下表是某典型日光伏高发、负荷晚高峰的测试结果Γ取值购电成本(元)储能循环次数光伏利用率求解时间(s)0确定性286503.298.5%0.830%293203.497.1%1.250%301803.695.2%1.5100%324504.191.3%2.1可以看到Γ从0到100%的变化过程中总成本增加了约13%光伏利用率从98.5%下降到91.3%。这说明完全最坏情况保护的代价是很高的实际项目中不建议取100%30%~50%是经济性和鲁棒性的较优平衡点。储能循环次数增加的原因也很直观为了应对光伏出力最坏情况下降需要提前腾出储能容量充放电频次自然上去了。这个指标对电池寿命影响很大如果项目运营周期是15年循环次数增加10%对全生命周期成本的影响不容忽视。4.2 不同不确定度区间的影响除了Γ不确定度区间的大小对结果影响也很大。我们用历史数据统计出的5%~95%分位数区间比固定±15%区间整体小一些但不同时段差异明显。中午12点到14点这个时段光伏不确定区间是±1.8MW而早晚时段只有±0.6MW。如果把全天统一设成±15%会过度保护——中午的偏差覆盖不足早晚的偏差保护过剩整体成本虚高。所以参数设定必须基于实际数据的分位数不能拍脑袋。多花一天时间整理历史预测误差分布效果比调一个月模型参数都好。4.3 与滚动时域控制的联动效果单时段静态鲁棒优化的意义有限我们把它作为滚动时域控制的优化核心。每15分钟滚动一次每次都把最新预测数据代入用鲁棒优化求解未来4小时的调度计划只执行第一个时段的结果。这样做的好处是鲁棒优化保证了每个滚动窗口内方案的可行性而滚动机制保证了方案能跟随最新的预测信息更新。实测下来在光伏波动剧烈的天气条件下每分钟光伏出力变化超过5%这种方案比确定性模型加事后修正的控制策略减少了约60%的功率越限事件。5. 工程化落地中的常见问题与排查技巧5.1 模型求解失败不可行问题诊断Gurobi报infeasible是最常见的事。我强烈建议用m.computeIIS()定位不可行约束的最小冲突集。析出的结果里如果看到功率平衡约束和SOC递推约束同时出现基本可以确定是初始SOC设置或储能容量约束导致的。if m.status GRB.INFEASIBLE: m.computeIIS() m.write(model.ilp) print(不可行约束已输出到 model.ilp请检查以下约束) for c in m.getConstrs(): if c.IISConstr: print(f - 约束 {c.ConstrName} 不可行)排查思路三步走先检查数据单位是否统一MW和MWh最容易混再检查约束方向符号最后检查SOC初值是否在可行域内。90%的问题都能被这三步覆盖。5.2 求解时间过长的处理策略鲁棒对等模型里有二元变量时求解时间可能从2秒跳到20秒。处理方案是给模型设置MIP gap容忍度。我们一般设MIPGap为1%这意味着求解器找到的最优解与理论最优解的差距在1%以内就停。调度场景下1%的成本差距对实际运行影响很小但求解时间可以大幅缩短。m.setParam(MIPGap, 0.01) m.setParam(TimeLimit, 30) m.setParam(Threads, 8)另一个实用技巧是根据调度时段截断模型未来4小时内的时段用精细模型4小时后的时段用聚合模型比如把储能SOC分辨率从15分钟降为1小时可以显著降低整数变量数量。5.3 数据驱动的参数校准不确定集参数的更新不要一成不变。我们项目里做了一个简易的月度校准机制每个月用上个月的实际预测误差数据重新计算分位数区间和Γ值动态更新不确定集参数。这个机制看起来土但对提升方案经济性帮助很大——不同季节的光伏波动特性差异显著固定参数容易在夏季过度保守或在冬季保护不足。5.4 数值稳定性问题大M法和宽泛边界值容易引发数值问题。我们踩过的一个典型坑是储能SOC被设为连续变量但取值范围跨越0~1而充放电功率变量取值范围是0~5MW两者量级差异导致求解器在判定可行性时出现数值误差。解决方案很简单统一量纲。要么把SOC按MWh为单位建模0~10要么把功率按相对值归一化。推荐用MWh为单位这样递推关系式中不需要额外乘容量系数也更直观。6. 项目经验总结与后续优化方向整个项目从前期的数据清洗、不确定集设计、模型搭建到最后的滚动控制集成前后花了两个月。最终交付的鲁棒优化调度模块在光伏和负荷存在±10%以内预测误差时可以保证功率平衡约束全天不违约同时成本比确定性方案平均高3%~8%这个代价在项目评审时被认为是可以接受的。后续还有几个优化方向我计划继续深入一是把光伏出力的时空相关性建模到不确定集里目前是逐时段独立的不确定参数实际上相邻时段的光伏出力存在明显的自相关性二是考虑把需求侧响应资源的调节成本也纳入优化模型现在可调负荷的调节成本是固定值实际运行中不同用户的响应意愿差异很大差异化定价后会有更好的经济性三是把鲁棒优化和模型预测控制做更深的融合用一个统一的框架处理多时间尺度协调问题而不是现在这样手动级联。这个项目做完最大的体会是鲁棒优化不是一个算法而是一整套工程方法论。不确定集怎么建、预算怎么取、结果怎么评估这些环节都需要结合项目实际数据反复迭代。算法本身解决的是数学问题而工程落地解决的是参数的合理性问题后者往往更花时间和精力。