新能源并网最让人头疼的是什么？不是设备本身，而是出力曲线那副"看心情"的劲儿。风一停、云一来，功率瞬间跳水，调度侧却必须拿出确定的方案应对不确定的发电。解决这个问题的标准技术路线，就是"场景生成与削减"。场景生成负责把随机性变成一组可计算的离散样本，场景削减则负责在保证精度的前提下把样本数量压缩到可操作的范围内。这套方法我在Matlab里完整跑过一遍，从风速分布建模、蒙特卡洛抽样，到同步回代削减和K-means聚类，每一步都踩过不少坑。这篇就把整条实现路径拆开讲清楚，适合刚接触新能源建模的研究生、做微电网或储能优化配置的工程师，以及任何想用Matlab把"不确定性"落到代码里的朋友参考。

1. 整体思路拆解：为什么要先生成、再削减

1.1 新能源出力的随机性本质

风电、光伏的出力本质上是一个随机过程。风速受气象条件、地形、湍流影响，光照强度则随云层移动、大气衰减变化，短时间内很难用确定性模型精确描述。但在做容量规划、经济调度、储能配置时，决策模型里的输入又必须是确定的数值。这就产生了矛盾：物理世界是随机的，优化模型却要求确定性输入。

场景法就是为了解决这个矛盾而存在的。它的核心逻辑是：把随机变量按其概率分布进行大量抽样，得到一组"可能发生但又各不相同"的出力曲线，每条曲线称为一个场景。场景集合并在一起，就能近似刻画原始随机分布的特征。简单说，就是用"抽样的离散集合"逼近"连续的随机分布"。

1.2 为什么必须做场景削减

蒙特卡洛抽样动辄生成几千上万条场景，直接把全部场景带入优化模型，计算量会大到无法接受。比如一个含风电的机组组合问题，每条场景对应一组约束，场景数上千时，求解时间会从分钟级恶化到小时级甚至无法求解。

但场景又不能随便删，删多了会丢失分布的关键信息，导致优化结果偏乐观或偏保守。场景削减的目标就是在"精度"与"规模"之间找一个平衡点：用尽量少的典型场景，最大程度保留原始场景集的概率分布特征。实际项目中，几千条场景削减到几十条甚至十几条，优化结果与真实情况的偏差可以控制在很小范围内。

我自己的经验是，场景数量每减少一个数量级，求解速度可能提升几十倍，而精度损失通常在几个百分点以内，完全在工程可接受范围内。所以这套"先生成、后削减"的流程，基本是所有含新能源不确定性优化的标准预处理步骤。

2. 场景生成的核心数学模型与Matlab实现

2.1 风速与光照的概率分布建模

场景生成的第一步，是确定随机变量服从什么分布。风电出力通常由风速驱动，工程上风速普遍用两参数Weibull分布描述：

f(v) = (k/c) * (v/c)^(k-1) * exp(-(v/c)^k)

其中k是形状参数，控制分布曲线的形态；c是尺度参数，控制风速整体量级。给定历史风速数据，可以用极大似然估计拟合成k和c。Matlab可以直接用wblrnd生成服从该分布的随机数，也可以用fitdist对历史数据做参数拟合。

光伏出力主要取决于光照强度，通常用Beta分布建模，其取值在[0,1]区间，正好对应归一化后的光照强度。Beta分布由α和β两个形状参数控制，Matlab中对应betarnd。同样可以用fitdist从历史辐照度数据中拟合参数。

这里要特别注意一个细节：风速和光照的时序相关性。直接在整体分布上抽样，生成的是独立同分布的序列，但真实的风速是连续变化的，前一时刻的风速会强烈影响后一时刻。忽略时序相关性生成的场景，一天内可能出现多次大起大落，明显不符合物理规律。

2.2 蒙特卡洛场景生成流程

蒙特卡洛生成场景的基本流程可以分为四步：

第一步，确定随机变量的概率分布模型以及参数。如果是做日前调度，通常会按小时划分时段，对每个时段单独拟合分布参数。

第二步，对每个时段的随机变量进行大量抽样。抽样次数决定了初始场景数。例如生成1000个场景，每个场景覆盖24小时，那么就需要抽24*1000个样本。Matlab中wblrnd(k, c, 1, 1000)一句就能完成一个时段的抽样。

第三步，根据风速-功率转换关系或光照-功率转换关系，把气象样本转换为出力数据。风机出力与风速的关系通常用分段函数描述：

P(v) = 0, v < v_in 或 v > v_out P(v) = P_rated * (v - v_in) / (v_rated - v_in), v_in ≤ v < v_rated P(v) = P_rated, v_rated ≤ v ≤ v_out

光伏出力则近似为光照强度乘以额定容量。这段转换逻辑不复杂，但容易出错的地方是分段点判断和单位归一化，建议单独封装成函数。

第四步，把逐时段的抽样结果组合成完整场景集。Matlab中一个常见做法是用三维数组存储：维度分别是"场景序号 × 时段数 × 电源类型"。

我补充一个实操建议：抽样完成后先做一次异常值检查。蒙特卡洛抽样虽然理论上分布正确，但有限样本中偶尔会出现极端值，比如Beta分布抽样得到接近0或接近1的异常光照强度。这些值虽然概率很小，但会大幅影响后续削减结果。我的做法是设定合理的物理边界，超出范围的样本直接剔除重抽，保证初始场景集的质量。

2.3 时序相关性的建模技巧

如果直接用独立抽样生成24小时场景，会发现相邻时段出力完全不连贯，看起来像白噪声。这在优化计算中虽然也能用，但结果可能偏保守，因为忽略了出力变化的惯性。处理时序相关性有几个常见办法。

最简单的是马尔可夫链法。把连续出力值离散化成若干个状态，统计历史数据的状态转移概率矩阵，然后按转移概率逐时段抽样。Matlab里可以用histcounts做状态划分，用累积概率矩阵配合rand实现状态转移抽样。这个方法的优点是实现简单、计算快，缺点是状态数太少会损失精度，状态数太多转移矩阵会变得稀疏。

更精细的做法是采用ARIMA时间序列模型。arima函数直接拟合历史出力序列，然后通过simulate函数生成大量未来路径。这个方法对单一电源的出力场景非常有效，尤其在数据质量好、时序规律明显的场合。但ARIMA也有局限性，它假设序列是平稳的，而风光出力受昼夜和季节影响，往往有明显的周期性和非平稳特征，建模前需要做差分处理。

我实际测试下来，如果只是想生成用于优化计算的日前场景，马尔可夫链配合适当的平滑处理已经够用。只有当研究重点在于长期时序特性，比如年尺度储能配置分析时，才值得上ARIMA或更复杂的生成模型。

3. 场景削减算法解析与Matlab代码实现

3.1 同步回代削减法：原理与实现

场景削减算法里最经典的当属同步回代削减法，也叫后向削减法。它的核心思想非常直观：每次迭代删除一个对场景集整体特征贡献最小的场景，并把被删除场景的概率累加到与它距离最近的保留场景上。

算法的具体步骤如下：

1. 计算所有场景两两之间的距离，通常用欧氏距离。假设每个场景是1行N列的时序向量，距离矩阵可以用pdist2一句算出。

2. 对每个场景，找到它与其他所有场景的最近距离以及对应的最近邻。

3. 找出所有"最近距离"中最小的那个，这个场景就是要删掉的场景。这一步的含义是：删除它造成的"信息损失"最小。

4. 把被删场景的概率加到其最近邻场景上，更新场景集。

5. 重复步骤2-4，直到场景数达到预设目标。

这个算法的好处是物理意义明确、实现简单，而且能保证保留下来的场景是原始场景中实际存在的样本，不会像K-means那样产生"虚拟场景"。但缺点也很明显：每删一个场景就要重算一次距离矩阵，时间复杂度高。初始场景1000个、削减到100个的迭代过程中，整体计算量存在压力。

Matlab中我建议用向量化操作来代替逐轮循环。比如用上三角矩阵去掉重复距离计算，用min函数一次找出最小距离。削减到目标场景数时的代码骨架大致是：

num_scenes = size(scenes, 1); probs = ones(num_scenes, 1) / num_scenes; target_num = 20; while num_scenes > target_num D = pdist2(scenes, scenes, 'euclidean'); D(1:num_scenes+1:end) = inf; % 对角元设无穷大 [min_dist, idx] = min(D, [], 2); [~, del_idx] = min(min_dist); neighbor = idx(del_idx); probs(neighbor) = probs(neighbor) + probs(del_idx); scenes(del_idx, :) = []; probs(del_idx) = []; num_scenes = num_scenes - 1; end

这段代码简洁但计算效率一般。如果要跑大场景集，建议把两两距离矩阵在外层预计算好，削减过程中只做局部更新，能省掉相当多重复计算。

3.2 K-means聚类削减及其改进方向

K-means聚类削减的思路与同步回代不同，它不直接删场景，而是把所有场景划分成K个簇，然后用每个簇的质心代表该簇的所有场景。质心场景的权重是该簇场景数量占总场景数量的比例。

K-means在Matlab里可以直接用kmeans函数，只需指定场景矩阵和聚类数K：

[idx, centroid] = kmeans(scenes, K, 'Distance', 'sqeuclidean', 'MaxIter', 500);

用完后统计每个簇的样本数，再按比例计算质心的概率权重。

与同步回代相比，K-means最大的优势是计算效率高，尤其场景数上万时，K-means仍然可以在秒级到分钟级完成。而且聚类后的"典型场景"往往具有明确的代表含义，便于分析场景的典型特征。

但K-means有一个被很多人忽略的问题：它生成的是质心场景，而不是真实场景。质心是簇内所有场景的算术平均，得到的曲线会趋向平滑，某些极端场景特征会被抹掉。对于优化问题，这种平滑化有时会导致结果偏保守，因为它低估了出力的波动幅度。

改进方向主要有两个。一是使用K-medoids算法，它选择的代表点是簇内离所有点最近的实际场景点，能保留真实曲线特征。Matlab中可以用kmedoids函数直接调用，代价是计算量略大。二是先用K-medoids聚类确定初始质心，再用K-means做精调，兼顾代表性和计算速度。我实测下来的经验是：如果后续优化模型对曲线波动敏感，优先用K-medoids；如果只关注总出力和期望值水平，K-means完全够用。

3.3 削减质量评估指标

削减做完，必须回答一个问题：削减后的场景集有没有"变质"。不看指标直接带入优化模型，是比较冒进的做法。我自己习惯用三个指标做评估。

第一个是场景总期望值误差。计算削减前后所有场景各时段的期望出力，对比最大偏差。这个指标直接反映优化结果是否会偏。如果偏差超过3%-5%，说明削减过度，需要增加场景数或换算法。

第二个是CDF对比。画出削减前后场景出力在每个时段的累积分布函数，观察两条曲线的贴合程度。这里用cdfplot就能快速完成。误差大的时段，说明削减算法在该时段丢失了分布细节。

第三个是相关性变化。计算削减前后场景序列的自相关系数或不同时段间的相关系数矩阵，对比差异。部分削减算法会对时序相关性造成破坏，尤其是K-means对每段时间独立聚类的时候。

我在项目中会用一张表记录不同场景数下的三个指标，对比后确定最终的场景数。没有单一准则适用所有情况，但有一个经验规律：场景数从1000削减到100的损失远小于从100削减到20的损失。也就是说越往后削减，边际信息损失越大，所以不要一味追求场景少。

4. 完整案例实操：风光联合场景生成与削减

4.1 数据准备与参数设置

接下来用一个完整案例走一遍流程。假设要为一个含风电场和光伏电站的微电网做日前调度准备场景输入，时间范围为24小时，时间步长1小时。

风速模型选用Weibull分布，假设基于历史数据拟合得到的形状参数k=2.3，尺度参数c=8.5。风机额定功率为1.5MW，切入风速3m/s，额定风速12m/s，切出风速25m/s。

光照模型选用Beta分布，按24个时段分别拟合。这里做一个简化假设：白天六个时段的光照Beta分布参数α=2.1、β=1.8，夜间时段不发电。光伏额定容量为1MW。

初始场景数设为1000，削减目标分别为50、30、20三种，通过指标对比确定最终选择。Matlab的随机数种子这里固定下来，方便复现结果：

rng(42);

4.2 生成与削减全流程代码实现

场景生成部分，先用wblrnd生成风速场景，再用betarnd生成光照场景，分别转换为出力曲线，按场景序号组织成矩阵。

% 参数定义 hours = 24; num_initial = 1000; k = 2.3; c = 8.5; v_in = 3; v_rated = 12; v_out = 25; P_rated_wind = 1.5; alpha = 2.1; beta = 1.8; P_rated_pv = 1.0; % 生成风速场景 wind_speed = wblrnd(k, c, num_initial, hours); % 风速转出力 wind_power = zeros(num_initial, hours); for i = 1:num_initial for t = 1:hours v = wind_speed(i, t); if v < v_in || v > v_out wind_power(i, t) = 0; elseif v >= v_rated wind_power(i, t) = P_rated_wind; else wind_power(i, t) = P_rated_wind * (v - v_in) / (v_rated - v_in); end end end % 生成光照场景 solar_power = zeros(num_initial, hours); for i = 1:num_initial for t = 6:17 % 假设白天时段 irrad = betarnd(alpha, beta); if irrad < 0.01 solar_power(i, t) = 0; else solar_power(i, t) = P_rated_pv * irrad; end end end % 组合出力场景并做归一化 total_power = wind_power + solar_power; total_power = total_power / max(total_power(:));

双层循环写起来直观，但效率不高。实际项目中我建议把风速转出力封装成向量化函数，用逻辑索引一次性处理大批量数据。上面是为了可读性保留了循环结构，读者可以自己改成向量版本。

削减阶段，分别跑同步回代和K-medoids聚类。以同步回为例：

target_nums = [50, 30, 20]; for tn = target_nums reduced_tn = backward_reduction(total_power, ones(num_initial,1)/num_initial, tn); % 存储结果用于后续指标评估 end

上面调用的backward_reduction函数就是第3.1节中那段循环代码的封装。

4.3 结果分析与误差评估

三种削减结果与初始场景集做对比后，我在实际运行中得到的典型数据是：

对于日前调度场景，期望值误差控制在2%以内通常可以接受，因此30个场景是性价比最高的选择。如果追求更精细的结果，就选50个场景，代价是求解时间明显上升。

可视化方面，我会画两张关键图。一张是削减前后的24小时期望出力对比曲线，看整体趋势是否一致。另一张是削减后的典型场景热力图，横轴是时段、纵轴是场景序号、颜色代表出力水平，能直观看出场景集是否覆盖了从低出力到高出力的各种典型情况。这两张图是检验场景质量最直观的方式，胜过看一堆数字指标。

5. 常见问题与调试经验

5.1 场景数如何选择

场景数没有标准答案，取决于下游模型的计算负担和精度阈值。有一个简单的实验方法：从100个场景开始，每次减少20个，观察期望出力误差的变化曲线。误差平缓的区域说明场景数还有压缩空间，误差开始明显上升的点就是临界场景数。

我个人在微电网优化中常用的场景数是20-50。机组组合这类大规模混合整数规划问题，场景数超过50求解时间往往很难接受；而储能容量配置这类线性规划问题，对场景数的容忍度更高，取100也没问题。总之要在求解器可承受范围内尽量保留更多场景。

5.2 计算性能优化技巧

场景生成和削减的计算瓶颈主要在两处。一是蒙特卡洛抽样本身，二是距离矩阵计算。

抽样环节可以用randraw工具包来生成各种自定义分布的随机数，比Matlab内置函数更灵活，但速度差别不大。真正影响性能的是转换逻辑里的循环，一定要向量化。比如风速转出力，用逻辑索引一次处理整个矩阵，可以把原来几秒的循环压缩到零点几秒。

距离矩阵计算在场景数超过5000时内存占用会变得很可观。一个5000×5000的双精度矩阵要占用200MB内存，迭代过程反复重算会非常慢。建议预先算一次距离矩阵，之后用索引更新，或者改用K-means这种不需要完整距离矩阵的算法。如果场景数量级上万，优先考虑K-means路线，不要硬跑同步回代。

5.3 我踩过的几个坑

Beta分布拟合在光照数据非零的情况下才有意义。如果历史数据里包含大量零值时段，整体拟合一个Beta分布会得到非常奇怪的参数。我踩过这个坑之后，改为对非零时段单独拟合，零值时段直接按0处理。

另一个坑是K-means聚类结果的随机性。K-means初始质心是随机选择的，同样的数据跑两次结果可能不同。解决办法是固定随机种子，或者使用多次运行取最优结果。kmeans函数可以通过'Replicates', 10参数让Matlab自动跑10次并返回最优聚类结果，强烈建议加上。

还有一点容易被忽略：场景生成时坐标和单位的归一化。风速用m/s、出力用MW、光照用百分比，直接混在一起算距离的话，量级差异会导致距离矩阵几乎由量级最大的变量主导。需要先把所有变量归一化到同一尺度，比如统一除以各自最大值或标准差。否则削减结果会偏向于照顾量级大的电源，另一个电源的场景特征可能被悄悄丢掉。

5.4 排查问题的基本流程

场景生成结果异常时，先不要急着调算法参数，按顺序排查这三件事。

第一，检查输入数据。历史风速和光照数据是否包含异常值、缺测值，时间序列是否有明显跳变。数据质量差是场景生成结果离谱的最常见原因。

第二，检查概率分布拟合结果。用plot把拟合曲线和实际数据直方图叠在一起看，直观判断拟合得好不好。仅靠拟合出的数值参数，很难发现分布形态不对。

第三，检查随机数生成过程。固定随机种子，重复生成两次初始场景，对比结果是否一致。如果不一致，说明代码某处存在未固定的随机源，比如K-means初始质心或抽样函数。固定rng后结果仍然不一致，就要检查是否有循环内重新设置随机种子的操作。

调试期间我习惯把中间变量逐一保存成mat文件，每一步结束后都画图检查。场景生成这类"随机性很强"的代码，一步错往往要到最后才发现结果异常，提前可视化能省下大量排查时间。

结尾分享

场景生成与削减这套流程，看起来是简单的统计抽样加聚类，真正跑通后才发现细节远比想象的多。分布参数拟合不好，后续场景质量全线崩塌；削减算法选错，精度和速度双双受损；连距离计算里量纲归一化这种小事都可能让削减结果偷偷跑偏。我在多次项目迭代中的一个体会是：场景削减不是一个"一步到位"的步骤，而是一个需要反复试算、不断与下游优化结果对照反馈的过程，留出足够的调参时间比追求一次跑通更重要。如果你也在用Matlab做新能源不确定性建模，建议先拿一套自己熟悉的历史数据，把生成、削减、评价这三步完整跑通，再逐步换成新数据和新场景。希望这篇分享能帮你少走几步弯路。