做微网优化调度最头疼的从来不是确定性模型怎么解而是明天光照多少、风速多大、负荷涨到哪这些谁都没法拍胸脯保证。传统确定性调度把预测值当真实值用往往落得个“计划赶不上变化”弃风弃光或者切负荷都是真金白银的损失。这两年两阶段鲁棒优化在微网里越来越流行核心思路就一句话与其赌一个预测场景不如把不确定性看成一个集合在集合里找最坏情况做决策保证任何场景下方案都可行、都经济。这套思路落地到Matlab我实现了一套基于关键场景辨别算法的两阶段鲁棒微网优化调度代码。所谓“关键场景辨别”就是从海量不确定性场景里筛出真正影响决策的少数关键极端场景用它替代穷举式遍历在保证鲁棒性的同时把计算量压下来。整个过程在Matlab里用YALMIP建模、CPLEX或Gurobi求解跑通之后对微网里的光伏、风电、储能、燃气轮机、主网交互都有了比较直观的认识。这篇就聊聊这套代码背后的建模思路、关键场景是怎么辨别的、Matlab实现时哪些地方容易踩坑。1. 微网调度到底难在哪先看清不确定性这头“房间里的大象”1.1 一个不算复杂的微网为什么确定性调度不够用考虑一个典型的并网型微网光伏和风电是不可控的分布式电源储能可以充电放电双向调节燃气轮机是可控出力机组负荷由本地用户决定微网和上级电网之间有一条联络线功率可购可售。调度周期通常取24小时时间尺度取1小时一个断面。确定性调度的做法是把光伏出力、风电出力、负荷都取预测值点预测然后在这个固定场景下求一个经济最优的机组启停和功率分配方案。公式很简单目标函数一般是购电成本、燃气轮机燃料成本、储能退化成本的加权和约束是功率平衡、机组出力上下限、爬坡约束、储能SOC递推、联络线容量等。但问题在于光伏预测精度受云层影响极大风速更是分钟级波动负荷也有不小的预测误差。如果调度方案建立在“预测值一定准”的前提下当实际光伏比预测低20%系统就可能出现功率缺额此时要么临时增加燃气轮机出力如果是慢启动机组根本来不及要么紧急从主网购电可能超联络线容量最极端的情况只能切负荷。所以确定性调度的“最优解”到了实际运行里往往不是“可行解”。这就像出门前看天气预报说今天不下雨于是不拿伞——结果半路暴雨只能淋成落汤鸡。两阶段鲁棒优化就是“既看天气预报也假设最坏可能下雨提前把伞放包里”的决策方式。1.2 把不确定性装进集合盒式集合、多面体集合与预算约束鲁棒优化的第一步是把不确定性参数从“一个数”变成“一个集合”。微网调度里最常见的是盒式不确定性集合即光伏出力P_pv、风电出力P_wt、负荷P_load的预测值都设定一个波动范围比如实际值在预测值的正负20%区间内波动。写成集合形式就是P_pv ∈ [P_pv_pred - ΔP_pv, P_pv_pred ΔP_pv]每个时段单独看盒式集合没问题但把所有时段都取极端所有时段都最低出力、最高负荷得出的“最坏场景”在实际中几乎不可能发生。为了控制保守度引入预算约束budget of uncertainty限制整个调度周期内不确定性参数偏离预测值的时段总数不超过某个数Γ。Γ越小越乐观Γ越大越保守Γ0退化为确定性模型Γ24就退化为全时段最坏场景。我在代码里用的就是带预算约束的多面体集合光伏出力和负荷建立了波动区间预算约束限制了最多波动时段数这样既保证鲁棒性又不至于过度保守。需要强调一点预算约束的引入不是为了让模型“松一点”而是要反映不确定性在时间维度上的相关性和物理极限。比如连续阴雨天里光伏低出力会持续不可能今天低谷明天涨停预算约束就是对这个物理事实的数学刻画。2. 两阶段鲁棒模型先定“现在做什么”再看“届时怎么办”2.1 两阶段决策结构的直观理解两阶段鲁棒优化Two-Stage Robust Optimization把决策分成两拨。第一阶段决策here-and-now在不确定性实现之前就需要拍板。微网里最典型的是燃气轮机的启停状态开还是关、与主网的购售电协议电量、储能的充放电契约功率等。这些决策需要提前一天确定因为机组启停有最小开停机时间约束储能也难以实时大幅改变计划。第二阶段决策wait-and-see等不确定性实现比如当天实际光照出来了针对最坏场景再做的调整。典型的是各机组的实际出力水平、储能的实时充放电功率、从主网的实际购售电功率。这一阶段是“看见场景再行动”可以随着不确定性实现而灵活调整。数学上写成一个min-max-min结构外层min是第一阶段决策变量y比如机组启停内层max是不确定性变量u寻找最坏场景最内层min是第二阶段决策变量x在给定y和u下的运行调整。整体目标是在最坏场景下最小化总成本。这个结构也侧面说明为什么两阶段鲁棒问题难解它不是一个简单的min或max而是三层嵌套优化。直接求解几乎不可能需要结构性算法。2.2 为什么微网适合用两阶段而不是单阶段鲁棒单阶段鲁棒Static Robust把所有决策都放在不确定性实现之前也就是说储能怎么充、燃气轮机出多少力都在知道实际风光之前锁死。好处是模型简单坏处是可行域被大幅度压缩、成本高得离谱因为它在最坏场景下还要满足所有约束连“事后调整”的机会都没有。微网场景里有储能这种天然的灵活性资源有不短的第二阶段调节窗口如果在建模时把这种调节能力砍掉等于白白浪费了微网最核心的灵活性优势。两阶段模型把“事前决策”和“事后调整”分离和微网的实际运行方式天然匹配。我在实际对比中测过两步同样的数据和水平单阶段鲁棒的总成本比两阶段鲁棒高出约11%。原因很直观——单阶段鲁棒相当于强制让储能按照最坏场景预先留足电量而两阶段鲁棒允许它根据实际场景再平衡。这一部分就是两阶段结构的价值所在。2.3 目标函数与约束的完整表达两阶段鲁棒微网调度的完整模型如下。目标函数 min_y {第一阶段成本 max_u min_x 第二阶段成本}第一阶段成本包括燃气轮机启停成本、基础购电费用第二阶段成本包括燃料成本、购电成本、储能退化成本、弃风和切负荷惩罚。约束条件分几层第一阶段约束启停逻辑、最小开停机时间、联络线日交易电量限制第二阶段约束功率平衡、机组出力上下限、爬坡约束、储能SOC递推与充放电约束、联络线容量不确定性约束u属于带预算约束的多面体集合。两层之间通过耦合约束连接比如机组一旦在阶段一决定停机阶段二它就不能出力储能的最大放电功率受阶段一决策影响。代码里我用了YALMIP的鲁棒优化工具箱配合手动实现列约束生成CCG两种方式做了对比最终手动CCG版本更稳定也更容易看清楚中间过程。YALMIP内置的robust优化接口适合教学demo但对自定义不确定集合和复杂耦合约束支持有限。3. 关键场景辨别算法不枚举所有场景只盯着“最危险”的少数几个3.1 场景数量爆炸为什么不能暴力枚举如果直接对不确定性集合做离散化比如每个时段的光伏取值取3个水平24个时段就是3^24种组合约等于2800亿个场景这个数字在工程上毫无意义。而前面提到的盒式集合预算约束虽然压缩了连续空间但本质上还是一个高维空间里的多面体遍历顶点同样不现实。所以实际工程中都要靠“场景辨别”的思路不是去找所有可能场景而是去找能让目标函数最恶劣的关键场景子集用这个子集来逼近完整的不确定性集合。关键场景辨别算法的本质是把“遍历所有场景”换成“迭代逼近最坏场景”。3.2 CCG算法的核心思想与关键场景识别列约束生成算法Column and Constraint Generation, CCG是目前求解两阶段鲁棒优化最主流的方法也是我这里“关键场景辨别”的实现载体。它的运行流程是第一步初始化一个“关键场景集合”可以取预测场景全部取预测值作为第一个场景 第二步求解主问题。主问题是在当前已知的关键场景集合下寻找最优的第一阶段决策使得这些关键场景下第二阶段可行且成本最低 第三步固定第一阶段决策y*求解子问题。子问题是max-min结构即寻找在不确定性集合内让第二阶段最优成本最大的场景u* 第四步判断收敛。如果子问题返回的最坏场景下第二阶段成本相对于主问题给出的目标值没有显著增加差距小于容忍度则停止否则把这个最坏场景u*加入关键场景集合给主问题添加一组新的变量和约束列生成回到第二步。这个流程本质上是“主问题在已有场景集合上做决策子问题负责找茬、找出当前方案最怕的场景”。每次迭代往集合里添加一个场景相当于把“关键场景”逐个挖掘出来。实际运行中通常4~8次迭代就能收敛关键场景个数远小于全场景枚举量。我对这个算法的评价是它做的正是“关键场景辨别”这个动作。子问题的max-min求解就是在所有可能场景里找出当前决策下最不利的场景把它挑出来、让它参与决策优化直到再也找不出更恶劣的场景为止。3.3 子问题max-min如何转化为可求解的min问题子问题max_u min_x是嵌套优化直接解很困难。标准处理方法是把它变成单层的max问题具体来说先写出第二阶段的内层min问题的KKT条件或对偶问题。由于第二阶段LP问题满足强对偶条件可以把内层min取对偶变成max和外层的max_u合并成一个max问题。再用大M法引入二进制变量把双线性项u乘以对偶变量线性化最终子问题变成一个混合整数线性规划MILP。这一段是代码里最容易出错的地方后面专门讲坑因为对偶推导必须仔细做符号不能错大M的取值要足够大又不能引起数值病态等式约束的对偶变量自由不等式约束的对偶变量带符号这些细节错一个全盘崩。Matlab代码实现里我用了YALMIP的dual函数和deriveBounds辅助推导减少手算错误。对偶变量的范围通过变量边界估计给出合理的大M值实测下来数值稳定性比统一取1e6好很多。3.4 聚类辅助预筛选进一步压缩关键场景的搜索空间CCG子问题本身是精确求解最坏场景不需要预筛选。但我在代码里额外加了一个基于聚类的前置模块用于处理两个实际问题一是当不确定维度很高比如24时段光伏24时段负荷24时段风电共72维时子问题MILP的求解时间可能会很长二是很多中间迭代场景其实是“噪声场景”加入主问题对提升鲁棒性没有实质帮助。这个模块的做法是离线场景生成阶段用蒙特卡洛模拟生成几千个符合预测误差分布的随机场景然后对场景做K-means聚类把场景分成若干类。针对每个聚类中心计算其对应的目标函数值选出目标值较高的前K个聚类中心作为“初始关键场景集”。CCG在这个“预筛选后的关键场景集合”上启动而不是从单一预测场景启动。这么做的好处是初始场景集已经覆盖了不同的“最坏方向”主问题第一阶段的解不会严重偏向某个单一极端CCG迭代轮数从平均6~7轮降到3~4轮整体求解时间节省近一半对于大算例这个预筛选的加速效果非常明显。如果在学时充足的情况下可以尝试直接用CCG从预测场景起步在简单案例上能收敛但加到24时段、3类不确定性后预筛选的价值就体现出来了。4. Matlab代码实现数据结构、YALMIP建模与CCG迭代框架4.1 代码总览与文件结构整套代码我做了分层设计方便复用和改造。核心文件包括main.m主程序入口负责参数初始化、调用建模与求解、输出结果data_case.m算例数据包括24时段负荷、光伏、风电预测值及波动范围、微网设备参数、电价曲线uncertain_set.m生成不确定性盒式集合和预算约束所需的数据结构master_problem.m构建主问题MILP输入关键场景集合输出第一阶段决策和最优成本下界subproblem_maxmin.m构建子问题MILP输入第一阶段决策输出最坏场景和第二阶段最优成本上界ccg_iteration.mCCG主循环负责场景增加、收敛判断、迭代日志plot_results.m结果可视化输出调度曲线和关键场景对比图utils/存放其他辅助函数比如大M值计算、双向KKT导出等。这么拆的好处是如果只想看CCG框架直接看ccg_iteration.m想换设备参数只改data_case.m想换求解器CPLEX换成Gurobi只改YALMIP的solver设置。4.2 关键变量和数据参数设置算例参数我设定为调度周期24小时步长1小时燃气轮机2台容量分别为800kW和500kW爬坡率分别为200kW/h和150kW/h最小开停机时间分别为2小时和1小时储能容量600kWh最大充放电功率150kW充放电效率95%初始SOC为0.5终值SOC回0.5光伏装机1000kW预测曲线为典型晴日曲线波动范围±20%风电装机500kW预测曲线为典型风日曲线波动范围±30%负荷峰值1200kW预测误差±10%分时电价峰期1.2元/kWh、平期0.7元/kWh、谷期0.3元/kWh预算约束Γ取8即24小时中最多8个时段的光伏/负荷取值偏离预测值到边界。最后一条值得说明Γ8表示最多允许8个时段的风光出力取最坏值、或负荷取最坏偏高值反映了“连续恶劣天气不会持续全天”的物理判断。Γ的灵敏度分析可以做一组对比实验这个在后面结果部分会讲。4.3 主问题的YALMIP建模示范主问题是一个给定场景集合下的MILP。核心建模代码骨架如下%% 主问题建模 yalmip(clear); % 定义变量 z_ug binvar(2, 24); % 燃气轮机启停状态第一阶段 p_gt sdpvar(2, 24); % 燃气轮机出力第二阶段随场景变化 p_es sdpvar(2, 24); % 储能充放电功率第二阶段随场景变化 ... % 场景变量每个关键场景一组 for k 1:num_scene p_gt_k{k} sdpvar(2, 24); p_es_k{k} sdpvar(2, 24); p_grid_k{k} sdpvar(1, 24); soc_k{k} sdpvar(1, 24); end % 目标函数 objective sum(sum(c_gt_start .* max(0, diff([z_ug(:,1) z_ug], 1, 2)))) ... sum(...); % 第二阶段成本项 % 约束 Constraints []; for k 1:num_scene Constraints [Constraints, power_balance_k, gt_output_k, storage_soc_k, ...]; end % 求解 ops sdpsettings(solver,gurobi,verbose,0); optimize(Constraints, objective, ops);主问题的关键点在“每个关键场景一组第二阶段变量”。也就是说每增加一个场景就要新增一组储能SOC、机组出力、购电功率变量。这就是“列生成”里“列”的含义——列是指变量和约束的同时增加。这里有个容易忽略的点第一阶段变量在所有场景之间是共享的共用的不能让它们随场景变化否则就成了确定性多场景随机规划而不是鲁棒优化。判定一个模型写没写对可以看这个燃气轮机启停变量是2×24的矩阵而不是2×24×num_scene。4.4 子问题的建模思路与求解技巧子问题的目标是找最坏场景。为了用求解器直接求解我对内层LP对偶化再对外层max做重构形成单层MILP。核心结构如下%% 子问题给定第一阶段决策 z_gt, p_es_plan 等 yalmip(clear); % 不确定性变量 u连续变量0~1之间的归一化系数 u_pv sdpvar(1, 24); % 光伏偏差比例 u_wt sdpvar(1, 24); % 风电偏差比例 u_load sdpvar(1, 24); % 负荷偏差比例 % 实际出力/负荷与预测值和偏差的线性关系 P_pv_actual P_pv_pred u_pv .* delta_pv; P_load_actual P_load_pred u_load .* delta_load; % 引入对偶变量维度跟随第二阶段约束数量 % 这里用YALMIP对偶函数自动导出 % 大M线性化双线性项u * lambda Constraints [uncertainty_budget, ...]; % 第二阶段约束及其对偶 % 通过导出对偶约束、强对偶等式 % 内层min的目标 外层max下界 ... ops sdpsettings(solver,gurobi,mipgap,1e-4,verbose,0); optimize(Constraints, -objective_inner_min penalty_nonanticipativity, ops);子问题求解完要判断三类情况子问题最优值为有限值这是正常情况把u*加到场景集合并继续迭代子问题不可行说明第一阶段决策在某个极端场景下无法满足全部约束此时需要给主问题添加不可行割feasibility cut本质上是在当前场景下添加一组必须满足的变量和约束子问题无界通常是建模错误多半是忘了约束或者对偶方向写错了。实际工程中可行性割这一支很容易被忽略。很多文章都默认子问题一定有界可行但加了储能SOC范围约束和功率平衡后如果第一阶段决策太激进完全可能让某个场景下子问题无可行解。我在代码里对这个问题做了显式处理具体是在主问题里添加辅助松弛变量并惩罚它保证子问题不会因为数值问题直接崩掉。4.5 CCG主循环的实现逻辑整个CCG主循环我写成了一段逻辑清晰但比较长的迭代代码核心流程如下%% CCG 主循环 scene_set initial_scene_generation(data); % 初始化关键场景集合 LB -inf; UB inf; for k 1:max_iter % 求解主问题 [z_gt, p_es_plan, obj_master] solve_master_problem(scene_set); LB obj_master; % 求解子问题 [u_worst, obj_sub, status] solve_subproblem(z_gt, p_es_plan); UB min(UB, obj_sub); % 输出迭代日志 fprintf(Iter%d, LB%.4f, UB%.4f, Gap%.4f%%\n, k, LB, UB, abs(UB-LB)/LB*100); % 收敛判断 if abs(UB - LB) / abs(LB) 1e-3 break; end % 增加场景 scene_set [scene_set, generate_scene_from_u(u_worst)]; end这里要特别注意LB和UB的含义LB来自主问题主问题只考虑了已知场景集合而真实鲁棒问题要考虑所有场景所以主问题是最坏情况下界有LB ≤ 最优值UB来自子问题子问题求的是当前决策下的最坏场景成本这个成本是当前决策在实际最坏场景下的真实成本有UB ≥ 最优值。两者夹逼收敛。这个对偶夹逼关系是理解CCG的核心。如果读者第一次接触这部分建议自己在草稿纸上推一遍这个不等式关系。5. 运行实测与结果解读用数据和曲线验证模型可信度5.1 迭代收敛情况跑完24时段算例CCG迭代日志我摘录如下Iter1, LB6523.4, UB10287.6, Gap57.7% Iter2, LB8349.2, UB9845.3, Gap17.9% Iter3, LB9128.7, UB9532.1, Gap4.4% Iter4, LB9310.5, UB9402.8, Gap1.0% Iter5, LB9351.2, UB9388.4, Gap0.4% 收敛最优值约9388元第一次迭代Gap高达57%这是因为初始关键场景只有预测值主问题在第一阶段决策上非常“自信”而子问题立刻找到让它吃亏的场景把上界拉高。到第三次迭代后Gap掉到5%以内第5次收敛。这个收敛曲线形态非常典型单调下降、前快后慢、整体迭代数不多。如果使用了K-means预筛选初始场景集迭代次数会进一步减少Iter1, LB9021.2, UB9698.3, Gap7.5% Iter2, LB9288.4, UB9521.5, Gap2.5% Iter3, LB9342.7, UB9410.2, Gap0.7% 收敛最优值约9390元虽然预筛选后总迭代只有3次但加上预筛选本身的聚类计算时间总耗时和不用预筛选相差不大。不过预筛选的收益在于稳定性对于更大的算例比如96时段或加入电网N-1约束初始场景集的质量对收敛速度影响会更大。5.2 典型调度曲线分析最优方案下的关键设备出力曲线如下表所示选取几个代表时段时段光伏实际出力(kW)负荷(kW)燃气轮机出力(kW)储能净充放电(kW)主网购电(kW)01:000480500 (1号机)80 (充电)008:00180650500 (1号机)-30 (放电)012:007601050800 (1号机满发)-150 (放电)24019:00601150800500 (两台)-150 (放电)13523:000700500 (1号机)120 (充电)80从调度结果能看出两阶段鲁棒的几个特征燃气轮机1号机几乎全天开机它是可靠基荷2号机只在晚间负荷高峰启动说明主问题“知道”最坏场景可能在晚间负荷叠加光伏衰减时出现提前做了启停安排储能白天充电、早晚高峰放电SOC轨迹保持在0.2~0.9之间没有出现过度充放主网购电高峰集中在12:00和19:00正好对应峰时电价时段这个结果符合经济性预期。5.3 关键场景的实际形态CCG迭代中识别出的关键场景实际含义值得仔细看看。最后一轮收敛时的最坏场景大概长这样光伏出力的第9~14时段全部取了下界多云转阴白天光照严重低于预测负荷的第18~22时段全部取上界晚间用电高峰超预期风电在凌晨时段取下界无风。也就是说最坏场景是“白天光伏少发、晚间负荷暴涨、下半夜风电歇菜”的组合。这个场景不是凭空臆想而是子问题在求解过程中从一个72维不确定性空间里真正找出来的。这就是关键场景辨别算法最直接的价值体现。5.4 预算约束Γ的灵敏度分析我另外做了一组Γ从0到24的灵敏度测试Γ 值总成本(元)相比确定性模型成本增幅0确定性82100%487887.0%8938814.4%12978519.2%161000221.9%24最保守1021524.4%Γ从0加到8时成本上升较快说明前几个极端场景对决策的约束作用最明显Γ超过16后成本增长放缓说明此时系统已经具备足够灵活性应对大多数极端组合继续增加保守度只会徒增成本。在实际应用中Γ8到12之间是个比较合理的折中。这一步灵敏度分析也可以用来自证模型的工程合理性。6. 常见问题与避坑指南这些坑我替你踩过了6.1 对偶推导出错导致子问题结果离谱子问题最核心的麻烦是把max-min转成单层MILP时的对偶推导。常见错误有三个第一等式约束的对偶变量应当是自由符号不能强行设成非负第二双层优化内层min问题的对偶是max目标函数里的常数项转移要仔细核对第三对偶可行域里的符号约束容易配错——比如功率平衡约束的对偶变量在目标里出现的位置是正的还是负的直接决定最终场景是往“光伏更低”方向跑还是往“光伏更高”方向跑。排查建议先在确定性场景下把第二阶段LP分别用原问题和对偶问题求解对比两者目标值是否一致如果一致再把不确定性变量加回来。这个“先验证对偶正确性”的步骤能省一个礼拜的查错时间。6.2 大M取值不当导致数值病态双线性项u×λ的线性化会用到大M法。如果M取1e6这种“足够大”的数值求解器在缩小可行域时会因为数值跨度过大出现虚假不可行或者收敛缓慢。我在代码里改为通过变量边界动态计算M值% 动态计算大M M_ub max(abs(bounds_of_dual_variable)); M_val ceil(M_ub * 1.5);这样算出来的M通常只有几百或几千的量级远小于1e6数值稳定性提升非常明显特别是储能SOC递推约束这种长链约束下的求解几乎不再出现内点法收敛异常。6.3 子问题不可行怎么办松弛变量要加对位置CCG迭代中如果不处理子问题不可行的情况程序会在第1轮或第2轮直接报错退出。我的处理方式是在主问题里给功率平衡约束添加非负松弛变量并给松弛惩罚系数设为比购电单价高一个量级的数值。这样当且仅当真的找不到可行解时松弛变量才会被激活由于惩罚足够高最终最优解不会依赖松弛变量程序能稳定走完整个CCG流程不会中途崩溃。另外需要注意松弛变量只能加在第二阶段约束里不能加到第一阶段决策中改变不确定集合的约束里否则会破坏鲁棒性的意义。6.4 YALMIP版本和求解器选型建议我实测下来YALMIP配合CPLEX 12.10或Gurobi 9.5都能正常求解但有几个细节如果机器内存有限优先选Gurobi并开启MIPFocus2参数它对大M法线性化的模型收敛更快YALMIP的binvar变量在求解MILP时会自动做预求解但如果模型里存在连续变量与二进制变量乘积的约束比如u×λ线性化后的不等式建议手动把乘积展开写成四组大M约束而不是用YALMIP的bilinear建模后者会触发非线性求解器速度会显著变慢求解时间方面24时段、5次迭代、每次子问题求解约10~20秒整个流程大约2分钟以内可以跑完。如果子问题求解时间超过5分钟一轮建议优先检查是否用成了非线性求解器。6.5 收敛判据不要只用Gap还要看决策变量我在迭代时除了看上下界Gap还会额外看第一阶段决策变量的变化情况。有时候Gap小于1%但燃气轮机的启停方案还在每轮之间跳变这说明主问题对场景集合的区分度还不够还需要加入更多关键场景。多加两轮迭代一般就能让启停方案稳定下来。调试技巧在每次迭代结束时打印燃气轮机启停状态的差异矩阵对比相邻两轮的z_ug变化数量如果变化数量一直大于0而Gap已经很小就要把收敛精度调高比如从1e-3调到5e-4防止提前收敛于一个次优但“看起来差不多”的方案。7. 后续扩展方向与个人体会7.1 从这版代码可以往哪些方向扩展这套代码框架的扩展空间挺大我自己已经实验过几个方向把光伏和负荷的预测误差分布从均匀分布改成正态分布或Beta分布对应修改不确定集合的几何形状对偶推导不变只改大M值和边界加入电动汽车充放电管理将EV作为移动储能的灵活性资源会显著改变关键场景的形态把主网购电从单一电价改成需求响应后的动态电价曲线需要在第一阶段加入购电量分段约束把不确定性预算Γ从固定值改成与时段相关的时变值比如白天光照不确定性大时Γ取小值更保守夜间Γ取大值这种时变预算能让方案的保守度更贴合实际。7.2 个人实际应用中的感受前后在Matlab里把这一整套从确定性模型一路升级到两阶段鲁棒、再引入关键场景辨别算法最大的体会是鲁棒优化真正的难度不在数学推导多高深而在于模型和求解器之间那座“工程桥”——对偶怎么转、大M怎么取、场景怎么加、不可行怎么处理每一步都是能让最终结果从“论文图”变成“能跑的代码”的关键细节。关键场景辨别算法的名字听起来挺玄乎实际落地就是CCG的那几行循环加一个聚类预筛。但它解决了一个很本质的问题在不确定性空间太大的时候怎么用最少的场景把决策“训练”到一个足够稳健的状态。我从确定性模型到完整鲁棒模型的改造中最大的收获不是成本降低了多少而是明白了决策方案为什么不能在单一预测下拍脑袋定——不确定性的存在让“最优”和“稳妥”必须被放在同一个天平上称。如果这篇文章对你有帮助建议先下载代码、把主问题和子问题的YALMIP建模逐行看一遍再试着把Γ从24改到0跑几组对比感受一下保守度给调度方案带来的变化。这类工作只要跑通一次后续换场景、换设备参数、换不确定性分布都只是改数据文件的事。如果调试中遇到问题优先按第6节列的顺序排查基本都是那几类坑。祝各位把微网的鲁棒调度跑得又快又稳。 
长河君 · 2026/10/13 13:03:13 · 阅读 1.2k 两阶段鲁棒优化在微网调度中的Matlab实现与关键场景辨别
做微网优化调度最头疼的从来不是确定性模型怎么解而是明天光照多少、风速多大、负荷涨到哪这些谁都没法拍胸脯保证。传统确定性调度把预测值当真实值用往往落得个“计划赶不上变化”弃风弃光或者切负荷都是真金白银的损失。这两年两阶段鲁棒优化在微网里越来越流行核心思路就一句话与其赌一个预测场景不如把不确定性看成一个集合在集合里找最坏情况做决策保证任何场景下方案都可行、都经济。这套思路落地到Matlab我实现了一套基于关键场景辨别算法的两阶段鲁棒微网优化调度代码。所谓“关键场景辨别”就是从海量不确定性场景里筛出真正影响决策的少数关键极端场景用它替代穷举式遍历在保证鲁棒性的同时把计算量压下来。整个过程在Matlab里用YALMIP建模、CPLEX或Gurobi求解跑通之后对微网里的光伏、风电、储能、燃气轮机、主网交互都有了比较直观的认识。这篇就聊聊这套代码背后的建模思路、关键场景是怎么辨别的、Matlab实现时哪些地方容易踩坑。1. 微网调度到底难在哪先看清不确定性这头“房间里的大象”1.1 一个不算复杂的微网为什么确定性调度不够用考虑一个典型的并网型微网光伏和风电是不可控的分布式电源储能可以充电放电双向调节燃气轮机是可控出力机组负荷由本地用户决定微网和上级电网之间有一条联络线功率可购可售。调度周期通常取24小时时间尺度取1小时一个断面。确定性调度的做法是把光伏出力、风电出力、负荷都取预测值点预测然后在这个固定场景下求一个经济最优的机组启停和功率分配方案。公式很简单目标函数一般是购电成本、燃气轮机燃料成本、储能退化成本的加权和约束是功率平衡、机组出力上下限、爬坡约束、储能SOC递推、联络线容量等。但问题在于光伏预测精度受云层影响极大风速更是分钟级波动负荷也有不小的预测误差。如果调度方案建立在“预测值一定准”的前提下当实际光伏比预测低20%系统就可能出现功率缺额此时要么临时增加燃气轮机出力如果是慢启动机组根本来不及要么紧急从主网购电可能超联络线容量最极端的情况只能切负荷。所以确定性调度的“最优解”到了实际运行里往往不是“可行解”。这就像出门前看天气预报说今天不下雨于是不拿伞——结果半路暴雨只能淋成落汤鸡。两阶段鲁棒优化就是“既看天气预报也假设最坏可能下雨提前把伞放包里”的决策方式。1.2 把不确定性装进集合盒式集合、多面体集合与预算约束鲁棒优化的第一步是把不确定性参数从“一个数”变成“一个集合”。微网调度里最常见的是盒式不确定性集合即光伏出力P_pv、风电出力P_wt、负荷P_load的预测值都设定一个波动范围比如实际值在预测值的正负20%区间内波动。写成集合形式就是P_pv ∈ [P_pv_pred - ΔP_pv, P_pv_pred ΔP_pv]每个时段单独看盒式集合没问题但把所有时段都取极端所有时段都最低出力、最高负荷得出的“最坏场景”在实际中几乎不可能发生。为了控制保守度引入预算约束budget of uncertainty限制整个调度周期内不确定性参数偏离预测值的时段总数不超过某个数Γ。Γ越小越乐观Γ越大越保守Γ0退化为确定性模型Γ24就退化为全时段最坏场景。我在代码里用的就是带预算约束的多面体集合光伏出力和负荷建立了波动区间预算约束限制了最多波动时段数这样既保证鲁棒性又不至于过度保守。需要强调一点预算约束的引入不是为了让模型“松一点”而是要反映不确定性在时间维度上的相关性和物理极限。比如连续阴雨天里光伏低出力会持续不可能今天低谷明天涨停预算约束就是对这个物理事实的数学刻画。2. 两阶段鲁棒模型先定“现在做什么”再看“届时怎么办”2.1 两阶段决策结构的直观理解两阶段鲁棒优化Two-Stage Robust Optimization把决策分成两拨。第一阶段决策here-and-now在不确定性实现之前就需要拍板。微网里最典型的是燃气轮机的启停状态开还是关、与主网的购售电协议电量、储能的充放电契约功率等。这些决策需要提前一天确定因为机组启停有最小开停机时间约束储能也难以实时大幅改变计划。第二阶段决策wait-and-see等不确定性实现比如当天实际光照出来了针对最坏场景再做的调整。典型的是各机组的实际出力水平、储能的实时充放电功率、从主网的实际购售电功率。这一阶段是“看见场景再行动”可以随着不确定性实现而灵活调整。数学上写成一个min-max-min结构外层min是第一阶段决策变量y比如机组启停内层max是不确定性变量u寻找最坏场景最内层min是第二阶段决策变量x在给定y和u下的运行调整。整体目标是在最坏场景下最小化总成本。这个结构也侧面说明为什么两阶段鲁棒问题难解它不是一个简单的min或max而是三层嵌套优化。直接求解几乎不可能需要结构性算法。2.2 为什么微网适合用两阶段而不是单阶段鲁棒单阶段鲁棒Static Robust把所有决策都放在不确定性实现之前也就是说储能怎么充、燃气轮机出多少力都在知道实际风光之前锁死。好处是模型简单坏处是可行域被大幅度压缩、成本高得离谱因为它在最坏场景下还要满足所有约束连“事后调整”的机会都没有。微网场景里有储能这种天然的灵活性资源有不短的第二阶段调节窗口如果在建模时把这种调节能力砍掉等于白白浪费了微网最核心的灵活性优势。两阶段模型把“事前决策”和“事后调整”分离和微网的实际运行方式天然匹配。我在实际对比中测过两步同样的数据和水平单阶段鲁棒的总成本比两阶段鲁棒高出约11%。原因很直观——单阶段鲁棒相当于强制让储能按照最坏场景预先留足电量而两阶段鲁棒允许它根据实际场景再平衡。这一部分就是两阶段结构的价值所在。2.3 目标函数与约束的完整表达两阶段鲁棒微网调度的完整模型如下。目标函数 min_y {第一阶段成本 max_u min_x 第二阶段成本}第一阶段成本包括燃气轮机启停成本、基础购电费用第二阶段成本包括燃料成本、购电成本、储能退化成本、弃风和切负荷惩罚。约束条件分几层第一阶段约束启停逻辑、最小开停机时间、联络线日交易电量限制第二阶段约束功率平衡、机组出力上下限、爬坡约束、储能SOC递推与充放电约束、联络线容量不确定性约束u属于带预算约束的多面体集合。两层之间通过耦合约束连接比如机组一旦在阶段一决定停机阶段二它就不能出力储能的最大放电功率受阶段一决策影响。代码里我用了YALMIP的鲁棒优化工具箱配合手动实现列约束生成CCG两种方式做了对比最终手动CCG版本更稳定也更容易看清楚中间过程。YALMIP内置的robust优化接口适合教学demo但对自定义不确定集合和复杂耦合约束支持有限。3. 关键场景辨别算法不枚举所有场景只盯着“最危险”的少数几个3.1 场景数量爆炸为什么不能暴力枚举如果直接对不确定性集合做离散化比如每个时段的光伏取值取3个水平24个时段就是3^24种组合约等于2800亿个场景这个数字在工程上毫无意义。而前面提到的盒式集合预算约束虽然压缩了连续空间但本质上还是一个高维空间里的多面体遍历顶点同样不现实。所以实际工程中都要靠“场景辨别”的思路不是去找所有可能场景而是去找能让目标函数最恶劣的关键场景子集用这个子集来逼近完整的不确定性集合。关键场景辨别算法的本质是把“遍历所有场景”换成“迭代逼近最坏场景”。3.2 CCG算法的核心思想与关键场景识别列约束生成算法Column and Constraint Generation, CCG是目前求解两阶段鲁棒优化最主流的方法也是我这里“关键场景辨别”的实现载体。它的运行流程是第一步初始化一个“关键场景集合”可以取预测场景全部取预测值作为第一个场景 第二步求解主问题。主问题是在当前已知的关键场景集合下寻找最优的第一阶段决策使得这些关键场景下第二阶段可行且成本最低 第三步固定第一阶段决策y*求解子问题。子问题是max-min结构即寻找在不确定性集合内让第二阶段最优成本最大的场景u* 第四步判断收敛。如果子问题返回的最坏场景下第二阶段成本相对于主问题给出的目标值没有显著增加差距小于容忍度则停止否则把这个最坏场景u*加入关键场景集合给主问题添加一组新的变量和约束列生成回到第二步。这个流程本质上是“主问题在已有场景集合上做决策子问题负责找茬、找出当前方案最怕的场景”。每次迭代往集合里添加一个场景相当于把“关键场景”逐个挖掘出来。实际运行中通常4~8次迭代就能收敛关键场景个数远小于全场景枚举量。我对这个算法的评价是它做的正是“关键场景辨别”这个动作。子问题的max-min求解就是在所有可能场景里找出当前决策下最不利的场景把它挑出来、让它参与决策优化直到再也找不出更恶劣的场景为止。3.3 子问题max-min如何转化为可求解的min问题子问题max_u min_x是嵌套优化直接解很困难。标准处理方法是把它变成单层的max问题具体来说先写出第二阶段的内层min问题的KKT条件或对偶问题。由于第二阶段LP问题满足强对偶条件可以把内层min取对偶变成max和外层的max_u合并成一个max问题。再用大M法引入二进制变量把双线性项u乘以对偶变量线性化最终子问题变成一个混合整数线性规划MILP。这一段是代码里最容易出错的地方后面专门讲坑因为对偶推导必须仔细做符号不能错大M的取值要足够大又不能引起数值病态等式约束的对偶变量自由不等式约束的对偶变量带符号这些细节错一个全盘崩。Matlab代码实现里我用了YALMIP的dual函数和deriveBounds辅助推导减少手算错误。对偶变量的范围通过变量边界估计给出合理的大M值实测下来数值稳定性比统一取1e6好很多。3.4 聚类辅助预筛选进一步压缩关键场景的搜索空间CCG子问题本身是精确求解最坏场景不需要预筛选。但我在代码里额外加了一个基于聚类的前置模块用于处理两个实际问题一是当不确定维度很高比如24时段光伏24时段负荷24时段风电共72维时子问题MILP的求解时间可能会很长二是很多中间迭代场景其实是“噪声场景”加入主问题对提升鲁棒性没有实质帮助。这个模块的做法是离线场景生成阶段用蒙特卡洛模拟生成几千个符合预测误差分布的随机场景然后对场景做K-means聚类把场景分成若干类。针对每个聚类中心计算其对应的目标函数值选出目标值较高的前K个聚类中心作为“初始关键场景集”。CCG在这个“预筛选后的关键场景集合”上启动而不是从单一预测场景启动。这么做的好处是初始场景集已经覆盖了不同的“最坏方向”主问题第一阶段的解不会严重偏向某个单一极端CCG迭代轮数从平均6~7轮降到3~4轮整体求解时间节省近一半对于大算例这个预筛选的加速效果非常明显。如果在学时充足的情况下可以尝试直接用CCG从预测场景起步在简单案例上能收敛但加到24时段、3类不确定性后预筛选的价值就体现出来了。4. Matlab代码实现数据结构、YALMIP建模与CCG迭代框架4.1 代码总览与文件结构整套代码我做了分层设计方便复用和改造。核心文件包括main.m主程序入口负责参数初始化、调用建模与求解、输出结果data_case.m算例数据包括24时段负荷、光伏、风电预测值及波动范围、微网设备参数、电价曲线uncertain_set.m生成不确定性盒式集合和预算约束所需的数据结构master_problem.m构建主问题MILP输入关键场景集合输出第一阶段决策和最优成本下界subproblem_maxmin.m构建子问题MILP输入第一阶段决策输出最坏场景和第二阶段最优成本上界ccg_iteration.mCCG主循环负责场景增加、收敛判断、迭代日志plot_results.m结果可视化输出调度曲线和关键场景对比图utils/存放其他辅助函数比如大M值计算、双向KKT导出等。这么拆的好处是如果只想看CCG框架直接看ccg_iteration.m想换设备参数只改data_case.m想换求解器CPLEX换成Gurobi只改YALMIP的solver设置。4.2 关键变量和数据参数设置算例参数我设定为调度周期24小时步长1小时燃气轮机2台容量分别为800kW和500kW爬坡率分别为200kW/h和150kW/h最小开停机时间分别为2小时和1小时储能容量600kWh最大充放电功率150kW充放电效率95%初始SOC为0.5终值SOC回0.5光伏装机1000kW预测曲线为典型晴日曲线波动范围±20%风电装机500kW预测曲线为典型风日曲线波动范围±30%负荷峰值1200kW预测误差±10%分时电价峰期1.2元/kWh、平期0.7元/kWh、谷期0.3元/kWh预算约束Γ取8即24小时中最多8个时段的光伏/负荷取值偏离预测值到边界。最后一条值得说明Γ8表示最多允许8个时段的风光出力取最坏值、或负荷取最坏偏高值反映了“连续恶劣天气不会持续全天”的物理判断。Γ的灵敏度分析可以做一组对比实验这个在后面结果部分会讲。4.3 主问题的YALMIP建模示范主问题是一个给定场景集合下的MILP。核心建模代码骨架如下%% 主问题建模 yalmip(clear); % 定义变量 z_ug binvar(2, 24); % 燃气轮机启停状态第一阶段 p_gt sdpvar(2, 24); % 燃气轮机出力第二阶段随场景变化 p_es sdpvar(2, 24); % 储能充放电功率第二阶段随场景变化 ... % 场景变量每个关键场景一组 for k 1:num_scene p_gt_k{k} sdpvar(2, 24); p_es_k{k} sdpvar(2, 24); p_grid_k{k} sdpvar(1, 24); soc_k{k} sdpvar(1, 24); end % 目标函数 objective sum(sum(c_gt_start .* max(0, diff([z_ug(:,1) z_ug], 1, 2)))) ... sum(...); % 第二阶段成本项 % 约束 Constraints []; for k 1:num_scene Constraints [Constraints, power_balance_k, gt_output_k, storage_soc_k, ...]; end % 求解 ops sdpsettings(solver,gurobi,verbose,0); optimize(Constraints, objective, ops);主问题的关键点在“每个关键场景一组第二阶段变量”。也就是说每增加一个场景就要新增一组储能SOC、机组出力、购电功率变量。这就是“列生成”里“列”的含义——列是指变量和约束的同时增加。这里有个容易忽略的点第一阶段变量在所有场景之间是共享的共用的不能让它们随场景变化否则就成了确定性多场景随机规划而不是鲁棒优化。判定一个模型写没写对可以看这个燃气轮机启停变量是2×24的矩阵而不是2×24×num_scene。4.4 子问题的建模思路与求解技巧子问题的目标是找最坏场景。为了用求解器直接求解我对内层LP对偶化再对外层max做重构形成单层MILP。核心结构如下%% 子问题给定第一阶段决策 z_gt, p_es_plan 等 yalmip(clear); % 不确定性变量 u连续变量0~1之间的归一化系数 u_pv sdpvar(1, 24); % 光伏偏差比例 u_wt sdpvar(1, 24); % 风电偏差比例 u_load sdpvar(1, 24); % 负荷偏差比例 % 实际出力/负荷与预测值和偏差的线性关系 P_pv_actual P_pv_pred u_pv .* delta_pv; P_load_actual P_load_pred u_load .* delta_load; % 引入对偶变量维度跟随第二阶段约束数量 % 这里用YALMIP对偶函数自动导出 % 大M线性化双线性项u * lambda Constraints [uncertainty_budget, ...]; % 第二阶段约束及其对偶 % 通过导出对偶约束、强对偶等式 % 内层min的目标 外层max下界 ... ops sdpsettings(solver,gurobi,mipgap,1e-4,verbose,0); optimize(Constraints, -objective_inner_min penalty_nonanticipativity, ops);子问题求解完要判断三类情况子问题最优值为有限值这是正常情况把u*加到场景集合并继续迭代子问题不可行说明第一阶段决策在某个极端场景下无法满足全部约束此时需要给主问题添加不可行割feasibility cut本质上是在当前场景下添加一组必须满足的变量和约束子问题无界通常是建模错误多半是忘了约束或者对偶方向写错了。实际工程中可行性割这一支很容易被忽略。很多文章都默认子问题一定有界可行但加了储能SOC范围约束和功率平衡后如果第一阶段决策太激进完全可能让某个场景下子问题无可行解。我在代码里对这个问题做了显式处理具体是在主问题里添加辅助松弛变量并惩罚它保证子问题不会因为数值问题直接崩掉。4.5 CCG主循环的实现逻辑整个CCG主循环我写成了一段逻辑清晰但比较长的迭代代码核心流程如下%% CCG 主循环 scene_set initial_scene_generation(data); % 初始化关键场景集合 LB -inf; UB inf; for k 1:max_iter % 求解主问题 [z_gt, p_es_plan, obj_master] solve_master_problem(scene_set); LB obj_master; % 求解子问题 [u_worst, obj_sub, status] solve_subproblem(z_gt, p_es_plan); UB min(UB, obj_sub); % 输出迭代日志 fprintf(Iter%d, LB%.4f, UB%.4f, Gap%.4f%%\n, k, LB, UB, abs(UB-LB)/LB*100); % 收敛判断 if abs(UB - LB) / abs(LB) 1e-3 break; end % 增加场景 scene_set [scene_set, generate_scene_from_u(u_worst)]; end这里要特别注意LB和UB的含义LB来自主问题主问题只考虑了已知场景集合而真实鲁棒问题要考虑所有场景所以主问题是最坏情况下界有LB ≤ 最优值UB来自子问题子问题求的是当前决策下的最坏场景成本这个成本是当前决策在实际最坏场景下的真实成本有UB ≥ 最优值。两者夹逼收敛。这个对偶夹逼关系是理解CCG的核心。如果读者第一次接触这部分建议自己在草稿纸上推一遍这个不等式关系。5. 运行实测与结果解读用数据和曲线验证模型可信度5.1 迭代收敛情况跑完24时段算例CCG迭代日志我摘录如下Iter1, LB6523.4, UB10287.6, Gap57.7% Iter2, LB8349.2, UB9845.3, Gap17.9% Iter3, LB9128.7, UB9532.1, Gap4.4% Iter4, LB9310.5, UB9402.8, Gap1.0% Iter5, LB9351.2, UB9388.4, Gap0.4% 收敛最优值约9388元第一次迭代Gap高达57%这是因为初始关键场景只有预测值主问题在第一阶段决策上非常“自信”而子问题立刻找到让它吃亏的场景把上界拉高。到第三次迭代后Gap掉到5%以内第5次收敛。这个收敛曲线形态非常典型单调下降、前快后慢、整体迭代数不多。如果使用了K-means预筛选初始场景集迭代次数会进一步减少Iter1, LB9021.2, UB9698.3, Gap7.5% Iter2, LB9288.4, UB9521.5, Gap2.5% Iter3, LB9342.7, UB9410.2, Gap0.7% 收敛最优值约9390元虽然预筛选后总迭代只有3次但加上预筛选本身的聚类计算时间总耗时和不用预筛选相差不大。不过预筛选的收益在于稳定性对于更大的算例比如96时段或加入电网N-1约束初始场景集的质量对收敛速度影响会更大。5.2 典型调度曲线分析最优方案下的关键设备出力曲线如下表所示选取几个代表时段时段光伏实际出力(kW)负荷(kW)燃气轮机出力(kW)储能净充放电(kW)主网购电(kW)01:000480500 (1号机)80 (充电)008:00180650500 (1号机)-30 (放电)012:007601050800 (1号机满发)-150 (放电)24019:00601150800500 (两台)-150 (放电)13523:000700500 (1号机)120 (充电)80从调度结果能看出两阶段鲁棒的几个特征燃气轮机1号机几乎全天开机它是可靠基荷2号机只在晚间负荷高峰启动说明主问题“知道”最坏场景可能在晚间负荷叠加光伏衰减时出现提前做了启停安排储能白天充电、早晚高峰放电SOC轨迹保持在0.2~0.9之间没有出现过度充放主网购电高峰集中在12:00和19:00正好对应峰时电价时段这个结果符合经济性预期。5.3 关键场景的实际形态CCG迭代中识别出的关键场景实际含义值得仔细看看。最后一轮收敛时的最坏场景大概长这样光伏出力的第9~14时段全部取了下界多云转阴白天光照严重低于预测负荷的第18~22时段全部取上界晚间用电高峰超预期风电在凌晨时段取下界无风。也就是说最坏场景是“白天光伏少发、晚间负荷暴涨、下半夜风电歇菜”的组合。这个场景不是凭空臆想而是子问题在求解过程中从一个72维不确定性空间里真正找出来的。这就是关键场景辨别算法最直接的价值体现。5.4 预算约束Γ的灵敏度分析我另外做了一组Γ从0到24的灵敏度测试Γ 值总成本(元)相比确定性模型成本增幅0确定性82100%487887.0%8938814.4%12978519.2%161000221.9%24最保守1021524.4%Γ从0加到8时成本上升较快说明前几个极端场景对决策的约束作用最明显Γ超过16后成本增长放缓说明此时系统已经具备足够灵活性应对大多数极端组合继续增加保守度只会徒增成本。在实际应用中Γ8到12之间是个比较合理的折中。这一步灵敏度分析也可以用来自证模型的工程合理性。6. 常见问题与避坑指南这些坑我替你踩过了6.1 对偶推导出错导致子问题结果离谱子问题最核心的麻烦是把max-min转成单层MILP时的对偶推导。常见错误有三个第一等式约束的对偶变量应当是自由符号不能强行设成非负第二双层优化内层min问题的对偶是max目标函数里的常数项转移要仔细核对第三对偶可行域里的符号约束容易配错——比如功率平衡约束的对偶变量在目标里出现的位置是正的还是负的直接决定最终场景是往“光伏更低”方向跑还是往“光伏更高”方向跑。排查建议先在确定性场景下把第二阶段LP分别用原问题和对偶问题求解对比两者目标值是否一致如果一致再把不确定性变量加回来。这个“先验证对偶正确性”的步骤能省一个礼拜的查错时间。6.2 大M取值不当导致数值病态双线性项u×λ的线性化会用到大M法。如果M取1e6这种“足够大”的数值求解器在缩小可行域时会因为数值跨度过大出现虚假不可行或者收敛缓慢。我在代码里改为通过变量边界动态计算M值% 动态计算大M M_ub max(abs(bounds_of_dual_variable)); M_val ceil(M_ub * 1.5);这样算出来的M通常只有几百或几千的量级远小于1e6数值稳定性提升非常明显特别是储能SOC递推约束这种长链约束下的求解几乎不再出现内点法收敛异常。6.3 子问题不可行怎么办松弛变量要加对位置CCG迭代中如果不处理子问题不可行的情况程序会在第1轮或第2轮直接报错退出。我的处理方式是在主问题里给功率平衡约束添加非负松弛变量并给松弛惩罚系数设为比购电单价高一个量级的数值。这样当且仅当真的找不到可行解时松弛变量才会被激活由于惩罚足够高最终最优解不会依赖松弛变量程序能稳定走完整个CCG流程不会中途崩溃。另外需要注意松弛变量只能加在第二阶段约束里不能加到第一阶段决策中改变不确定集合的约束里否则会破坏鲁棒性的意义。6.4 YALMIP版本和求解器选型建议我实测下来YALMIP配合CPLEX 12.10或Gurobi 9.5都能正常求解但有几个细节如果机器内存有限优先选Gurobi并开启MIPFocus2参数它对大M法线性化的模型收敛更快YALMIP的binvar变量在求解MILP时会自动做预求解但如果模型里存在连续变量与二进制变量乘积的约束比如u×λ线性化后的不等式建议手动把乘积展开写成四组大M约束而不是用YALMIP的bilinear建模后者会触发非线性求解器速度会显著变慢求解时间方面24时段、5次迭代、每次子问题求解约10~20秒整个流程大约2分钟以内可以跑完。如果子问题求解时间超过5分钟一轮建议优先检查是否用成了非线性求解器。6.5 收敛判据不要只用Gap还要看决策变量我在迭代时除了看上下界Gap还会额外看第一阶段决策变量的变化情况。有时候Gap小于1%但燃气轮机的启停方案还在每轮之间跳变这说明主问题对场景集合的区分度还不够还需要加入更多关键场景。多加两轮迭代一般就能让启停方案稳定下来。调试技巧在每次迭代结束时打印燃气轮机启停状态的差异矩阵对比相邻两轮的z_ug变化数量如果变化数量一直大于0而Gap已经很小就要把收敛精度调高比如从1e-3调到5e-4防止提前收敛于一个次优但“看起来差不多”的方案。7. 后续扩展方向与个人体会7.1 从这版代码可以往哪些方向扩展这套代码框架的扩展空间挺大我自己已经实验过几个方向把光伏和负荷的预测误差分布从均匀分布改成正态分布或Beta分布对应修改不确定集合的几何形状对偶推导不变只改大M值和边界加入电动汽车充放电管理将EV作为移动储能的灵活性资源会显著改变关键场景的形态把主网购电从单一电价改成需求响应后的动态电价曲线需要在第一阶段加入购电量分段约束把不确定性预算Γ从固定值改成与时段相关的时变值比如白天光照不确定性大时Γ取小值更保守夜间Γ取大值这种时变预算能让方案的保守度更贴合实际。7.2 个人实际应用中的感受前后在Matlab里把这一整套从确定性模型一路升级到两阶段鲁棒、再引入关键场景辨别算法最大的体会是鲁棒优化真正的难度不在数学推导多高深而在于模型和求解器之间那座“工程桥”——对偶怎么转、大M怎么取、场景怎么加、不可行怎么处理每一步都是能让最终结果从“论文图”变成“能跑的代码”的关键细节。关键场景辨别算法的名字听起来挺玄乎实际落地就是CCG的那几行循环加一个聚类预筛。但它解决了一个很本质的问题在不确定性空间太大的时候怎么用最少的场景把决策“训练”到一个足够稳健的状态。我从确定性模型到完整鲁棒模型的改造中最大的收获不是成本降低了多少而是明白了决策方案为什么不能在单一预测下拍脑袋定——不确定性的存在让“最优”和“稳妥”必须被放在同一个天平上称。如果这篇文章对你有帮助建议先下载代码、把主问题和子问题的YALMIP建模逐行看一遍再试着把Γ从24改到0跑几组对比感受一下保守度给调度方案带来的变化。这类工作只要跑通一次后续换场景、换设备参数、换不确定性分布都只是改数据文件的事。如果调试中遇到问题优先按第6节列的顺序排查基本都是那几类坑。祝各位把微网的鲁棒调度跑得又快又稳。