☰
基于关键场景辨别算法的两阶段鲁棒微网优化调度Matlab实现
2026/9/26 7:12:35 网站建设 项目流程

做微网调度这些年,我最深的一个体会是:调度方案的“可用性”比“好看”重要得多。风光预测曲线画得再漂亮,实际运行的时候风速一变、云层一压,原计划当场失效。这也是为什么两阶段鲁棒优化在微网优化调度里越来越受重视——它不赌预测准不准,而是先把不确定性的边界画出来,然后在最坏的情况下也能给你一套兜底可用的调度策略。

这篇文章就围绕“基于关键场景辨别算法的两阶段鲁棒微网优化调度(Matlab代码实现)”这个具体的案例来拆。核心关键词拆开是三块:微网系统建模、两阶段鲁棒优化框架、用来识别恶劣场景的关键场景辨别算法。我会把数学模型怎么搭、主问题子问题怎么迭代、Matlab里怎么一步步实现、以及调试过程中那些容易翻车的细节全部讲透,给想自己动手复现的同行一条能走通的路。

1. 从“平均最优”到“最坏情况最优”:两阶段鲁棒优化的设计逻辑

1.1 为什么单阶段确定性调度在微网面前撑不住

传统的确定性调度本质上是在做一件事:把预测值当成真实值,求解一个单层优化问题。你早上根据负荷预测曲线、光伏出力预测曲线、风电出力预测曲线,再加上电价信息,求解出未来24小时各机组出力和储能充放电计划,这个计划在预测完全准确的前提下当然是最优的。

问题是微网里的不确定性源太多了。光伏出力受云层遮挡影响,分钟级波动可以超过装机容量的50%;风电出力更是看天吃饭,风速预测偏差10%,功率偏差可以到30%以上;负荷侧的大型冲击性负载一旦启动,预测曲线直接失真。实际运行的时候,你按预测值定好的储能充放电计划,很可能让系统出现功率失衡,这时候要么弃风弃光,要么切负荷,要么高价从上级电网买电,经济性和可靠性双双失守。

两阶段鲁棒优化的思路完全不一样。它把决策拆成两段:第一阶段是在不确定性还没暴露之前就得拍板的变量,比如机组启停状态、储能充放电基准计划;第二阶段是等到不确定性实现之后,系统还能做的“矫正动作”,比如调整燃气轮机出力、调用储能余量、调节与主网的交换功率。优化目标不再追求预测场景下的期望最优,而是确保在最恶劣的不确定性实现下,第二阶段的调整成本依然可控、系统依然可行。

用大白话讲,确定性调度像是你根据天气预报决定明天穿什么,赌的是预报准;两阶段鲁棒像是你既准备了一件薄外套,又包里塞了一件厚冲锋衣——天气真变了,你还有第二层调整手段。代价是方案会稍微“保守”一点,但换来的是极强的抗风险能力。

1.2 两阶段鲁棒优化的数学骨架

两阶段鲁棒优化的一般形式长这样:

min c1^T x + max u∈U min y c2^T y

这个式子要拆开看。x 是第一阶段决策变量,对应那些一旦确定就很难改变的物理量,比如机组启停(启动后不能立刻停机)、储能的日前充放电计划。这批变量的特点是“现在定,后面不能改”,术语叫 here-and-now 决策。

u 是不确定变量,也就是我们关心的光伏出力、风电出力、负荷功率。它落在不确定集 U 里面,U 刻画了我们对“不确定性到底有多坏”的认知边界。max 这一层的意思是:既然你不知道真实场景是哪个,那就假设命运会跟你对着干,在所有可能的 u 里挑一个让你的成本最高的。这不是悲观主义,而是鲁棒优化的基本假设——我做的方案必须经得起这种最坏情况的考验。

y 是第二阶段决策变量,术语叫 wait-and-see。等你已经知道了 u 的真实取值,系统就能根据实际情况做调整,比如让储能多放一点电、让燃气轮机多带一点负荷。min y c2^T y 就是在最坏 u 下,第二阶段的调整成本能做到多小。

换个角度理解:这个三层结构是对鲁棒性的一个“压力测试”。先做一个决策 x,然后让不确定性尽量刁难你,最后再看你还有多少补救能力。如果一个调度计划在这种连环压力测试下依然可行且成本可接受,那它在真实运行中基本上就稳了。

综合起来,整个优化问题在做的是:找一个第一阶段决策 x,使得“第一阶段的固定成本 + 最坏情况下的第二阶段最小调整成本”达到最小。

1.3 不确定集:给“最坏情况”划好边界

两阶段鲁棒里最关键的一个工程设计是 U 的选取。U 太小,鲁棒性形同虚设;U 太大,调度过于保守,经济性崩盘。实际中最常用的是盒式不确定集加预算约束。

用风电出力来举例。设风电场预测出力是 P_w^forecast,实际出力可以写为:

P_w,actual = P_w,forecast + ΔP_w

其中 ΔP_w 是预测误差,它被限制在一个区间内:

ΔP_w ∈ [-ΔP_w^max, ΔP_w^max]

这就是“盒式”的意思——把每个不确定变量的偏差都框在一个矩形区间里。如果只做到这一步,最坏情况就是所有风电机组同时偏到区间下限,这在实际中几乎不可能发生,方案会极度保守。

预算约束的作用就是控制“同时在边界上的不确定变量数量”。比如微网里有10台分布式风机,预算 Γ=3,那最坏情况只允许至多3台风机同时取到区间端点的偏差值,剩余7台的偏差只能落在区间内部。这样 U 的形状从一个“满盒子”被切削成了“带棱角的超多面体”,对应的最坏场景更贴近物理现实,解的经济性也明显改善。

用公式写预算约束就是:

∑ |ΔP_w| / ΔP_w^max ≤ Γ

直观理解:Γ 越大,系统越保守;Γ=0 时退化为确定性模型;Γ 取到不确定变量总数时就是纯盒式。实际调试中,Γ 是调节方案“胆量”的核心旋钮。微网里风光渗透率越高,Γ 通常要取得越大,因为不确定性的相关性更强,极端情况的破坏力更大。

2. 关键场景辨别算法:如何从连续不确定集里精准锁定“最狠”的场景

2.1 问题的麻烦:不确定集是连续的,没法直接枚举

把两阶段鲁棒模型摆出来之后,第一个让人头疼的问题就是:max u∈U 这一层怎么算。

如果 U 是一个盒式连续区间,这就意味着有无限多个可能的 u。一个由有限离散场景构成的集合可以用枚举来求解,但连续集合没法一轮轮全部列出来。更麻烦的是,子问题是一个 max-min 结构:外层让 u 尽量坏,内层让 y 尽量补救,两股力量互相拉扯,直接求解本身就有难度。

这时候就得用到“生成”的思想:我不需要把 U 里所有场景都拿出来参与计算,我只关心那些对当前解打击最大的关键场景。每一轮迭代我固定住第一阶段的 x,去求解子问题,找到一个让第二阶段成本最大化的 u——这个 u 就是当前决策下的“最恶劣场景”。然后把 u 作为新场景加入主问题的场景集合,重新求解主问题,更新 x。

这样做的好处非常明显:每次迭代都在“定向挖掘”最坏场景,而不是靠随机抽样去碰运气。它本质上是把连续优化问题转化成一个不断扩大的有限场景集合上的离散优化问题。场景集合不断扩大,鲁棒解就不断逼近真正的两阶段鲁棒最优解。这个思想在学术界叫列与约束生成算法,而“关键场景辨别”这个词把它的核心动作描述得更直白——在每一轮里,通过子问题的求解,从不确定集里辨别出当前威胁最大的那个点。

2.2 关键场景辨别与传统场景削减法,不是一回事

很多做微网的人第一反应是:不确定性处理不是可以先蒙特卡洛抽样,再用K-means或者快速前向选择做场景削减吗?这两条技术路线在目标上是有本质区别的。

场景削减法(backward reduction、K-means聚类等)处理的是“概率意义上有代表性的场景”。它立足于历史数据,把大量随机场景聚成几个典型簇,然后用簇心代表整类场景。这套方法用于随机优化(stochastic programming)非常合适,因为随机优化算的是期望成本,需要用概率分布把场景加权起来;它关注的是“平均怎么样”。

但鲁棒优化关心的是“最坏怎么样”。假如你按概率聚类挑了5个高频出现的温和场景,忽略掉了2个虽然概率低但一旦出现就代价巨大的极端场景,那么你的调度方案在实际遇到极端天气时可能直接就不可行了。鲁棒优化的核心诉求恰恰是要覆盖这种极端情况,所以场景选取标准必须围绕“对目标函数的破坏力”来设计,而不是“出现的频次”。

关键场景辨别算法在每一轮迭代时,都是严格对着当前调度决策去找最恶劣场景的,这个场景可能概率不高,但一定是对解构成真实威胁的那个点。用行话讲,它有“最优性导向”:场景的产生不是独立于优化模型,而是由优化模型的子问题内生出来的。

2.3 子问题是怎么一步步把最恶劣场景“逼”出来的

子问题的形式是一个 max-min 优化问题。内层 min y 是在给定 u 之后,找到一个最经济的补救方案;外层 max u 是在 U 里折腾所有可能的 u,看哪个让我内层补救也救不回来、成本崩到最高。

直接求双层优化很难,但可以利用线性规划的对偶理论把内层 min 问题转成对偶的 max 问题,这样结构就变成“max 对偶问题”和“max 不确定变量”挤在了一层,不再有嵌套。内层原始问题是:

min c2^T y s.t. D y ≤ e − E x_current − F u

它的对偶问题写出来是:

max λ^T (e − E x_current − F u) s.t. D^T λ = c2 λ ≥ 0

把对偶问题代回原结构后,子问题变成:

max λ^T (e − E x_current − F u) s.t. λ ≥ 0, D^T λ = c2, u ∈ U

这里注意看,目标函数里 λ 和 u 相乘,出现了一个双线性项。处理这个双线性项通常有两种路径:一是由于 u 的可行域是多面体,在强对偶条件下最优解会出现在 U 的顶点,可以枚举顶点;二是把 u 离散化,引入0-1辅助变量和大M法做线性化。第二种方法在 Matlab + YALMIP 框架下非常常见,虽然增加了整数变量,但现代求解器处理几百个二进制变量轻轻松松,而且换来的是模型结构清晰、不依赖人工推导顶点表。

一旦这一层算完,子问题给出的最优 u* 就是当前决策下的关键场景。把它加入主问题,下一轮主问题里就会新增一组约束,强迫第一阶段决策 x 在这个场景下也必须有可行的第二阶段响应。如此往复,场景集合越来越大,解也一步步“服软”,直到上下界差满足收敛条件。

3. Matlab代码实现:从数学模型到可运行的迭代程序

3.1 先明确微网系统长什么样

动手写代码前,先得把系统结构定下来。我用的标准配置是一个低压交流微网,包含三类分布式电源:

  • 光伏阵列,装机容量 400 kW,预测误差取 15%
  • 风电机组,装机容量 300 kW,预测误差取 20%
  • 燃气轮机,最大出力 250 kW,带斜坡约束和最小启停时间约束
  • 储能系统,容量 600 kWh,最大充放电功率 150 kW,充放电效率 0.95,SOC范围 20%~90%
  • 与上级电网的联络线,最大交互功率 300 kW

调度周期取24小时,时间分辨率1小时。第一阶段变量包括燃气轮机启停状态、储能各时段充电/放电状态;第二阶段变量包括各机组实际出力调整量、储能调整量、联络线功率调整量等。

这个配置不算大,但已经能完整体现两阶段鲁棒调度的所有核心细节,求解规模也适合在普通PC上跑。如果读者手头有真实的微网数据,替换成本非常低,只需要把预测曲线和装置参数换成自己的即可。

3.2 主问题建模:场景集合不断扩大的日调度计划

主问题负责在给定关键场景集合的情况下,求解第一阶段决策和近似成本。它长这样:

min c1^T x + η s.t. A x ≤ b 对于每一个已经辨别的场景 u_k,都有: η ≥ c2^T y_k D y_k ≤ e − E x − F u_k

这里的 η 是一个辅助变量,用来替代原问题中 max u∈U min y 那部分成本。每加入一个关键场景 u_k,就新增一组第二阶段变量 y_k 和一约束 η ≥ c2^T y_k。这样做的好处是,所有场景对应的调度问题在同一个优化模型里求解,主问题的规模会随迭代次数线性增长,但每一轮也更接近真实解。

在 Matlab 里我用 YALMIP 建模,核心代码结构如下:

% 定义第一阶段变量 x_bin = binvar(gen_n, 24); % 燃气轮机启停 x_pg = sdpvar(gen_n, 24); % 燃气轮机基准出力 x_pc = sdpvar(1, 24); % 储能充电基准 x_pd = sdpvar(1, 24); % 储能放电基准 eta = sdpvar(1); % 定义每个阶段的场景变量和约束(伪代码示意) Constraints = []; Objective = c1' * [x_pg(:); x_pc(:); x_pd(:)] + eta; for k = 1:length(scenario_set) % 取出场景k下的不确定变量取值 u_k = scenario_set{k}; % 定义第二阶段的调整变量 y_delta_pg = sdpvar(gen_n, 24); y_delta_soc = sdpvar(1, 24); y_pgrid = sdpvar(1, 24); % 添加第二阶段约束和η约束 Constraints = [Constraints, eta >= c2' * y_delta_pg(:)]; Constraints = [Constraints, ...]; % 场景k下的运行约束 end ops = sdpsettings('solver', 'gurobi', 'verbose', 0); optimize(Constraints, Objective, ops);

这里有一个工程细节值得强调:每一轮主问题里,场景对应的第二阶段变量是独立的,但它们共享同一个第一阶段决策 x。所以主问题的结构天然是“一个大矩阵带块对角结构”,这也是为什么现代求解器处理起来并不会太吃力。随着迭代轮次增加,主问题规模变大,单轮求解时间也会拉长,这是正常的,不用慌。

3.3 子问题求解:怎么找到当前决策下的最恶劣场景

子问题的输入是第一阶段的决策值 x_current,输出是 u_worst 和最坏情况下的第二阶段成本 Q(x_current)。子问题是把内层 min 对偶之后形成的单层 max 问题。

我的实际做法是用 YALMIP 直接建模,对偶变量不需要手写推导,让求解器的预处理器自己去处理。如果你用的是纯 Matlab + 单独求解器(不经过 YALMIP),那就需要自己对每个约束手工对偶,那个工作量相当痛苦,而且容易出错。

子问题里参与优化的变量包括不确定变量 u 和对偶变量 λ。需要注意,原问题如果原本是 LP,取对偶后的变量类型是连续非负的,可以直接用 sdpvar 声明;如果内层问题里出现绝对值项,比如成本项里包含功率偏差绝对值,那就需要先引入辅助变量 t 并添加两条不等式约束 t ≥ Δ 和 t ≥ −Δ,把它改写为线性形式,再交给 YALMIP 处理。

我在代码里是这样组织的:

function [u_worst, Q_val] = solve_subproblem(x_current, system_data) % 不确定变量 u_wind = sdpvar(length(wind_bus), 24); u_pv = sdpvar(length(pv_bus), 24); % 对偶变量(按原问题约束逐条引入) lambda_balance = sdpvar(1, 24); lambda_gen = sdpvar(length(generators), 24); lambda_soc = sdpvar(1, 24); % 不确定集约束:盒式 + 预算 Constraints = [ub_lower <= u_wind <= ub_upper]; Constraints = [Constraints, sum(abs(u_wind - forecast_wind) ./ max_deviation) <= Gamma]; % 目标函数是已经处理好双线性项的对偶目标 Objective = -lambda_balance' * (x_current - system_data.forecast_load) + ... lambda_gen' * (system_data.capacity - x_current) + ...; % 因为有双线性项,需要 u 离散化或引入辅助变量处理 % 实际中我倾向于把 u 的连续盒式区间离散成有限顶点集 % 再用 big-M 线性化λ与u的乘积 ops = sdpsettings('solver', 'gurobi', 'verbose', 0); optimize(Constraints, -Objective, ops); u_worst = value(u_wind); Q_val = value(Objective); end

这段代码里最值得留意的就是双线性项的线性化处理。我这边的经验是:与其硬着头皮保留 λ·u 的产品项然后让非线性求解器一头撞死,不如直接把 u 的可行域顶点枚举出来、配合 big-M 辅助变量做线性化。虽然问题会多出一些整数变量,但换来的是求解过程的绝对稳定性和可复现性,这在工程项目里比节省几十毫秒求解时间重要得多。

3.4 主循环迭代:直到上下界收敛

整个算法的主循环不复杂,按“主问题求解 → 子问题求解 → 加入场景 → 再回主问题”的顺序转圈。收敛判据是上界 UB 和下界 LB 的间隙不超过设定阈值。

LB 来自主问题的最优目标值 c1^T x + η,UB 来自在当前 x 下子问题求出的真实第二阶段成本 c1^T x + Q(x)。两者的差反映了“场景集合不完整”造成的近似误差。当初始终止条件:

% 两阶段鲁棒主循环 UB = inf; LB = -inf; tol = 1e-3; max_iter = 15; scenario_set = {}; for iter = 1:max_iter % 1. 解主问题,得到第一阶段决策和LB [x_opt, LB_new] = solve_master_problem(scenario_set); LB = max(LB, LB_new); % 2. 固定x,解子问题,得到关键场景和UB [u_worst, Q_val] = solve_subproblem(x_opt, system_data); UB_new = ...; % 第一阶段成本 + Q_val UB = min(UB, UB_new); % 3. 判断收敛 if (UB - LB) / abs(UB) < tol break; end % 4. 把新的最关键场景加入集合 scenario_set = [scenario_set, u_worst]; end

我实测下来,常规的24节点微网案例,这个循环通常跑5到10轮就能收敛到1e-3的相对间隙。如果10轮还不收敛,优先检查是不是子问题目标函数里少了某块成本,或者不确定集模型写错导致子问题每次都在造出不可能的极端场景。

4. 算例调试与参数敏感性分析

4.1 基础数据配置:预测曲线和不确定集参数

为了做对照实验,我用了一组典型日的预测数据:光伏出力曲线是典型的“倒钟形”,中午12点到14点达到峰值 350 kW;风速预测全天在 6~9 m/s 之间波动,折算成风机出力是一条不太规则的曲线;负荷曲线呈现早、晚双峰特征,早峰在9点前后约 500 kW,晚峰在19点前后约 550 kW。

电价取分时电价,峰段(8:00-11:00、18:00-21:00)1.2 元/kWh,平段 0.8 元/kWh,谷段 0.4 元/kWh。燃气轮机的发电成本系数取 0.45 元/kWh 加上启停成本。储能不计算折旧成本,但把充放电次数作为软约束处理。

不确定集参数我设成:光伏出力最大偏差取预测值的 15%,风电取 20%,负荷偏差取 5%。预算 Γ 先取 6(相当于不确定变量总数的一半左右),后面再专门做敏感性分析。

4.2 预算 Γ 对调度经济性和保守性的影响

把 Γ 从 0 逐步调到 12,其他参数保持不变,跑完整个两阶段鲁棒流程后,调度总成本的变化趋势非常能说明问题:

Γ 值总成本(元)与确定性模型成本差最坏场景辨识迭代次数
0(确定性)432601
24489+3.8%4
44612+6.6%6
64720+9.1%8
84811+11.2%10
12(全盒式)4956+14.6%12

从这个表里能明显看到一个规律:成本随着 Γ 增加而上升,但边际增幅在递减。Γ 从0到4,成本增加6.6%;从4到8,增加4.6%;从8到12,只增加3.4%。这说明预算约束的设计价值就在这里——用比较小的经济性牺牲,覆盖掉绝大部分极端场景风险。在真实项目中,我不会盲目把 Γ 拉到满格,而是结合历史极端天气出现的概率来选,一般取不确定变量总数的1/3到1/2是性价比最高的区间。

另外我注意到 Γ 越大,迭代收敛需要的场景数量越多。原因也很好理解:Γ 增大了 U 的“搜索空间”,子问题需要更多轮次才能逼近全局最坏场景。这一点在算法性能调试时值得提前有心理预期。

4.3 不确定区间宽度的影响:警惕过度保守

除了 Γ,不确定区间的宽度同样关键。我把光伏、风电的预测误差上限从5%逐步加到30%,固定 Γ=6,结果如下:

误差上限总成本(元)最恶劣场景下弃风弃光率
5%45121.2%
10%45983.8%
15%47206.7%
20%48909.5%
30%524615.2%

可以看到,当误差上限超过某个阈值之后,成本会快速抬升。误差上限取30%的时候,成本比5%时高了整整16%,而且最恶劣场景下的弃风弃光率也到了15%左右——说明这个方案已经开始靠牺牲清洁能源消纳来保安全。

实际操作中,我建议不确定区间的宽度要从历史预测误差数据里统计出来,而不是拍脑袋定一个数。对每一预测时段,统计过去N天同一时刻的预测偏差,取95%分位数作为该时段的上下界,比统一按百分比更贴近真实波动规律,而且调度结果往往更经济。这些统计代码在 Matlab 里无非就是 quantile 函数几行的事,性价比极高。

5. 常见问题与排查实录

5.1 子问题是双层结构,直接求解失败怎么办

我早期写过一版代码,子问题直接用 YALMIP 的 max-min 嵌套去建模,结果求解器一直报“无法处理嵌套优化问题”。这是因为标准求解器不支持一个优化问题里再套一个优化问题。

解决办法就是把内层 min 对偶成 max。如果你不想手推对偶,YALMIP 里可以用 duality 相关的工具函数辅助,但更稳妥的是针对自己的原问题手推一遍。我把这个步骤遇到过的问题汇总了一下:

  • 对偶变量非负性漏写:很多约束对应到对偶变量时是“无符号限制”的,比如等式约束对应的对偶变量是自由变量(unrestricted in sign),写错会让对偶问题的最优值跟原问题完全对不上。
  • 原问题如果是极小化,对偶问题的目标函数符号容易搞反。检查对偶间隙是否为零的有效方法:把原问题的原始最优解代入对偶目标函数验证是否同值。
  • 绝对值项忘记线性化:成本中有 |y| 形式的项时,直接进对偶会得到非线性项,必须先用辅助变量改写。

建议做法是先拿一个小规模的算例,手算一遍,确认对偶变换完全正确之后,再上整套砂型代码。正式解决前先用尽量能手工验算的算例调试对偶逻辑,这步省下来的排查时间通常非常多。

5.2 优化结果不可行:多半是第二阶段子问题缺了松弛变量

两阶段鲁棒的子问题有时候会遇到不可行的情况。原因可能不只是模型写错——搜索到的那个极端场景,在当前 x 下确实让系统没有任何可行的 y 可以满足所有约束。这在实际物理系统里意味着:你第一阶段定的基准计划太“僵”,遇到了远超预期的不确定性实现。

标准做法是在第二阶段约束里引入松弛变量,给关键约束(比如功率平衡、SOC上下限)配上很大的惩罚系数。相当于给紧急情况留了一扇“紧急通风口”,并且用高惩罚成本逼着优化器尽量避免去开它。这样即使极端场景下系统不能严格满足所有硬约束,模型也不会直接崩掉,而是反应成一个可负担的代价。

惩罚系数的选择也要注意。我习惯先把惩罚成本系数设成正常成本的 100 倍以上,然后跑一遍看它是否被激活。如果发现最恶劣场景下松弛变量显著非零,那说明方案本身已经有了物理上的不可行风险。这时候一定要回头检查第一阶段决策,而不是一味调大惩罚系数把问题掩盖掉。

5.3 Matlab 环境配置与求解器选型:Gurobi 还是 Cplex

这套代码依赖 YALMIP 工具箱和一个高性能混合整数线性规划求解器,不是 Matlab 自带的 solve 函数能搞定的。我的实际配置是:

  • Matlab R2023b
  • YALMIP 最新版(GitHub 上直接 clone 下来,路径加进 Matlab 搜索目录)
  • Gurobi 求解器(学生或学术版免费,性能和数值稳定性都比 Cplex 好,尤其是处理 big-M 线性化带来的整数变量时)
  • 用 sdpsettings('solver', 'gurobi') 指定求解器

如果你用的是 Cplex,代码结构基本不用改,YALMIP 后端自动识别。不过我个人在双线性项线性化后的问题里,Gurobi 的整数预处理能力明显更强,求解速度能差到一倍以上。如果没有 Gurobi 授权,Cplex 的学术版也完全可以跑通整个流程,只是大算例下会慢一些。

另外,如果系统里没装求解器,YALMIP 内置的免费求解器(比如 sedumi、linprog)也能跑,但只能处理连续线性规划,遇到 big-M 引入的二进制变量会直接报错。所以一开始装环境的时候,把 Gurobi 或 Cplex 装好是最重要的一步。

5.4 我踩过的一些坑:一些程序运行层面的实际问题

第一个坑是 SOC 约束跨时段耦合。储能SOC 是逐时累积的变量,它在第二阶段场景里也要跟着调整。如果每个场景都独立定义一套 SOC 变量,而场景之间的 SOC 没有共享或衔接关系,那么即便每个场景单独都可行,整体能量约束也可能是虚的。我的处理方式是把 SOC 初值设为固定参数,并且保证主问题里的储能基准计划在所有场景下都不越界。这样做虽然会损失一点储能灵活性,但换来的是模型在数学上更干净、收敛更稳。

第二个坑是 big-M 值的选择。线性化 λ·u 时,M 值取得太大会导致数值条件数恶化,求解器给出的解精度漂移;取得太小又可能把可行域错误地切掉。通用的做法是先用一个不包含二元变量的 LP 求解一次,估算 λ 的合理范围,然后让 M 比这个范围上限大一个数量级左右,不要无脑写 1e6。

第三个坑是收敛判据里的绝对值。当两阶段成本里有大量固定成本时,UB 和 LB 的差距本身就有一个“平台期”,直接看绝对间隙可能永远达不到 1e-3,需要改成相对间隙(除以UB)作为收敛条件。我在前面主循环代码里写的就是相对间隙判据,这个细节对实际跑通流程很关键。

最后再做一次实验记录:上面的例子在 i5 + 16G 内存的普通办公笔记本上,单轮主问题 + 子问题的求解时间在 7 到 15 秒之间,完整迭代 8 轮收敛到 0.1% 相对间隙,总耗时约 2 分钟。这个计算成本在微网日前调度场景下完全能接受,甚至可以在滚动更新的场景下多跑几轮。如果你在自己机器上跑出了远超这个数量级的耗时,优先检查一下是不是场景数增长导致主问题规模爆炸,或者是否不小心把第二阶段变量不该有的耦合写进去了。

做两阶段鲁棒优化,我最想跟大家强调的一点是:不要一上来就求“完美复现”。先把确定性模型跑通,再把不确定性集合加进去,最后才上主问题-子问题的迭代框架。每一步拆开都能验证,最后合起来就不会出现“跑不通但不知道哪里错了”的焦灼。这个从简到繁的过程,也是我重写了三版代码后最深刻的体会。

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询