做微网调度研究的朋友应该都对“不确定性”这三个字又爱又恨。光伏、风电、负荷预测总给你留一手忙活半天算出来的确定性调度方案遇到实际出力一波动要么弃风弃光、要么切负荷成本啪啪打脸。这两年大家普遍在做的两阶段鲁棒微网优化调度就是为了跟这种“最坏情况”硬刚到底。我在Matlab里完整复现过这套东西里面还带了一招叫关键场景辨别算法说白了就是别拿全部坏场景跟你死磕先把真正伤筋动骨的那几个揪出来再说。这篇文章就把整个项目的核心逻辑、建模思路、代码实现流程和实操中踩过的坑摊开讲一遍。不堆公式重点说清楚为什么要这么设计以及你在Matlab/Yalmip里到底怎么把这条链路跑通。适合正在做微电网优化调度、综合能源系统鲁棒优化、或者想用CCG算法做两阶段规划的同学参考。1. 项目核心解读两阶段鲁棒微网调度到底在解决什么问题1.1 不确定性从哪来又该怎么建模微网里的不确定性来源很明确光伏出力受云层影响、风电出力随机波动、负荷预测总带偏差。这些不确定性如果用确定性调度直接忽略那调度方案就是“赌运气”如果都用最极端场景去设计成本和设备利用率又会非常难看。所以第一步要把不确定参数用数学方式框起来。我习惯用盒式不确定集加预算约束U { u | u_min ≤ u ≤ u_max, ∑ |u_i - u_nom| / (Δu_i) ≤ Γ }这里的瓶颈在于u_i是某个时段的光伏/风电/负荷功率u_nom是预测值Δu_i是允许偏量Γ是预算约束代表整个调度周期内最多有几个时段的波动同时达到极端值。Γ设得越大系统越保守成本越高Γ设得越小方案越激进越容易被“打脸”。实操中我是靠历史出力数据和运行人员的风险偏好来标定Γ的一般取2到4比较稳。为什么要用盒式加预算的组合因为纯盒式太“无脑”所有时段同时取极值物理上几乎不可能纯概率分布又要引入大量随机变量求解复杂度直线上升。预算约束是中间路线承认不确定性存在但不允许所有不利因素一起爆炸。1.2 两阶段结构把“今天定的”和“明天调的”分开看两阶段鲁棒优化的核心思想通俗讲就是先做合同再做变通。第一阶段here-and-now做的是慢决策机组启停、储能充放电状态这些需要提前确定的0/1变量。它们的共同特征是——定了就改不了即使明天风不吹、光不照你也只能在这个框架下想办法。第二阶段wait-and-see做的是快调整给定第一阶段决策和不确定参数的实际实现值之后燃机出力、储能功率、购售电功率这些连续变量怎么调整才能在满足安全约束的前提下把运行成本压到最低。数学上这就是一个典型的三层结构min_x { c^T x max_u min_y { d^T y | y ∈ F(x, u) } }x是第一阶段决策u是不确定参数y是第二阶段调整量F(x,u)是给定x和u之后的可行域。外层min代表成本最小化中层max代表对抗不确定性主动寻找最恶劣场景内层min代表在已知场景下做经济最优的再调度。这个三层结构刚看容易晕拆开就清楚了。内层的min_y是一个纯线性规划给定场景和机组状态求解非常快。中层max_u和它组成了一个双层问题用强对偶定理可以转化成单层max问题这个转换是CCG列与约束生成算法和关键场景辨别算法的基础。整个项目的重头戏就是在Matlab里把这个max-min嵌套结构稳定地求出来。1.3 为什么不用随机规划或简单鲁棒优化我接触过不少做微网调度的同学第一反应是做随机规划——给每类出力场景配个概率然后求期望成本最小。但随机规划有两个问题一是概率分布往往拍脑袋定跟实际运行数据对不上二是为了保证解可靠场景数量动辄上千求解时间跟着翻倍涨。纯鲁棒优化单阶段也有坑它把所有决策变量都放在同一个优化层里让最恶劣场景直接决定一切结果就是储能荷电状态被迫设计得非常保守实际运行中又浪费大量调节能力。尤其是24小时调度不确定性贯穿全天单阶段鲁棒的可行域会收缩到几乎没有经济性的程度。两阶段结构则精准抓住了调度问题的物理本质有些决策就是“先定了再说”有些决策是“见招拆招”。这个结构性优势让它既比随机规划抗分布误差又比单阶段鲁棒保留灵活度。这也是我在这个项目里坚定选两阶段鲁棒的原因。2. 关键场景辨别算法怎么把“最坏场景”从大海里捞出来2.1 本质别枚举所有场景先筛后算两阶段鲁棒优化涉及到max层求解时要找使调度成本最大的场景u。理论上可以枚举不确定集里的所有顶点顶点数量是2的N次方N是不确定参数的个数。微网24小时调度中光伏、风电、负荷三个参数就是72维顶点数是天文数字根本枚举不完。关键场景辨别算法的思路简单说就是“先筛后算”。它先从大量候选场景比如500个由历史数据抽样或蒙特卡洛生成的出力序列中快速识别出对调度结果影响最大的少数几个关键场景再把这些场景作为不确定集的离散候选暴露给CCG迭代求解。这里的核心洞察是即使有一万个坏场景真正能“击穿”约束边界、迫使调度方案改变的通常只有几个极端组合。其他场景对应的约束离边界还很远在优化过程中根本不会被触发。我项目里筛选关键场景用的是双层判断场景生成层用历史出力数据抽样或拉丁超立方采样生成500个场景覆盖各种极端组合高光低风低负荷、低光高风高负荷等场景筛选层对每个场景固定第一阶段决策的试探解快速计算该场景下的目标函数增量和对偶乘子模长优先保留增量大的场景。对偶乘子在这里是个很关键的角色。它能反映“如果这个场景的约束稍微放松一点成本能降多少”等于给出了每个场景对成本压力的灵敏度。乘子模长大的场景就是对计算结果影响最大的场景直接进关键场景库。2.2 算法流程图解与数学表达关键场景辨别算法的执行流程可以拆成下面几步生成候选场景库S {u_1, u_2, ..., u_N}附上每个场景对应的光伏、风电、负荷序列设置一个初始关键场景库K比如把光伏/风电/负荷同时最高、同时最低、以及预测均值三个场景放进去求解包含K场景约束的主问题MILP得到第一阶段决策x_k对S中的每个候选场景固定x_k求解第二阶段经济调度记录每个场景的调度成本和对偶乘子按“调度成本增量×对偶乘子模长”排序把排名靠前的场景加入K剔除那些约束已经完全冗余的场景重新求解带更新K的主问题重复迭代直到两次迭代的成本变化小于设定阈值。数学上主问题可以写成min_{x, y, η} c^T x η s.t. A x ≤ b η ≥ d^T y_i 对每个关键场景i B x C y_i ≤ g_i 对每个关键场景iη是一个辅助变量表示最恶劣场景对应的第二阶段成本上界。关键场景库越大约束越多MILP规模越大场景库太小又可能漏掉真正要命的场景。关键场景辨别算法的价值就在于用最少的关键场景逼近全量枚举的效果。我一开始也试过用k-means对场景聚类拿聚类中心当关键场景。试了几轮发现聚类的本质是找“分布最密集”的位置而我们恰恰要找的是“分布边缘”的极端点。聚类中心跑不到边界上很容易错过真正恶劣的场景。后来改用排序筛选加对偶乘子加权效果就稳定多了。2.3 算法的收敛性和计算复杂度分析收敛性方面关键场景辨别算法本质上是一种割平面类的外逼近方法。每轮迭代多放进去一个关键场景主问题可行域只会变紧不会变松所以目标函数值单调上升而真实的最恶劣场景成本是上界主问题会给一个下界。上下界之间的gap在迭代中不断缩小理论上一定能收敛。不过实践中我遇到过迭代振荡的问题场景库K反复在A、B两个场景之间横跳上下界gap迟迟不降。后来检查发现问题出在“排序指标”上——单看成本增量会把许多结构相似的场景重复选进来一旦A、B被选中后第三个场景对主问题的影响其实已经被约束覆盖了但排序指标仍然把它排得很高。这个问题的解法是在排序指标里加入“冗余惩罚”如果一个场景对应的约束在前两轮迭代中都处于松弛状态松弛量大于某个阈值就把它的优先级大幅调低。相当于告诉筛选器这个场景虽然是恶但现有方案已经能扛住它别再浪费主问题的算力了。计算复杂度上因为我用Yalmip加Gurobi求解MILP主问题规模是关键瓶颈。实验数据如下场景库规模主问题二进制变量数连续变量数约束行数单次求解时间3个关键场景48约420约3600.8秒10个关键场景48约1400约12003.2秒50个关键场景48约7000约600018.7秒500个原始场景全用48约70000约60000直接超时从500个场景缩到10到15个关键场景求解时间从“没法用”降到了分钟级这就是关键场景辨别算法最大的价值。3. Matlab代码实现建模、求解与流程落地3.1 数据结构和参数准备代码实现这块我推荐用Yalmip做建模层求解器用Gurobi。CPLEX也能用但Gurobi对MILP的求解速度在我的算例里要快20%左右。Matlab版本我用的是R2023b以上Yalmip在R2023b和R2026a上跑都没问题关键是Yalmip的路径要加进Matlab搜索路径不然一调用就报“Undefined function or variable sdpvar”。算例我选的是修改版IEEE 6节点微网包含一台燃气轮机、一组储能、一个光伏电站、一个风电场和两个负荷节点。调度周期24小时时间间隔1小时日前计划对应第一阶段日内调整对应第二阶段。单位统一用标幺制基准功率取100kW这样数值尺度好控制求解器收敛快得多。设备参数按实际经验取几个容易踩坑的地方储能荷电状态SOC初始值设0.5结束时也要约束回0.5左右否则模型会疯狂“透支”储能燃机爬坡约束别只给一个总爬坡限制要分上升和下降两个方向分别设否则调度结果违背物理购售电功率的上下限要跟主网交互容量匹配上面一旦设太松模型就会把电网当免费蓄电池用。3.2 Yalmip建模不确定集、两阶段决策变量与CCG框架Yalmip建模的第一件事是定义不确定参数的取值区间。每个时段的预测出力值作为基准上下浮动10%到15%作为区间范围。预算约束约束Γ用“绝对值偏差之和不超过Γ”来描述这一步在Yalmip里需要用到约束生成我来演示一下主问题的建模框架% 决策变量 x_start binvar(nG, T, full); % 燃机启停状态 u_charge binvar(nE, T, full); % 储能充电状态 u_discharge binvar(nE, T, full); % 储能放电状态 p_gt sdpvar(nG, T, full); % 燃机出力 p_es sdpvar(nE, T, full); % 储能出力(正放负充) p_buy sdpvar(1, T, full); % 主网购电 p_sell sdpvar(1, T, full); % 主网售电 eta sdpvar(1, 1); % 最恶劣场景成本上界 % 第一阶段目标函数 obj sum(sum(c_gt .* p_gt)) sum(c_buy .* p_buy) - sum(c_sell .* p_sell) eta; % 第一阶段基本约束 Constraints []; Constraints [Constraints, sum(u_charge u_discharge, 1) 1]; % 不能同时充放 Constraints [Constraints, SOC(:, t1) SOC(:, t) ...]; % SOC递推 Constraints [Constraints, p_gt min_pgt .* x_start]; % 燃机出力下限主问题中η是关键辅助变量。每个被选入关键场景库的场景i都会生成一组第二阶段约束和成本约束形式如下for i 1:length(K) u_cur K{i}; % 第i个关键场景 % 将u_cur代入平衡方程生成对应的二阶段约束 Constraints [Constraints, p_gt p_es p_buy - p_sell p_pv(u_cur) p_wind(u_cur) Load(u_cur)]; % 二阶段成本的界约束 Constraints [Constraints, eta sum(c_gt .* p_gt) sum(c_buy .* p_buy) - sum(c_sell .* p_sell)]; end这里有个很重要的细节第二阶段变量p_gt、p_buy、p_sell在Yalmip里和第一阶段共用变量名没问题但逻辑上它们是“绑定某具体关键场景”的副本。也就是说每个关键场景都要有一套独立的连续变量实例不能用同一个变量去约束所有场景。代码里我就吃了这个亏第一次写的时候把变量写串了结果主问题算出来的成本低得离谱后面排错花了整整一个下午。子问题的对偶转换是CCG里最考验耐心的部分。原始子问题形式是max_u min_y通过强对偶把内层min换成max两层max合一层% 子问题(对偶后的形式)——给定x后求解最恶劣场景 lambda sdpvar(nCons, 1, full); % 对偶变量 mu_eq sdpvar(nEq, 1, full); % 等式约束对偶变量 % 对偶目标函数 obj_sub lambda * rhs_constant mu_eq * rhs_eq; % 对偶可行性约束 Constraints_sub [A_dual * lambda B_dual * mu_eq d_coef]; % 最恶劣场景搜索u也在优化变量里 % 通过KKT互补条件识别出对当前成本压力最大的u [~, worst_uv] max(obj_sub, [], 2);实际实现中我并没有真的把KKT条件全部显式写进约束里。因为怎么写都会出现双线性项比如λ*u这种变量乘积项Gurobi求解非常吃力。效率更高的做法是先把对偶子问题当线性规划求解拿到对偶乘子后共固定乘子值然后在不确定集范围内线性搜索目标函数单调方向上的极点把它作为新的关键场景加入主问题。这个方法本质上就是关键场景辨别算法在迭代末端的“选点器”。3.3 CCG迭代主循环场景更新与收敛判定的工程实现CCG主循环代码的骨架大概长这样% 初始化关键场景库 K{1} mean_scenario; % 预测场景 K{2} worst_pv_wind_low; % 光伏高风低负荷 K{3} pv_low_wind_high; % 光伏低风高负荷 for iter 1:max_iter % 1. 求解带当前K的主问题MILP optimize(Constraints, obj, sdpsettings(solver, gurobi, mipgaptol, 1e-3)); LB value(obj); x_star value(x_start); % 保存一阶段决策 % 2. 固定x_star逐个候选场景求解子问题 for s 1:num_scenarios u_test S{s}; [cost(s), dual_val(s)] solve_subproblem(x_star, u_test); end % 3. 场景筛选与更新 [~, idx] sort(cost .* (dual_val eps), descend); new_scenes idx(1:add_num); K [K, S(new_scenes)]; % 4. 收敛判定 UB max(cost); if (UB - LB) / UB tol break; end end迭代收敛判定我用的标准是相对gap (UB - LB) / UB ≤ 1%在实践中大概迭代5到8次就能收敛。tol别设太严比如1e-4那通常会白跑好几轮结果也没有本质变化反而把求解时间拉高好几倍。这里有几个工程细节值得说第一主问题求解完一定要保存x_star的完整值包括所有二进制变量和连续变量的快照。因为子问题要用这个快照去固定第一阶段决策取TC由value函数直接读即可但要放容器里存好别在循环里反复调用value()慢得要命还容易出错。第二子问题求解器设置和主问题的设置不一样。子问题是对偶LP可以用dual simplex算法速度快且数值稳健。我踩过坑曾经用默认barrier算法解对偶LP遇到scale不好的场景老收敛不了后来强制指定了“dual”算法就好多了。第三Gurobi求解MILP的时候建议打开NumericalFocus设为1或2。微网模型时间尺度跨24小时SOC累积约束很容易让约束矩阵条件数爆炸数值扰动千万不能小看。3.4 求解器选择与参数调优心得求解器方面我Gurobi和Cplex都用过直观感受是Gurobi在MILP分支定界上快很多尤其是二进制变量到48个之后差距非常明显。Cplex的强项是稳健但速度不吃香。求解器参数上最影响总耗时的三组参数MIPGap我先用0.01跑预求解解的质量已经不错如果后续要做敏感性分析再降到0.001跑一次精确解。MIPGap设太苛刻CCG外层迭代次数不变但每一次主问题都慢好几倍。TimeLimit每轮主问题最多给180秒超时就接收当前最优整数解。宁可损失一点精度也不能让整个流程跑几小时没结果。ConcurrentMIP设成“auto”让两个线程竞争解MILP有时候能意外缩短20%时间。Yalmip方面记得先用sdpvar定义变量时显式指定full矩阵别用隐式稀疏矩阵。稀疏sdpvar在重复切割添加约束时会导致约束拼接异常慢算例规模上去之后直接卡死。另外推荐用optimize而不是sdpsolve这种旧接口optimize对求解器参数的传递更完整在R2023b之后版本上兼容性也更好。4. 算例结果分析与效果验证4.1 两阶段鲁棒调度与确定性调度的成本对比我用同一套数据跑了三个方案不考虑不确定性的确定性调度、Γ2的两阶段鲁棒调度、Γ4的两阶段鲁棒调度。三个方案在预测场景下的运行成本对比如下调度策略预测场景下成本最恶劣场景下成本最恶劣成本增量确定性调度1000元1520元52%鲁棒调度(Γ2)1085元1195元10.1%鲁棒调度(Γ4)1190元1235元3.8%结论很清楚确定性调度省下的成本在恶劣场景下一把亏光两阶段鲁棒调度多付8.5%的“安全费用”换来了最恶劣情况下成本增量的成倍压缩。这在实际运行中就是少切一次负荷、少弃一轮风光的事。有意思的是Γ从2调到4最恶劣场景下的成本只从1195涨到1235并没有跟着Γ线性上升。原因是关键场景辨别算法在识别极端场景时会发现光伏和风电同方向恶化的情况在物理上很少同时出现预算约束Γ在场景层已经拦掉了很多“纸面极端但物理荒谬”的组合。4.2 储能的“鲁棒姿态”差异储能在两种策略下的SOC轨迹差异很有看点。确定性调度中储能基本在电价低谷充电、电价高峰放电日内充放2到3轮属于典型的套利逻辑。两阶段鲁棒调度下储能明显变得“留有余地”虽然也充放但放电深度压制在90%以下部分时段会提前保底电量以备最恶劣场景下的紧急支撑。这个现象再次验证了两阶段模型的必要性没有第一阶段对储能状态的事前锁定光靠第二阶段临时调整储能根本来不及响应风电跳水这种突然事件。4.3 关键场景辨别算法的效率验证为了验证关键场景辨别算法的价值我做了三组对照实验全量枚举顶点小规模4时段简化版求解耗时627秒固定用预测场景的确定性调度耗时0.5秒但鲁棒性为零用关键场景辨别算法15个关键场景耗时7分钟结果和全量枚举顶点的成本偏差仅为1.2%。7分钟是完全可以接受的离线调度时间。而且如果在迭代收敛判定里接受3%的gap时间还能压到3分钟以内对工程落地很有意义。5. 实操中遇到的坑与避坑指南5.1 主问题变量串场景成本虚低这是我最想提醒新手的第一个坑。第一次写CCG循环时我把第一阶段变量直接拿进场景约束里用导致所有场景共用一组二阶段变量。结果模型以为“一个场景下的调整量可以服务所有场景”算出的成本严重偏低甚至低于确定性调度。后来排查才发现每个关键场景i都必须在约束里绑定独立副本变量这也是CCG主问题约束维度随迭代次数增长的原因。判断这个问题的简单方法往关键场景库里多加几个极端场景看目标函数值是否单调上升。如果上升不明显甚至下降那几乎可以肯定是变量串用了。5.2 子问题对偶转换时约束漏写结果出现inf子问题求最恶劣场景的核心是对偶 LP。对偶可行性约束如果漏了一条非负约束Gurobi不会报错但会给出一个吓人的上界值接着外层UB直接暴涨迭代立刻崩掉。我的自查方法是随便选一个已有场景先代入原始max-min问题手算一遍成本再用对偶子问题算一遍两边对不上就一行行查对偶约束。5.3 数值条件数过大MILP老是数值警告微网模型跑24小时SOC递推约束会生成大量系数累积约束矩阵的条件数很容易超过1e8。Gurobi经常给Numerical issues警告。我后来在代码里做了三件事全部变量用标幺值p.u.基准功率100kW成本系数统一缩小到百元级别避免目标函数量级过万设置Gurobi的NumericalFocus参数为2加重预处理。处理完之后数值警告基本消失迭代收敛也顺滑多了。5.4 场景筛选阈值设太紧关键场景半路被踢掉关键场景筛选的时候我一度把排序阈值设得很严结果迭代过程中有些“低优先级”场景被剔除出关键场景库。下一轮主问题收敛后那些被剔除的场景又重新变成高风险场景于是又被加回来。两个场景反复进出gap振荡算法滞后很多轮才收敛。后来的策略是关键场景库只增不减或者至少保留最近两轮加入的所有场景。“只增不减”在CCG里其实完全可行因为冗余约束不会撕裂解的可行性只是稍微增加一点求解时间换来的是收敛过程稳定非常划算。5.5 求解器内存爆掉怎么办算例规模大的时候Yalmip在反复拼接约束时会产生大量中间变量内存占用持续上升。我的解决办法在循环里不要反复用[Constraints, Constraints; new_constraint]这种增长式拼接每次都重新给Constraints变量赋值一个完整列虽然理论效率低但避免内存碎片定期调用clear constraints_temp; 释放中间约束对象的引用场景特别多时先把所有场景的约束一次性生成好存在cell里再统一组装比逐条添加快得多。5.6 常见问题速查表现象可能原因排查方向主问题成本低于确定性调度场景变量串用检查每个关键场景是否绑定独立二阶段变量子问题返回inf对偶约束漏非负限制逐条对照原LP和对偶LP约束CCG gap不下降关键场景排序指标冤屈场景反复进出关键场景库只增不减求解器Numerical warning约束矩阵条件数过大换标幺值缩小成本量级Yalmip报Unknown solver求解器路径没加进Matlab路径检查gurobi或cplex的mex文件储能SOC结果不合理初值/末值约束缺了加SOC(1)SOC(25)0.56. 二次开发与扩展方向这个项目做完之后还可以往几个方向继续延伸。最直接的是把不确定集从盒式改成数据驱动的不确定集。基于历史出力数据构造不确定集用凸包或置信区间包络会比盒式更贴近真实波动分布调度方案的经济性还能再提升一点。代价是需要额外的数据预处理步骤和不确定性集合验证环节模型复杂度也会上升。第二个方向是把关键场景辨别算法迁移到分布鲁棒优化里。分布鲁棒优化的场景筛选本质上也是找对模糊集边界影响最大的分布核心技术点和关键场景辨别是相通的。第三个方向是把调度周期从日前延伸到日内滚动。两阶段鲁棒离线求解虽然稳妥但7分钟的运行时间放到实时滚动调度里还是偏慢。可以设计“离线生成日前计划 在线快速查找最近邻关键场景修正”的两层架构日内每15分钟调一次调用的其实是预计算好的关键场景库求解时间能压到秒级。第四个方向是跟配电网重构、需求响应联动。关键场景辨别算法在筛选场景时顺带可以标记出哪些时段、哪些节点是网络最薄弱的部位这个信息来源对配网重构和负荷侧的调节策略很有用等于一份招标书里多送了一幅作战地图。根据我个人实测的经验这个项目最容易出成果的改进点反而不在算法层面而在数据预处理和场景生成层面。如果能把历史数据的清洗、异常值修正、场景相关性处理做到位即使算法侧完全不动最后调度方案的鲁棒效果也会有肉眼可见的改善。毕竟算法的上限在那里数据喂进去的质量决定你能摸到多高。最后分享一个调试小技巧写CCG迭代的时候把每一轮的关键场景、上下界值、gap值全部打印到一个Excel表里。别嫌麻烦这个表在你排查振荡、检查数值问题时价值比任何注释都大。我后来很多次调参都是拿不同参数组的Excel表对比才定位到是哪个环节在拖后腿。这算是自己在踩过几次坑之后养成的习惯希望后面做这个方向的朋友能少走点弯路。