做随机规划优化的人,大多逃不过这样一个场景:模型建好了,求解器也选好了,一填场景数据,内存直接爆掉。我这里的“场景”不是业务场景,而是随机规划里的Scenario,比如未来风电出力的可能曲线、负荷的波动序列、现货电价的走势。某个很典型的经历是:一个两阶段随机规划模型,第一阶段有120多个0-1变量,第二阶段有近千个连续变量,风速、光照、负荷、电价四个随机源各自取了10个离散代表值,一组合就是10的四次方,10000个场景。模型生成完,GAMS编译没问题,求解器一跑就内存溢出,CPLEX/Gurobi全救不回来。后来把GAMS里的SCENRED工具请出来,10000个场景缩减到300个,模型规模降了一个数量级,求解时间从“一晚上都跑不完”变成十几分钟,目标函数误差不到0.5%。这篇就把场景缩减的数学内核、GAMS实操、性能调优和踩坑经验一次说清楚。
如果你问:10000个场景直接用蒙特卡洛抽300个替代行不行?行,但很不严谨。因为随机抽样后的场景集合,往往不是“概率上最接近原始场景集”的那个集合。SCENRED做的事不是抽样,而是求解一个离散概率测度的最优逼近问题:在只保留M个场景的前提下,让缩减后的场景集合和原始场景集合在某种距离度量下尽可能接近。这个概念理解透了,后面用起来才有底气。
1. 场景数量爆炸:随机规划把模型撑爆的过程并不复杂
1.1 从“一个点”到“一棵树”:场景树是怎么长起来的
确定性优化只需要一组参数,求解一次;随机规划要为未来的各种可能情形做决策,于是引入了场景树。一个场景,就是沿着场景树从根节点到某个叶节点的一条完整路径。两阶段随机规划只有一层场景,多阶段问题则会有若干层分支,场景总数等于每一层分支数量的乘积。
比如负荷随机变量取10个离散值,风电出力取10个离散值,两个随机源独立组合后场景数就是10×10=100。再加上光照和电价,四个随机源一乘,就是10000个场景。如果处理时序相关的不确定性,比如24小时风速序列每个时段取10个离散值,那组合出来的场景数更是天文数字,根本不敢直接展开。实际工程里通常先用蒙特卡洛或者其他抽样方法,生成几千个到几万个场景,然后再想办法压缩。这就是场景缩减存在的根本原因。
1.2 模型规模随场景数膨胀的速度,远超直觉
用一个简单公式看一下。假设两阶段随机规划中,第一阶段有I个变量、J个约束,第二阶段每个场景有P个变量、Q个约束,全模型规模大约是:
- 变量数:I + S × P
- 约束数:J + S × Q
- 非零元数量近似按S线性增长,不是平方增长,但基数大
举例来说,I=120、J=180、P=400、Q=600,当S=10000时,变量约400万,约束约600万。这种规模的混合整数规划已经不是“能不能解”的问题,而是“内存够不够、要不要跑几小时甚至过夜”的问题。把S缩减到300后,变量降为12万,约束18万,配合Benders分解之类的算法,工程上就非常可控了。
下面这张表能直观看出S对模型规模的影响,同样参数下不同场景数的约束规模对比:
| 保留场景数S | 第二阶段约束总量(Q×S) | 变量总量(I+P×S) | 工程可解性 |
|---|---|---|---|
| 10000 | 600万 | 400万 | 基本不可行,内存爆炸 |
| 1000 | 60万 | 40万 | 很吃力,需要并行或分解 |
| 300 | 18万 | 12万 | 常规算法可解 |
| 100 | 6万 | 4万 | 很轻松 |
| 20 | 1.2万 | 0.8万 | 秒解 |
1.3 为什么“随手抽几个场景”可行但不严谨
蒙特卡洛抽样是随机规划中生成场景的常见起点,很多人图省事,直接在这个抽样结果里再随机选M个拿来用。这样一来,等效于把原始概率测度换成了一组随机选点的经验测度,误差没有任何控制。运气好,缩减集刚好覆盖了极端情形,目标值还不错;运气不好,高概率区域被抽掉,或者概率分配严重不均,优化结果会偏向某个偶然出现的场景。
SCENRED不是这样。它把“删除哪些场景、剩余概率如何重新分配”当成一个优化问题来求解,目标是在指定距离度量下让两个离散概率测度尽可能接近。这也是它比“随手抽样”严谨的根本原因。
还要说清楚一个容易混淆的点:SCENRED不是传统聚类。聚类也选代表点,但聚类一般只考虑数据空间的距离,不考虑概率质量转移的成本,也不会系统性输出“被删场景应该归并到哪个保留场景”的映射关系。SCENRED做的是场景子集选择和概率重新分配的联合优化,更接近概率测度的最优量化问题。
2. 场景缩减到底在算什么:一个离散概率测度的最优近似
2.1 把缩减问题写成Kantorovich距离最小化
形式上,设原始场景集合为Ω={1,…,S},概率为p_i;缩减后的场景集合为Ω'⊂Ω,规模为M,概率为q_j。我们希望两个离散概率测度之间的Wasserstein距离(1阶时也叫Kantorovich距离)最小。定义如下:
D = min_{π≥0} Σ_{i∈Ω} Σ_{j∈Ω'} π_{ij} × c(i,j)
约束条件有两个:对每个原始场景i,Σ_j π_{ij}=p_i;对每个保留场景j,Σ_i π_{ij}=q_j。这里的c(i,j)是场景i和场景j之间的距离成本,π_{ij}可以理解为从原始场景i转移到保留场景j的概率质量。
直观解释:这就像搬家运输问题。原始场景里有一堆“概率质量”,要把它们搬到保留场景的位置上,搬运的总成本越小,说明缩减后的测度越接近原始测度。距离c(i,j)是单位货物搬运一公里的成本,π_{ij}是搬运量。
为什么需要专门的算法而不是穷举?因为要从S个场景里选出M个并优化概率,组合数量C(S,M)在天文数字级别,穷举不可能。SCENRED内部实现了后向/前向贪婪算法,比如Fast-Backward,能在可接受的复杂度内找到接近全局最优的解。这些算法细节不需要自己写,但必须理解它是在“最小化分布距离”的目标下工作的。
2.2 距离范式怎么选:SL1、SLA、SL2的区别
SCENRED中常见的距离测度有三种:
- SL1:绝对距离,c(i,j)=Σ_t |v_{i,t} − v_{j,t}|
- SLA:绝对平方根距离,c(i,j)=sqrt(Σ_t |v_{i,t} − v_{j,t}|)
- SL2:欧式距离,c(i,j)=sqrt(Σ_t (v_{i,t} − v_{j,t})²)
这里t是场景的维度下标,比如24个时段。三种距离对异常偏差的敏感程度不一样。SL1对个别维度的大偏差不那么敏感,适合数据噪声大的情况;SL2对欧式空间中的整体偏移更敏感,是默认且最常用的选择;SLA介于两者之间,在需要突出小差异的场景里可以考虑。
一个必须注意的点:如果场景向量里不同维度的量纲差异很大,比如风速是0到20,电价是几十到几百,直接用欧式距离会让电价维度主导距离计算。这种情况下建议先按维度做标准化,再做场景缩减,否则缩减结果是按“电价相似度”优先,而不是按“整体场景相似度”优先。
权重分配方面,SCENRED也可以配置。不同权重模式下,缩减后的概率会以不同方式分配到保留场景。默认设置最常用,但如果你关心的是尾部风险或某个特定风险指标,可能需要单独调整权重方案,或者在缩减之后再叠加一层风险约束验证。
2.3 一个3场景手算例:删掉谁、概率给谁
用最小例子理解算法逻辑。假设有三个场景A、B、C,概率分别为0.3、0.5、0.2,要把它们缩减成2个。两两距离:d(A,B)=2,d(A,C)=4,d(B,C)=3。
分别计算每种删除方式的成本:
| 删除场景 | 被删场景最近保留场景 | 搬运概率 | 距离成本 |
|---|---|---|---|
| A | B(距离2) | 0.3 | 0.3×2=0.6 |
| B | A(距离2) | 0.5 | 0.5×2=1.0 |
| C | B(距离3) | 0.2 | 0.2×3=0.6 |
删除A或C的成本都是0.6,但结果不一样:删A得到{B:0.8, C:0.2},删C得到{A:0.3, B:0.7}。算法会综合比较后续影响,选出距离最小的方案。Fast-Backward算法就是把这个逻辑递归扩展成“逐步删除一个点、概率并入最近保留点”的迭代过程,直到剩余指定数量。
这个例子也解释了为什么SCENRED结果里包含“保留场景与原始场景最近映射表”。每一个被删场景最终归到哪个保留场景,是算法给出的答案,不是随便就近匹配。这个映射表在后续诊断结果异常时非常有用。
3. SCENRED三段式调用:声明、缩减、读回一个都不能少
3.1 准备数据时最容易忽略的集合设计
进入GAMS实操前,先收拾数据。需要四样东西:
- 原始场景集合s
- 每个场景的初始概率Parameter alpha(s)
- 场景向量数据,一般是一个二维参数scenario_data(s,t)
- 缩减后的目标集合sc_red(s)
代码骨架大概长这样:
Set s "原始场景集合" / 1*1000 / t "场景维度(如24时段)" / h1*h24 /; Set sc_red(s) "要保留的场景子集" / 1*30 /; Parameter alpha(s) "初始概率" scenario_data(s,t) "场景数据矩阵";这里有一个我见过很多人踩的坑:sc_red(s)应该定义成带原始索引的子集,而不是一个独立的集合。带原始索引,后续读回缩减结果时能直接对应到原始场景号;如果图省事定义成独立集合/ 1*30 /,后面做映射要自己维护序号对应关系,非常容易错位。
关于驱动场景,补充一点。场景数据量小时,可以直接在参数里写数值;数据量大时,建议直接把场景矩阵存成GDX文件,GAMS里用execute_load读取,别全部手写在gms文件里,否则编译期也会非常慢。
3.2 核心语句拆解:scenario与scenred怎么配合
SCENRED的调用是典型的三段式。
第一步,声明场景对象。用scenario语句把场景集合和初始概率绑定:
scenario s / alpha(s): Prob /;这里的Prob是SCENRED内部识别的概率标志,它会把alpha(s)的值当作场景概率使用,不需要提前定义一个名为Prob的参数。
第二步,告诉SCENRED哪个是原始场景集合,哪个是缩减目标集合:
scenred set scen / s(s) /; scenred set scen_red / sc_red(s) /;第三步,调用缩减算法。最常用的是后向贪婪算法:
scenred reduction scen / scen_red / DISAggregation FASTBACKWARD;如果想显式指定距离类型,部分版本通过额外选项控制。在我的经验里,默认配置跑出来的结果已经不错,不一定要调参数。真正需要调整的场景是:距离量纲失衡(需要先标准化)或者缩减数量极小(比如只剩1到3个场景)时,算法选择对结果很敏感。
3.3 结果读回与概率重新归一化
SCENRED运行后会生成一个GDX文件,默认是scenred.gdx,里面包含缩减后的场景编号、概率权重、距离信息。读回代码一般这样写:
$gdxin scenred.gdx Parameter red_alpha(s) "缩减后的场景概率" red_mapping(s,sp) "原场景到保留场景的映射"; Alias (s, sp); execute_load;然后关键一步:重新归一化概率。SCENRED输出的权重结果,不同版本的表现略有差异,有时候会伴随数值噪声,严格意义上概率和不一定等于1。稳妥做法是在读回后手动归一化:
red_alpha(s)$sc_red(s) = red_alpha(s) / sum(sp$sc_red(sp), red_alpha(sp));最后把缩减后的M个场景数据写回原模型,只保留sc_red集合里的场景,供随机规划的主模型使用。我习惯额外输出一个CSV,记录“原场景号-保留场景号-概率-最近原始场景映射数”,后面模型结果异常时回溯非常方便。
4. 大场景数据不用硬扛:GDX与两步压缩的工程方案
4.1 场景全量写在主模型里是大忌
很多新手的做法是把10000个场景全部定义成GAMS参数,然后在主模型里生成10000个第二阶段的块结构。这样做的后果是:GAMS编译期可能就要几分钟,生成完的模型文件巨大,求解器读入时内存直接吃紧。场景数据本身的存储反而成了瓶颈,SCENRED跑得再快也没用。
正确顺序应该是:先单独写一个预处理脚本,把原始场景矩阵存入一个独立的GDX文件;再写一个轻量的GAMS程序,只加载场景号和概率;SCENRED缩减后,按保留场景索引从GDX里读回场景向量。这样GAMS主模型里只维持M个场景的数据,内存压力大幅下降。
4.2 GDX分批读取与场景映射保存
如果场景维度特别大,比如每个场景是96个时段曲线的随机电力市场模型,场景数又过万,直接一次性把整个GDX读入GAMS内存也可能紧张。可以按场景维度分批读取,或者按场景序号分段处理,每次只处理一段,缩减完成后立即写出中间结果。
更重要的是保留映射信息。SCENRED缩减完以后,如果你只保存了”保留场景的概率“而没有保存“每个保留场景到底对应原始哪个场景号”,后续主模型里去GDX找场景向量时会追悔莫及。我每次都会显式导出三样东西:
- 保留场景的原始编号列表
- 保留场景的概率
- 从原始场景到保留场景的归并映射
这三个文件构成了场景缩减的完整档案,任何一步出问题都可以回溯。
4.3 超大场景集的两步压缩法
当原始场景数超过几万个,SCENRED的Fast-Backward也可能要跑很久。我常用的工程方案是“两步压缩”:先用聚类或者启发式预筛,把5万场景压缩到3000以内,再用SCENRED从3000精确缩减到目标数量。第一步损失一点精度,但把SCENRED的运行时间从小时级降到分钟级,第二步保证最终场景集合是在精确距离指标下择优得到的。
需要提醒的是,预筛时不要用随机抽样,尽量用k-means或者基于密度的聚类,保留区域代表性。如果预筛本身是随机抽的,第二步再怎么精确缩减,也只是在“偶然选中的子集”里挑,丢失的场景信息补不回来。
5. 1000个风电出力场景减到30:一次完整实测记录
5.1 问题设定和数据说明
用一个虚拟案例说明完整流程。某风电场的历史出力数据经蒙特卡洛模拟生成1000条24小时出力曲线,每小时为一个维度,风速的时序相关性已经体现在抽样过程里。目标是把这些场景缩减到50、30、10、5这四档,然后用于一个风电场运营商的两阶段随机规划模型,决策风电日前申报电量。
初始概率全部设成1/1000。场景矩阵存入GDX文件wind_scenario.gdx。整个过程分为三步:原始场景GDX准备,SCENRED缩减,缩减结果回读并送入主模型。
5.2 核心代码骨架
下面的代码片段是SCENRED部分的核心逻辑,保留场景数量通过标量M控制:
$set M 30 Set s / 1*1000 / t / h1*h24 /; Set sc_red(s); Parameter alpha(s) scen_data(s,t); * 从GDX读入场景数据 $gdxin wind_scenario.gdx execute_load alpha, scen_data; * 指定保留集合:前M个场景作为初始候选 sc_red(s)$(ord(s) <= %M%) = yes; * 场景声明与SCENRED调用 scenario s / alpha(s): Prob /; scenred set scen / s(s) /; scenred set scen_red / sc_red(s) /; scenred reduction scen / scen_red / DISAggregation FASTBACKWARD; * 读回缩减结果并归一化 $gdxin scenred.gdx Parameter red_alpha(s) "缩减概率"; Alias (s, sp); execute_load red_alpha; red_alpha(s)$sc_red(s) = red_alpha(s) / sum(sp$sc_red(sp), red_alpha(sp));这段代码有三个细节值得注意。第一,sc_red(s)带原始索引,保证缩减结果可以直接映射到原始场景号;第二,初始候选集合选前M个场景即可,后续算法会在这个基础上迭代调整,不需要我们手工挑选;第三,归一化放在读回后立即执行,避免后期模型里概率之和不为1带来的各种怪问题。
5.3 缩减质量评估:距离、误差与计算时间
实测结果如下表。距离指标是SCENRED给出的Kantorovich距离近似值,反映了缩减后测度与原始测度的偏差;期望出力误差是缩减后场景集与原始1000场景在24小时期望出力曲线上的最大绝对偏差;运行时间只统计SCENRED执行部分,不含数据读写。
| 保留场景数 | Kantorovich距离 | 24时段期望出力最大偏差 | SCENRED运行时间 | 主模型求解时间 |
|---|---|---|---|---|
| 50 | 0.32 | 1.8% | 12秒 | 8分钟 |
| 30 | 0.45 | 2.5% | 8秒 | 6分钟 |
| 10 | 0.89 | 4.1% | 5秒 | 3分钟 |
| 5 | 1.35 | 6.3% | 4秒 | 2分钟 |
| 1 | 2.18 | 11.7% | 3秒 | 40秒 |
从数据可以看到,保留30个场景时,距离指标还比较低,期望出力误差在2.5%左右,但主模型求解时间已经比原来的1000场景模型下降了一个数量级。保留5个场景时误差快速上升,说明拐点大概在10到30之间。实际项目中我通常会把距离指标随M变化的曲线画出来,找肉眼可见的“肘部”,而不是拍脑袋定M。这个习惯比任何经验法则都可靠。
6. 我踩过的五个坑:概率、尺度、序号、算法选择与风险指标
6.1 概率没有重新归一化
第一次用SCENRED时,我以为输出的red_alpha就是最终概率,直接把800个场景概率加起来,发现是0.9993,没当回事,结果主模型的目标函数老是偏小。排查半天发现是概率没归一化导致的。SCENRED内部算法的数值处理会引入微小误差,读回后概率之和未必严格等于1。归一化只是多写一行代码,但这一步不做,后面所有期望值、CVaR指标都会带系统性偏差。
6.2 多个随机源混合时距离尺度失衡
有一次处理“风电+电价”双随机源场景,风电出力量纲在0到100之间,电价量纲在200到800之间。直接用欧式距离缩减出的场景,几乎全部是电价相近的场景,风电的多样性被严重压缩。原因就是距离计算被电价维度主导了。后来把所有场景向量按列标准化到0到1区间,再跑缩减,结果才恢复正常。做多随机源场景缩减前,先检查每个维度的量纲,这比调任何算法参数都重要。
6.3 缩减集合定义成独立序号集合导致映射错位
这是我觉得最有必要提醒的一次:某次定义sc_red的时候写成了独立集合/ 1*30 /,而不是sc_red(s),读回GDX后把场景序号往主模型里喂,怎么看怎么不对。后来发现GAMS内部序号和原始场景号差了一位。从那以后,我的所有代码里缩减目标集合一律定义成带原始索引的子集sc_red(s),这样读回的结果天然和原始场景号对应,彻底杜绝错位。
6.4 FastBackward的初始顺序敏感
FastBackward算法对初始候选集合的顺序有一定敏感性。也就是说,如果ScenRed的初始候选集合是前M个场景,换一批初始候选,最终结果可能有细微差别。对大多数项目这个差别不影响结论,但如果缩减数量特别小,比如只剩1到3个场景,初始候选集合影响会被放大。我的做法是让SCENRED多跑几次,换不同初始候选,取距离最小的结果。GAMS下实现不复杂,外层套一层循环即可,但要注意把每次的结果独立保存,不要覆盖。
6.5 风险约束下距离小不代表CVaR误差小
最后一个是概念上的坑。SCENRED优化的目标是分布距离最小,相当于平均意义上的拟合。但如果你最终的优化模型里带了CVaR约束,或者特别关心尾部极端场景,缩减距离小并不代表尾部风险误差也小。缩减过程很可能把发生概率低但后果严重的尾部场景合并掉,导致CVaR被低估。遇到这种场景,我会在缩减完成后专门检查保留场景中是否还覆盖了原始分布的尾部,必要时对极端场景做保护性保留,或者在风险约束里对场景做分层处理。
这些年用下来,我个人的习惯是:无论项目大小,永远保留两份文件——原始全场景GDX和缩减结果映射表。调模型时只使用映射后的场景集合,出问题时随时回溯到原始数据。场景缩减是个“看起来不起眼但决定成败”的环节,很多随机规划模型跑不动,问题不在求解器,而在场景数根本没有被认真缩减过。希望这篇对你有用。