做新能源电力系统仿真的朋友十有八九都遇到过这类问题风光出力数据拿来直接用模型规模稍微大一点就卡得不行不处理数据误差又大得没法交代。场景生成与场景削减就是专门解决这个矛盾的。我手头这套MATLAB代码基于概率距离快速削减法能把初始生成的大量风光场景精简成少量有代表性的典型场景同时还保留每个场景出现的大致概率。它适合做微电网规划、含可再生能源的电力系统经济调度、可靠性评估这类工作是新能源领域很实用的一个数据处理环节。这套方法的核心逻辑并不复杂先用概率分布模型生成成百上千条可能的风光出力曲线再用概率距离作为衡量标准把相似的曲线合并、叠概率最终得到几十条既能覆盖典型情况、又不至于让优化模型算爆的场景。接下来我把整个方案的思路、生成过程、削减算法和MATLAB实现细节完整拆开尤其是那些代码里看不出来的坑一次讲清楚。1. 项目背景与核心设计思路1.1 为什么不能直接使用单条时序曲线在风电、光伏并网分析里出力本身是随机过程受风速、光照、温度、云量等因素影响今天中午的光伏出力曲线和明天中午的几乎不可能完全一样。如果我们只用一条“典型日曲线”去做调度或容量配置等于默认了未来只有一种天气情况这对决策来说风险太大。反过来如果保留一整年8760小时的历史数据优化模型里变量太多求解时间会指数级上升严重拖慢工作节奏。所以主流做法分两步走。第一步叫场景生成根据风速、光照的概率分布蒙特卡洛采样出大量可能的出力时序曲线每条曲线描述一种可能出现的情况。第二步叫场景削减在保留概率特性的前提下把几百上千条曲线压缩到几十条让后续计算量降到可接受范围。场景削减不是简单的“删除多余的线”而是要保证削减后的场景集合和原始场景集合在统计意义上足够接近比如均值、方差、峰值分布都要对得上。1.2 概率距离快速削减法的选型逻辑市面上常用场景削减算法有K-means聚类、层次聚类、后向削减等。先说K-means它确实快但有个问题聚类中心是平均值会把场景的尖峰、陡坡这类细节磨平而且聚类结果不保留概率语义每个簇到底代表多大的发生概率需要额外计算。层次聚类也有类似问题场景的代表性更多依赖聚类半径的设置主观性比较强。概率距离快速削减法属于经典的最优场景削减近似算法。它的数学背景是在给定场景数量上限的条件下寻找一个场景子集使得削减前后场景集合的Kantorovich距离最小。精确求解这个问题是一个大规模组合优化问题计算量极大不适合在线计算。概率距离快速削减法用贪心策略逼近最优解每次迭代只合并“距离最近的一对场景”并把被合并场景的概率累加到保留场景上循环往复直到达到目标场景数量。这种方法既保留了概率信息又避免了大规模优化求解实测下来结果非常接近精确最优解。1.3 这个方法能解决什么实际问题举个例子我在做微电网容量配置时初始生成1500个风光场景优化模型一跑就是几个小时。用概率距离削减到30个场景后求解时间缩短到几分钟而配置结果和原始1500个场景跑出来的结果相差不到5%。这说明削减过程没有丢掉关键的不确定性信息。对做调度、做可靠性评估、做风电消纳分析的人来说这套方法能把“算不动”的问题变成“算得动”的问题同时保证决策质量不滑坡。2. 场景生成从概率分布到时序场景2.1 风速与光照的概率分布建模场景生成的前提是给风速和光照建立概率模型。风速在工程上通常用两参数威布尔分布描述概率密度函数是[ f(v) \frac{k}{c}\left(\frac{v}{c}\right)^{k-1} \exp\left(-\left(\frac{v}{c}\right)^k\right) ]其中 ( k ) 是形状参数控制分布形态( c ) 是尺度参数控制整体风速水平。这组参数可以直接根据历史风速数据用极大似然估计求得没有历史数据时也可以参考经验值沿海风场 ( k ) 一般在2.0到2.5之间( c ) 在6到9之间。光照强度通常用Beta分布描述定义在0到1之间通过均值和方差来确定两个形状参数 ( \alpha ) 和 ( \beta )[ \alpha \mu \left( \frac{\mu(1-\mu)}{\sigma^2} - 1 \right) ][ \beta (1-\mu) \left( \frac{\mu(1-\mu)}{\sigma^2} - 1 \right) ]这里的 ( \mu ) 和 ( \sigma^2 ) 是太阳辐照度归一化后的均值和方差。实际项目里如果只有日最大辐照度的历史数据也可以先拟合日最大辐照度的分布再叠加日内变化曲线来生成时序场景。2.2 风电功率与光伏出力转换模型有了风速场景和光照场景还需要把它们换算成电功率场景。风电功率转换是分段函数典型参数包括切入风速 ( v_{cin} )、切出风速 ( v_{cout} )、额定风速 ( v_r )、额定功率 ( P_r )。常用的线性化模型是当 ( v v_{cin} ) 或 ( v v_{cout} ) 时( P_w 0 )当 ( v_{cin} \le v \le v_r ) 时( P_w P_r \cdot \frac{v - v_{cin}}{v_r - v_{cin}} )当 ( v_r \le v \le v_{cout} ) 时( P_w P_r )光伏出力简化模型通常写成[ P_{pv} \eta \cdot S \cdot A ]其中 ( \eta ) 是光电转换效率( S ) 是太阳辐照度( A ) 是光伏阵列总面积。如果考虑温度影响可以再加一个温度修正系数但大多数场景生成场景下简化模型已经够用。2.3 生成带时序相关性的初始场景这里有一个容易忽略的细节随机生成的风速时序如果每个时段独立采样得到的曲线会非常“毛糙”跟真实风速的慢变特性完全不符。为了解决这个问题需要在采样时引入时间相关性。我用的是AR(1)模型。对每个场景先采样第一个时段的风速值然后按下面的递推式生成后续时段[ v(t) \alpha v(t-1) (1-\alpha)\mu \varepsilon(t) ]其中 ( \alpha ) 是自相关系数一般取0.7到0.9( \mu ) 是风速均值( \varepsilon ) 是服从威布尔分布残差的随机项。这样生成的每条风速曲线都有自然的“惯性”不会出现前一小时风速12米/秒、后一小时骤降到2米/秒这种不合理的跳跃。采样方式上我强烈建议用拉丁超立方采样代替纯蒙特卡洛采样。拉丁超立方会把概率空间分层保证每个区间都有样本点避免随机采样常见的“挤成一团”问题。用同样数量的样本拉丁超立方覆盖的分布空间更完整生成出来的场景集质量更高。3. 场景削减概率距离快速削减法算法剖析3.1 场景间距离的定义概率距离快速削减法的第一步是定义“两个场景有多像”。最常用的距离是欧氏距离对两个维度相同的时间序列 ( x_i ) 和 ( x_j )[ d(i,j) \sqrt{\sum_{t1}^{T} \left( x_i(t) - x_j(t) \right)^2 } ]其中 ( T ) 是时段数比如一天24小时。如果场景规模特别大也可以省略开方运算直接用平方欧氏距离做比较结果排序完全一致还能省一点计算时间。这里要注意如果风光联合场景是由风电场景和光伏场景拼接而成的必须对两个部分分别计算距离再求和而不能直接把数值混在一起算否则量纲不一致会导致光伏部分被风电部分“带偏”。3.2 迭代合并流程与概率更新规则概率距离快速削减法的核心循环非常清晰用伪代码表示输入初始场景集 ( S {s_1, s_2, ..., s_N} )对应概率 ( P {p_1, p_2, ..., p_N} ) 目标场景数 ( N_{target} )。计算所有场景对之间的距离矩阵 ( D )。找出距离最小的一对场景 ( (i^, j^) )。将场景 ( j^* ) 的概率累加到场景 ( i^* ) 上( p_{i^} p_{i^} p_{j^*} )。从场景集中删除场景 ( j^* )。更新距离矩阵中与 ( i^* ) 相关的距离值。如果当前场景数量大于 ( N_{target} )回到步骤3否则输出削减后的场景集和概率向量。合并时是保留还是删除某个场景有几种变体。经典做法是直接删除其中一个保留另一个的原始数据概率叠加。另一种做法是生成两个场景的加权平均作为新场景。我测试过两种方案加权平均会让峰值被削弱不利于分析极端场景所以推荐直接删除法。3.3 关键参数初始场景数与目标场景数这两个参数对结果影响很大。初始场景数太少分布覆盖不完整削减结果偏差大初始场景数太多距离矩阵计算时间长内存占用也高。我常用的配置是初始场景1000到2000个目标场景30到50个。1000个初始场景生成30个典型场景基本能保留原始场景集的95%以上统计特征如果只留10个场景虽然计算更快但极端出力场景容易被合并掉可靠性评估会偏乐观。另外要特别注意概率向量的归一化。每次合并后概率总量保持不变但最后最好再统一除以总和避免浮点误差累计导致概率之和不为1。3.4 概率距离削减法与其他方法的场景对比为了直观感受差异我做过一个对比测试同一组1000个场景分别用概率距离快速削减法和K-means聚类削减到20个然后计算削减后场景集与原始场景集的均值误差和方差误差。结果显示概率距离法的均值误差在2%以内方差误差在5%以内K-means的均值误差接近3%方差误差却超过10%。原因还是我之前说的K-means聚类中心作为代表场景会把极端时段的出力拉平导致方差被系统性低估在处理风光这类波动性很强的数据时特别明显。4. MATLAB代码实现与性能优化4.1 整体函数结构与数据流这套MATLAB代码按功能拆成场景生成和场景削减两个独立模块中间用矩阵接口衔接方便替换成自己的模型。代码的大致结构是% 主脚本入口 % 1. 设置参数 params.weibull_k 2.1; % 风速形状参数 params.weibull_c 7.5; % 风速尺度参数 params.rho 0.85; % 风速时序自相关系数 params.P_r 1500; % 风机额定功率 kW params.v_cin 3; % 切入风速 params.v_r 12; % 额定风速 params.v_cout 25; % 切出风速 params.N_scenes 1000; % 初始场景数 params.N_target 30; % 削减后场景数 params.horizon 24; % 时段数 % 2. 生成初始风电场景 scenes_wind generate_wind_scenes(params); % 3. 生成初始光伏场景 scenes_pv generate_pv_scenes(params); % 4. 拼接成风光联合场景 scenes_all [scenes_wind, scenes_pv]; % 5. 等概率初始化 prob ones(params.N_scenes, 1) / params.N_scenes; % 6. 场景削减 [scenes_red, prob_red] reduce_scenes(scenes_all, prob, params.N_target);4.2 场景生成模块核心代码风速场景生成的核心是AR(1)递推加拉丁超立方采样。下面这个函数实现了单个风速场景的时序生成function scene gen_wind_scene_ar1(params) % 生成一个风电出力场景24小时 horizon params.horizon; k params.weibull_k; c params.weibull_c; rho params.rho; % 生成均匀随机数并转换为威布尔分布采样 % 使用威布尔分布逆变换采样 u rand(horizon, 1); % 第一时段从威布尔分布采样 v zeros(horizon, 1); v(1) c * (-log(1 - u(1)))^(1/k); % 后续时段用AR(1)递推 mu c * gamma(1 1/k); for t 2:horizon % 条件期望 噪声项 % 这里简化为新采样值围绕条件期望波动 target_mean rho * v(t-1) (1 - rho) * mu; % 从威布尔分布重新采样并调整到目标均值附近 v_new c * (-log(1 - rand()))^(1/k); % 通过残差项调整相关性 v(t) target_mean 0.3 * (v_new - mu) (1 - 0.3) * (v(t-1) - target_mean); v(t) max(v(t), 0); end % 转换为风电出力 P zeros(horizon, 1); for t 1:horizon if v(t) params.v_cin || v(t) params.v_cout P(t) 0; elseif v(t) params.v_r P(t) params.P_r * (v(t) - params.v_cin) / (params.v_r - params.v_cin); else P(t) params.P_r; end end scene P; end注意真实项目中通常先用拉丁超立方采样生成所有的初始风速值再通过AR模型重排时序这样既能保证分布覆盖又能保留时间相关性。上面这段代码为了便于理解做了一定简化但核心逻辑——分布采样 自相关处理 功率转换——是一致的。光伏场景生成也是类似逻辑先对日最大辐照度做Beta分布采样再用典型日辐照度曲线按比例缩放最后乘上转换效率和面积得到功率时序。4.3 场景削减模块核心代码概率距离快速削减法的核心实现如下function [scenes_red, prob_red] reduce_scenes(scenes, prob, N_target) % scenes: N * T 矩阵每行是一个场景 % prob: N * 1 向量每个场景的概率 % N_target: 目标场景数量 N size(scenes, 1); sc scenes; pr prob(:); while N N_target % 计算距离矩阵使用向量化方法避免双层循环 D pdist2(sc, sc, squaredeuclidean); % 平方欧氏距离 % 对角线设置为无穷大排除场景自身的距离 D(1:N1:end) inf; % 找到最小距离对 [~, idx] min(D(:)); [i, j] ind2sub([N, N], idx); % 合并场景将j的概率累加到i删除j pr(i) pr(i) pr(j); % 删除场景j sc(j, :) []; pr(j) []; N N - 1; end % 概率归一化 pr pr / sum(pr); scenes_red sc; prob_red pr; end这段代码的逻辑和算法原理完全对应每次迭代只找全局最近的场景对概率叠加后删除。由于每次删除后都要重算距离矩阵当初始场景数超过2000时循环次数多、矩阵规模大计算时间会显著上升。如果遇到这种情况可以用局部更新代替全局重算只更新与新场景有关的那一行一列距离其余保持不变。4.4 性能优化技巧我曾经用1500个初始场景跑削减到30最初版本的代码跑了将近4分钟优化后只用了40秒。关键优化点有三个。第一距离计算用向量化。不要写双层for循环去计算场景间距离直接用pdist2函数MATLAB底层对矩阵运算是并行优化的速度能快一个数量级。第二删除场景时用布尔索引批量操作不要用循环逐个删除。上面代码里用sc(j, :) []就是一次到位。第三如果场景数量特别大比如超过5000可以用KD树做最近邻搜索避免每次都计算全量距离矩阵。MATLAB的knnsearch支持KD树对这种场景削减问题能大幅提速。还有一个细节保存削减过程中的中间结果。每次迭代其实对应一个场景数量后面如果要分析“削减到50个场景时效果如何”“削减到30个场景时效果如何”直接从中间结果取即可不用重新跑。5. 常见问题与调参实战5.1 削减后的场景太集中失去了波动性出现这种情况通常有三个原因。一是初始场景数太少分布覆盖不全耦合削减后自然就聚到少数几个模式上。二是场景间距离计算没有加权重光伏部分数值大、风电部分数值小距离被光伏主导导致风电的波动特征被忽略。三是目标场景数设得太少比如从1000个直接砍到5个再好的算法也会丢掉大量信息。我的建议是目标场景数不要低于初始场景数的2%到3%。如果削减后的场景集合不能满足可靠性评估需求可以分两步走先整体削减出常规场景再单独做极端天气场景手动加入削减结果里保证典型场景和极端场景都有覆盖。5.2 削减前后统计特征误差怎么验证拿到削减结果后不能只看曲线好不好看要做量化验证。我常用的验证指标有三个。第一个是概率距离值削减后场景集合和原始场景集合的Kantorovich距离越小越好一般控制在原始场景集合均值的5%以内。第二个是均值与方差误差分别对比每个时段削减前后出力均值、方差的相对误差。第三个是用KS检验对比单时段出力的累积分布函数看削减后的分布是否和原始分布一致。这三个指标都能写进MATLAB代码里自动输出一份验证报告。5.3 风光场景的相关性怎么处理风速和光照本质上不是独立的很多地区存在“风小光强”或“互补”的特性。如果忽略这种相关性生成的联合场景会出现一些概率极低的组合比如大风天光照极强这种实际很少发生的情况导致调度结果偏保守或偏乐观。处理办法是在采样时对风速和光照的随机数施加相关性。常见方式是生成相关系数矩阵做Cholesky分解得到相关的标准正态随机变量再通过逆变换映射到威布尔分布和Beta分布。具体来说先采样一组相关的正态变量再用正态CDF反变换得到均匀分布的随机数最后用均匀随机数反变换采样风速和光照。这个过程在MATLAB里只需要几行代码就能实现但对场景质量的提升非常明显。风速和光照的相关系数可以根据历史数据估算没有数据时可以取-0.2到0.3之间的经验值。5.4 运行内存不足怎么办初始场景2000个每个24小时存成矩阵是 2000×24双精度占内存约1.5MB并不大。但如果把场景维度扩展到逐小时8760小时矩阵就是 2000×8760瞬间涨到140MB距离矩阵更是达到 2000×2000×8字节等于32MB加起来还能接受但再大就危险了。建议用单精度存储场景矩阵scenes single(scenes)内存直接减半。距离矩阵也可以分块计算避免一次性生成完整矩阵。如果场景数上万建议改用采样加聚类的方式先粗减到1000个再用概率距离快速削减法精减到目标数量两阶段削减的速度远快过一阶段硬算。5.5 调参经验速查参数推荐范围影响初始场景数500-2000太少覆盖不足太多计算慢目标场景数20-50太少损失波动性太多削减意义小风速自相关系数0.7-0.9控制时序平滑程度风光相关系数-0.2-0.3控制联合场景的合理性场景距离权重按额定功率归一化防止量纲不一致导致距离偏差最后一列是重点场景距离计算前最好先把风电和光伏出力都归一化到各自的额定功率再拼接场景。这样距离度量的物理意义变成了“归一化出力曲线的相似度”风电和光伏在距离计算中权重均衡削减结果更合理。结尾一点个人体会这套代码我在多个项目里反复用过最深的体会是场景削减不是“数据处理完事”它直接决定后续优化模型能不能算出靠谱的结果。概率距离快速削减法胜在原理清晰、可解释性强削减后的场景保留了原始概率语义这对调度方案的风险分析特别有用。后来我还试过用Wasserstein距离替代Kantorovich距离做削减效果接近数值上更稳定如果你有精力可以往这个方向扩展。再分享一个小技巧削减后的场景建议顺便输出每个场景的累积概率分布图和几个典型场景的时序可视化贴在报告里对方哪怕不懂算法一眼也能看出“这30条曲线确实覆盖了各种天气情况”这一点在项目汇报时价值极大。 
长河君 · 2026/10/13 15:12:23 · 阅读 1.2k 风光场景生成与削减：概率距离快速削减法MATLAB实现
做新能源电力系统仿真的朋友十有八九都遇到过这类问题风光出力数据拿来直接用模型规模稍微大一点就卡得不行不处理数据误差又大得没法交代。场景生成与场景削减就是专门解决这个矛盾的。我手头这套MATLAB代码基于概率距离快速削减法能把初始生成的大量风光场景精简成少量有代表性的典型场景同时还保留每个场景出现的大致概率。它适合做微电网规划、含可再生能源的电力系统经济调度、可靠性评估这类工作是新能源领域很实用的一个数据处理环节。这套方法的核心逻辑并不复杂先用概率分布模型生成成百上千条可能的风光出力曲线再用概率距离作为衡量标准把相似的曲线合并、叠概率最终得到几十条既能覆盖典型情况、又不至于让优化模型算爆的场景。接下来我把整个方案的思路、生成过程、削减算法和MATLAB实现细节完整拆开尤其是那些代码里看不出来的坑一次讲清楚。1. 项目背景与核心设计思路1.1 为什么不能直接使用单条时序曲线在风电、光伏并网分析里出力本身是随机过程受风速、光照、温度、云量等因素影响今天中午的光伏出力曲线和明天中午的几乎不可能完全一样。如果我们只用一条“典型日曲线”去做调度或容量配置等于默认了未来只有一种天气情况这对决策来说风险太大。反过来如果保留一整年8760小时的历史数据优化模型里变量太多求解时间会指数级上升严重拖慢工作节奏。所以主流做法分两步走。第一步叫场景生成根据风速、光照的概率分布蒙特卡洛采样出大量可能的出力时序曲线每条曲线描述一种可能出现的情况。第二步叫场景削减在保留概率特性的前提下把几百上千条曲线压缩到几十条让后续计算量降到可接受范围。场景削减不是简单的“删除多余的线”而是要保证削减后的场景集合和原始场景集合在统计意义上足够接近比如均值、方差、峰值分布都要对得上。1.2 概率距离快速削减法的选型逻辑市面上常用场景削减算法有K-means聚类、层次聚类、后向削减等。先说K-means它确实快但有个问题聚类中心是平均值会把场景的尖峰、陡坡这类细节磨平而且聚类结果不保留概率语义每个簇到底代表多大的发生概率需要额外计算。层次聚类也有类似问题场景的代表性更多依赖聚类半径的设置主观性比较强。概率距离快速削减法属于经典的最优场景削减近似算法。它的数学背景是在给定场景数量上限的条件下寻找一个场景子集使得削减前后场景集合的Kantorovich距离最小。精确求解这个问题是一个大规模组合优化问题计算量极大不适合在线计算。概率距离快速削减法用贪心策略逼近最优解每次迭代只合并“距离最近的一对场景”并把被合并场景的概率累加到保留场景上循环往复直到达到目标场景数量。这种方法既保留了概率信息又避免了大规模优化求解实测下来结果非常接近精确最优解。1.3 这个方法能解决什么实际问题举个例子我在做微电网容量配置时初始生成1500个风光场景优化模型一跑就是几个小时。用概率距离削减到30个场景后求解时间缩短到几分钟而配置结果和原始1500个场景跑出来的结果相差不到5%。这说明削减过程没有丢掉关键的不确定性信息。对做调度、做可靠性评估、做风电消纳分析的人来说这套方法能把“算不动”的问题变成“算得动”的问题同时保证决策质量不滑坡。2. 场景生成从概率分布到时序场景2.1 风速与光照的概率分布建模场景生成的前提是给风速和光照建立概率模型。风速在工程上通常用两参数威布尔分布描述概率密度函数是[ f(v) \frac{k}{c}\left(\frac{v}{c}\right)^{k-1} \exp\left(-\left(\frac{v}{c}\right)^k\right) ]其中 ( k ) 是形状参数控制分布形态( c ) 是尺度参数控制整体风速水平。这组参数可以直接根据历史风速数据用极大似然估计求得没有历史数据时也可以参考经验值沿海风场 ( k ) 一般在2.0到2.5之间( c ) 在6到9之间。光照强度通常用Beta分布描述定义在0到1之间通过均值和方差来确定两个形状参数 ( \alpha ) 和 ( \beta )[ \alpha \mu \left( \frac{\mu(1-\mu)}{\sigma^2} - 1 \right) ][ \beta (1-\mu) \left( \frac{\mu(1-\mu)}{\sigma^2} - 1 \right) ]这里的 ( \mu ) 和 ( \sigma^2 ) 是太阳辐照度归一化后的均值和方差。实际项目里如果只有日最大辐照度的历史数据也可以先拟合日最大辐照度的分布再叠加日内变化曲线来生成时序场景。2.2 风电功率与光伏出力转换模型有了风速场景和光照场景还需要把它们换算成电功率场景。风电功率转换是分段函数典型参数包括切入风速 ( v_{cin} )、切出风速 ( v_{cout} )、额定风速 ( v_r )、额定功率 ( P_r )。常用的线性化模型是当 ( v v_{cin} ) 或 ( v v_{cout} ) 时( P_w 0 )当 ( v_{cin} \le v \le v_r ) 时( P_w P_r \cdot \frac{v - v_{cin}}{v_r - v_{cin}} )当 ( v_r \le v \le v_{cout} ) 时( P_w P_r )光伏出力简化模型通常写成[ P_{pv} \eta \cdot S \cdot A ]其中 ( \eta ) 是光电转换效率( S ) 是太阳辐照度( A ) 是光伏阵列总面积。如果考虑温度影响可以再加一个温度修正系数但大多数场景生成场景下简化模型已经够用。2.3 生成带时序相关性的初始场景这里有一个容易忽略的细节随机生成的风速时序如果每个时段独立采样得到的曲线会非常“毛糙”跟真实风速的慢变特性完全不符。为了解决这个问题需要在采样时引入时间相关性。我用的是AR(1)模型。对每个场景先采样第一个时段的风速值然后按下面的递推式生成后续时段[ v(t) \alpha v(t-1) (1-\alpha)\mu \varepsilon(t) ]其中 ( \alpha ) 是自相关系数一般取0.7到0.9( \mu ) 是风速均值( \varepsilon ) 是服从威布尔分布残差的随机项。这样生成的每条风速曲线都有自然的“惯性”不会出现前一小时风速12米/秒、后一小时骤降到2米/秒这种不合理的跳跃。采样方式上我强烈建议用拉丁超立方采样代替纯蒙特卡洛采样。拉丁超立方会把概率空间分层保证每个区间都有样本点避免随机采样常见的“挤成一团”问题。用同样数量的样本拉丁超立方覆盖的分布空间更完整生成出来的场景集质量更高。3. 场景削减概率距离快速削减法算法剖析3.1 场景间距离的定义概率距离快速削减法的第一步是定义“两个场景有多像”。最常用的距离是欧氏距离对两个维度相同的时间序列 ( x_i ) 和 ( x_j )[ d(i,j) \sqrt{\sum_{t1}^{T} \left( x_i(t) - x_j(t) \right)^2 } ]其中 ( T ) 是时段数比如一天24小时。如果场景规模特别大也可以省略开方运算直接用平方欧氏距离做比较结果排序完全一致还能省一点计算时间。这里要注意如果风光联合场景是由风电场景和光伏场景拼接而成的必须对两个部分分别计算距离再求和而不能直接把数值混在一起算否则量纲不一致会导致光伏部分被风电部分“带偏”。3.2 迭代合并流程与概率更新规则概率距离快速削减法的核心循环非常清晰用伪代码表示输入初始场景集 ( S {s_1, s_2, ..., s_N} )对应概率 ( P {p_1, p_2, ..., p_N} ) 目标场景数 ( N_{target} )。计算所有场景对之间的距离矩阵 ( D )。找出距离最小的一对场景 ( (i^, j^) )。将场景 ( j^* ) 的概率累加到场景 ( i^* ) 上( p_{i^} p_{i^} p_{j^*} )。从场景集中删除场景 ( j^* )。更新距离矩阵中与 ( i^* ) 相关的距离值。如果当前场景数量大于 ( N_{target} )回到步骤3否则输出削减后的场景集和概率向量。合并时是保留还是删除某个场景有几种变体。经典做法是直接删除其中一个保留另一个的原始数据概率叠加。另一种做法是生成两个场景的加权平均作为新场景。我测试过两种方案加权平均会让峰值被削弱不利于分析极端场景所以推荐直接删除法。3.3 关键参数初始场景数与目标场景数这两个参数对结果影响很大。初始场景数太少分布覆盖不完整削减结果偏差大初始场景数太多距离矩阵计算时间长内存占用也高。我常用的配置是初始场景1000到2000个目标场景30到50个。1000个初始场景生成30个典型场景基本能保留原始场景集的95%以上统计特征如果只留10个场景虽然计算更快但极端出力场景容易被合并掉可靠性评估会偏乐观。另外要特别注意概率向量的归一化。每次合并后概率总量保持不变但最后最好再统一除以总和避免浮点误差累计导致概率之和不为1。3.4 概率距离削减法与其他方法的场景对比为了直观感受差异我做过一个对比测试同一组1000个场景分别用概率距离快速削减法和K-means聚类削减到20个然后计算削减后场景集与原始场景集的均值误差和方差误差。结果显示概率距离法的均值误差在2%以内方差误差在5%以内K-means的均值误差接近3%方差误差却超过10%。原因还是我之前说的K-means聚类中心作为代表场景会把极端时段的出力拉平导致方差被系统性低估在处理风光这类波动性很强的数据时特别明显。4. MATLAB代码实现与性能优化4.1 整体函数结构与数据流这套MATLAB代码按功能拆成场景生成和场景削减两个独立模块中间用矩阵接口衔接方便替换成自己的模型。代码的大致结构是% 主脚本入口 % 1. 设置参数 params.weibull_k 2.1; % 风速形状参数 params.weibull_c 7.5; % 风速尺度参数 params.rho 0.85; % 风速时序自相关系数 params.P_r 1500; % 风机额定功率 kW params.v_cin 3; % 切入风速 params.v_r 12; % 额定风速 params.v_cout 25; % 切出风速 params.N_scenes 1000; % 初始场景数 params.N_target 30; % 削减后场景数 params.horizon 24; % 时段数 % 2. 生成初始风电场景 scenes_wind generate_wind_scenes(params); % 3. 生成初始光伏场景 scenes_pv generate_pv_scenes(params); % 4. 拼接成风光联合场景 scenes_all [scenes_wind, scenes_pv]; % 5. 等概率初始化 prob ones(params.N_scenes, 1) / params.N_scenes; % 6. 场景削减 [scenes_red, prob_red] reduce_scenes(scenes_all, prob, params.N_target);4.2 场景生成模块核心代码风速场景生成的核心是AR(1)递推加拉丁超立方采样。下面这个函数实现了单个风速场景的时序生成function scene gen_wind_scene_ar1(params) % 生成一个风电出力场景24小时 horizon params.horizon; k params.weibull_k; c params.weibull_c; rho params.rho; % 生成均匀随机数并转换为威布尔分布采样 % 使用威布尔分布逆变换采样 u rand(horizon, 1); % 第一时段从威布尔分布采样 v zeros(horizon, 1); v(1) c * (-log(1 - u(1)))^(1/k); % 后续时段用AR(1)递推 mu c * gamma(1 1/k); for t 2:horizon % 条件期望 噪声项 % 这里简化为新采样值围绕条件期望波动 target_mean rho * v(t-1) (1 - rho) * mu; % 从威布尔分布重新采样并调整到目标均值附近 v_new c * (-log(1 - rand()))^(1/k); % 通过残差项调整相关性 v(t) target_mean 0.3 * (v_new - mu) (1 - 0.3) * (v(t-1) - target_mean); v(t) max(v(t), 0); end % 转换为风电出力 P zeros(horizon, 1); for t 1:horizon if v(t) params.v_cin || v(t) params.v_cout P(t) 0; elseif v(t) params.v_r P(t) params.P_r * (v(t) - params.v_cin) / (params.v_r - params.v_cin); else P(t) params.P_r; end end scene P; end注意真实项目中通常先用拉丁超立方采样生成所有的初始风速值再通过AR模型重排时序这样既能保证分布覆盖又能保留时间相关性。上面这段代码为了便于理解做了一定简化但核心逻辑——分布采样 自相关处理 功率转换——是一致的。光伏场景生成也是类似逻辑先对日最大辐照度做Beta分布采样再用典型日辐照度曲线按比例缩放最后乘上转换效率和面积得到功率时序。4.3 场景削减模块核心代码概率距离快速削减法的核心实现如下function [scenes_red, prob_red] reduce_scenes(scenes, prob, N_target) % scenes: N * T 矩阵每行是一个场景 % prob: N * 1 向量每个场景的概率 % N_target: 目标场景数量 N size(scenes, 1); sc scenes; pr prob(:); while N N_target % 计算距离矩阵使用向量化方法避免双层循环 D pdist2(sc, sc, squaredeuclidean); % 平方欧氏距离 % 对角线设置为无穷大排除场景自身的距离 D(1:N1:end) inf; % 找到最小距离对 [~, idx] min(D(:)); [i, j] ind2sub([N, N], idx); % 合并场景将j的概率累加到i删除j pr(i) pr(i) pr(j); % 删除场景j sc(j, :) []; pr(j) []; N N - 1; end % 概率归一化 pr pr / sum(pr); scenes_red sc; prob_red pr; end这段代码的逻辑和算法原理完全对应每次迭代只找全局最近的场景对概率叠加后删除。由于每次删除后都要重算距离矩阵当初始场景数超过2000时循环次数多、矩阵规模大计算时间会显著上升。如果遇到这种情况可以用局部更新代替全局重算只更新与新场景有关的那一行一列距离其余保持不变。4.4 性能优化技巧我曾经用1500个初始场景跑削减到30最初版本的代码跑了将近4分钟优化后只用了40秒。关键优化点有三个。第一距离计算用向量化。不要写双层for循环去计算场景间距离直接用pdist2函数MATLAB底层对矩阵运算是并行优化的速度能快一个数量级。第二删除场景时用布尔索引批量操作不要用循环逐个删除。上面代码里用sc(j, :) []就是一次到位。第三如果场景数量特别大比如超过5000可以用KD树做最近邻搜索避免每次都计算全量距离矩阵。MATLAB的knnsearch支持KD树对这种场景削减问题能大幅提速。还有一个细节保存削减过程中的中间结果。每次迭代其实对应一个场景数量后面如果要分析“削减到50个场景时效果如何”“削减到30个场景时效果如何”直接从中间结果取即可不用重新跑。5. 常见问题与调参实战5.1 削减后的场景太集中失去了波动性出现这种情况通常有三个原因。一是初始场景数太少分布覆盖不全耦合削减后自然就聚到少数几个模式上。二是场景间距离计算没有加权重光伏部分数值大、风电部分数值小距离被光伏主导导致风电的波动特征被忽略。三是目标场景数设得太少比如从1000个直接砍到5个再好的算法也会丢掉大量信息。我的建议是目标场景数不要低于初始场景数的2%到3%。如果削减后的场景集合不能满足可靠性评估需求可以分两步走先整体削减出常规场景再单独做极端天气场景手动加入削减结果里保证典型场景和极端场景都有覆盖。5.2 削减前后统计特征误差怎么验证拿到削减结果后不能只看曲线好不好看要做量化验证。我常用的验证指标有三个。第一个是概率距离值削减后场景集合和原始场景集合的Kantorovich距离越小越好一般控制在原始场景集合均值的5%以内。第二个是均值与方差误差分别对比每个时段削减前后出力均值、方差的相对误差。第三个是用KS检验对比单时段出力的累积分布函数看削减后的分布是否和原始分布一致。这三个指标都能写进MATLAB代码里自动输出一份验证报告。5.3 风光场景的相关性怎么处理风速和光照本质上不是独立的很多地区存在“风小光强”或“互补”的特性。如果忽略这种相关性生成的联合场景会出现一些概率极低的组合比如大风天光照极强这种实际很少发生的情况导致调度结果偏保守或偏乐观。处理办法是在采样时对风速和光照的随机数施加相关性。常见方式是生成相关系数矩阵做Cholesky分解得到相关的标准正态随机变量再通过逆变换映射到威布尔分布和Beta分布。具体来说先采样一组相关的正态变量再用正态CDF反变换得到均匀分布的随机数最后用均匀随机数反变换采样风速和光照。这个过程在MATLAB里只需要几行代码就能实现但对场景质量的提升非常明显。风速和光照的相关系数可以根据历史数据估算没有数据时可以取-0.2到0.3之间的经验值。5.4 运行内存不足怎么办初始场景2000个每个24小时存成矩阵是 2000×24双精度占内存约1.5MB并不大。但如果把场景维度扩展到逐小时8760小时矩阵就是 2000×8760瞬间涨到140MB距离矩阵更是达到 2000×2000×8字节等于32MB加起来还能接受但再大就危险了。建议用单精度存储场景矩阵scenes single(scenes)内存直接减半。距离矩阵也可以分块计算避免一次性生成完整矩阵。如果场景数上万建议改用采样加聚类的方式先粗减到1000个再用概率距离快速削减法精减到目标数量两阶段削减的速度远快过一阶段硬算。5.5 调参经验速查参数推荐范围影响初始场景数500-2000太少覆盖不足太多计算慢目标场景数20-50太少损失波动性太多削减意义小风速自相关系数0.7-0.9控制时序平滑程度风光相关系数-0.2-0.3控制联合场景的合理性场景距离权重按额定功率归一化防止量纲不一致导致距离偏差最后一列是重点场景距离计算前最好先把风电和光伏出力都归一化到各自的额定功率再拼接场景。这样距离度量的物理意义变成了“归一化出力曲线的相似度”风电和光伏在距离计算中权重均衡削减结果更合理。结尾一点个人体会这套代码我在多个项目里反复用过最深的体会是场景削减不是“数据处理完事”它直接决定后续优化模型能不能算出靠谱的结果。概率距离快速削减法胜在原理清晰、可解释性强削减后的场景保留了原始概率语义这对调度方案的风险分析特别有用。后来我还试过用Wasserstein距离替代Kantorovich距离做削减效果接近数值上更稳定如果你有精力可以往这个方向扩展。再分享一个小技巧削减后的场景建议顺便输出每个场景的累积概率分布图和几个典型场景的时序可视化贴在报告里对方哪怕不懂算法一眼也能看出“这30条曲线确实覆盖了各种天气情况”这一点在项目汇报时价值极大。