1. 为什么风光场景生成非得做削减这一刀做新能源电力系统规划或调度优化的朋友对场景法应该都不陌生。风电和光伏出力受天气影响具有强烈的随机性和波动性我们没办法用一个确定的曲线去描述明天风多大、太阳多强只能通过历史数据或概率模型生成大量可能出现的出力场景。理论上场景数量越多对不确定性刻画越精确但随之而来的问题是计算量呈指数级膨胀。举个例子一台风机的出力场景如果生成1000个那么一个含10台风机、5个光伏电站的系统组合起来的场景总数就是天文数字。调度模型里每个场景都对应一组约束条件求解器面对成千上万的场景约束轻则求解时间从分钟级变成小时级重则直接内存爆炸无法求解。所以业界普遍的做法是先生成大量场景再用某种方法削减到少量有代表性的场景用这几十个场景去逼近原始上千场景的统计特性。这个思路听起来简单但削减方法选得不好结果会非常难看。我见过有人直接用K-means聚类做场景削减聚类中心倒是选出来了但聚类算法本质上是把场景当成点来划分它压根不考虑场景之间的转移概率和时序相关性削减后的典型场景集在概率分布上往往失真后续优化结果偏差很大。而概率距离快速削减法也叫基于Kantorovich距离的场景削减法简称SD场景削减是目前学术界和工业界都认可的主流方法。它的核心逻辑不是找聚类中心而是通过计算场景之间的概率距离矩阵逐步剔除对整体概率分布影响最小的场景并把被剔除场景的概率累加到距离它最近的保留场景上。每一步削减都保证原始场景集与削减后场景集的概率距离增量最小因此能在给定场景数条件下最大程度保留原始分布特性。MATLAB实现这套方法代码量不大但有几个坑非常隐蔽。本文我用一个完整的风电光伏联合场景生成与削减案例把原理、代码、调参经验和验证方法一次讲透。适合正在做电力系统随机优化、微电网容量配置、风电光伏消纳分析的研究生和工程师参考。2. 风电光伏联合场景生成先弄懂我们削减的对象是什么2.1 风力和光伏出力的数学模型怎么建要削减场景第一步是得有场景。场景生成的方法很多常见的有基于历史数据的直接采样、基于ARIMA时间序列模型、基于马尔可夫链、基于Copula函数、基于概率密度函数的蒙特卡洛采样等。不同的方法适合不同的数据条件。在实际项目中我通常采用概率分布模型时序相关性还原的方式。风电出力近似服从Weibull分布光伏出力近似服从Beta分布但直接按这两个分布独立采样出来的场景毫无时序逻辑——风不会上一秒满发下一秒突然归零光伏也不会中午出力呈白噪声式跳动。所以完整的场景生成链路是对风速历史数据拟合Weibull分布通过风速-功率转换曲线得到风电出力序列。对光照强度历史数据拟合Beta分布通过光照-功率转换模型得到光伏出力序列。引入自相关系数和互相关性让生成的场景在时间维度上平滑、在风-光之间具有合理的相关性。用蒙特卡洛采样生成N个初始场景每个场景包含24个时段的双维出力风电光伏。需要说明的是场景削减算法本身并不关心场景是怎么生成的它只负责在给定初始场景集和概率权重通常假设等概率的前提下找到最优的保留子集。所以哪怕你用别的方法生成场景只要格式是场景数×时段数×变量数后面的削减代码可以直接复用。2.2 初始场景规模如何确定初始场景规模N的选取直接影响生成精度和削减效果。太少了削减后的典型场景缺乏代表性太多了计算概率距离矩阵时的时间和内存开销会显著增加。我的经验法则是单变量场景至少1000个起步风光联合场景建议2000到5000个。1000个以上场景经过削减到10-20个统计特性保留效果较好如果初始只有200个场景削减到10个之后误差会明显偏大。这里贴一段生成初始风电场景的MATLAB核心代码采用Weibull分布采样并加入时序自相关处理% 参数设置 nScenarios 2000; % 初始场景数 nPeriods 24; % 调度时段数小时 k 2.1; % Weibull形状参数根据历史风速拟合 c 8.5; % Weibull尺度参数m/s vRated 12; % 额定风速 vCutIn 3; % 切入风速 vCutOut 25; % 切出风速 Prated 1.5; % 额定功率 MW % 生成独立Weibull风速样本 windSpeed wblrnd(c, k, nScenarios, nPeriods); % 加入时序平滑简单一阶自回归让相邻时段风速平滑过渡 rho 0.85; % 自相关系数 for t 2:nPeriods windSpeed(:, t) rho * windSpeed(:, t-1) ... (1-rho) * windSpeed(:, t) randn(nScenarios, 1) * 0.3; windSpeed(windSpeed 0) 0; end % 风速→功率转换 windPower zeros(nScenarios, nPeriods); for s 1:nScenarios for t 1:nPeriods v windSpeed(s, t); if v vCutIn || v vCutOut windPower(s, t) 0; elseif v vRated windPower(s, t) Prated; else % 线性区间的近似功率曲线 windPower(s, t) Prated * (v - vCutIn) / (vRated - vCutIn); end end end光伏场景的生成逻辑类似区别在于Beta分布拟合的是光照强度的归一化值并且时序上要体现白天出力、晚上归零的强周期特征。联合场景就是把两个矩阵横向拼接每个场景行向量维度变成4824时段风电24时段光伏。3. 概率距离快速削减法原理每一步剔除都让分布损失最小3.1 Kantorovich距离到底是什么很多教程一上来就堆公式把人看懵。我用一句话概括概率距离削减法的本质它在反复问一个问题——如果把某个场景删掉把它的概率转移给离它最近的另一个场景整体概率分布变化有多大如果这个变化是所有候补删除项中最小的那就删它。衡量变化有多大的标尺就是Kantorovich距离也叫Wasserstein距离。两个场景集合之间的Kantorovich距离可以理解为把一个概率分布搬到另一个概率分布所需的最小代价这里的代价就是场景差异乘以转移概率。对于离散场景集Kantorovich距离的计算公式可以转化为一个最优运输问题。但在场景削减这个特定情境下我们不需要求解完整的最优运输问题因为削减过程是逐步进行的每一步从一个场景集中删除一个场景并把其概率累加到最近邻场景上。这种迭代贪婪策略得到的削减结果在每一步都是局部最优的实际应用中已经足够逼近全局最优。3.2 算法的完整流程拆解概率距离快速削减法的具体步骤如下初始化有N个场景每个场景概率权重p_i 1/N也可根据历史频率赋权。计算所有场景两两之间的欧氏距离d(i, j)得到N×N对称距离矩阵。距离定义可以是逐时段绝对误差之和L1范数也可以是欧氏距离L2范数常用的是L2范数。对于每个场景i找到与其距离最近的场景jj ≠ i记录距离d_min(i) min d(i, j)。计算每个场景的可删除代价cost(i) p_i × d_min(i)。这个值的物理含义是删除场景i后概率损失p_i乘以它不得不转移给最近邻居的距离代价d_min(i)。选出cost最小的场景k把它删除并将其概率p_k累加到最近邻居j上p_j p_k p_j。更新距离矩阵删除第k行第k列重复步骤3-5直到保留的场景数量达到预设值M。用伪代码表达就是while nScenarios nTarget % 对每个场景找最近邻居距离 minDist zeros(nScenarios, 1); minIdx zeros(nScenarios, 1); for i 1:nScenarios dist_i distMatrix(i, :); dist_i(i) inf; % 排除自身 [minDist(i), minIdx(i)] min(dist_i); end % 计算删除代价 deleteCost probWeights .* minDist; % 找到代价最小的场景 [~, delIdx] min(deleteCost); % 概率转移给最近邻居 probWeights(minIdx(delIdx)) probWeights(minIdx(delIdx)) probWeights(delIdx); % 删除场景和对应行列 distMatrix(delIdx, :) []; distMatrix(:, delIdx) []; probWeights(delIdx) []; nScenarios nScenarios - 1; end注意这里有个细节距离矩阵每删除一个场景就要更新一次因为删除后某些场景的最近邻居可能发生变化必须重新计算最近邻距离。这也是该算法时间复杂度为O(N³)的原因——每轮删除都要重新扫描一遍距离矩阵。不过对于几千个场景来说MATLAB跑起来也就几秒到几十秒完全可以接受。3.3 为什么不是一次删一批有一种快速近似做法是预先计算所有场景的距离矩阵然后一次性选出多个待删除对象再统一做概率转移。这样做速度更快但存在一个问题第一批被判定为可删的场景中可能会有两个场景互为最近邻居。如果同时删除它们各自的概率都要转移到对方身上就出现了概率无处安放的矛盾。更合理的做法是每删除一个场景立即更新距离矩阵并重新计算最近邻保证每一步操作都是基于最新场景集的。虽然计算量更大但结果更可靠。对于超大规模场景上万个可以考虑用KD-Tree等近似最近邻搜索加速但对于电力系统场景削减的典型规模暴力计算完全够用。4. MATLAB完整实现从距离矩阵到削减主循环4.1 场景切削函数封装实际项目中我不会把主脚本堆成一坨而是封装成一个通用函数输入初始场景矩阵和期望保留数输出削减后的场景、对应概率和削减过程中的误差记录。以下是完整的MATLAB函数代码function [reducedScenarios, reducedProb, distHistory] ... scenarioReduction(scenarios, probs, targetCount) % 基于概率距离快速削减法的场景削减函数 % 输入 % scenarios - 初始场景矩阵维度 [N, T]N为场景数T为时段数含多变量拼接 % probs - 场景概率向量维度 [N, 1]缺省时默认等概率 % targetCount - 目标保留场景数 % 输出 % reducedScenarios - 削减后场景矩阵维度 [targetCount, T] % reducedProb - 削减后场景概率向量维度 [targetCount, 1] % distHistory - 每次删减后的Kantorovich距离记录用于观察削减误差 N size(scenarios, 1); if nargin 2 || isempty(probs) probs ones(N, 1) / N; end if nargin 3 error(必须指定目标场景数); end if targetCount N error(目标场景数必须小于初始场景数); end % 计算初始距离矩阵欧氏距离 distMatrix squareform(pdist(scenarios, euclidean)); % 避免对角线为0置为无穷大 distMatrix(1:N1:end) inf; distHistory zeros(N - targetCount, 1); iter 0; while size(distMatrix, 1) targetCount iter iter 1; currentN size(distMatrix, 1); % 找每个场景的最近邻居 [minDist, minIdx] min(distMatrix, [], 2); % 计算每个场景被删除时造成的概率距离增量 deleteCost probs(:) .* minDist(:); % 选择代价最小的场景删除 [~, delIdx] min(deleteCost); % 概率累加到最近邻居 neighborIdx minIdx(delIdx); probs(neighborIdx) probs(neighborIdx) probs(delIdx); % 记录当前的累积削减误差用当前最小代价近似单步损失 distHistory(iter) deleteCost(delIdx); % 从距离矩阵和概率向量中删除该场景 distMatrix(delIdx, :) []; distMatrix(:, delIdx) []; probs(delIdx) []; end % 最后提取保留的场景需要知道原场景索引 % 由于不断删行距离矩阵的行序和场景矩阵的行序在删除后不一致 % 这里需要用逻辑索引追踪原位置简化方式是通过在削减过程中同步删除场景矩阵的行来实现。 % 由于函数输入未返回索引这里采用同步删除方式重建场景矩阵上面代码有个需要注意的地方我在函数开头提到需要追踪原场景索引。实际上最稳妥的做法是在删除距离矩阵行列时同步删除场景矩阵中的对应行否则最后不知道哪些原始场景被保留了。修正后的核心循环如下reducedScenarios scenarios; idxMap (1:N); while size(reducedScenarios, 1) targetCount currentN size(reducedScenarios, 1); distMat squareform(pdist(reducedScenarios, euclidean)); distMat(1:currentN1:end) inf; [minDist, minIdx] min(distMat, [], 2); deleteCost probs(:) .* minDist(:); [~, delIdx] min(deleteCost); neighborIdx minIdx(delIdx); probs(neighborIdx) probs(neighborIdx) probs(delIdx); % 同步删除场景行和概率 reducedScenarios(delIdx, :) []; probs(delIdx) []; end每轮循环重新调用pdist计算距离矩阵比维护一个动态更新的距离矩阵更直观。虽然理论上多次调用pdist会增加计算量但对几千个场景、48维的场景矩阵来说MATLAB性能完全扛得住代码可读性却高了很多。4.2 主脚本生成场景→削减→可视化主脚本的逻辑很清晰生成2000个风光联合场景削减到10个然后画出削减前后对比图。代码如下clear; clc; close all; rng(2024); % 固定随机种子保证结果可复现 % 生成初始场景以2.2节的windPower为基础增加光伏部分 % 此处直接载入预生成的场景矩阵实际使用时替换为自己的场景生成代码 % scenarios: [2000, 48]前24列风电后24列光伏 load(initialScenarios.mat, scenarios); N size(scenarios, 1); probs ones(N, 1) / N; % 设置削减目标数 M 10; % 调用削减函数 [redScen, redProb, histErr] scenarioReduction(scenarios, probs, M); % 输出结果 fprintf(初始场景数: %d削减后场景数: %d\n, N, M); fprintf(削减后概率分布: \n); disp(redProb(:)); % 可视化削减前后某几个场景的对比 figure; subplot(2,1,1); plot(1:24, scenarios(1:50, 1:24), Color, [0.7 0.7 0.7]); hold on; plot(1:24, redScen(:, 1:24), LineWidth, 1.5); xlabel(时段/h); ylabel(风电出力/MW); title(风电场景削减前后对比); legend(初始场景, 削减后典型场景); grid on; subplot(2,1,2); plot(1:24, scenarios(1:50, 25:48), Color, [0.9 0.7 0.7]); hold on; plot(1:24, redScen(:, 25:48), LineWidth, 1.5); xlabel(时段/h); ylabel(光伏出力/MW); title(光伏场景削减前后对比); legend(初始场景, 削减后典型场景); grid on;运行后削减后的10个典型场景会在2000个原始场景的云图中清晰勾勒出主要分布模式。我截取一张典型输出效果描述风电场景削减后的曲线基本覆盖了从低出力到高出力的各个典型水平光伏场景则呈现出明显的午间高峰形态概率最大的场景是整个分布中最普通的那种——白天光伏高发、晚间风电接力。5. 削减效果怎么评估三种Kantorovich距离验证方法5.1 削减前后概率距离的直接计算最直观的评估方法是计算削减前场景集与削减后场景集之间的Kantorovich距离。这个距离越小说明削减后的场景集对原始分布的还原度越高。MATLAB中可以通过调用Wasserstein距离的相关实现或直接按定义用线性规划求解最优运输问题。对于中等规模问题用transport函数或者CVX求解即可。但实际中我们更关心的是随着目标场景数M从1增加到50削减误差怎样变化把这个变化曲线画出来就能直观看到用多少个场景能在精度和计算量之间取得平衡。我的经验是对于风力发电场景M10时距离通常能降到初始距离的5%以内继续增加M边际收益开始大幅递减。M_values 1:10:50; errorValues zeros(size(M_values)); for i 1:length(M_values) [~, ~, ~] deal([]); [redScen, redProb, histErr] scenarioReduction(scenarios, probs, M_values(i)); errorValues(i) sum(histErr); % 累计削减误差 end figure; plot(M_values, errorValues, o-, LineWidth, 1.5); xlabel(保留场景数); ylabel(累计削减误差); title(削减误差随目标场景数的变化); grid on;5.2 与K-means和K-medoids聚类的对比很多场景削减文献会对比概率距离法与聚类法的差异。我建议你在自己的算例上也做一次对比因为在某些数据结构下K-medoids和概率距离法的结果很接近但概率距离法在概率分配上天然更合理。K-means聚类削减的问题在于聚类得到的簇中心是一个虚拟场景它可能不代表任何真实存在的出力模式而且它无法保留场景的出现概率概念。K-medoids虽然选取的是真实场景作为簇中心但每个簇内场景被视为等权重忽略了簇内场景的密度信息。概率距离法则通过概率转移机制让概率密度高的区域保留更多场景概率密度低的区域被合并。直观理解如果某个区域有400个相似场景削减后这个区域可能保留4个场景每个概率约0.1如果另一区域只有50个场景削减后可能只保留1个场景。这个特性对后续优化计算特别友好因为概率权重可以被直接用作随机规划中的场景概率。5.3 进一步检验用削减后场景求解优化问题对比这是最实际也最有说服力的验证方式。假设你正在做微电网容量配置先用2000个场景求解一版结果再用削减后的10个场景求解对比两类结果的目标函数值偏差。如果偏差在可接受范围内比如5%以内说明削减方案有效可以放心在更大规模求解中使用。具体操作上你需要一个两阶段随机规划模型如投资决策运行模拟场景不确定性体现在风电光伏出力上。用削减后场景求解得到的投资方案再拿到全场景下做运行模拟回代检验计算期望运行成本。这个回代过程基本涵盖了随机优化中最关键的一步——场景削减与模型求解的闭环验证。6. 实操中的关键坑与调参经验6.1 归一化处理风电和光伏的量纲差异是隐形杀手这是我最想强调的一点。风电出力单位是MW数值可能在0到1500之间光伏出力如果做归一化处理数值在0到1之间。如果你把这两个量直接拼接到一个48维向量里计算欧氏距离距离计算会被风电的数值尺度主导光伏的差异在距离中几乎不起作用。最终削减出来的场景风电部分很精确光伏部分则非常粗糙。正确做法是先归一化。把风电和光伏都标幺到各自的装机容量下即除以额定功率然后再拼接计算距离。这样每个维度的量纲一致距离度量才能均衡反映两类电源的差异。削减完成后如果需要还原到实际功率再乘以各自的装机容量即可。% 归一化后再削减 windNormalized windPower / windCapacity; pvNormalized pvPower / pvCapacity; scenariosNormalized [windNormalized, pvNormalized]; % 削减... % 还原 windReduced reducedScenarios(:, 1:24) * windCapacity; pvReduced reducedScenarios(:, 25:48) * pvCapacity;6.2 距离度量选择L1还是L2默认的欧氏距离L2对异常场景更敏感会让削减结果更偏向于方差较大的区域曼哈顿距离L1则更稳健在一些风光场景数据中表现更好。我建议两种都跑一遍对比削减后的累积距离误差。如果数据中有明显的离群场景L1通常更稳。另外有一种进阶做法是在距离计算中引入时段权重。比如在调度问题中峰时段的出力偏差对系统运行的影响远大于谷时段可以给峰时段赋予更高权重。具体做法是把原始场景向量乘上一个权重向量再计算欧氏距离。这个技巧在工程实践中很实用能让削减结果更贴合约束的关键区域。6.3 目标场景数M怎么定M的取值没有绝对标准但有几个参考维度随机规划求解能力如果用的是Cplex/Gurobi求解混合整数规划M在10-20个时求解速度通常可接受如果模型简单如线性规划M可以放到50个。分布还原精度需求用5.1节的误差曲线观察误差下降的拐点取拐点处的M值。对比文献常用值多数期刊论文中M取10、20、50三档做敏感性分析。我个人的习惯是先取M20作为基准然后做M10和M50的敏感性分析把三条结果曲线放在同一张图里审稿人或评审专家看到这种分析会觉得你考虑得很周全。6.4 运行效率优化技巧当初始场景数达到5000以上时每次循环都重算距离矩阵会变得有些吃力。两个优化方向第一用向量化操作替代循环。MATLAB的pdist本身是C语言实现的效率已经很高关键是尽量避免在循环里做额外的矩阵复制。第二可以在每轮削减前先检查最近邻矩阵是否可以增量更新只有当被删除场景的邻居关系受影响时才重算对应行。不过这个优化实现复杂度较高一般场景规模下收益不大。如果场景数量过万建议先用K-means粗削减到500个再用概率距离法精确削减到目标数。这种两级削减策略兼顾速度和精度是工程中很实用的做法。7. 案例把风光场景削减用到微电网容量规划中7.1 问题设定我做一个简化版的案例演示。假设一个独立微电网含风力发电、光伏发电和储能系统需要确定风电、光伏和储能的配置容量目标是在满足负荷需求的前提下年综合成本最小。风速和光照的不确定性用场景描述。初始生成2000个风光联合场景用本文方法削减到20个典型场景。7.2 求解结果对比用削减后的20个场景求解容量配置再将配置结果放回2000个全场景中做运行模拟验证实际期望成本。结果是这样的数据仅供演示不同算例会浮动配置方案风电容量/MW光伏容量/MW储能容量/MWh全场景回代成本/万元全场景求解2000场景15.28.640.51236.8削减后求解20场景15.18.941.01241.3K-means削减后求解16.87.938.21295.7从结果看概率距离法削减后求解的决策变量与全场景求解结果非常接近回代成本偏差仅0.36%而K-means削减后求解的偏差达到4.8%。这个案例很直观地说明了削减方法选择对优化决策的真实影响。7.3 从案例中得出的经验第一场景削减的价值不仅在于压缩计算量更在于保留概率分布的关键特征。选对方法20个场景能扛住2000个场景的统计精度选错方法50个场景也照样失真。第二评估削减质量不能只看削减前后的概率距离必须结合下游具体问题如容量配置结果做闭环检验。有时候概率距离根小但关键尾部分布被削没了导致规划方案在极端场景下失去鲁棒性。第三风光联合削减不是简单的分别削减再把结果拼起来。联合削减能保留风-光之间的时序相关性和互补性这对含高比例新能源的微电网规划至关重要。如果你分别对风电和光伏做削减得到的典型场景组合可能根本不存在于原始联合分布中。8. MATLAB实现中的代码细节与扩展思路8.1 函数化封装的设计心得我在写scenarioReduction函数时刻意把输入输出设计成纯矩阵操作不依赖全局变量。这样好处很明显可以直接放进parfor循环做批量实验比如同时测试不同的M值、不同的距离度量、不同的初始场景集互不干扰。推荐你把函数文件保存为scenarioReduction.m放在工作目录下然后写一个runExperiments.m主脚本集中管理参数设置、调用和结果输出。我所有的MATLAB项目都遵循这个习惯数据生成、核心算法、实验分析三个层次分离改一个模块不影响其他模块。8.2 可视化场景质量的小技巧除了简单的曲线图我还会画一种削减效果热力图横轴是24小时时段纵轴是不同的削减后场景颜色表示出力水平。这种图可以快速看出削减后的场景是否覆盖了全时段的各种出力水平。另外把削减后场景的概率权重画成一个饼图或条形图能直观看到概率分布是否过于集中——如果出现某个场景权重超过0.3的情况说明削减后代表性可能不足建议增加场景数。8.3 扩展到多变量和时序耦合场景上面的例子只有风、光两个变量实际中可能还需要考虑负荷的不确定性、电价的不确定性等。实现上只需把更多变量的时序数据拼接到场景矩阵中维度变为[T × 变量数]。归一化时注意每个变量单独处理避免某一变量主导距离计算。对于需要考虑季节差异的场景集比如四季各有不同特性可以分季节分别生成和削减再合并成包含季节标签的完整场景集。在随机规划中这种做法还能同时表达季节概率和场景条件概率的层次结构。8.4 与Matlab工具箱的配合如果你的MATLAB安装了Optimization Toolbox可以把削减后的场景直接接入随机规划求解框架。场景概率向量对应随机规划中的场景权重削减后的场景矩阵作为不确定参数输入到约束中。另外Statistics and Machine Learning Toolbox中的pdist函数是计算距离矩阵的基础工具如果你要做K-means对比实验kmeans函数也是现成的。整套代码不需要任何额外的第三方工具箱基础的MATLAB环境加Statistics工具箱就完全够用。这点对很多用学校正版License或者教育版的朋友很友好不用担心依赖问题。9. 写在最后场景削减的实际项目经验小结概率距离快速削减法是我在多个新能源规划项目中验证过的最可靠场景削减方案。它的优势不是某个特定算例的巧合而是有严格的最优运输理论支撑——每一步削减都保证了概率距离损失最小整个削减过程可解释、可验证、可复现。分享几条我在实际项目中总结的经验初始场景质量决定削减效果上限。生成场景时务必保证概率分布拟合的准确性时序相关性处理不到位削减再精确也没用。归一化是所有场景削减的第一步不要跳过去。不同量纲的变量拼在一起距离计算就失真了。削减后的验证必须闭环。算一下削减前后场景的均值、方差、相关系数等统计指标不够最好能代入你的实际优化模型看决策差异。MATLAB代码注意随机种子设置。场景生成有随机性不固定rng种子前后两次结果可能明显不同论文复现性会受质疑。如果读者把场景数从10加到50计算时间翻倍但误差改善甚微说明关键信息已经被10个场景捕捉到了没必要继续增加。最后再分享一个小技巧在削减过程中可以记录每一步的单步削减误差即每次被删除场景的p_i × d_min(i)画成误差曲线观察衰减趋势。如果发现中间某一步误差突然跳升说明当前场景集中存在一个孤立场景——它与所有其他场景距离都很远删它会造成很大分布损失。这种情况下保留它比多保留两个相似场景更有价值你甚至可以把这种现象作为自适应选择目标场景数的依据。这套方法我已经在多个项目里跑顺了代码量不大原理也不复杂但每次用都能体会到用数学优化替代拍脑袋削减的踏实感。如果你正在为风光场景数量发愁不妨把文中的代码跑一遍应该能帮你少走不少弯路。 基于Kantorovich距离的风光场景削减方法及MATLAB实现
1. 为什么风光场景生成非得做削减这一刀做新能源电力系统规划或调度优化的朋友对场景法应该都不陌生。风电和光伏出力受天气影响具有强烈的随机性和波动性我们没办法用一个确定的曲线去描述明天风多大、太阳多强只能通过历史数据或概率模型生成大量可能出现的出力场景。理论上场景数量越多对不确定性刻画越精确但随之而来的问题是计算量呈指数级膨胀。举个例子一台风机的出力场景如果生成1000个那么一个含10台风机、5个光伏电站的系统组合起来的场景总数就是天文数字。调度模型里每个场景都对应一组约束条件求解器面对成千上万的场景约束轻则求解时间从分钟级变成小时级重则直接内存爆炸无法求解。所以业界普遍的做法是先生成大量场景再用某种方法削减到少量有代表性的场景用这几十个场景去逼近原始上千场景的统计特性。这个思路听起来简单但削减方法选得不好结果会非常难看。我见过有人直接用K-means聚类做场景削减聚类中心倒是选出来了但聚类算法本质上是把场景当成点来划分它压根不考虑场景之间的转移概率和时序相关性削减后的典型场景集在概率分布上往往失真后续优化结果偏差很大。而概率距离快速削减法也叫基于Kantorovich距离的场景削减法简称SD场景削减是目前学术界和工业界都认可的主流方法。它的核心逻辑不是找聚类中心而是通过计算场景之间的概率距离矩阵逐步剔除对整体概率分布影响最小的场景并把被剔除场景的概率累加到距离它最近的保留场景上。每一步削减都保证原始场景集与削减后场景集的概率距离增量最小因此能在给定场景数条件下最大程度保留原始分布特性。MATLAB实现这套方法代码量不大但有几个坑非常隐蔽。本文我用一个完整的风电光伏联合场景生成与削减案例把原理、代码、调参经验和验证方法一次讲透。适合正在做电力系统随机优化、微电网容量配置、风电光伏消纳分析的研究生和工程师参考。2. 风电光伏联合场景生成先弄懂我们削减的对象是什么2.1 风力和光伏出力的数学模型怎么建要削减场景第一步是得有场景。场景生成的方法很多常见的有基于历史数据的直接采样、基于ARIMA时间序列模型、基于马尔可夫链、基于Copula函数、基于概率密度函数的蒙特卡洛采样等。不同的方法适合不同的数据条件。在实际项目中我通常采用概率分布模型时序相关性还原的方式。风电出力近似服从Weibull分布光伏出力近似服从Beta分布但直接按这两个分布独立采样出来的场景毫无时序逻辑——风不会上一秒满发下一秒突然归零光伏也不会中午出力呈白噪声式跳动。所以完整的场景生成链路是对风速历史数据拟合Weibull分布通过风速-功率转换曲线得到风电出力序列。对光照强度历史数据拟合Beta分布通过光照-功率转换模型得到光伏出力序列。引入自相关系数和互相关性让生成的场景在时间维度上平滑、在风-光之间具有合理的相关性。用蒙特卡洛采样生成N个初始场景每个场景包含24个时段的双维出力风电光伏。需要说明的是场景削减算法本身并不关心场景是怎么生成的它只负责在给定初始场景集和概率权重通常假设等概率的前提下找到最优的保留子集。所以哪怕你用别的方法生成场景只要格式是场景数×时段数×变量数后面的削减代码可以直接复用。2.2 初始场景规模如何确定初始场景规模N的选取直接影响生成精度和削减效果。太少了削减后的典型场景缺乏代表性太多了计算概率距离矩阵时的时间和内存开销会显著增加。我的经验法则是单变量场景至少1000个起步风光联合场景建议2000到5000个。1000个以上场景经过削减到10-20个统计特性保留效果较好如果初始只有200个场景削减到10个之后误差会明显偏大。这里贴一段生成初始风电场景的MATLAB核心代码采用Weibull分布采样并加入时序自相关处理% 参数设置 nScenarios 2000; % 初始场景数 nPeriods 24; % 调度时段数小时 k 2.1; % Weibull形状参数根据历史风速拟合 c 8.5; % Weibull尺度参数m/s vRated 12; % 额定风速 vCutIn 3; % 切入风速 vCutOut 25; % 切出风速 Prated 1.5; % 额定功率 MW % 生成独立Weibull风速样本 windSpeed wblrnd(c, k, nScenarios, nPeriods); % 加入时序平滑简单一阶自回归让相邻时段风速平滑过渡 rho 0.85; % 自相关系数 for t 2:nPeriods windSpeed(:, t) rho * windSpeed(:, t-1) ... (1-rho) * windSpeed(:, t) randn(nScenarios, 1) * 0.3; windSpeed(windSpeed 0) 0; end % 风速→功率转换 windPower zeros(nScenarios, nPeriods); for s 1:nScenarios for t 1:nPeriods v windSpeed(s, t); if v vCutIn || v vCutOut windPower(s, t) 0; elseif v vRated windPower(s, t) Prated; else % 线性区间的近似功率曲线 windPower(s, t) Prated * (v - vCutIn) / (vRated - vCutIn); end end end光伏场景的生成逻辑类似区别在于Beta分布拟合的是光照强度的归一化值并且时序上要体现白天出力、晚上归零的强周期特征。联合场景就是把两个矩阵横向拼接每个场景行向量维度变成4824时段风电24时段光伏。3. 概率距离快速削减法原理每一步剔除都让分布损失最小3.1 Kantorovich距离到底是什么很多教程一上来就堆公式把人看懵。我用一句话概括概率距离削减法的本质它在反复问一个问题——如果把某个场景删掉把它的概率转移给离它最近的另一个场景整体概率分布变化有多大如果这个变化是所有候补删除项中最小的那就删它。衡量变化有多大的标尺就是Kantorovich距离也叫Wasserstein距离。两个场景集合之间的Kantorovich距离可以理解为把一个概率分布搬到另一个概率分布所需的最小代价这里的代价就是场景差异乘以转移概率。对于离散场景集Kantorovich距离的计算公式可以转化为一个最优运输问题。但在场景削减这个特定情境下我们不需要求解完整的最优运输问题因为削减过程是逐步进行的每一步从一个场景集中删除一个场景并把其概率累加到最近邻场景上。这种迭代贪婪策略得到的削减结果在每一步都是局部最优的实际应用中已经足够逼近全局最优。3.2 算法的完整流程拆解概率距离快速削减法的具体步骤如下初始化有N个场景每个场景概率权重p_i 1/N也可根据历史频率赋权。计算所有场景两两之间的欧氏距离d(i, j)得到N×N对称距离矩阵。距离定义可以是逐时段绝对误差之和L1范数也可以是欧氏距离L2范数常用的是L2范数。对于每个场景i找到与其距离最近的场景jj ≠ i记录距离d_min(i) min d(i, j)。计算每个场景的可删除代价cost(i) p_i × d_min(i)。这个值的物理含义是删除场景i后概率损失p_i乘以它不得不转移给最近邻居的距离代价d_min(i)。选出cost最小的场景k把它删除并将其概率p_k累加到最近邻居j上p_j p_k p_j。更新距离矩阵删除第k行第k列重复步骤3-5直到保留的场景数量达到预设值M。用伪代码表达就是while nScenarios nTarget % 对每个场景找最近邻居距离 minDist zeros(nScenarios, 1); minIdx zeros(nScenarios, 1); for i 1:nScenarios dist_i distMatrix(i, :); dist_i(i) inf; % 排除自身 [minDist(i), minIdx(i)] min(dist_i); end % 计算删除代价 deleteCost probWeights .* minDist; % 找到代价最小的场景 [~, delIdx] min(deleteCost); % 概率转移给最近邻居 probWeights(minIdx(delIdx)) probWeights(minIdx(delIdx)) probWeights(delIdx); % 删除场景和对应行列 distMatrix(delIdx, :) []; distMatrix(:, delIdx) []; probWeights(delIdx) []; nScenarios nScenarios - 1; end注意这里有个细节距离矩阵每删除一个场景就要更新一次因为删除后某些场景的最近邻居可能发生变化必须重新计算最近邻距离。这也是该算法时间复杂度为O(N³)的原因——每轮删除都要重新扫描一遍距离矩阵。不过对于几千个场景来说MATLAB跑起来也就几秒到几十秒完全可以接受。3.3 为什么不是一次删一批有一种快速近似做法是预先计算所有场景的距离矩阵然后一次性选出多个待删除对象再统一做概率转移。这样做速度更快但存在一个问题第一批被判定为可删的场景中可能会有两个场景互为最近邻居。如果同时删除它们各自的概率都要转移到对方身上就出现了概率无处安放的矛盾。更合理的做法是每删除一个场景立即更新距离矩阵并重新计算最近邻保证每一步操作都是基于最新场景集的。虽然计算量更大但结果更可靠。对于超大规模场景上万个可以考虑用KD-Tree等近似最近邻搜索加速但对于电力系统场景削减的典型规模暴力计算完全够用。4. MATLAB完整实现从距离矩阵到削减主循环4.1 场景切削函数封装实际项目中我不会把主脚本堆成一坨而是封装成一个通用函数输入初始场景矩阵和期望保留数输出削减后的场景、对应概率和削减过程中的误差记录。以下是完整的MATLAB函数代码function [reducedScenarios, reducedProb, distHistory] ... scenarioReduction(scenarios, probs, targetCount) % 基于概率距离快速削减法的场景削减函数 % 输入 % scenarios - 初始场景矩阵维度 [N, T]N为场景数T为时段数含多变量拼接 % probs - 场景概率向量维度 [N, 1]缺省时默认等概率 % targetCount - 目标保留场景数 % 输出 % reducedScenarios - 削减后场景矩阵维度 [targetCount, T] % reducedProb - 削减后场景概率向量维度 [targetCount, 1] % distHistory - 每次删减后的Kantorovich距离记录用于观察削减误差 N size(scenarios, 1); if nargin 2 || isempty(probs) probs ones(N, 1) / N; end if nargin 3 error(必须指定目标场景数); end if targetCount N error(目标场景数必须小于初始场景数); end % 计算初始距离矩阵欧氏距离 distMatrix squareform(pdist(scenarios, euclidean)); % 避免对角线为0置为无穷大 distMatrix(1:N1:end) inf; distHistory zeros(N - targetCount, 1); iter 0; while size(distMatrix, 1) targetCount iter iter 1; currentN size(distMatrix, 1); % 找每个场景的最近邻居 [minDist, minIdx] min(distMatrix, [], 2); % 计算每个场景被删除时造成的概率距离增量 deleteCost probs(:) .* minDist(:); % 选择代价最小的场景删除 [~, delIdx] min(deleteCost); % 概率累加到最近邻居 neighborIdx minIdx(delIdx); probs(neighborIdx) probs(neighborIdx) probs(delIdx); % 记录当前的累积削减误差用当前最小代价近似单步损失 distHistory(iter) deleteCost(delIdx); % 从距离矩阵和概率向量中删除该场景 distMatrix(delIdx, :) []; distMatrix(:, delIdx) []; probs(delIdx) []; end % 最后提取保留的场景需要知道原场景索引 % 由于不断删行距离矩阵的行序和场景矩阵的行序在删除后不一致 % 这里需要用逻辑索引追踪原位置简化方式是通过在削减过程中同步删除场景矩阵的行来实现。 % 由于函数输入未返回索引这里采用同步删除方式重建场景矩阵上面代码有个需要注意的地方我在函数开头提到需要追踪原场景索引。实际上最稳妥的做法是在删除距离矩阵行列时同步删除场景矩阵中的对应行否则最后不知道哪些原始场景被保留了。修正后的核心循环如下reducedScenarios scenarios; idxMap (1:N); while size(reducedScenarios, 1) targetCount currentN size(reducedScenarios, 1); distMat squareform(pdist(reducedScenarios, euclidean)); distMat(1:currentN1:end) inf; [minDist, minIdx] min(distMat, [], 2); deleteCost probs(:) .* minDist(:); [~, delIdx] min(deleteCost); neighborIdx minIdx(delIdx); probs(neighborIdx) probs(neighborIdx) probs(delIdx); % 同步删除场景行和概率 reducedScenarios(delIdx, :) []; probs(delIdx) []; end每轮循环重新调用pdist计算距离矩阵比维护一个动态更新的距离矩阵更直观。虽然理论上多次调用pdist会增加计算量但对几千个场景、48维的场景矩阵来说MATLAB性能完全扛得住代码可读性却高了很多。4.2 主脚本生成场景→削减→可视化主脚本的逻辑很清晰生成2000个风光联合场景削减到10个然后画出削减前后对比图。代码如下clear; clc; close all; rng(2024); % 固定随机种子保证结果可复现 % 生成初始场景以2.2节的windPower为基础增加光伏部分 % 此处直接载入预生成的场景矩阵实际使用时替换为自己的场景生成代码 % scenarios: [2000, 48]前24列风电后24列光伏 load(initialScenarios.mat, scenarios); N size(scenarios, 1); probs ones(N, 1) / N; % 设置削减目标数 M 10; % 调用削减函数 [redScen, redProb, histErr] scenarioReduction(scenarios, probs, M); % 输出结果 fprintf(初始场景数: %d削减后场景数: %d\n, N, M); fprintf(削减后概率分布: \n); disp(redProb(:)); % 可视化削减前后某几个场景的对比 figure; subplot(2,1,1); plot(1:24, scenarios(1:50, 1:24), Color, [0.7 0.7 0.7]); hold on; plot(1:24, redScen(:, 1:24), LineWidth, 1.5); xlabel(时段/h); ylabel(风电出力/MW); title(风电场景削减前后对比); legend(初始场景, 削减后典型场景); grid on; subplot(2,1,2); plot(1:24, scenarios(1:50, 25:48), Color, [0.9 0.7 0.7]); hold on; plot(1:24, redScen(:, 25:48), LineWidth, 1.5); xlabel(时段/h); ylabel(光伏出力/MW); title(光伏场景削减前后对比); legend(初始场景, 削减后典型场景); grid on;运行后削减后的10个典型场景会在2000个原始场景的云图中清晰勾勒出主要分布模式。我截取一张典型输出效果描述风电场景削减后的曲线基本覆盖了从低出力到高出力的各个典型水平光伏场景则呈现出明显的午间高峰形态概率最大的场景是整个分布中最普通的那种——白天光伏高发、晚间风电接力。5. 削减效果怎么评估三种Kantorovich距离验证方法5.1 削减前后概率距离的直接计算最直观的评估方法是计算削减前场景集与削减后场景集之间的Kantorovich距离。这个距离越小说明削减后的场景集对原始分布的还原度越高。MATLAB中可以通过调用Wasserstein距离的相关实现或直接按定义用线性规划求解最优运输问题。对于中等规模问题用transport函数或者CVX求解即可。但实际中我们更关心的是随着目标场景数M从1增加到50削减误差怎样变化把这个变化曲线画出来就能直观看到用多少个场景能在精度和计算量之间取得平衡。我的经验是对于风力发电场景M10时距离通常能降到初始距离的5%以内继续增加M边际收益开始大幅递减。M_values 1:10:50; errorValues zeros(size(M_values)); for i 1:length(M_values) [~, ~, ~] deal([]); [redScen, redProb, histErr] scenarioReduction(scenarios, probs, M_values(i)); errorValues(i) sum(histErr); % 累计削减误差 end figure; plot(M_values, errorValues, o-, LineWidth, 1.5); xlabel(保留场景数); ylabel(累计削减误差); title(削减误差随目标场景数的变化); grid on;5.2 与K-means和K-medoids聚类的对比很多场景削减文献会对比概率距离法与聚类法的差异。我建议你在自己的算例上也做一次对比因为在某些数据结构下K-medoids和概率距离法的结果很接近但概率距离法在概率分配上天然更合理。K-means聚类削减的问题在于聚类得到的簇中心是一个虚拟场景它可能不代表任何真实存在的出力模式而且它无法保留场景的出现概率概念。K-medoids虽然选取的是真实场景作为簇中心但每个簇内场景被视为等权重忽略了簇内场景的密度信息。概率距离法则通过概率转移机制让概率密度高的区域保留更多场景概率密度低的区域被合并。直观理解如果某个区域有400个相似场景削减后这个区域可能保留4个场景每个概率约0.1如果另一区域只有50个场景削减后可能只保留1个场景。这个特性对后续优化计算特别友好因为概率权重可以被直接用作随机规划中的场景概率。5.3 进一步检验用削减后场景求解优化问题对比这是最实际也最有说服力的验证方式。假设你正在做微电网容量配置先用2000个场景求解一版结果再用削减后的10个场景求解对比两类结果的目标函数值偏差。如果偏差在可接受范围内比如5%以内说明削减方案有效可以放心在更大规模求解中使用。具体操作上你需要一个两阶段随机规划模型如投资决策运行模拟场景不确定性体现在风电光伏出力上。用削减后场景求解得到的投资方案再拿到全场景下做运行模拟回代检验计算期望运行成本。这个回代过程基本涵盖了随机优化中最关键的一步——场景削减与模型求解的闭环验证。6. 实操中的关键坑与调参经验6.1 归一化处理风电和光伏的量纲差异是隐形杀手这是我最想强调的一点。风电出力单位是MW数值可能在0到1500之间光伏出力如果做归一化处理数值在0到1之间。如果你把这两个量直接拼接到一个48维向量里计算欧氏距离距离计算会被风电的数值尺度主导光伏的差异在距离中几乎不起作用。最终削减出来的场景风电部分很精确光伏部分则非常粗糙。正确做法是先归一化。把风电和光伏都标幺到各自的装机容量下即除以额定功率然后再拼接计算距离。这样每个维度的量纲一致距离度量才能均衡反映两类电源的差异。削减完成后如果需要还原到实际功率再乘以各自的装机容量即可。% 归一化后再削减 windNormalized windPower / windCapacity; pvNormalized pvPower / pvCapacity; scenariosNormalized [windNormalized, pvNormalized]; % 削减... % 还原 windReduced reducedScenarios(:, 1:24) * windCapacity; pvReduced reducedScenarios(:, 25:48) * pvCapacity;6.2 距离度量选择L1还是L2默认的欧氏距离L2对异常场景更敏感会让削减结果更偏向于方差较大的区域曼哈顿距离L1则更稳健在一些风光场景数据中表现更好。我建议两种都跑一遍对比削减后的累积距离误差。如果数据中有明显的离群场景L1通常更稳。另外有一种进阶做法是在距离计算中引入时段权重。比如在调度问题中峰时段的出力偏差对系统运行的影响远大于谷时段可以给峰时段赋予更高权重。具体做法是把原始场景向量乘上一个权重向量再计算欧氏距离。这个技巧在工程实践中很实用能让削减结果更贴合约束的关键区域。6.3 目标场景数M怎么定M的取值没有绝对标准但有几个参考维度随机规划求解能力如果用的是Cplex/Gurobi求解混合整数规划M在10-20个时求解速度通常可接受如果模型简单如线性规划M可以放到50个。分布还原精度需求用5.1节的误差曲线观察误差下降的拐点取拐点处的M值。对比文献常用值多数期刊论文中M取10、20、50三档做敏感性分析。我个人的习惯是先取M20作为基准然后做M10和M50的敏感性分析把三条结果曲线放在同一张图里审稿人或评审专家看到这种分析会觉得你考虑得很周全。6.4 运行效率优化技巧当初始场景数达到5000以上时每次循环都重算距离矩阵会变得有些吃力。两个优化方向第一用向量化操作替代循环。MATLAB的pdist本身是C语言实现的效率已经很高关键是尽量避免在循环里做额外的矩阵复制。第二可以在每轮削减前先检查最近邻矩阵是否可以增量更新只有当被删除场景的邻居关系受影响时才重算对应行。不过这个优化实现复杂度较高一般场景规模下收益不大。如果场景数量过万建议先用K-means粗削减到500个再用概率距离法精确削减到目标数。这种两级削减策略兼顾速度和精度是工程中很实用的做法。7. 案例把风光场景削减用到微电网容量规划中7.1 问题设定我做一个简化版的案例演示。假设一个独立微电网含风力发电、光伏发电和储能系统需要确定风电、光伏和储能的配置容量目标是在满足负荷需求的前提下年综合成本最小。风速和光照的不确定性用场景描述。初始生成2000个风光联合场景用本文方法削减到20个典型场景。7.2 求解结果对比用削减后的20个场景求解容量配置再将配置结果放回2000个全场景中做运行模拟验证实际期望成本。结果是这样的数据仅供演示不同算例会浮动配置方案风电容量/MW光伏容量/MW储能容量/MWh全场景回代成本/万元全场景求解2000场景15.28.640.51236.8削减后求解20场景15.18.941.01241.3K-means削减后求解16.87.938.21295.7从结果看概率距离法削减后求解的决策变量与全场景求解结果非常接近回代成本偏差仅0.36%而K-means削减后求解的偏差达到4.8%。这个案例很直观地说明了削减方法选择对优化决策的真实影响。7.3 从案例中得出的经验第一场景削减的价值不仅在于压缩计算量更在于保留概率分布的关键特征。选对方法20个场景能扛住2000个场景的统计精度选错方法50个场景也照样失真。第二评估削减质量不能只看削减前后的概率距离必须结合下游具体问题如容量配置结果做闭环检验。有时候概率距离根小但关键尾部分布被削没了导致规划方案在极端场景下失去鲁棒性。第三风光联合削减不是简单的分别削减再把结果拼起来。联合削减能保留风-光之间的时序相关性和互补性这对含高比例新能源的微电网规划至关重要。如果你分别对风电和光伏做削减得到的典型场景组合可能根本不存在于原始联合分布中。8. MATLAB实现中的代码细节与扩展思路8.1 函数化封装的设计心得我在写scenarioReduction函数时刻意把输入输出设计成纯矩阵操作不依赖全局变量。这样好处很明显可以直接放进parfor循环做批量实验比如同时测试不同的M值、不同的距离度量、不同的初始场景集互不干扰。推荐你把函数文件保存为scenarioReduction.m放在工作目录下然后写一个runExperiments.m主脚本集中管理参数设置、调用和结果输出。我所有的MATLAB项目都遵循这个习惯数据生成、核心算法、实验分析三个层次分离改一个模块不影响其他模块。8.2 可视化场景质量的小技巧除了简单的曲线图我还会画一种削减效果热力图横轴是24小时时段纵轴是不同的削减后场景颜色表示出力水平。这种图可以快速看出削减后的场景是否覆盖了全时段的各种出力水平。另外把削减后场景的概率权重画成一个饼图或条形图能直观看到概率分布是否过于集中——如果出现某个场景权重超过0.3的情况说明削减后代表性可能不足建议增加场景数。8.3 扩展到多变量和时序耦合场景上面的例子只有风、光两个变量实际中可能还需要考虑负荷的不确定性、电价的不确定性等。实现上只需把更多变量的时序数据拼接到场景矩阵中维度变为[T × 变量数]。归一化时注意每个变量单独处理避免某一变量主导距离计算。对于需要考虑季节差异的场景集比如四季各有不同特性可以分季节分别生成和削减再合并成包含季节标签的完整场景集。在随机规划中这种做法还能同时表达季节概率和场景条件概率的层次结构。8.4 与Matlab工具箱的配合如果你的MATLAB安装了Optimization Toolbox可以把削减后的场景直接接入随机规划求解框架。场景概率向量对应随机规划中的场景权重削减后的场景矩阵作为不确定参数输入到约束中。另外Statistics and Machine Learning Toolbox中的pdist函数是计算距离矩阵的基础工具如果你要做K-means对比实验kmeans函数也是现成的。整套代码不需要任何额外的第三方工具箱基础的MATLAB环境加Statistics工具箱就完全够用。这点对很多用学校正版License或者教育版的朋友很友好不用担心依赖问题。9. 写在最后场景削减的实际项目经验小结概率距离快速削减法是我在多个新能源规划项目中验证过的最可靠场景削减方案。它的优势不是某个特定算例的巧合而是有严格的最优运输理论支撑——每一步削减都保证了概率距离损失最小整个削减过程可解释、可验证、可复现。分享几条我在实际项目中总结的经验初始场景质量决定削减效果上限。生成场景时务必保证概率分布拟合的准确性时序相关性处理不到位削减再精确也没用。归一化是所有场景削减的第一步不要跳过去。不同量纲的变量拼在一起距离计算就失真了。削减后的验证必须闭环。算一下削减前后场景的均值、方差、相关系数等统计指标不够最好能代入你的实际优化模型看决策差异。MATLAB代码注意随机种子设置。场景生成有随机性不固定rng种子前后两次结果可能明显不同论文复现性会受质疑。如果读者把场景数从10加到50计算时间翻倍但误差改善甚微说明关键信息已经被10个场景捕捉到了没必要继续增加。最后再分享一个小技巧在削减过程中可以记录每一步的单步削减误差即每次被删除场景的p_i × d_min(i)画成误差曲线观察衰减趋势。如果发现中间某一步误差突然跳升说明当前场景集中存在一个孤立场景——它与所有其他场景距离都很远删它会造成很大分布损失。这种情况下保留它比多保留两个相似场景更有价值你甚至可以把这种现象作为自适应选择目标场景数的依据。这套方法我已经在多个项目里跑顺了代码量不大原理也不复杂但每次用都能体会到用数学优化替代拍脑袋削减的踏实感。如果你正在为风光场景数量发愁不妨把文中的代码跑一遍应该能帮你少走不少弯路。