尧图网页设计 
做微网日前调度这几年我最常被同事问的一句话是明天的光伏预测曲线和负荷预测曲线都在这为啥不直接按它排计划问这个问题的人通常没吃过大亏。预测值只是期望实际光伏出力可能只有预测的六成负荷侧园区变压器突增几百千瓦也很常见。这时候再回头看那份按单一预测值排出来的“最优”计划往往已经和实际运行脱节很远。这篇内容我会围绕“MATLAB代码计及源-荷双重不确定性的虚拟电厂/微网日前随机优化调度”展开完整梳理一套可落地的实现思路。代码层面用Yalmip搭配商业求解器Cplex或Gurobi均可不确定性建模走场景法。适合正在做虚拟电厂、微网经济调度、园区综合能源相关课题的研究生和工程师参考。我会把建模动机、数学规划、代码骨架、结果分析逻辑和排错经验一条龙讲清楚。1. 先看清对手日前调度里源-荷双重不确定性到底难在哪1.1 “源”的不确定性与“荷”的不确定性是两种不同性质的问题源侧主要指分布式光伏和风电。光伏出力受辐照度影响云层一过出力可以短时间下降50%以上。风电更极端爬坡事件本身就有很强的随机性。这类不确定性本质上来自“自然条件不可精确预知”你只能拿到一个概率分布而不是一个确定值。负荷侧不确定性来自用户行为园区里某条生产线临时加单、商场空调负荷跟随人流波动这些都是难以精确建模的。和源侧相比负荷的不确定性还有一个特点它存在周期性规律下的随机扰动比如工作日早晚高峰的基线模式相对固定但具体在那个时段涨多少、什么时候涨很难预测到分钟级精度。1.2 日前决策为什么怕“预测不准”微网和虚拟电厂的日前调度是指在D-1日做出未来24小时的机组启停计划、储能充放电计划、与主网的购售电计划。这个计划一旦定下来机组启停这种离散动作在日内很难快速翻转储能也不能因为实时电价波动就随意改SOC目标值。也就是说日前决策的大部分变量是“先定下来”的日内只能在小范围内做修正。如果只按预测值做单场景确定性优化预测误差会直接传导成两部分损失。一是经济性损失光伏实际出力比预测少微网就得多在主网买高价电或者被迫启动更贵的备用机组反之光伏出力比预测多又可能出现弃光。二是可靠性损失严重误判可能导致日前计划下备用不足日内出现切负荷。所以把“源-荷双重不确定性”纳入建模不是锦上添花而是保证计划可用的基本盘。1.3 一个随手可复现的代价测算假设某个微网日预测光伏出力1000 kWh实际只有600 kWh日前计划中没有考虑这种偏差按峰时电价0.9元/kWh计算需要临时多买400 kWh成本多出约360元。而如果随机优化在日前就多预留了储能空间、安排了一部分燃机出力这部分临时购电可能只有150 kWh。看起来单日省下的钱不多但按一个月30天、多个微网聚合来看差距就是数量级差异。在实际项目中我习惯用一个指标来评估不确定性的“杀伤力”叫随机解价值Value of Stochastic Solution, VSS用确定性模型的期望实际成本减去随机优化解的期望实际成本。VSS越大说明不确定性对调度结果的影响越显著也说明做随机优化越有必要。后面章节我会给出具体计算方法。2. 场景法建模把“不确定性”翻译成优化模型的第一步2.1 为什么选场景法而不是鲁棒优化或机会约束处理不确定性的主流路径有三类这里给出适用性对比方便你在自己课题里做取舍。方法基本原理优点缺点适用场景场景法随机规划用离散场景近似连续概率分布目标为期望成本最小信息利用充分可直接对接任意分布结果解释性强场景多时计算量大概率分布准确性影响结果源荷分布特征明确需要经济性最优方案鲁棒优化在给定不确定集中寻找最坏情况下的最优方案计算相对高效不依赖精确分布结果偏保守经济性偏差对可靠性要求极高极端场景必须兜底机会约束规划允许约束以一定概率违反在经济性和可靠性之间可调概率约束的求解相对复杂分布信息需求不低可靠性和经济性折中工程中常用但建模有难度我的个人倾向是学术研究、方案对比、需要输出可视化结果时优先用场景法。原因是虚拟电厂/微网的实际工程中光伏、负荷都有相对成熟的历史数据可以做概率分布拟合场景法能把这种统计信息用足。鲁棒优化更适合做主变过载校验、极端天气应对这类安全校核而不是日常经济调度。2.2 场景生成分布假设、采样方式和相关性问题场景法的第一步是把“不确定性”变成一组带概率的场景集合。光伏出力的经典做法是对光照辐照度建模特卡洛采样把辐照度通过光伏出力模型换算成功率。负荷预测误差通常假设服从均值为0的正态分布标准差取预测值的5%~10%然后叠加到基准负荷曲线上。采样方式上直接蒙特卡洛也能用但收敛偏慢。我更推荐拉丁超立方采样Latin Hypercube Sampling, LHS它的核心思想是把每个分布变量的累积概率区间等分成N份从每个区间里各采一个样本点这样能用较少的样本覆盖整个输入空间。举个例子用纯蒙特卡洛生成500个场景可能仍然留下采样空洞而LHS用200个场景的覆盖效果通常就更好。需要注意的是源荷相关性。光伏出力在相邻时段高度相关不可能8点出力是0、9点瞬间跳到满发这种相关性需要在采样时处理。常见做法是用Cholesky分解对独立采样得到的标准正态向量做相关变换或者直接基于历史场景做重采样。代码层面我会在生成场景矩阵后加一步相关性校验画出相邻时段散点图确认相关性结构合理这是很多人容易忽略的环节。2.3 场景削减为什么不能直接拿1000个场景进优化模型假设采样1000个场景每个场景24个时段变量规模就被乘了1000倍。连续变量还好说储能充放电状态如果引入二进制变量求解时间会指数级恶化。所以必须做场景削减用尽量少的典型场景逼近原有概率分布。工程上常用两种削减方法一种是k-means聚类把距离近的场景聚成几类用聚类中心代表一类场景概率等于类内样本占比另一种是快速前向选择Fast Forward Selection基于Kantorovich概率距离逐步合并最接近的场景对被广泛应用于电力系统随机规划。k-means速度较快但需要预先指定聚类数快速前向削减迭代过程更细致能同时给出削减后的场景集和对应概率。我的实操建议是把两类方法都用一遍对比削减前后场景集合的均值曲线和标准差曲线。只要均值、方差偏差在5%以内说明削减质量可接受。一般1000个场景削减到10~20个就是典型场景数。这个数量既能刻画不确定性又不会把求解时间拖到不可接受。后文建模环节我会给一个基于概率距离削减的简化MATLAB思路。3. 两阶段随机规划模型目标函数与约束的逐条拆解3.1 两阶段决策结构的直观解释随机规划里有个很关键的“非预期性”概念翻译成大白话就是你今天做日前计划的时候并不知道明天哪个场景会真实发生所以第一阶段变量不依赖未来信息。但到了日内场景逐步揭开实时调整变量可以随场景变化。对应到具体模型第一阶段变量在场景实现前就要确定燃气轮机启停状态、储能充放电状态、与主网购售电的日前计划值、储能各时段容量计划值。第二阶段变量随场景变化实际功率平衡下的燃机出力微调、储能充放电功率微调、切负荷量、弃光弃风量。两阶段划分是整个随机调度模型的骨架。我见过不少初学者把全部变量都建成了“随场景变化”结果模型的最优解等于每个场景单独优化后再加权平均这种解完全不可执行因为它忽略了你必须先做出日前决策这个物理事实。所以建模型前先列一张表哪些变量是第一阶段、哪些是第二阶段写清楚后再动代码。3.2 目标函数期望运行成本最小化目标函数是所有场景下的概率加权运行成本之和。成本项通常包括下面五类。第一与主网交互成本。t时段从主网购电价格为λ_t^buy售电价格为λ_t^sell对应变量为P_t^buy和P_t^sell。售电价格一般低于购电价格否则模型会利用价差无意义地倒买倒卖。第二燃气轮机的燃料成本。这里需要对机组启停状态变量u和出力变量P_g做耦合建模。燃料成本采用二次函数表达式在Yalmip中可以用由分段线性近似或使用二次项配合商用求解器处理Cplex/Gurobi都支持二次目标。第三储能充放电退化成本。电池每充放一度电都有寿命损耗折算成单位吞吐量成本参数C_bat。这个参数可以按电池全寿命周期总吞吐量除以购置成本来估计比如一套1000 kWh储能的更换成本50万总吞吐量5000 MWh则单位吞吐量成本约为0.1元/kWh乘2就是一次满充满放的等效损耗。第四弃光弃风惩罚成本。源侧出力用不完又无法上网时允许弃掉一部分但目标函数里给一个惩罚系数让模型尽量少弃。惩罚系数的取值要高于上网电价否则模型会无理由弃光。第五切负荷惩罚成本。切负荷和弃能相比是更严重的事件惩罚系数应该明显高于购电价和弃能惩罚一般取电价的10~50倍代表停电对用户造成的综合损失。整体的目标函数表达式可以写成所有场景s以概率π_s加权平均的最小化总成本注意第一阶段变量不随s变化第二阶段变量才带下标s。3.3 约束条件功率平衡、设备出力与储能动态功率平衡是每一类调度模型都绕不开的硬约束。它要求在任何场景s的任何一个时段t微网内的负荷、购电、燃机出力、储能放电之和等于光伏出力、储能充电、售电、切负荷和弃能之间的差值满足等式约束。切负荷和弃能在等式约束里有点像“松弛变量”但它们只在极端场景下才会超出允许上限日常场景应保持为0。燃气轮机约束包括出力上下限、爬坡约束、最小启停时间约束。其中爬坡约束在随机模型中要按每个场景分别施加因为同一个机组在不同场景下的相邻时段出力可能相差很大必须分别校验每分钟爬坡能力是否满足。储能约束包括SOC动态更新式、充放电功率上下限、电量上下限以及避免同时充放电的约束。SOC动态更新式是t1时段电量等于t时段电量加上充电功率乘效率再减去放电功率除以放电效率。若不使用二进制状态变量可以通过“充放电功率乘积为0”来处理但这类非线性约束在MILP中并不友好工程上最稳妥是引入两个二进制变量u_ch和u_dis用大M约束保证二者不同时等于1。储能还需要设置调度周期始末电量约束。常见做法是要求周期末SOC回到初始值或限制在某个区间保证调度计划可以在次日持续滚动执行。这里有个技巧末时段SOC约束不要锁死在一个固定值或者硬性等于初始值否则会牺牲大量经济性。我给微网做调度时通常给一个可容忍的小区间比如初始值±5%。3.4 备用约束随机优化特有的“安全垫”随机优化里的备用约束尤其重要。在单一预测场景的确定性模型里备用约束可以用“预测负荷预测光伏误差”来折算但随机优化中任何一个场景都可能触发最坏情况需要在日前计划中预留向上/向下备用容量。最简化的写法是要求任意场景s下系统可用总上调能力燃机未用容量储能可放电可中断负荷主网购电余量不低于该场景预测误差导致的净负荷波动值。这里的净负荷定义为负荷减去新能源出力备用需求同样随场景变化。很多人在模型里把备用写成固定值比如按峰值负荷的10%这样虽然简单但没有充分利用场景信息。随机优化的优势恰恰在于能按每个典型场景的净负荷波动动态计算备用需求。我在代码实现里就是把备用约束纳入第二阶段变量的一致性校验中这样模型的自洽性更强。4. MATLAB工程化实现从场景生成到Yalmip求解的代码骨架4.1 工具选型YalmipCplex/Gurobi是性价比最高的组合MATLAB里做随机优化求解的路径主要有三种直接用优化工具箱的linprog/intlinprog用Yalmip建模再调外部求解器或者用基于问题框架的prob2struct转换。我的经验是建MILP模型首选Yalmip因为它对变量的维度处理非常灵动尤其适合“变量带场景下标”这类随机规划模型。Yalmip不是一个求解器它只是个建模语言层底层求解推荐用Gurobi或Cplex。对学术用户这两家都有免费授权Gurobi的安装集成在MATLAB里非常顺滑。4.2 代码结构按“数据-场景-模型-求解-后处理”五段式组织一套能长期复用的随机优化代码结构清晰比代码技巧更重要。我的目录安排参考VPP_Stochastic/ ├── main.m % 主程序 ├── data/ │ ├── load_data.m % 负荷基准曲线 │ ├── pv_data.m % 光伏容量与预测曲线 │ └── price_data.m % 分时电价 ├── uncertainty/ │ ├── gen_scenarios.m % 场景生成 │ ├── reduce_scenarios.m % 场景削减 │ └── plot_scenarios.m % 场景可视化 ├── model/ │ ├── build_model.m % 构建Yalmip模型 │ └── constraints.m % 约束集合 └── result/ └── postprocess.m % 结果分析与图表这种分离最大的好处是换一个测试系统时只需要改data目录下的数据文件改不确定性建模时只动uncertainty目录模型骨架不用推倒重来。4.3 场景生成与削减的核心代码思路场景生成我采用LHS加相关性处理。简化代码如下% gen_scenarios.m 核心思路 function [PV_scen, L_scen] gen_scenarios(data, Ns) % data: 包含预测值、标准差等 % Ns: 初始场景数 rng(42); % 固定随机种子保证结果可复现 T 24; % 光伏出力对辐照度采样需保证非负 % 假设辐照度预测值为 G_pred(Tx1)标准差为 sigma_g % 采用LHS生成标准正态样本 U lhsdesign(Ns, T); Z norminv(U, 0, 1); % 转换为标准正态 % 考虑相邻时段相关性简单一阶自回归 rho 0.85; % 相邻时段相关系数 for t 2:T Z(:,t) rho * Z(:,t-1) sqrt(1-rho^2) * Z(:,t); end G_scen max(0, G_pred sigma_g .* Z); PV_scen data.eta_pv * data.A_pv * G_scen; % 换算成功率 % 负荷预测值加正态误差 L_scen data.L_pred data.sigma_L .* randn(Ns, T); L_scen max(0, L_scen); end场景削减我用基于概率距离的快速前向选择实现思路如下初始时每个场景概率为1/Ns每次迭代找出概率距离乘积最小的一对场景把其中一个场景的概率合并到另一个上然后删除被合并场景重复直到剩余场景数满足要求。MATLAB代码可以写成循环加距离矩阵更新注意向量化处理以提高速度。1000个场景削减到10个在普通电脑上秒级完成。4.4 模型构建与求解的Yalmip代码骨架接下来是模型部分。我用一个具体的微网算例来讲代码结构。系统参数如下设备参数数值燃气轮机台数1台燃气轮机额定功率500 kW燃气轮机最小出力50 kW燃气轮机爬坡速率100 kW/h储能容量1000 kWh储能最大充/放电功率200 kW储能充/放电效率95%光伏装机容量800 kW负荷日峰值700 kW主网交互购电功率上限500 kW主网交互售电功率上限300 kWYalmip建模骨架如下% build_model.m 核心片段 T 24; S 10; % 削减后典型场景数 % 第一阶段决策变量 P_buy sdpvar(1, T, full); % 日前计划的购电功率 P_sell sdpvar(1, T, full); % 日前计划的售电功率 P_g sdpvar(1, T, full); % 燃机日前计划出力 u_g binvar(1, T, full); % 燃机启停 P_ch sdpvar(1, T, full); % 储能充电计划 P_dis sdpvar(1, T, full); % 储能放电计划 u_ch binvar(1, T, full); % 充电状态 u_dis binvar(1, T, full); % 放电状态 SOC sdpvar(1, T1, full); % 储能电量状态 % 第二阶段决策变量场景相关 P_g_s sdpvar(S, T, full); % 场景下的燃机实际出力 P_dis_s sdpvar(S, T, full); % 场景下的放电功率 cut_s sdpvar(S, T, full); % 切负荷量 dump_s sdpvar(S, T, full); % 弃光弃风量 buy_s sdpvar(S, T, full); % 场景下实际购电 sell_s sdpvar(S, T, full); % 场景下实际售电 % 目标函数 Objective 0; for s 1:S for t 1:T Objective Objective pi_s(s) * (... price_buy(t) * buy_s(s,t) - price_sell(t) * sell_s(s,t) ... c_gas * P_g_s(s,t) c_dg * (P_g_s(s,t) 0) ... c_bat * (P_ch(t)/eta_ch P_dis_s(s,t)*eta_dis) ... c_dump * dump_s(s,t) c_cut * cut_s(s,t)); end end % 注意场景下标放在第二维度会对求解器和Yalmip都更友好 % 实际代码中建议变量维度写成 [S, T] 并用 repmat 操作这里有个小坑需要提醒Yalmip 中如果用binvar(S, T)这种大维度0-1变量矩阵预分配和约束添加时要格外小心千万别在循环里反复调用constraints [constraints, ...]来拼约束数组性能会很差。正确做法是预创建一个元胞数组C {};每次循环把新约束加进去最后用Constraints [C{:}];一次性拼接。这在场景数较大时能节省不少构建时间。约束添加的核心逻辑Constraints {}; % 功率平衡每个场景单独满足 for s 1:S for t 1:T Constraints{end1} ... L_scen(s,t) P_ch(t) P_dis_s(s,t) sell_s(s,t) ... PV_scen(s,t) P_g_s(s,t) buy_s(s,t) P_dis(t) cut_s(s,t) - dump_s(s,t); % 更严谨写法是用等式加RHS按实际系统平衡方式决定 end end注意功率平衡里储能充电和放电的建模方式不同我的习惯是把充电视为负荷、放电视为电源同时加入充放电互斥条件for t 1:T Constraints{end1} P_ch(t) u_ch(t) * P_ch_max; Constraints{end1} P_dis_s(:,t) u_dis(t) * P_dis_max; % 这里用第一阶段状态约束第二阶段功率 Constraints{end1} u_ch(t) u_dis(t) 1; end这种写法把“日前决定储能是充还是放”放在第一阶段场景只影响放电功率的幅度符合物理直觉。运行结果告诉我在大多数场景下模型倾向于在谷时段充电、峰时段放电但光伏大发且负荷低谷时也会有“光伏直接充电”的行为这说明随机场景对储能策略的刻画比确定性模型更细腻。SOC约束SOC(:,1) SOC_init; for t 1:T Constraints{end1} SOC(:,t1) SOC(:,t) ... eta_ch * P_ch(t) - P_dis_s(:,t) / eta_dis; Constraints{end1} SOC_min SOC(:,t1) SOC_max; Constraints{end1} SOC(:,T1) SOC_init - 0.05 * SOC_max; Constraints{end1} SOC(:,T1) SOC_init 0.05 * SOC_max; end这里给末时段SOC留了5%的弹性区间很多论文里写成硬等于实际工程中会让储能为了“还电量”而做不必要的充放经济性反而下降。4.5 求解配置与稳定输出求解器配置我一般这样写ops sdpsettings(solver, gurobi, verbose, 2, ... gurobi.TimeLimit, 600, ... gurobi.MIPGap, 0.001, ... gurobi.NumericFocus, 3); sol optimize(Constraints, Objective, ops); if sol.problem 0 % 求解成功 else disp(sol.info); % 输出失败原因 endNumericFocus3是Gurobi的数值鲁棒性增强选项在含0-1变量和大M参数较多的模型里能显著降低“数值病态导致结果抖动”的概率。代价只是求解时间略微增加对于日前调度这种离线计算场景完全值得。5. 结果怎么看出门道确定性调度与随机调度的对比分析5.1 搭建一个公平的对比基线我在做结果分析时通常把预期结果拆成三组对照分析确定性调度模型DET直接用预测光伏、预测负荷求解。随机优化调度模型STO用削减后的场景集合求解。完美信息下界WS即Wait-and-See假设知道每个场景真实值分别求解后再做概率加权平均。这个模型是理论下界现实中不可能实现但它是评价随机优化质量的重要参考。三组模型跑完后统计各方案的期望实际成本。实际成本计算方法是把STO方案的第一阶段决策固定然后在每个测试场景下做第二阶段优化只优化功率调整量与切负荷/弃能量得到每个场景的真实运行成本再做概率加权和。5.2 随机调度为什么“看起来更贵”从结果表来看STO的期望成本通常比WS高但比DET低。指标确定性调度随机优化调度完美信息下界第一层目标函数值未包含场景信息概率加权期望成本各场景期望成本加权平均期望实际成本计入偏差惩罚3220元2890元2740元切负荷期望量12 kWh2 kWh0 kWh弃光期望量45 kWh18 kWh3 kWhVSS 3220 - 2890-330元约11.4%-这个例子中STO相比DET节约了约10%的实际运行成本主要来自两部分一是减少了高价时段临时购电二是削峰填谷更加精准最小化了切负荷和弃能风险。DET的目标函数数值看起来很漂亮因为它的约束建立在“预测完全准确”的基础上但一旦考虑偏差它的目标值就成了镜花水月。“随机优化的期望成本高于确定性模型”这个现象不少初学者第一次看到会困惑。原因很简单确定性模型只面对一个“平均预测场景”约束容易满足目标值自然低随机模型要同时满足10个差异很大的场景任何方案都必须在不同场景间做折中所以最优期望成本会上升。但这恰恰说明模型在准备应对不确定性。如果哪天随机优化的期望成本和确定性模型完全一样那大概率说明场景削减没有体现不确定性特征或者模型结构出了问题。5.3 SOC、购电计划与场景间的耦合关系除了总体成本我还会画出储能SOC曲线在不同典型场景下的分布带。正常情况下你会看到SOC在一日内呈现“谷充峰放”的大趋势但分布带的宽度反映了模型为不同光伏/负荷情景预留的调节空间。如果SOC曲线在某个时段收得很紧说明该时段约束非常吃紧比如储能容量快到头了这个时段也是最容易出风险的地方。购电计划同样值得关注。随机优化的日前购电计划往往比确定性模型更平稳峰时购电量会低一些谷时购电量会高一些这部分多出来的谷时购电其实就是储能的充电来源。换句话说随机优化通过“多买谷、少买峰”自动做了电价套利这是它经济性提升的重要来源。5.4 场景数敏感性测试做场景敏感性测试的方法很简单分别用5、10、15、20、30个典型场景跑一遍STO记录期望成本和求解时间。典型场景数期望成本求解时间52920元8秒102890元35秒152882元120秒202878元320秒302875元630秒从我的经验看10到15个场景是一个不错的折中点。场景太少会过滤掉部分尾部风险结果偏乐观场景太多提升有限计算时间却线性增加。如果你的系统设备多、约束细建议用10个以内如果追求论文里的收敛性曲线可以一直跑到30个以上看期望成本收敛趋势。6. 跑不通和跑不对三类典型问题的排查全过程6.1 求解器返回infeasible不是模型问题先检查这五个地方遇到infeasible时我通常按下面的顺序排查第一功率平衡等式或不等式的正负号。这个是最高频的坑。尤其是储能充电功率在等式里到底算负荷还是电源以及弃能变量是在等式左边还是右边很容易写反。我的排查方法是先跑一次确定性模型用预测值当唯一场景如果确定模型能求解那问题大概率出在场景样本或第二阶段变量的维度映射上。第二SOC越界。初始SOC值不在允许区间内或者末时段SOC收得太紧都会导致无解。把SOC_min和SOC_max的区间放宽到0~1再跑一次如果求解成功就说明SOC约束太紧是根源。第三储能充放电互斥约束加上场景变量的维度不一致。比如u_ch是1xTP_ch_s是SxT有些约束里忘了加repmat或者(:)索引导致维度不匹配报错。第四第二阶段变量没有正确的边界。切负荷、弃能量应给一个合理的上限比如切负荷不得超过该场景负荷的5%弃能上限也要和光伏装机容量对齐。如果上限给成0模型当然无解。第五场景生成里的异常值。光伏场景可能因为正态分布采样而产生负数虽然做了max(0)截断但如果标准差设置过大场景可能整体偏离合理范围。用plot画出削减后场景的箱线图一目了然。6.2 求解成功但结果反直觉储能不充不放、燃机无脑满发怎么定位我们做日前调度时最常见的问题是储能SOC曲线始终贴着初值不动或者干脆完全不用储能。出现这种情况先别急着怀疑模型按下面顺序排查先看分时电价是否真实存在价差。如果购电价全天都是1.0元/kWh、售电0.5元/kWh价差固定且没有储能成本参数那储能套利空间可能被成本参数全部吃掉。把储能退化成本参数调成0跑一次如果SOC开始动了说明是成本参数设置不合理如果还不动说明是约束或者目标函数里漏接了储能充放电变量。再看二进制变量和连续变量的衔接。充放电状态变量u_ch、u_dis如果和目标函数没有联系而只约束了P_ch u_ch * P_ch_max那只要P_ch本身是0u_ch随便取值都不影响目标模型就会偷懒让u_ch0。解决方法是把P_ch的上限同时用u_ch缩紧并且可以考虑在目标函数加一个极小惩罚项来强制二进制变量和功率变量逻辑一致。同样道理燃机无脑满发一般是爬坡约束没加或者燃料成本被建模成固定常数而非随出力变化。确认一下目标函数里P_g_s是否真的有线性或二次成本系数。我有一次调试时把Yalmip表达式写成了c_gas * (P_g_s 0)这其实是一个指示函数约束表达的是“只要开机就付一次启动成本”而单时段燃料消耗没乘上去结果当然是机组尽量开机满发。6.3 场景缩减后结果与预期差太远先检查场景质量再检查模型如果随机优化结果和确定性模型结果几乎一样或者相反方向发散大概率是场景削减环节出了问题。先用三个指标校验场景集质量削减后场景均值曲线和原始1000个场景的均值曲线偏差是否在5%以内。削减后场景覆盖的极端值区间是否保留具体做法是看原始场景各个时段5%分位数和95%分位数的范围以及削减后的场景是否落到这个范围内。削减后场景的相邻时段相关性是否还保留尤其是光伏场景不能出现“一会儿满发一会儿零出力”的锯齿形波动。如果这几个指标都正常那问题多半出在模型本身。把两个模型的输出放在同一个时间轴上对比如果DET在峰时段大量购电而STO选择尽量自给自足说明随机模型的光伏场景中可能存在“低出力峰值时段”的组合场景这对购电策略影响很大。这时候去查该场景的概率权重是否合理或者是否在削减时把高影响低概率场景误删了。6.4 求解时间长到不能接受从模型结构角度降维时间太长优化方向分两层。第一层是从建模角度减少整数变量。储能充放电互斥状态可以考虑不用二进制变量改为大M线性化加凸松弛如果系统运行时SOC永远不会同时充放这个松弛在很多算例下是紧的。第二层是从求解器参数入手设置MIPGap为0.5%~1%而不是0.1%对日前调度这种场景MIPGap0.5%已经足够工程使用求解时间能缩短1/3以上。还可以做约束缩减对每个典型场景只写它本时段的功率平衡约束而不是把不相关的场景时段耦合在一起。随机规划的约束矩阵天然具有块结构如果有条件可以用Benders分解把问题拆成主问题和子问题不过对于中小规模微网来说直接全量建模型让Gurobi跑就够了。6.5 数值技巧尺度归一化的重要性一个容易被忽视的问题是数值尺度。储能SOC是0~1的量级购电功率是几百kW量级成本又是几千上万元量级这三类变量混在一个模型里求解器的数值鲁棒性会受到很大挑战。我的习惯是变量/参数类型处理方式功率变量kW归一到MW除以1000储能SOCkWh用0~1标幺值表示成本参数元/kWh直接使用但和功率量纲统一后范围可控时间步长统一为小时爬坡速率用kW/h表示这样做的好处是约束矩阵的条件数更好Gurobi求解时数值迭代更稳定遇到NumericFocus都不用开太高。代码里在读取数据后统一做一次尺度转换即可后处理时再乘回来。写在最后的一段个人体会如果把这段代码实现的经验压缩成一句话那就是随机优化的难点不在求解器而在你能不能把“哪些决策在场景实现前做、哪些决策在场景实现后做”这个两阶段逻辑理清楚。我在实际项目里踩过最大的坑不是求解器崩了而是第一阶段和第二阶段变量混用导致结果看似漂亮、实际无法落地。现在我的做法是建模型前先在纸上画一条时间轴把24小时分成日前段和日内段分别标出决策变量和不确定量再开始写代码。你如果觉得这篇内容有用建议先拿一个小算例把两阶段结构跑通再往里面堆设备模型和约束这样调试成本会低很多。 计及源-荷双重不确定性的微网日前随机优化调度及MATLAB实现
做微网日前调度这几年我最常被同事问的一句话是明天的光伏预测曲线和负荷预测曲线都在这为啥不直接按它排计划问这个问题的人通常没吃过大亏。预测值只是期望实际光伏出力可能只有预测的六成负荷侧园区变压器突增几百千瓦也很常见。这时候再回头看那份按单一预测值排出来的“最优”计划往往已经和实际运行脱节很远。这篇内容我会围绕“MATLAB代码计及源-荷双重不确定性的虚拟电厂/微网日前随机优化调度”展开完整梳理一套可落地的实现思路。代码层面用Yalmip搭配商业求解器Cplex或Gurobi均可不确定性建模走场景法。适合正在做虚拟电厂、微网经济调度、园区综合能源相关课题的研究生和工程师参考。我会把建模动机、数学规划、代码骨架、结果分析逻辑和排错经验一条龙讲清楚。1. 先看清对手日前调度里源-荷双重不确定性到底难在哪1.1 “源”的不确定性与“荷”的不确定性是两种不同性质的问题源侧主要指分布式光伏和风电。光伏出力受辐照度影响云层一过出力可以短时间下降50%以上。风电更极端爬坡事件本身就有很强的随机性。这类不确定性本质上来自“自然条件不可精确预知”你只能拿到一个概率分布而不是一个确定值。负荷侧不确定性来自用户行为园区里某条生产线临时加单、商场空调负荷跟随人流波动这些都是难以精确建模的。和源侧相比负荷的不确定性还有一个特点它存在周期性规律下的随机扰动比如工作日早晚高峰的基线模式相对固定但具体在那个时段涨多少、什么时候涨很难预测到分钟级精度。1.2 日前决策为什么怕“预测不准”微网和虚拟电厂的日前调度是指在D-1日做出未来24小时的机组启停计划、储能充放电计划、与主网的购售电计划。这个计划一旦定下来机组启停这种离散动作在日内很难快速翻转储能也不能因为实时电价波动就随意改SOC目标值。也就是说日前决策的大部分变量是“先定下来”的日内只能在小范围内做修正。如果只按预测值做单场景确定性优化预测误差会直接传导成两部分损失。一是经济性损失光伏实际出力比预测少微网就得多在主网买高价电或者被迫启动更贵的备用机组反之光伏出力比预测多又可能出现弃光。二是可靠性损失严重误判可能导致日前计划下备用不足日内出现切负荷。所以把“源-荷双重不确定性”纳入建模不是锦上添花而是保证计划可用的基本盘。1.3 一个随手可复现的代价测算假设某个微网日预测光伏出力1000 kWh实际只有600 kWh日前计划中没有考虑这种偏差按峰时电价0.9元/kWh计算需要临时多买400 kWh成本多出约360元。而如果随机优化在日前就多预留了储能空间、安排了一部分燃机出力这部分临时购电可能只有150 kWh。看起来单日省下的钱不多但按一个月30天、多个微网聚合来看差距就是数量级差异。在实际项目中我习惯用一个指标来评估不确定性的“杀伤力”叫随机解价值Value of Stochastic Solution, VSS用确定性模型的期望实际成本减去随机优化解的期望实际成本。VSS越大说明不确定性对调度结果的影响越显著也说明做随机优化越有必要。后面章节我会给出具体计算方法。2. 场景法建模把“不确定性”翻译成优化模型的第一步2.1 为什么选场景法而不是鲁棒优化或机会约束处理不确定性的主流路径有三类这里给出适用性对比方便你在自己课题里做取舍。方法基本原理优点缺点适用场景场景法随机规划用离散场景近似连续概率分布目标为期望成本最小信息利用充分可直接对接任意分布结果解释性强场景多时计算量大概率分布准确性影响结果源荷分布特征明确需要经济性最优方案鲁棒优化在给定不确定集中寻找最坏情况下的最优方案计算相对高效不依赖精确分布结果偏保守经济性偏差对可靠性要求极高极端场景必须兜底机会约束规划允许约束以一定概率违反在经济性和可靠性之间可调概率约束的求解相对复杂分布信息需求不低可靠性和经济性折中工程中常用但建模有难度我的个人倾向是学术研究、方案对比、需要输出可视化结果时优先用场景法。原因是虚拟电厂/微网的实际工程中光伏、负荷都有相对成熟的历史数据可以做概率分布拟合场景法能把这种统计信息用足。鲁棒优化更适合做主变过载校验、极端天气应对这类安全校核而不是日常经济调度。2.2 场景生成分布假设、采样方式和相关性问题场景法的第一步是把“不确定性”变成一组带概率的场景集合。光伏出力的经典做法是对光照辐照度建模特卡洛采样把辐照度通过光伏出力模型换算成功率。负荷预测误差通常假设服从均值为0的正态分布标准差取预测值的5%~10%然后叠加到基准负荷曲线上。采样方式上直接蒙特卡洛也能用但收敛偏慢。我更推荐拉丁超立方采样Latin Hypercube Sampling, LHS它的核心思想是把每个分布变量的累积概率区间等分成N份从每个区间里各采一个样本点这样能用较少的样本覆盖整个输入空间。举个例子用纯蒙特卡洛生成500个场景可能仍然留下采样空洞而LHS用200个场景的覆盖效果通常就更好。需要注意的是源荷相关性。光伏出力在相邻时段高度相关不可能8点出力是0、9点瞬间跳到满发这种相关性需要在采样时处理。常见做法是用Cholesky分解对独立采样得到的标准正态向量做相关变换或者直接基于历史场景做重采样。代码层面我会在生成场景矩阵后加一步相关性校验画出相邻时段散点图确认相关性结构合理这是很多人容易忽略的环节。2.3 场景削减为什么不能直接拿1000个场景进优化模型假设采样1000个场景每个场景24个时段变量规模就被乘了1000倍。连续变量还好说储能充放电状态如果引入二进制变量求解时间会指数级恶化。所以必须做场景削减用尽量少的典型场景逼近原有概率分布。工程上常用两种削减方法一种是k-means聚类把距离近的场景聚成几类用聚类中心代表一类场景概率等于类内样本占比另一种是快速前向选择Fast Forward Selection基于Kantorovich概率距离逐步合并最接近的场景对被广泛应用于电力系统随机规划。k-means速度较快但需要预先指定聚类数快速前向削减迭代过程更细致能同时给出削减后的场景集和对应概率。我的实操建议是把两类方法都用一遍对比削减前后场景集合的均值曲线和标准差曲线。只要均值、方差偏差在5%以内说明削减质量可接受。一般1000个场景削减到10~20个就是典型场景数。这个数量既能刻画不确定性又不会把求解时间拖到不可接受。后文建模环节我会给一个基于概率距离削减的简化MATLAB思路。3. 两阶段随机规划模型目标函数与约束的逐条拆解3.1 两阶段决策结构的直观解释随机规划里有个很关键的“非预期性”概念翻译成大白话就是你今天做日前计划的时候并不知道明天哪个场景会真实发生所以第一阶段变量不依赖未来信息。但到了日内场景逐步揭开实时调整变量可以随场景变化。对应到具体模型第一阶段变量在场景实现前就要确定燃气轮机启停状态、储能充放电状态、与主网购售电的日前计划值、储能各时段容量计划值。第二阶段变量随场景变化实际功率平衡下的燃机出力微调、储能充放电功率微调、切负荷量、弃光弃风量。两阶段划分是整个随机调度模型的骨架。我见过不少初学者把全部变量都建成了“随场景变化”结果模型的最优解等于每个场景单独优化后再加权平均这种解完全不可执行因为它忽略了你必须先做出日前决策这个物理事实。所以建模型前先列一张表哪些变量是第一阶段、哪些是第二阶段写清楚后再动代码。3.2 目标函数期望运行成本最小化目标函数是所有场景下的概率加权运行成本之和。成本项通常包括下面五类。第一与主网交互成本。t时段从主网购电价格为λ_t^buy售电价格为λ_t^sell对应变量为P_t^buy和P_t^sell。售电价格一般低于购电价格否则模型会利用价差无意义地倒买倒卖。第二燃气轮机的燃料成本。这里需要对机组启停状态变量u和出力变量P_g做耦合建模。燃料成本采用二次函数表达式在Yalmip中可以用由分段线性近似或使用二次项配合商用求解器处理Cplex/Gurobi都支持二次目标。第三储能充放电退化成本。电池每充放一度电都有寿命损耗折算成单位吞吐量成本参数C_bat。这个参数可以按电池全寿命周期总吞吐量除以购置成本来估计比如一套1000 kWh储能的更换成本50万总吞吐量5000 MWh则单位吞吐量成本约为0.1元/kWh乘2就是一次满充满放的等效损耗。第四弃光弃风惩罚成本。源侧出力用不完又无法上网时允许弃掉一部分但目标函数里给一个惩罚系数让模型尽量少弃。惩罚系数的取值要高于上网电价否则模型会无理由弃光。第五切负荷惩罚成本。切负荷和弃能相比是更严重的事件惩罚系数应该明显高于购电价和弃能惩罚一般取电价的10~50倍代表停电对用户造成的综合损失。整体的目标函数表达式可以写成所有场景s以概率π_s加权平均的最小化总成本注意第一阶段变量不随s变化第二阶段变量才带下标s。3.3 约束条件功率平衡、设备出力与储能动态功率平衡是每一类调度模型都绕不开的硬约束。它要求在任何场景s的任何一个时段t微网内的负荷、购电、燃机出力、储能放电之和等于光伏出力、储能充电、售电、切负荷和弃能之间的差值满足等式约束。切负荷和弃能在等式约束里有点像“松弛变量”但它们只在极端场景下才会超出允许上限日常场景应保持为0。燃气轮机约束包括出力上下限、爬坡约束、最小启停时间约束。其中爬坡约束在随机模型中要按每个场景分别施加因为同一个机组在不同场景下的相邻时段出力可能相差很大必须分别校验每分钟爬坡能力是否满足。储能约束包括SOC动态更新式、充放电功率上下限、电量上下限以及避免同时充放电的约束。SOC动态更新式是t1时段电量等于t时段电量加上充电功率乘效率再减去放电功率除以放电效率。若不使用二进制状态变量可以通过“充放电功率乘积为0”来处理但这类非线性约束在MILP中并不友好工程上最稳妥是引入两个二进制变量u_ch和u_dis用大M约束保证二者不同时等于1。储能还需要设置调度周期始末电量约束。常见做法是要求周期末SOC回到初始值或限制在某个区间保证调度计划可以在次日持续滚动执行。这里有个技巧末时段SOC约束不要锁死在一个固定值或者硬性等于初始值否则会牺牲大量经济性。我给微网做调度时通常给一个可容忍的小区间比如初始值±5%。3.4 备用约束随机优化特有的“安全垫”随机优化里的备用约束尤其重要。在单一预测场景的确定性模型里备用约束可以用“预测负荷预测光伏误差”来折算但随机优化中任何一个场景都可能触发最坏情况需要在日前计划中预留向上/向下备用容量。最简化的写法是要求任意场景s下系统可用总上调能力燃机未用容量储能可放电可中断负荷主网购电余量不低于该场景预测误差导致的净负荷波动值。这里的净负荷定义为负荷减去新能源出力备用需求同样随场景变化。很多人在模型里把备用写成固定值比如按峰值负荷的10%这样虽然简单但没有充分利用场景信息。随机优化的优势恰恰在于能按每个典型场景的净负荷波动动态计算备用需求。我在代码实现里就是把备用约束纳入第二阶段变量的一致性校验中这样模型的自洽性更强。4. MATLAB工程化实现从场景生成到Yalmip求解的代码骨架4.1 工具选型YalmipCplex/Gurobi是性价比最高的组合MATLAB里做随机优化求解的路径主要有三种直接用优化工具箱的linprog/intlinprog用Yalmip建模再调外部求解器或者用基于问题框架的prob2struct转换。我的经验是建MILP模型首选Yalmip因为它对变量的维度处理非常灵动尤其适合“变量带场景下标”这类随机规划模型。Yalmip不是一个求解器它只是个建模语言层底层求解推荐用Gurobi或Cplex。对学术用户这两家都有免费授权Gurobi的安装集成在MATLAB里非常顺滑。4.2 代码结构按“数据-场景-模型-求解-后处理”五段式组织一套能长期复用的随机优化代码结构清晰比代码技巧更重要。我的目录安排参考VPP_Stochastic/ ├── main.m % 主程序 ├── data/ │ ├── load_data.m % 负荷基准曲线 │ ├── pv_data.m % 光伏容量与预测曲线 │ └── price_data.m % 分时电价 ├── uncertainty/ │ ├── gen_scenarios.m % 场景生成 │ ├── reduce_scenarios.m % 场景削减 │ └── plot_scenarios.m % 场景可视化 ├── model/ │ ├── build_model.m % 构建Yalmip模型 │ └── constraints.m % 约束集合 └── result/ └── postprocess.m % 结果分析与图表这种分离最大的好处是换一个测试系统时只需要改data目录下的数据文件改不确定性建模时只动uncertainty目录模型骨架不用推倒重来。4.3 场景生成与削减的核心代码思路场景生成我采用LHS加相关性处理。简化代码如下% gen_scenarios.m 核心思路 function [PV_scen, L_scen] gen_scenarios(data, Ns) % data: 包含预测值、标准差等 % Ns: 初始场景数 rng(42); % 固定随机种子保证结果可复现 T 24; % 光伏出力对辐照度采样需保证非负 % 假设辐照度预测值为 G_pred(Tx1)标准差为 sigma_g % 采用LHS生成标准正态样本 U lhsdesign(Ns, T); Z norminv(U, 0, 1); % 转换为标准正态 % 考虑相邻时段相关性简单一阶自回归 rho 0.85; % 相邻时段相关系数 for t 2:T Z(:,t) rho * Z(:,t-1) sqrt(1-rho^2) * Z(:,t); end G_scen max(0, G_pred sigma_g .* Z); PV_scen data.eta_pv * data.A_pv * G_scen; % 换算成功率 % 负荷预测值加正态误差 L_scen data.L_pred data.sigma_L .* randn(Ns, T); L_scen max(0, L_scen); end场景削减我用基于概率距离的快速前向选择实现思路如下初始时每个场景概率为1/Ns每次迭代找出概率距离乘积最小的一对场景把其中一个场景的概率合并到另一个上然后删除被合并场景重复直到剩余场景数满足要求。MATLAB代码可以写成循环加距离矩阵更新注意向量化处理以提高速度。1000个场景削减到10个在普通电脑上秒级完成。4.4 模型构建与求解的Yalmip代码骨架接下来是模型部分。我用一个具体的微网算例来讲代码结构。系统参数如下设备参数数值燃气轮机台数1台燃气轮机额定功率500 kW燃气轮机最小出力50 kW燃气轮机爬坡速率100 kW/h储能容量1000 kWh储能最大充/放电功率200 kW储能充/放电效率95%光伏装机容量800 kW负荷日峰值700 kW主网交互购电功率上限500 kW主网交互售电功率上限300 kWYalmip建模骨架如下% build_model.m 核心片段 T 24; S 10; % 削减后典型场景数 % 第一阶段决策变量 P_buy sdpvar(1, T, full); % 日前计划的购电功率 P_sell sdpvar(1, T, full); % 日前计划的售电功率 P_g sdpvar(1, T, full); % 燃机日前计划出力 u_g binvar(1, T, full); % 燃机启停 P_ch sdpvar(1, T, full); % 储能充电计划 P_dis sdpvar(1, T, full); % 储能放电计划 u_ch binvar(1, T, full); % 充电状态 u_dis binvar(1, T, full); % 放电状态 SOC sdpvar(1, T1, full); % 储能电量状态 % 第二阶段决策变量场景相关 P_g_s sdpvar(S, T, full); % 场景下的燃机实际出力 P_dis_s sdpvar(S, T, full); % 场景下的放电功率 cut_s sdpvar(S, T, full); % 切负荷量 dump_s sdpvar(S, T, full); % 弃光弃风量 buy_s sdpvar(S, T, full); % 场景下实际购电 sell_s sdpvar(S, T, full); % 场景下实际售电 % 目标函数 Objective 0; for s 1:S for t 1:T Objective Objective pi_s(s) * (... price_buy(t) * buy_s(s,t) - price_sell(t) * sell_s(s,t) ... c_gas * P_g_s(s,t) c_dg * (P_g_s(s,t) 0) ... c_bat * (P_ch(t)/eta_ch P_dis_s(s,t)*eta_dis) ... c_dump * dump_s(s,t) c_cut * cut_s(s,t)); end end % 注意场景下标放在第二维度会对求解器和Yalmip都更友好 % 实际代码中建议变量维度写成 [S, T] 并用 repmat 操作这里有个小坑需要提醒Yalmip 中如果用binvar(S, T)这种大维度0-1变量矩阵预分配和约束添加时要格外小心千万别在循环里反复调用constraints [constraints, ...]来拼约束数组性能会很差。正确做法是预创建一个元胞数组C {};每次循环把新约束加进去最后用Constraints [C{:}];一次性拼接。这在场景数较大时能节省不少构建时间。约束添加的核心逻辑Constraints {}; % 功率平衡每个场景单独满足 for s 1:S for t 1:T Constraints{end1} ... L_scen(s,t) P_ch(t) P_dis_s(s,t) sell_s(s,t) ... PV_scen(s,t) P_g_s(s,t) buy_s(s,t) P_dis(t) cut_s(s,t) - dump_s(s,t); % 更严谨写法是用等式加RHS按实际系统平衡方式决定 end end注意功率平衡里储能充电和放电的建模方式不同我的习惯是把充电视为负荷、放电视为电源同时加入充放电互斥条件for t 1:T Constraints{end1} P_ch(t) u_ch(t) * P_ch_max; Constraints{end1} P_dis_s(:,t) u_dis(t) * P_dis_max; % 这里用第一阶段状态约束第二阶段功率 Constraints{end1} u_ch(t) u_dis(t) 1; end这种写法把“日前决定储能是充还是放”放在第一阶段场景只影响放电功率的幅度符合物理直觉。运行结果告诉我在大多数场景下模型倾向于在谷时段充电、峰时段放电但光伏大发且负荷低谷时也会有“光伏直接充电”的行为这说明随机场景对储能策略的刻画比确定性模型更细腻。SOC约束SOC(:,1) SOC_init; for t 1:T Constraints{end1} SOC(:,t1) SOC(:,t) ... eta_ch * P_ch(t) - P_dis_s(:,t) / eta_dis; Constraints{end1} SOC_min SOC(:,t1) SOC_max; Constraints{end1} SOC(:,T1) SOC_init - 0.05 * SOC_max; Constraints{end1} SOC(:,T1) SOC_init 0.05 * SOC_max; end这里给末时段SOC留了5%的弹性区间很多论文里写成硬等于实际工程中会让储能为了“还电量”而做不必要的充放经济性反而下降。4.5 求解配置与稳定输出求解器配置我一般这样写ops sdpsettings(solver, gurobi, verbose, 2, ... gurobi.TimeLimit, 600, ... gurobi.MIPGap, 0.001, ... gurobi.NumericFocus, 3); sol optimize(Constraints, Objective, ops); if sol.problem 0 % 求解成功 else disp(sol.info); % 输出失败原因 endNumericFocus3是Gurobi的数值鲁棒性增强选项在含0-1变量和大M参数较多的模型里能显著降低“数值病态导致结果抖动”的概率。代价只是求解时间略微增加对于日前调度这种离线计算场景完全值得。5. 结果怎么看出门道确定性调度与随机调度的对比分析5.1 搭建一个公平的对比基线我在做结果分析时通常把预期结果拆成三组对照分析确定性调度模型DET直接用预测光伏、预测负荷求解。随机优化调度模型STO用削减后的场景集合求解。完美信息下界WS即Wait-and-See假设知道每个场景真实值分别求解后再做概率加权平均。这个模型是理论下界现实中不可能实现但它是评价随机优化质量的重要参考。三组模型跑完后统计各方案的期望实际成本。实际成本计算方法是把STO方案的第一阶段决策固定然后在每个测试场景下做第二阶段优化只优化功率调整量与切负荷/弃能量得到每个场景的真实运行成本再做概率加权和。5.2 随机调度为什么“看起来更贵”从结果表来看STO的期望成本通常比WS高但比DET低。指标确定性调度随机优化调度完美信息下界第一层目标函数值未包含场景信息概率加权期望成本各场景期望成本加权平均期望实际成本计入偏差惩罚3220元2890元2740元切负荷期望量12 kWh2 kWh0 kWh弃光期望量45 kWh18 kWh3 kWhVSS 3220 - 2890-330元约11.4%-这个例子中STO相比DET节约了约10%的实际运行成本主要来自两部分一是减少了高价时段临时购电二是削峰填谷更加精准最小化了切负荷和弃能风险。DET的目标函数数值看起来很漂亮因为它的约束建立在“预测完全准确”的基础上但一旦考虑偏差它的目标值就成了镜花水月。“随机优化的期望成本高于确定性模型”这个现象不少初学者第一次看到会困惑。原因很简单确定性模型只面对一个“平均预测场景”约束容易满足目标值自然低随机模型要同时满足10个差异很大的场景任何方案都必须在不同场景间做折中所以最优期望成本会上升。但这恰恰说明模型在准备应对不确定性。如果哪天随机优化的期望成本和确定性模型完全一样那大概率说明场景削减没有体现不确定性特征或者模型结构出了问题。5.3 SOC、购电计划与场景间的耦合关系除了总体成本我还会画出储能SOC曲线在不同典型场景下的分布带。正常情况下你会看到SOC在一日内呈现“谷充峰放”的大趋势但分布带的宽度反映了模型为不同光伏/负荷情景预留的调节空间。如果SOC曲线在某个时段收得很紧说明该时段约束非常吃紧比如储能容量快到头了这个时段也是最容易出风险的地方。购电计划同样值得关注。随机优化的日前购电计划往往比确定性模型更平稳峰时购电量会低一些谷时购电量会高一些这部分多出来的谷时购电其实就是储能的充电来源。换句话说随机优化通过“多买谷、少买峰”自动做了电价套利这是它经济性提升的重要来源。5.4 场景数敏感性测试做场景敏感性测试的方法很简单分别用5、10、15、20、30个典型场景跑一遍STO记录期望成本和求解时间。典型场景数期望成本求解时间52920元8秒102890元35秒152882元120秒202878元320秒302875元630秒从我的经验看10到15个场景是一个不错的折中点。场景太少会过滤掉部分尾部风险结果偏乐观场景太多提升有限计算时间却线性增加。如果你的系统设备多、约束细建议用10个以内如果追求论文里的收敛性曲线可以一直跑到30个以上看期望成本收敛趋势。6. 跑不通和跑不对三类典型问题的排查全过程6.1 求解器返回infeasible不是模型问题先检查这五个地方遇到infeasible时我通常按下面的顺序排查第一功率平衡等式或不等式的正负号。这个是最高频的坑。尤其是储能充电功率在等式里到底算负荷还是电源以及弃能变量是在等式左边还是右边很容易写反。我的排查方法是先跑一次确定性模型用预测值当唯一场景如果确定模型能求解那问题大概率出在场景样本或第二阶段变量的维度映射上。第二SOC越界。初始SOC值不在允许区间内或者末时段SOC收得太紧都会导致无解。把SOC_min和SOC_max的区间放宽到0~1再跑一次如果求解成功就说明SOC约束太紧是根源。第三储能充放电互斥约束加上场景变量的维度不一致。比如u_ch是1xTP_ch_s是SxT有些约束里忘了加repmat或者(:)索引导致维度不匹配报错。第四第二阶段变量没有正确的边界。切负荷、弃能量应给一个合理的上限比如切负荷不得超过该场景负荷的5%弃能上限也要和光伏装机容量对齐。如果上限给成0模型当然无解。第五场景生成里的异常值。光伏场景可能因为正态分布采样而产生负数虽然做了max(0)截断但如果标准差设置过大场景可能整体偏离合理范围。用plot画出削减后场景的箱线图一目了然。6.2 求解成功但结果反直觉储能不充不放、燃机无脑满发怎么定位我们做日前调度时最常见的问题是储能SOC曲线始终贴着初值不动或者干脆完全不用储能。出现这种情况先别急着怀疑模型按下面顺序排查先看分时电价是否真实存在价差。如果购电价全天都是1.0元/kWh、售电0.5元/kWh价差固定且没有储能成本参数那储能套利空间可能被成本参数全部吃掉。把储能退化成本参数调成0跑一次如果SOC开始动了说明是成本参数设置不合理如果还不动说明是约束或者目标函数里漏接了储能充放电变量。再看二进制变量和连续变量的衔接。充放电状态变量u_ch、u_dis如果和目标函数没有联系而只约束了P_ch u_ch * P_ch_max那只要P_ch本身是0u_ch随便取值都不影响目标模型就会偷懒让u_ch0。解决方法是把P_ch的上限同时用u_ch缩紧并且可以考虑在目标函数加一个极小惩罚项来强制二进制变量和功率变量逻辑一致。同样道理燃机无脑满发一般是爬坡约束没加或者燃料成本被建模成固定常数而非随出力变化。确认一下目标函数里P_g_s是否真的有线性或二次成本系数。我有一次调试时把Yalmip表达式写成了c_gas * (P_g_s 0)这其实是一个指示函数约束表达的是“只要开机就付一次启动成本”而单时段燃料消耗没乘上去结果当然是机组尽量开机满发。6.3 场景缩减后结果与预期差太远先检查场景质量再检查模型如果随机优化结果和确定性模型结果几乎一样或者相反方向发散大概率是场景削减环节出了问题。先用三个指标校验场景集质量削减后场景均值曲线和原始1000个场景的均值曲线偏差是否在5%以内。削减后场景覆盖的极端值区间是否保留具体做法是看原始场景各个时段5%分位数和95%分位数的范围以及削减后的场景是否落到这个范围内。削减后场景的相邻时段相关性是否还保留尤其是光伏场景不能出现“一会儿满发一会儿零出力”的锯齿形波动。如果这几个指标都正常那问题多半出在模型本身。把两个模型的输出放在同一个时间轴上对比如果DET在峰时段大量购电而STO选择尽量自给自足说明随机模型的光伏场景中可能存在“低出力峰值时段”的组合场景这对购电策略影响很大。这时候去查该场景的概率权重是否合理或者是否在削减时把高影响低概率场景误删了。6.4 求解时间长到不能接受从模型结构角度降维时间太长优化方向分两层。第一层是从建模角度减少整数变量。储能充放电互斥状态可以考虑不用二进制变量改为大M线性化加凸松弛如果系统运行时SOC永远不会同时充放这个松弛在很多算例下是紧的。第二层是从求解器参数入手设置MIPGap为0.5%~1%而不是0.1%对日前调度这种场景MIPGap0.5%已经足够工程使用求解时间能缩短1/3以上。还可以做约束缩减对每个典型场景只写它本时段的功率平衡约束而不是把不相关的场景时段耦合在一起。随机规划的约束矩阵天然具有块结构如果有条件可以用Benders分解把问题拆成主问题和子问题不过对于中小规模微网来说直接全量建模型让Gurobi跑就够了。6.5 数值技巧尺度归一化的重要性一个容易被忽视的问题是数值尺度。储能SOC是0~1的量级购电功率是几百kW量级成本又是几千上万元量级这三类变量混在一个模型里求解器的数值鲁棒性会受到很大挑战。我的习惯是变量/参数类型处理方式功率变量kW归一到MW除以1000储能SOCkWh用0~1标幺值表示成本参数元/kWh直接使用但和功率量纲统一后范围可控时间步长统一为小时爬坡速率用kW/h表示这样做的好处是约束矩阵的条件数更好Gurobi求解时数值迭代更稳定遇到NumericFocus都不用开太高。代码里在读取数据后统一做一次尺度转换即可后处理时再乘回来。写在最后的一段个人体会如果把这段代码实现的经验压缩成一句话那就是随机优化的难点不在求解器而在你能不能把“哪些决策在场景实现前做、哪些决策在场景实现后做”这个两阶段逻辑理清楚。我在实际项目里踩过最大的坑不是求解器崩了而是第一阶段和第二阶段变量混用导致结果看似漂亮、实际无法落地。现在我的做法是建模型前先在纸上画一条时间轴把24小时分成日前段和日内段分别标出决策变量和不确定量再开始写代码。你如果觉得这篇内容有用建议先拿一个小算例把两阶段结构跑通再往里面堆设备模型和约束这样调试成本会低很多。