这个标题是我在电力系统优化方向绕了很久之后回头再看依然觉得含金量很高的一个题目。很多人第一次看到“基于线性准则的考虑风力发电不确定性的分布鲁棒优化机组组合”会被这一串术语吓住觉得门槛太高但实际上把它拆开看它解决的是一个非常接地气的工程问题风电出力说不准我到底该让火电机组怎么启停、怎么出力才能在成本和安全之间找到那个最优的平衡点。这篇文章我会把这套方法从原理到Matlab实现完整梳理一遍包括为什么用分布鲁棒而不是传统随机规划、线性准则里的线性决策规则到底怎么用、两阶段模型里那些约束是怎么写出来的以及我实际调代码时踩过的坑。内容偏硬核但会尽量讲人话适合正在做机组组合、风电并网或者鲁棒优化的研究生和工程师参考。1. 问题本质机组组合为什么要跟风电不确定性较劲1.1 从确定性模型说起为什么不能直接用点预测传统的机组组合(Unit Commitment, UC)是电力系统调度里最经典的问题。目标函数很直观给定负荷曲线和机组参数决定每台火电机组在哪些时段开机、哪些时段停机、每个时段出多少功率使得总发电成本最低同时满足负荷平衡、备用容量、机组爬坡、最小启停时间等一系列约束。如果没有风电这个问题就是一个确定性的混合整数规划(MIP)用商业求解器可以可靠求解。但风电一进来事情就变复杂了。风电出力不归调度员管它由天气决定预测误差随预测提前期增长而显著增大。比如提前24小时预测一个风电场的出力误差可能在20%到30%以上遇到天气过程转变时误差更大。如果我们直接把风电预测值当成确定值代入模型得到的机组组合方案往往在实时运行中出问题风电预测偏高时系统被迫紧急下调火电出力甚至弃风预测偏低时可能面临旋转备用不足极端情况下需要切负荷。我见过一个实际案例某地区风电场预测出力150MW实际出力只有80MW当天正好是晚高峰火电机组爬坡速度跟不上最后是紧急启动了调峰机组才稳住频率。这就是确定性UC方案的脆弱性。1.2 随机规划、鲁棒优化各自的问题为了应对不确定性学术界和工业界探索了两条技术路线。第一条是随机规划(Stochastic Programming)。它的思路是假设风电出力服从某个已知的概率分布一般是基于历史数据的经验分布通过蒙特卡洛采样生成大量场景然后求解一个期望成本最小化的优化问题。理论上很完美但实际用起来有几个尴尬之处。第一真实的概率分布我们并不知道只能从历史数据里估计。风速概率分布通常被拟合成Weibull分布但实际风电出力受风机控制策略、尾流效应、风机故障等多重因素影响跟理论分布的偏差并不小。如果估计出的分布不准随机规划的最优解在真实场景下表现可能很差。第二场景数量要足够多才能逼近真实分布但场景越多计算量越大。机组组合本身已经是含大量整数变量的MIP问题了再叠加几百上千个场景模型规模急剧膨胀求解时间让人难以接受。第三条路线是鲁棒优化(Robust Optimization)。它不假设分布已知而是构造一个不确定集合U要求优化方案在U中最坏情况下依然可行。这思路很稳但过度保守。因为不确定集合基本是按照极端情况设计的比如出力上下界区间而真实运行中大多数情况下风电出力都在区间中部不会总是逼近边界。鲁棒优化等于每天都在按最坏情况做打算调度成本会明显偏高实际运行中它的经济性往往过不了关。随机规划怕分布不准鲁棒优化怕过度保守。前者需要精确的概率信息而不可得后者不需要概率信息但代价太高。分布鲁棒优化(Distributionally Robust Optimization, DRO)就是在两者之间找第三条路。1.3 分布鲁棒优化的核心思想不猜分布圈出分布的范围分布鲁棒优化的想法很直接既然精确概率分布拿不到那我们就构造一个分布的集合叫模糊集(Ambiguity Set)把真实分布可能所处的位置圈定在这个集合里。优化的目标是在模糊集内最坏情况下的期望成本最小化。换句话说传统鲁棒优化是在跟“不确定的风电出力”博弈分布鲁棒优化是在跟“不确定的风电出力分布”博弈层次更高、信息利用更充分。它既用了分布信息通过模糊集的结构又保持了对分布估计误差的抵御能力。模糊集的构造方式有很多种。经典的是矩模糊集约束分布的均值、协方差落在给定范围内更现代的是Wasserstein距离模糊集以经验分布为中心限定真实分布到经验分布的Wasserstein距离不超过某个半径。后者在理论上有一个很漂亮的性质随着历史数据增多模糊集半径可以收缩到零分布鲁棒解收敛到真实随机规划解。这也是为什么Wasserstein模糊集近年来成为热门方向。而标题里的“线性准则”指的是在第二阶段决策中引入线性决策规则(Linear Decision Rule)把原本复杂的两阶段min-max-min问题转化为一个规模可处理的优化问题。这个思路后面我会详细展开。2. 线性准则把复杂的min-max-min问题变成可解的优化模型2.1 两阶段框架里的min-max-min结构和难点分布鲁棒机组组合本质上是一个三层的min-max-min问题。最外层min是机组组合决策启停、开机、出力的日前计划中间层max是在模糊集内寻找最坏情况分布内层min是在最坏分布下进行实时经济调度。这个结构的理解方式可以这样想调度员先做决定然后大自然选择对我们最不利的风电分布最后在这个条件下我们以最低成本调整运行方式。三层嵌套导致直接求解在计算上几乎不可行必须做结构上的松弛或改写。标题中的“基于线性准则”正是为了处理这个困难而引入的思路。没有这个改造的话内层的调度优化在每个不确定场景下都要重新求解一个优化问题在模糊集内要对无数个分布求期望这在数学上是极度困难的。2.2 线性决策规则和仿射策略为什么二阶段决策要设计成这种形式线性决策规则的思路是从鲁棒控制里的仿射策略( Affine Policy)借鉴来的。它的核心想法是第二阶段的可调变量比如火电机组的实时出力调整量不确定参数风电预测误差的显式线性函数。用数学式子表述就是实时调整量 基准值 系数矩阵 × 风电预测误差。这里的预测误差是随机变量通过这种线性关系调度变量被显式表达成了不确定参数的函数。引入线性决策规则后min-max-min问题中的内层min函数变成了一个关于不确定参数的线性函数期望运算和求最大值就可以通过对偶理论交换顺序原来的复杂结构被大幅简化。最终问题变成形式为最小化锥规划或半定规划的三层目标可交给商业求解器处理。2.3 线性准则与Wasserstein模糊集的搭配逻辑使用Wasserstein模糊集时对偶转化的过程会特别顺利。Wasserstein距离的Wasserstein模糊集允许通过Lagrange对偶和强对偶定理把内部的最大化问题转化为关于对偶变量的一些约束而这些约束往往是线性矩阵不等式(LMI)或二阶锥约束(SOC)能够被现代数值求解器有效处理。这里有一个比较重要的理论结果当第二阶段决策采用线性决策规则时Wasserstein模糊集下的分布鲁棒机组组合模型可以通过有限维线性矩阵不等式精确描述。这意味着原本无限维的问题被转化成了有限维的凸优化问题求解的可靠性大幅提升。我在实验中发现线性决策规则的决策质量与完全自由的第二阶段决策相比一般只有几个百分点的差距但计算时间可能相差一个数量级以上。在工程实践中这个取舍是非常值得的。2.4 为什么不能跳过线性准则直接求解有人可能会问既然线性决策规则有近似误差为什么不用完全自由的第二阶段决策?从理论上讲在一般情况下内层min问题与中间max问题无法简单交换顺序必须通过强对偶、次梯度分解或Benders分解处理而每次评估最坏情况分布都需要求解一个大规模优化问题嵌套在机组组合的主循环里计算量极其庞大。我在一篇文献里看到过对比数据使用线性决策规则的模型求解时间大约在分钟级而完全自由的二阶锥重构需要小时级甚至天级。对于日前机组组合这种需要每天滚动求解的应用场景线性决策规则几乎是唯一工程上可行的选择。这也是标题中“基于线性准则”是这个方法能落地的原因。3. Matlab代码实现从模糊集到机组组合约束的完整搭建3.1 工具选型YalmipGurobi搭建优化框架我实现这套模型用的工具链是MatlabYalmipGurobi。Yalmip是一个Matlab的建模工具箱它最大的价值在于把优化问题描述为高级语法然后自动调用底层的求解器。对于这种带有线性矩阵不等式和混合整数变量的复杂模型Yalmip可以大幅减少建模的时间和出错概率。Gurobi是目前工业界广泛使用的求解器处理混合整数规划、二阶锥规划的速度和稳定性都相当出色。对学术研究和工程项目来说Gurobi有免费学术许可申请过程比较简单。安装好之后核心调用语句就两行ops sdpsettings(solver, gurobi, verbose, 2); optimize(constraints, objective, ops);实际场景建模的核心是把模糊集、决策变量和约束条件准确表达成Yalmip的语法。后面我会逐步展开。3.2 风电场景库的生成与处理分布鲁棒优化的数据基础是风电出力场景库。我的做法分四步第一步是获取原始数据。可以使用风电场公开数据集比如NREL美国国家可再生能源实验室的风电数据或者国内风电场SCADA系统导出的实际出力数据。时间分辨率通常选15分钟或1小时建议积累至少一年以上以保证覆盖不同季节和天气类型。第二步是数据清洗。剔除负出力、超过额定容量的异常记录以及长时间恒值的数据点。清洗时要注意保留极端出力事件因为极端场景对模糊集的边界形状影响很大。我后面会专门讲这个坑。第三步是场景划分和降维。按季节和典型时段把数据分成多个子集对每个子集用K-means聚类或同步回代消除法进行场景约简。比如原始数据有8760个小时点通过聚类压缩成50个代表性场景每个场景带一个权重代表该类场景在历史数据中的出现频率。第四步是构造经验分布并计算Wasserstein球的半径。在Matlab中可以用wasserstein函数或自己实现最优传输算法来计算两个分布之间的距离。半径的选择与样本数量有关理论上有O(N^(-1/d))的收缩速率实际中我会通过交叉验证选择一个平衡鲁棒性和经济性的半径值。3.3 两阶段机组组合模型的目标函数与约束逐条拆解下面我把模型核心部分逐条列出给出Matlab代码和解释。决策变量分为两个阶段。第一阶段是常规机组的启停状态u(i,t)、开机动作v(i,t)、停机动作w(i,t)以及基准出力p0(i,t)。第二阶段是实时调整量r(i,t,xi)在采用线性决策规则后r是关于随机变量xi的仿射函数。% 常规机组参数 % n_gen: 机组数量, T: 时段数, S: 场景数量 % 第一阶段决策变量 u binvar(n_gen, T, full); % 启停状态 p0 sdpvar(n_gen, T, full); % 基准出力 % 第二阶段决策变量线性决策规则调整量 a b * xi a sdpvar(n_gen, T, full); % 仿射系数 b sdpvar(n_gen, T, full); % 仿射系数 % 风电出力不确定量 xi 在 [-1, 1] 区间波动 % 实时出力 p0 a b .* xi目标函数是最小化总的期望成本包含燃料成本、启停成本和最坏情况下的调整成本。% 成本系数: c2 二次项, c1 一次项, c0 常数项; startup_cost 启动成本 objective 0; for t 1:T for i 1:n_gen % 燃料成本用分段线性近似 objective objective c1(i) * p0(i,t) c0(i) * u(i,t); % 启动成本 objective objective startup_cost(i) * v(i,t); end end约束条件主要分几类。这里我列最重要的几个。功率平衡约束是每个时段必须满足的系统总出力等于总负荷加网损的基本条件。基准出力加期望调整量等于净负荷% 负荷平衡在标称场景下 constraints []; for t 1:T constraints [constraints, sum(p0(:,t)) sum(baseline_wind(:,t)) load(t)]; end不确定性下的功率平衡要求任意xi在最坏情况下也要满足平衡for t 1:T % a和b需要满足对任意xi都成立 constraints [constraints, sum(a(:,t)) sum(b(:,t) .* xi) sum(wind_forecast(:,t)) sum(xi) 0]; end机组出力上下限约束。这个约束在引入仿射策略后要写成对所有xi都成立的形式利用xi的取值区间转化为对系数a和b的线性约束for i 1:n_gen for t 1:T % Pmin * u p0 a b*xi Pmax * u constraints [constraints, Pmin(i) * u(i,t) p0(i,t) a(i,t) - norm(b(i,t), 1) 0]; constraints [constraints, p0(i,t) a(i,t) norm(b(i,t), 1) Pmax(i) * u(i,t)]; end end这里对任意xi∈[-1,1]成立的区间约束自然等价于对a的约束和对b的1-范数约束。这个转化是线性决策规则建模里最经典的一个trick。爬坡约束。相邻时段的出力变化量必须限制在爬坡速率以内for i 1:n_gen for t 2:T constraints [constraints, p0(i,t) a(i,t) - (p0(i,t-1) a(i,t-1)) - (norm(b(i,t),1) norm(b(i,t-1),1)) -ramp_down(i)]; constraints [constraints, p0(i,t) a(i,t) - (p0(i,t-1) a(i,t-1)) (norm(b(i,t),1) norm(b(i,t-1),1)) ramp_up(i)]; end end最小启停时间约束可以通过经典的三个不等式组实现for i 1:n_gen for t 1:T if t min_up(i) constraints [constraints, sum(u(i, t-min_up(i)1:t)) min_up(i) * v(i,t)]; end end end备用容量约束。系统必须预留足够的旋转备用。这里需要把不确定量的边界考虑进去for t 1:T constraints [constraints, sum(Pmax(i) .* u(:,t)) load(t) reserve(t) max_wind_deviation(t)]; end3.4 Wasserstein模糊集与分布鲁棒约束的对偶转化模糊集的建模是分布鲁棒优化的核心。在Yalmip里可以这样建立Wasserstein模糊集% 经验分布场景 % scenario_set: 场景矩阵, S行N列 % wasserstein_radius: Wasserstein半径 % 引入对偶变量 lambda 0 lambda sdpvar(1, 1); constraints [constraints, lambda 0]; % 对偶转化后的约束 % 这里将max_{P in ball} E[Q(xi)] 转化为有限维约束 for s 1:S constraints [constraints, ...]; end对偶转化的核心推导思路是最坏情况分布下的期望成本可以通过引入对偶变量等价为对偶函数在某些约束下的最大值而这个对偶函数是有限维且凸的。完整的推导涉及Wasserstein距离的定义、最优传输的对偶形式、Lagrange乘子法以及强对偶条件。我在附录里会给出一个详细的推导示意。3.5 完整求解流程和参数设置整个Matlab代码的执行流程如下% 主程序: DRO_UC_main.m % 1. 加载场景数据和机组参数 % 2. 构造模糊集参数 % 3. 定义决策变量 % 4. 构造目标函数和约束条件 % 5. 调用求解器求解 % 6. 输出机组组合结果 % 7. 后处理绘图和对比分析 % 关键参数 T 24; % 调度时段数小时 S 100; % 场景数 n_gen 10; % 机组数 wasserstein_radius 0.1; % Wasserstein半径需调参求解结束后输出几类结果各机组的启停状态表、各时段出力计划、总成本、最坏情况分布下的期望成本以及与传统随机规划、确定性优化的成本对比。这里要提醒一下Gurobi在求解混合整数二阶锥规划时有一些参数可以调节比如MIPGap可以设置容忍度LazyConstraints用于添加惰性约束。我在大批量实验时会把MIPGap设置到1e-3速度提升明显且精度足够。4. 常见问题与排查技巧实录4.1 求解器报错“Infeasible”或“Numerical issues”这是我被问得最多的问题。模型不可行绝大多数情况不是数学推导错了而是约束写得过于严格或者参数量级差异太大导致数值不稳定。一套系统性的排查思路第一步是定位冲突约束。先用Yalmip的check函数逐个约束组检查可行性找出是哪组约束在作怪。% 检查约束可行性 check(constraints(1:10));第二步是检查量级一致性。比如机组容量是MW级别而风电波动是kW级别或者成本是$/MWh而备用约束是MW数值差很多个数量级求解器很容易出现数值问题。我的习惯是把所有参数统一归一到相同量纲比如都用pu标幺值。第三步是检查线性决策规则的约束转化是否遗漏了某一项。特别是爬坡约束如果只约束了基准值而忘了处理仿射项就会出现看似合理但实际不可行的情况。第四步是尝试放宽备用约束或增大Wasserstein半径判断模型对参数的敏感性。如果模型在任何一个参数组合下都不可行那大概率是建模本身的问题如果只是某些参数组合不可行那就是参数取值不合理。4.2 Wasserstein半径怎么选才合适Wasserstein半径是一个超参数直接影响鲁棒性和经济性的平衡。半径越大模糊集越大解越保守但越稳健半径越小解越接近经验分布下的随机规划解但鲁棒性变弱。一个实际可行的调参方法是用验证集做交叉验证。把历史数据切成训练集和验证集用训练集构造模糊集并求解然后在验证集上评估实际运行成本和约束违反率。随着半径增大违反率会下降但成本会上升选取违反率降到可接受阈值的最小半径。我在实验中还观察到半径和场景数量的关系也有规律。理论上随着样本数N增大半径按O(N^(-1/(d1)))的速率缩小。如果数据量很大半径可以取得很小数据量少时则需要更大的半径来覆盖分布估计误差。4.3 求解时间爆炸为什么加了线性准则还是慢即便用了线性决策规则分布鲁棒机组组合依然可能很慢尤其是场景数多、机组多的情况下。有几个优化方向值得尝试。第一是减少场景数。如果场景数从500减到50计算时间可能缩短一个数量级。可以用K-means或者快速前向选择法挑选具有代表性的场景尽量保留分布形状。第二是减少整数变量的对称性。多个同类型机组会产生大量对称解相应地给求解器增加负担。Yalmip中可以添加对称破缺约束比如相同类型的机组按照编号顺序开机大幅提升分支定界效率。第三是设置合理的求解器参数。Gurobi的MIPFocus参数可以调整搜索策略对于我们的模型设置MIPFocus1侧重于寻找可行解在工程中往往更高效。4.4 结果异常机组启停频繁震荡怎么处理我遇到过不止一次这样的情况求解出来的机组组合方案中某些机组的启停状态在相邻时段之间来回切换站在工程角度说是完全不能接受的。原因往往是目标函数中启动成本设置不合理或者爬坡约束太松导致求解器倾向于用频繁启停来平衡出力。解决办法是在目标函数中对启停动作加惩罚系数或者直接设置最小运行/停机时间约束这两种方法在工程上都很有效。比如在某台机组刚启动后设置至少运行4小时才能停机这可以用一个简单的整数不等式实现。4.5 数据问题场景库缺失或异常怎么处理这类问题最容易出现在风电场数据缺失或通讯中断的时候尤其是新建风电场。先判断数据缺失类型连续缺失一般是通讯或SCADA故障零星缺失可能是传感器偶尔失灵或数据传输丢包。如果缺失比例低于5%简单插值就够了如果缺失比例较高建议用同一风电场内相邻风机或功率特性相似的机组做参照通过功率曲线的相关性进行估算。补全数据后需要重新验证场景库的统计特征确认机组组合方案没有因为数据变化出现剧烈变动。这里有一个容易被忽视的坑很多文献把极端出力场景当作坏数据直接剔除但对于分布鲁棒优化极端场景恰恰决定了模糊集的边界。清洗掉极端风电场景后系统更容易出现备用不足的告警。处理风电数据时要谨慎尽量保留真实的极端事件只剔除确认的异常记录。4.6 另一个教训Matlab版本差异与随机数种子我在两个不同版本Matlab上跑过同一份代码结果调度成本差了1%左右。最初以为是逻辑错误排查很久后发现是随机数流机制变了导致场景生成结果不同。这个问题在论文复现时特别容易遇到。解决方案是显式设置随机数种子并固定随机数生成器类型% 固定随机数种子,保证实验可复现 rng(2024, twister); % 如果涉及并行计算,还需要固定worker上的随机流 spmd rng(labindex * 1000 2024, twister); end另外如果使用并行工具箱跑场景生成每个worker的随机数流默认是不一致的可以根据worker编号设置不同但固定的种子既能保证可复现又能避免所有worker生成相同场景。5. 数据驱动鲁棒优化的场景化落地路径代码跑通只是第一步要让分布鲁棒优化机组组合真正发挥作用需要把场景库构建、模糊集设定、求解策略和实际业务数据打通。这里分享我梳理的一条落地路径供团队或研究组参考。5.1 场景库构建从数据清洗到场景约简整个流程的第一步是构建高质量的场景库。数据来源建议从风电场SCADA系统导出至少一年的历史出力数据采样间隔选15分钟或1小时。如果只有几周数据场景库的代表性会明显不足。数据清洗时剔除异常记录如负出力、超过额定容量、长时间恒值的故障数据但只剔除确认异常不要把极端场景一起丢掉。然后按季节和时段划分场景子集因为风电出力在夏季和冬季、白天和夜间的统计特性差异很大混合建模会模糊这些差异。如果场景数量过多比如超过5000个建议用K-means或同步回代消除法将场景压缩到50~100个代表性场景并记录每个代表性场景的权重。这样可以在不显著损失分布特征的前提下大幅降低计算量。5.2 模糊集构建两种实用方案的对比模糊集的选择直接决定优化结果和计算复杂度。矩模糊集实现简单、求解效率高但可能过于保守因为实际风速分布往往不是简单的高斯分布。Wasserstein距离模糊集以经验分布为中心用Wasserstein距离圈定半径对尾部风险的刻画更准更贴合风电出力分布特征但求解时要额外处理对偶变换复杂度会高一些。从我的实际测试看如果目标是快速验证算法思路先用矩模糊集把流程跑通再切换到Wasserstein模糊集做精细分析是性价比比较高的路径。5.3 两阶段模型的业务落地平抑偏差与现货市场在实际电力市场中分布鲁棒机组组合的价值主要体现在两处。一是日前计划阶段用分布鲁棒优化生成常规机组启停计划相比传统的确定性优化可以提前预留更多灵活调节空间而不是事后被动应对风电偏差。二是实时调度阶段在日内滚动修正阶段根据最新风电预测调整机组出力但启停状态通常保持日前计划的决策避免机组频繁启停带来的机械损耗和成本上升。我参与过的一个风电场并网项目中风电渗透率大约20%。确定性优化方案在个别强风日会频繁触发弃风和紧急调峰而分布鲁棒方案通过提前在模糊集边界内预留调节容量显著减少了这类事件。虽然日常运行成本略有上升大约2%~3%但在大风季和极端天气频发的时段系统的安全裕度明显更好综合来看更划算。6. 复盘与延展从机组组合到更广的鲁棒决策分布鲁棒优化最吸引我的地方是它提供了一种“中庸而不保守”的决策框架。传统随机优化需要知道精确的概率分布但现实中我们拿到的都是历史数据的经验分布直接假设分布已知本质上是在用一个不确定的输入换一个确定的输出。而传统鲁棒优化又太悲观所有不确定性都按最坏情况处理付出的成本代价往往高得离谱。分布鲁棒优化的定位正好在两者之间——它承认分布本身也是不确定的但通过模糊集把这种不确定圈定在可控范围内。用一句话总结就是与其精确地错不如在合理的范围内稳妥地对。这套方法完全可以迁移到其他含高比例可再生能源的决策场景。微电网能量管理中分布式光伏和负荷的不确定性同样适合用分布鲁棒框架建模储能容量配置与调度中可以将电量不确定性与储能SOC约束耦合优化储能充放电策略虚拟电厂聚合调度中多个分布式资源的不确定性汇聚后模糊集的构建方式会更有意思多区电力系统联络线计划中把区间耦合和风电不确定性同时放进鲁棒优化模型里是工程价值很高的方向。如果后续你有兴趣我也可以单独写一篇关于MatlabYalmip环境下分布鲁棒机组组合的完整代码解析把每个约束的实现细节逐个过一遍包括对偶变换后约束形式怎么写、KKT条件如何处理互补松弛等。毕竟这类带Wasserstein模糊集的模型第一次写对偶代码确实容易劝退。实践出真知建议你先找一份开源的风电数据集把场景库和模糊集的流程跑通再对照机组组合结果去理解每一个约束和参数的含义。等整条链路熟了再根据自己的业务场景调整目标函数和约束条件会比直接套用网上现成代码去调参有收获得多。