我做微网调度这行也有几年了最早用确定性模型做日前计划遇到过一次挺狼狈的事傍晚光伏出力突然掉到预测值的70%储能又在早高峰把电放完了结果只能按实时电价从主网高价补电那天的运行成本比平时多出将近20%。从那之后我开始认真对待不确定性建模试过随机优化最后落在了两阶段鲁棒这条路上。后来为了缓解鲁棒优化“过于保守”的毛病引入了关键场景辨别算法效果比较理想。今天这篇就来好好拆解这套“基于关键场景辨别算法的两阶段鲁棒微网优化调度”方案——包括建模思路、求解流程和Matlab实现经验把能落地的细节都写出来。适合正在做微电网日前调度、或者刚接触两阶段鲁棒优化的研究生和工程师参考。1. 微网日前调度的不确定性为什么确定性模型不够用1.1 一个让我下定决心改模型的真实算例我负责的微网不算大包含一台2MW的微型燃气轮机、1.2MWh储能、2MW风电、1MW光伏以及一条与主网交互上限2MW的联络线。调度周期是24小时时间分辨率1小时。最早版本就是教科书式的确定性经济调度把风电、光伏、负荷都取预测值然后求解一个线性规划目标函数是燃料成本加购电成本最小。平时这事没什么问题直到有一次晚高峰光伏骤降。当天的光伏预测显示傍晚还能有0.6MW出力实际只有0.18MW少了0.42MW。储能按照确定性方案在凌晨低谷给充满、早高峰放空了傍晚根本没有余量。燃气轮机已经压到上限联络线也满了最后只能靠临时可中断负荷协议硬切了0.25MW负荷代价是按合同赔付加上从主网高价补电的费用当日总成本比正常高出近20%。这个案例让我意识到一件事微网日前调度的瓶颈往往不是设备模型有多复杂而是你对手里预测数据的“信任程度”有问题。确定性模型把预测值当作确定值相当于默认“风光预测量多少实际就发多少”。但实际上风电出力在日前时间尺度上的预测误差通常在±15%~±25%之间光伏在阴雨天误差更大负荷也不是刚性不变的。误差一旦在某个时段集中暴露后续所有机组出力和储能安排都会跟着失衡。1.2 不确定性从哪里来以及传统方法各自的短板微网的不确定性来源主要有三个可再生能源出力波动主要是风电和光伏。风电的随机性来自风速的间歇性光伏来自云层遮挡和天气突变。负荷需求波动尤其是有居民和商业负荷混合的微网日常峰谷变化之外还有突发性的负荷启停。市场价格波动对并网型微网来说从主网购电的实时市场价格本身也是不确定的。针对这些不确定性传统做法有两条路线。一条是随机规划需要假设不确定参数的概率分布然后通过蒙特卡洛抽样或场景树来描述未来可能发生的情况。随机规划的优点是经济性较好因为它按概率加权平均求期望成本最低缺点是它对概率分布的准确性极其敏感——如果你假设的正态分布和实际偏态误差对不上算出来的方案会很“脆”。而且场景数量一多模型规模呈指数膨胀求解时间让人抓狂。另一条就是鲁棒优化。鲁棒优化不依赖概率分布只需要知道不确定参数的取值范围然后去求“在最坏的情况下怎么让自己损失最小”。这个思路从工程直觉上完全说得通也正好符合调度人员“先把最困难的情况兜住再谈优化”的心态。但它的经典版本有个被很多人诟病的问题过度保守。因为最坏情况往往意味着风最小、光最小、负荷最大同时发生而这种组合实际出现的概率极低。为了这个几乎不出现的场景你会牺牲大量日常运行经济性。所以正确的做法不是“要不要用鲁棒优化”而是“用多保守的鲁棒优化”。这就引出了今天文章的核心关键场景辨别算法。它要做的事就是帮你从庞大的不确定集合里找出那些“真正配得上被当作最坏场景”的边界点然后只针对这些关键场景做鲁棒保护。2. 两阶段鲁棒模型的数学骨架min-max-min怎么落到微网上2.1 第一阶段和第二阶段决策的划分依据两阶段鲁棒的核心思想是把决策按照“看到不确定性之前”和“看到不确定性之后”分为两类。拿微网日前调度举例。第一阶段决策是必须在不确定性实现之前就要拍板的业界常称为 here-and-now 决策。典型的包括燃气轮机的开机/停机状态、联络线购售电状态、储能的充放电状态逻辑什么时候允许充、什么时候允许放。这些决策的特点是调整代价极高比如燃气轮机启动需要时间、启停次数还有寿命损耗不可能根据光伏每五分钟波动一次就频繁切换。所以它们要在“不知道明天天气怎么样”的时候定下来。第二阶段决策是在不确定性实现之后可以快速响应调整的也就是 wait-and-see 决策。典型包括燃气轮机的实际出力增量、储能的实际充放电功率、弃风弃光量、切负荷量。这些变量在一天之内可以平滑调节所以它们可以等“看清了风光实际出力”再决定。这种“先定大框架、再随机应变”的思路比起单阶段鲁棒优化更贴近微网的实际运行逻辑。在数学上它天然形成一个 min-max-min 结构外层 min最小化第一阶段固定成本 最坏场景下的第二阶段调整成本中间 max在最坏的自然对手——即不确定参数向不利于我们的方向取值时让第二阶段成本尽量高内层 min在给定的最坏场景下调度员通过调节第二阶段变量把运行成本压到最低。这个结构就是两阶段鲁棒优化最基本的样子。2.2 模型方程与box不确定集合用公式写出来会更清晰。我这里用的是一个常见的微网抽象模型变量含义如下x第一阶段决策变量包括燃气轮机启停状态、购售电状态等y第二阶段决策变量包括燃气轮机出力调整、储能充放电功率、弃风弃光量、切负荷量u不确定参数具体为风电出力、光伏出力、负荷的任意组合目标函数min [ c1^T x max_{u∈U} min_{y∈Ω(x,u)} ( c2^T y 惩罚项 ) ]约束条件分成三个部分A x ≥ d表示只有第一阶段变量的约束比如燃气轮机最小运行状态约束、联络线购售电状态互斥约束B y C u ≤ h - D x表示第二阶段变量和不确定参数共同满足的运行约束比如功率平衡、储能SOC递推、线路潮流限制y ∈ Y(x,u)表示第二阶段变量的可行域包括出力上下限、爬坡约束等。这里的不确定集合U最常用的就是box不确定集。以风电为例U_w { P_w | P_w^pre - Δ_w^- ≤ P_w ≤ P_w^pre Δ_w^ }其中 P_w^pre 是预测值Δ 是允许的最大偏差。每个时段都可以定义各自的上下界。box集合的好处是结构简单、对偶变换干净CCG算法写起来不复杂。但刚才说过真正的难点不在box本身而在box怎么定——是全取理论边界还是用数据识别出来的关键场景来压缩边界。这正是下一节重点展开的内容。需要强调的是鲁棒优化并不排斥第二阶段变量在多个时段之间的耦合。比如储能SOC是跨时段的递推变量在最坏场景下子问题里的min-y本质上是一个跨24小时的动态优化而不是t个独立的单时段问题。这一点在Matlab实现里一定要留心否则很容易把储能的时序关系写错导致结果毫无物理意义。3. 关键场景辨别算法它解决的是鲁棒优化“过度保守”的问题3.1 全集合鲁棒为什么保守最恶劣场景几乎不会出现理解关键场景辨别算法首先要接受一个反直觉的事实最坏情况不是“一个”场景而是“一堆”边界点组成的集合。在box不确定集上风电、光伏、负荷各有上下界夹在一起最坏的组合方式不是唯一的。经典两阶段鲁棒在处理的时候内层max会让不确定参数向使成本最大化的方向靠最后算出来的“最坏场景”往往是所有不利偏差同时取极值的那个点。问题恰恰出在这里。举例来说如果风电预测是12MWbox给到±20%那就是风电最低9.6MW光伏预测在白天有6MWbox给到±15%最低5.1MW负荷预测峰值3.2MWbox给到±5%最高3.36MW。三者同时取最不利值意味着一个“大风力骤减同时阴天遮光同时负荷暴涨”的世界。这种组合不是说绝对不可能但概率微乎其微。标准鲁棒优化的缺点就是它不为这个概率操心只要这个场景属于不确定集合它就强迫你为它预留资源。结果就是燃气轮机多开、储能提前预留、联络线少售电日常运行成本被白白抬高。这也是很多工程师对鲁棒优化敬而远之的原因。3.2 关键场景辨别的三类常用实现路径所谓关键场景辨别本质上是给鲁棒优化加一个“数据筛选”的前置步骤从不确定集合中挑出那些对调度决策真正构成压力、真正处于约束边界附近的场景忽略那些虽然存在但影响不大的场景。工程上我见过三类做法第一类按偏差幅度排序筛选。把历史N天的风电、光伏、负荷预测误差按绝对值大小排序取某个置信水平比如95%以内的偏差作为box边界超过的分位数直接砍掉。这是最简单的做法优点是快缺点是它只关注了单变量的幅度没考虑变量之间的相关性和时间耦合。比如风电在15点的偏差大但那个时段负荷本身很低系统根本不在乎而光伏在18点的偏差只有10%却正好撞上晚高峰风险反而更大。第二类按运行压力指标筛选。先固定一个初始调度策略然后逐一对每个历史场景计算“运行压力”比如净负荷缺口、储能SOC越限程度、备用不足量等。压力超过阈值的场景被标记为关键场景之后再在这些关键场景上重新优化调度策略。这个做法的优点是与系统实际运行高度相关挑出来的场景确实是会让调度员头疼的场景缺点是要设计一个合理的压力指标指标选得不好筛选效果就打折扣。第三类结合CCG迭代轨迹进行关键场景辨识。在两阶段鲁棒求解过程中子问题每次都会给出一个当前最恶劣场景;随着迭代次数增加这些最恶劣场景会逐渐稳定在少数几个边界点上。我们把迭代过程中出现过的场景记录下来统计哪些场景反复出现、哪些场景导致过约束触发的对偶乘子显著不为零这些就是系统中的“天然关键场景”。这种做法的好处是不需要额外的压力指标直接用优化算法本身的反馈来判断关键性是我个人最推荐的一种。3.3 我用Matlab实现时选择的指标与筛选流程在这个项目里我用的是第二类加第三类的结合原因很简单第二类做法的物理意义直观第三类做法的自动化程度高。具体流程如下第一步从历史数据中提取预测误差样本。我取了过去90天的风电、光伏、负荷预测值和实际值按小时计算每个时段的偏差形成90×24×3的误差矩阵。然后用核密度估计把这些误差拟合成联合分布再抽样生成2000个初始场景。第二步计算每个场景的严重度。我选了一个物理意义很直接的压力指标固定一组参考调度策略后每个场景下的净负荷缺口总量与储能被迫越限量之和。净负荷缺口越大说明这个场景越容易造成供需失衡储能被迫越限说明系统已经动用了最后的调节手段。对2000个场景计算这个指标后按降序排列。第三步抽取TOP-N场景作为关键场景。我取了前50个然后检查这些场景覆盖的逐时段不确定性区间用它们重新构造压缩后的box不确定集。压缩后的box不是简单取50个场景中每个时段的min/max而是在此基础上保留5%的余量防止箱型边界贴得太紧导致预测误差稍微超界就出问题。第四步将压缩后的不确定集交给两阶段鲁棒模型在CCG迭代中求解。迭代结束后再回头比较子问题实际识别出的最恶劣场景是否都在50个关键场景集合内。如果出现了集合外的场景说明TOP-N取少了适当扩大N再算一次。这套流程我跑下来效果比较理想但有一点要提醒指标设计直接影响筛选质量没有一种指标是普适的。如果你的微网是孤岛型切负荷代价极高压力指标里就该加大切负荷惩罚如果你是并网型且电价波动大压力指标里就该突出高价购电风险。别指望一个指标打天下。4. CCG迭代求解主问题、子问题与关键场景如何衔接4.1 主问题MP带上已发现场景的扩展模型两阶段鲁棒模型不能直接喂给求解器因为max-min的内层结构不是商业求解器能直接处理的形式。标准的解法是CCG列与约束生成Column-and-Constraint Generation它把原问题拆成一个主问题和一个子问题通过迭代不断逼近最优解。主问题的形式是在已知一系列“候选最恶劣场景”u^(1), u^(2), ..., u^(k) 的前提下求解一个包含第一阶段变量x和第二阶段变量y^(1), y^(2), ..., y^(k) 的扩展模型。每个候选场景对应一组完整的第二阶段变量和约束目标函数是min c1^T x η其中η是一个辅助变量用来表示最坏场景下的第二阶段成本上界。约束除了第一阶段约束外还要对每个已知场景u^(k)加上η ≥ c2^T y^(k) B y^(k) C u^(k) ≤ h - D x y^(k) ∈ Y(x, u^(k))这里要注意每个y^(k)都带有场景索引但它们共享同一组x。共享x正是两阶段结构的体现第一阶段决策只能有一个不能随着场景改变第二阶段决策可以随场景变化而变化。主问题是一个标准的MILP如果含0-1变量或LP如果不含0-1变量可以直接用Gurobi或Cplex求解。求解主问题得到的x^(k)就是当前迭代下的最优第一阶段方案。4.2 子问题SP最恶劣场景识别与对偶变换子问题做的事情是固定第一阶段决策x x^(k)然后求解max_{u∈U} min_{y∈Ω(x^(k),u)} c2^T y这个max-min结构同样不能直接求解。但好消息是内层min是一个关于y的线性规划LP只要y的可行域是有界的就可以用强对偶定理把min问题转成max问题。对偶之后原来的max-min就变成了一个带双线性项的单层max问题max_{u, λ} ( λ^T (h - D x - C u) )s.t. λ^T B ≥ c2u ∈ U, λ ∈ Λ其中λ是对偶变量。问题在于约束里出现了λ和u相乘的项 λ^T C u这是一个双线性项非凸不能直接交给求解器。处理双线性项有两条常见路径。第一条是大M线性化。如果u是box连续变量把u离散成顶点枚举用二进制变量α表示是否取某个顶点然后用大M法把 λ·u 乘积线性化。第二条是枚举顶点。对于box不确定集u的最优解一定在box的顶点上。每个顶点对应一个确定的u值把顶点提前枚举出来分别求解只关于λ的LP取其中目标值最大的顶点作为最坏场景。第二条路径在不确定参数维度不高时非常有效我的算例里风光负荷三个不确定源24时段顶点数虽然理论上有2^(72)个但结合关键场景辨别的压缩实际需要枚举的候选顶点很少。4.3 迭代流程与收敛判据CCG的完整迭代流程我写成下面这样方便直接照着写代码初始化下界LB-∞上界UB∞迭代计数k1初始场景集合为空或放入一个保守的初始场景比如所有不确定参数取中间值。求解主问题MP得到最优目标值obj_MP和第一阶段解x^(k)。更新下界LB obj_MP。固定x x^(k)求解子问题SP得到最恶劣场景u^(k**)和最优值obj_SP。更新上界UB min(UB, c1^T x^(k) obj_SP)。检查收敛条件如果UB - LB ≤ ε停止迭代输出x^(k) 和UB作为最优调度方案。否则把新发现的场景u^(k**) 以及对应的第二阶段变量y^(k**) 加入主问题的约束集合k k 1返回步骤2。关于收敛判据工程上我一般设0.5%的相对间隙或者绝对间隙50元人民币。千万别设成0那会让迭代次数陡增而且意义不大。实际中CCG在这个问题上的收敛速度通常很快——标准两阶段鲁棒大概在5~9次迭代内收敛关键场景辨别后不确定集更紧凑往往3~6次就够了。关键场景辨别在这里的衔接点是初始场景集合不空手启动而是用关键场景集合中压力最大的几个场景作为初始u。这样做的好处至少有两个一是主问题第一轮就能给出一个比较合理的x后续子问题识别出的新场景和初始场景高度重合收敛自然快二是避免子问题在最开始因为x太离谱而给出一些极端但无关紧要的场景把迭代带偏。5. Matlab代码实现框架与求解器配置5.1 变量与数据结构设计文章标题里提到了Matlab代码实现这里我必须多说几句。Matlab里做这类问题我不推荐直接用linprog或intlinprog拼模型因为两阶段鲁棒模型的约束和变量会随着CCG迭代不断动态增加原生工具箱处理起来很别扭。我用的是YALMIP Gurobi/Cplex的组合。YALMIP负责建模底层求解器负责算动态加约束用YALMIP的约束列表非常简单。变量定义上我建议把一阶段、二阶段、不确定参数的索引分开管理否则迭代几次之后自己都容易绕晕。我的习惯是T 24; % 调度时段数按小时 % 一阶段变量 u_MT_on binvar(1, T, full); % 燃气轮机开停状态 v_buy binvar(1, T, full); % 购电状态 v_sell binvar(1, T, full); % 售电状态 % 二阶段变量每个候选场景一组用元胞数组存 y_MT cell(1, Kmax); % 燃气轮机出力调整量 y_P_cha cell(1, Kmax); % 储能充电功率 y_P_dis cell(1, Kmax); % 储能放电功率 y_curtail_w cell(1, Kmax); % 弃风 y_curtail_pv cell(1, Kmax); % 弃光 y_shed cell(1, Kmax); % 切负荷量不确定参数本身在子问题里是决策变量在主问题里是已知常数。所以我不会把u定义成全局变量而是在构造子问题时单独用sdpvar定义并施加box约束。5.2 主问题和子问题的代码骨架主问题的构造思路把每一轮发现的新场景u_found存成一个矩阵U_found每一列是一个场景的24时段风光出力负荷取值。然后在主问题里对每个场景生成一组二阶段变量和约束最后把所有约束拼起来。核心代码骨架大致是% 主问题约束集合 MP_cons []; MP_cons [MP_cons, ...]; % 第一阶段约束、启停逻辑约束等 for k 1:size(U_found, 2) u_k U_found(:, k); % 提取第k个场景 % 为第k个场景定义二阶段变量 y_MT{k} sdpvar(1, T, full); y_P_dis{k} sdpvar(1, T, full); ... % 添加功率平衡约束 MP_cons [MP_cons, ... y_MT{k} y_P_dis{k}*eff_dis - y_P_cha{k} u_k(1,:) u_k(2,:) ... P_grid{k} u_k(3,:) y_curtail_w{k} y_curtail_pv{k}]; % 添加储能时序约束 MP_cons [MP_cons, SOC{k}(t1) SOC{k}(t) y_P_cha{k}(t)*eff_cha - y_P_dis{k}(t)/eff_dis]; ... % 添加η ≥ 二阶段成本约束 MP_cons [MP_cons, eta sum( ... )]; end % 目标函数一阶段成本 eta MP_obj sum(MT_cost * u_MT_on) sum(buy_price .* v_buy .* ...) eta; ops sdpsettings(solver, gurobi, verbose, 0); sol_MP optimize(MP_cons, MP_obj, ops);子问题的构造就要小心了。固定x之后需要把一阶段变量的值代入约束然后把内层min对偶化。实现上有两种方式手动对偶把内层LP写成对偶形式再把外层max和对偶后的max合并处理双线性项后求解。直接用YALMIP的dualize函数自动生成对偶问题这个省事但有时候数值稳定性不如手动写。我实际用的是手动对偶因为可以完全控制双线性项的线性化方式。关键代码片段如下% 固定x后把一阶段参数代入 x_fixed value(x); % 定义不确定参数u作为sdpvar u_w sdpvar(1, T, full); u_pv sdpvar(1, T, full); u_load sdpvar(1, T, full); % box约束关键场景辨别后压缩的边界 SP_cons [u_w U_w_lb; u_w U_w_ub; ...]; % 内层y的对偶问题变量是对偶乘子lambda lambda sdpvar(size(..., 1), 1, full); SP_obj lambda * (h_vec - D * x_fixed - C * u_vector); % 双线性项线性化大M法 M_val 5 * max(max(abs(U_ub - U_lb)), 1); w sdpvar(length(lambda), 3*T, full); % lambda与u的乘积项 z binvar(3*T, 1, full); % 顶点选择变量 for i 1:length(lambda) for j 1:3*T SP_cons [SP_cons, w(i,j) U_lb(j)*z(j)]; SP_cons [SP_cons, w(i,j) U_ub(j)*z(j)]; SP_cons [SP_cons, w(i,j) lambda(i) - U_ub(j)*(1-z(j))]; SP_cons [SP_cons, w(i,j) lambda(i) - U_lb(j)*(1-z(j))]; end end % 子问题目标max对偶目标 SP_obj -sum(sum(w .* ... )) ...; sol_SP optimize(SP_cons, -SP_obj, ops);代码里有个细节值得说明大M线性化处理lambda连续变量和u二进制变量的乘积时我把u的连续变量离散到了box顶点上。这个离散化对于box不确定集是精确的因为线性目标函数在box的极值点取得最优。如果x固定后内层y的约束不是纯线性比如含有二次项那就要另行处理但微网调度线性化程度已经足够高这个限制基本不影响工程实践。5.3 求解器选择与参数调优我用的是Gurobi原因很简单主问题是MILPGurobi的MIP求解器在同类问题上的速度比Cplex快不少而且许可相对容易拿到。如果没有Gurobi用Cplex也行YALMIP的建模代码几乎不用改只需要把solver参数从gurobi改成cplex。有几个求解器参数我基本上是固定设置的ops sdpsettings(solver, gurobi, verbose, 0); ops.gurobi.MIPGap 1e-3; % MIP相对间隙 ops.gurobi.TimeLimit 300; % 主问题单次求解时间上限300秒 ops.gurobi.NumericFocus 1; % 数值稳定性优先这里特别想强调TimeLimit。CCG每次迭代都要重求一次主问题如果不设时间限制某个迭代轮次可能因为MIP分支定界太深而卡死10分钟整体算下来完全不可接受。我设了300秒上限配合0.5%的目标间隙整体算例通常在25分钟内稳定收敛。如果你只是小规模测试直接把TimeLimit设成60秒就够了。6. 算例对比关键场景辨别在成本与可靠性之间找到了平衡6.1 算例系统与参数为了把三种方案的差别说清楚我设计了一个可以直接对比的算例。系统参数如下设备参数数值燃气轮机额定出力2 MW燃气轮机最小出力0.4 MW燃气轮机单位燃料成本42 元/MWh储能容量1.2 MWh储能最大充/放电功率0.4 MW储能充放电效率0.95风电装机容量2 MW光伏装机容量1 MW联络线购售电上限2 MW切负荷惩罚成本8000 元/MWh日前风电预测误差box偏差±20%光伏预测误差box偏差±15%负荷预测误差box偏差±5%分时电价峰时08:00-11:00、18:00-21:001.2元/kWh平时0.8元/kWh谷时23:00-06:000.4元/kWh。关键场景辨别环节我按第3.3节的流程从2000个抽样场景中筛出TOP-50作为关键场景并据此把风电box压缩到约±14%光伏压缩到约±9%负荷不变±5%因为负荷预测精度还好。压缩后的不确定集合覆盖了原本约80%的联合分布范围这个覆盖率数据和分位数设定匹配。6.2 三种方案的运行结果与成本分析三种方案分别是方案A确定性优化使用预测值不考虑不确定性。方案B标准两阶段鲁棒优化box不确定集取理论全边界不压缩。方案C关键场景辨别 两阶段鲁棒优化不确定集来自关键场景统计。我从三个维度对比结果日均运行成本、最恶劣场景下的成本、典型极端场景下的切负荷量。方案日均运行成本元最坏场景成本元约束违反/切负荷情况A 确定性2340031600典型极端场景下切负荷0.25 MWhB 标准两阶段鲁棒2890027800所有box场景内均不切负荷C 关键场景辨别鲁棒2540027100在压缩box内基本不切负荷解释一下这张表方案A在预测值下成本最低只有23400元但一旦遇到“光伏低出力撞晚高峰”的典型极端场景实际成本会飙到31600元多出来的8200元就是高价购电和切负荷赔偿。方案B把所有理论极端场景都兜住了代价是日常成本抬升到28900元比确定性高23.5%。方案C把不切负荷的保证范围压缩到关键场景覆盖区间日常成本只比确定性高8.5%但在关键场景覆盖范围内依然足够可靠。我自己的分析是方案B看似最安全实际是在为一个几乎不会发生的多变量同时极端事件每年多付上百万元成本。方案C通过关键场景辨别把保护范围从“理论所有可能的世界”缩到“数据上有意义的世界”在可靠性和经济性之间找到了一个漂亮的平衡点。迭代收敛情况也值得记录方法迭代次数总耗时秒B 标准两阶段鲁棒81520C 关键场景辨别鲁棒4830方案C不仅成本更低求解还更快。这是因为压缩后的不确定集让子问题需要枚举的候选顶点大幅减少主问题每轮新增的场景也更少整体计算负担明显下降。在实际工程里这意味着当日可以更快拿到优化结果给调度员留出更多复核时间。这个算例也印证了我一直以来的观点鲁棒优化的保守性不是算法带来的而是不确定集合边界定义带来的。解同样的鲁棒模型你把box画得更大方案就必然更保守你把box画得合理方案就务实得多。关键场景辨别算法的本质就是亲手把这个“box”画得更贴近现实。7. 实际运行踩坑记录与效率调优建议7.1 子问题对偶不可行与无界问题CCG第一次跑通之前我被两个问题折磨了很久这里单独拎出来说。第一个是对偶不可行。子问题内层LP的对偶成立需要满足一个前提原始问题可行且对偶可行。但如果第一阶段决策x给得太离谱比如储能SOC在某个时段已经被推到了上限之外内层min问题本身就不可行对偶问题就是不可行或无界状态。求解器报infeasible或unbounded整个CCG就不知道该怎么继续。我踩坑之后的做法是在模型里给所有关键运行约束加上一个松弛变量并给松弛量配上极高的惩罚系数。比如储能SOC上下限约束改成SOC_min - s_{soc,t} ≤ SOC(t) ≤ SOC_max s_{soc,t}其中 s_{soc,t} ≥ 0惩罚系数取10000元/MWh。这个松弛项不影响正常情况下的解——因为惩罚够高正常时没人会去用它但一旦子问题遇到x导致的越限它就不会直接返回infeasible而是返回一个带惩罚成本的结果CCG可以继续迭代。这招在工程稳定性上非常管用。第二个是双线性项的大M取值问题。M取小了会错误地把可行域“削掉”一块导致解出来的所谓最坏场景不是真正的最坏M取大了数值精度下降求解器在分支定界时会不断报数值警告。我的经验是M值按物理意义来定比如不确定参数的偏差上限本身再加一个10%的裕度比随便给一个1e6强得多。7.2 大M取值、求解器数值稳定性与代码容错关于大M具体再展开一点。我在第5.2节的代码里写的是M_val 5 * max(...)但实际项目里我应该把M取得和物理范围严格对应。比如lambda是某个平衡约束的对偶乘子它的量纲是元/MWh实际值通常在几十元上下u是风电出力范围在0~2MW。lambda和u的乘积最坏量级是几十元×2MW所以大M取100就完全够用而不是随便写个1000。给每个lambda配各自的M_val比全部共用一个万能M稳定得多。数值稳定性方面我建议把所有目标函数和约束里的系数统一化到同一量级。燃料成本如果按元/MWh计电价也按元/kWh会差1000倍混合起来求解器可能出现数值问题。我的习惯是电价统一换成元/MWh这样整个目标函数的系数都落在几十到几百这个量级Gurobi内部preSolve就不容易出幺蛾子。代码容错还有一个细节每次求解后检查sol.info。如果求解器的输出是infeasible别急着看结果先把当前迭代的主问题/子问题约束导出成LP文件或Matlab的debug信息快速定位是哪条约束把可行域挖空的。我一般用如下代码if sol_MP.problem ~ 0 warning(主问题求解失败%s, sol_MP.info); save(dump_MP.mat, MP_cons, MP_obj, ops); end这段代码看起来简单但真的能帮你省下大量排查时间。7.3 计算效率的优化思路效率调优这块我总结几条实际有效的经验初始场景集合别空着。CCG空集启动第一轮主问题往往很差导致子问题发现一些“过渡性最坏场景”。这些场景在后面迭代里会被淘汰但白白增加了主问题的变量数。我用关键场景中压力最大的5个场景做初始集合迭代次数明显下降。场景入库去重。迭代过程中子问题可能反复给出同一个场景新场景加入主问题前先判断一下是否已经存在。我按参数向量做绝对差判重阈值取1e-6几百次迭代下来能少加不少重复约束。主问题提供初始可行解。每次迭代主问题求解前直接把上一轮的一阶段解作为初始解传入Gurobi的x0字段MIP求解耗时能再降20%左右。这个操作对大规模微网尤其划算。子问题按并行拆分时段。如果局部不确定参数之间的耦合很弱可以把24小时切成几个独立时段分别求子问题再取最大值。不过储能SOC跨时段约束会让这种拆分变复杂我的建议是先把串行版本跑稳定再考虑并行别一上来就搞复杂化。最后再分享一个我在实际运行中的体会关键场景辨别算法并不是要把鲁棒优化变得很复杂相反它追求的是“少算点没用的场景多保护点该保护的运行边界”。如果你正在做微网调度也不一定要一上来就上全套两阶段鲁棒。可以先拿历史数据算一算看看你系统里真正贵的不确定性是什么——是风光波动大还是负荷峰谷差大然后针对最贵的那个不确定性来源设计关键场景筛选指标。这个思路比堆满一整页数学公式更有工程价值。