风电不确定性下多目标优化调度：场景生成、NSGA-II与滚动优化

简介：这是一篇发表于《黑龙江电力》2014年第5期的学术论文PDF，面向电力系统运行分析、调度及风电并网研究领域的工程师和硕博研究生。论文针对风电机组并网规模扩大带来的火电机组频繁启停、运行效率低等问题，构建了考虑风电不确定性的电力系统多目标优化调度模型，引入maxmin函数和s支配的多目标粒子群算法进行求解，并通过算例验证了模型在经济性与平稳性上的表现。资源包为单个PDF文件，约367KB，内容包含调度目标函数、机组启停排序策略、梯形隶属度函数建模及算例分析，可帮助读者快速了解2014年前后风电不确定性处理与多目标调度建模的主流思路。目前已有153人学习浏览，适合用作电力系统优化调度课程论文、科研立项或技术报告撰写时的参考文献。

1. 风电不确定性下的多目标优化调度，问题不只是“多一个随机源”

风电出力是电力系统调度里最难伺候的电源。天气预报给出的日前风速曲线，到了实时运行阶段可能整体平移，也可能突然剧烈波动；一个容量百兆瓦的风电场，出力偏差动辄二三十兆瓦。如果调度模型把风电当成确定值来用，所有备用容量和线路潮流校验都建立在“偏差不会发生”的幻觉上。与此同时，降碳指标要求调度员必须同时压低发电成本和碳排放量，两个目标互相拉扯，问题就从“怎么把发电计划排出来”变成了“在风电不确定性场景集合下，怎么找一组成本、排放、可靠性都说得过去的帕累托解”。这篇博客会把从不确定性建模、多目标优化算法到滚动调度的完整链条拆开，给出可以直接动手改的代码和参数。

2. 风电不确定性建模：场景生成、削减与鲁棒区间

任何多目标优化调度模型，第一步不是写目标函数，而是决定不确定性用什么结构进入模型。常见做法有两种：场景法和区间鲁棒法。场景法生成一组可能的风电出力时序，给每条时序配一个发生概率；区间鲁棒法则只保留预测值的置信区间。场景法信息量多，但会放大求解规模；区间鲁棒计算快，但容易过度保守。工程上两者经常搭配使用：日前计划用场景法做成本与排放推演，安全校验用区间鲁棒检查极端偏差。

2.1 用拉丁超立方抽样生成初始风速场景

风速的概率分布通常用双参数 Weibull 拟合，形状参数 k 和尺度参数 c 可以通过场站历史风速的最大似然估计得到。风电功率与风速不是线性关系，典型风机在切入风速 3m/s 左右开始发电，在额定风速 12m/s 附近进入限功率状态，超过切出风速则停机。生成大量初始场景最直接的是蒙特卡洛采样，但那样会浪费大量样本在概率密度高的中间区域；拉丁超立方抽样（LHS）通过分层采样让场景均匀覆盖整个概率空间，适合先铺开、再削减。

import numpy as np from scipy.stats import weibull_min # 风速拟合参数：k=2.0, c=8.5，对应年均风速约 7.5m/s 的场站 n_scenarios = 500 n_horizon = 24 # 拉丁超立方抽样：把 [0,1] 分成 500 个等距区间，每个区间随机取一个分位数 quantiles = (np.arange(n_scenarios) + np.random.rand(n_scenarios, n_horizon)) / n_scenarios wind_speed = weibull_min.ppf(quantiles, 2.0, scale=8.5) # 分段简化功率曲线：切入 3m/s、额定 12m/s、切出 25m/s def speed_to_power(v, rated=2.0): power = np.zeros_like(v) linear_zone = (v >= 3) & (v < 12) rated_zone = (v >= 12) & (v < 25) power[linear_zone] = rated * (v[linear_zone] - 3) / (12 - 3) power[rated_zone] = rated return power wind_power = speed_to_power(wind_speed)

weibull_min.ppf是 Weibull 分布的分位数函数，位置参数 2.0 是形状参数 k，scale=8.5是尺度参数 c。LHS 把每个时段分成 500 个等概率区间，在每列内只取一个随机分位点，因此 500 条曲线在单时段上都服从同样的边缘分布，又避免了蒙特卡洛重复采样造成的空隙。speed_to_power里的rated是风电场额定总容量，单位 MW；如果场站额定容量是 50MW，就把 rated 改成 50。这段简化功率曲线没有考虑尾流效应和风机控制策略，用于调度原型验证足够；真实工程应换用场内功率预测系统给出的功率曲线。

生成 500 个场景后还有一个实际问题：相邻时段的采样是独立的，曲线上的波动可能比真实风速平滑曲线更毛糙。可以在采样后做一次一阶指数平滑，或者用带时间相关性的 ARIMA 生成风速序列，再做功率变换。更省事的做法是直接把历史预测误差叠加到预测曲线上，这样得到的场景天然带有时间相关性。

2.2 基于 K-means 的场景削减与概率重标定

直接从 2.1 得到的 500 条 24 维曲线大多相似，都塞进优化模型会让机组组合问题的规模爆炸。常见做法是用聚类算法把相似曲线归类，然后取每类的中心曲线作为代表性场景，并把该类的样本数占比作为场景概率。聚类目标可以是全部 24 小时出力，也可以按关键时段加权。

from sklearn.cluster import KMeans # 用 24 小时出力向量做 K-means 聚类，聚成 10 类 n_clusters = 10 kmeans = KMeans(n_clusters=n_clusters, random_state=42).fit(wind_power) centers = kmeans.cluster_centers_ # shape (10, 24) # 统计每个类的样本数，归一化为概率权重 counts = np.bincount(kmeans.labels_, minlength=n_clusters) weights = counts / counts.sum() # 检查削减前后平均曲线是否接近，避免把关键形态聚丢 print(np.dot(weights, centers).mean(axis=0)) print(wind_power.mean(axis=0))

centers是 10 条代表性风电出力曲线，weights是它们对应的概率，两者相加会近似等于原始场景总体的均值，可以用来快速验证削减过程没有系统性偏移。random_state=42是为了结果可复现；实际工程中应多试几个随机种子，看概率权重是否稳定。n_clusters 的选择会影响调度结果保守度：10~30 个场景对带旋转备用的机组组合比较常见，少于 5 个场景容易把晚高峰风电出力高的形态合并掉，多于 30 个场景会让每个时段的备用约束核算显著变慢。

还有一个容易被忽略的细节：K-means 聚类中心是均值曲线，不是原始场景的一员，所以它天然会把极端的“零出力”和“满发”场景平滑掉。如果调度结果用于安全校核，应该用每个类内最恶劣的原始场景替代聚类中心，或者保留每个类的 90% 分位曲线。这样得出的经济成本会更保守，但线路潮流和备用校验会可靠很多。

场景削减不是越精细越好。电力系统优化调度模型里，每多一个场景，功率平衡约束、备用约束和网络安全约束都会多一组；求解时间随场景数近似线性增长，而多目标优化算法的种群评估又要把所有场景全部算一遍，所以最终场景数要量力而行。

2.3 区间鲁棒与置信度选择

场景法适合回答“如果风电这样波动，成本和排放是多少”，区间鲁棒适合回答“极端情况下，备用和线路是否还能扛住”。区间鲁棒把每个时段的风电出力限制成预测均值加减一个带宽，再在约束最坏情况下求保护。为了避免过度保守，可以引入不确定预算：每个时段内最多只有少数风电场同时取到极端偏差。

表格里的“多一组对偶变量”是指在鲁棒优化处理中使用强对偶把最小化问题改写成可求解形式。实际工程中很少有人每次重新推导对偶；对于线性约束，可以直接用盒子鲁棒加不确定预算，把每个时段的风电上下偏差乘上一个 0/1 整数变量，再引入预算上限，这样模型仍然是混合整数线性规划。调度人员可以通过历史预测误差的分位数确定带宽：例如取过去 30 天同一风电场逐小时预测误差，统计到 90% 分位数，让覆盖率达到 90% 以上。

实际操作中，我一般会把场景法和区间鲁棒配合：用场景法算日前多目标曲线，用区间鲁棒对“系统净负荷”的极端上下界做一次安全校核，校核不过的时段再给备用约束附加条件。这个组合既能保留概率信息，又不会让优化模型被场景数量淹死。

3. 多目标优化调度模型的变量设计、目标与约束

不确定性场景准备好以后，下一步是把调度问题写成一个可以用算法求解的多目标优化模型。写模型之前先想清楚三件事：决策变量是什么，目标函数有几个，约束条件哪些会被场景数放大。下面按这个顺序展开。

3.1 决策变量和目标函数怎么设

在火电占比仍然较高的电力系统里，调度决策变量通常包括各台火电机组在每一个时段的有功出力、机组启停状态，以及风电场的弃风量。如果系统带储能，还要加充电、放电以及充放电状态变量；如果是微网，还要外加与主网的购售电功率和分时电价下的交互变量。决策变量太多会让多目标进化算法的搜索空间指数增长，所以能用连续变量表达的不要都用整数变量。机组启停在日前调度里必须用 0/1 变量，但弃风量可以设置成连续变量并限制在 0 到预测出力之间。

目标函数最常见的两元组合是“系统总发电成本”和“碳排放量”。发电成本包含火电机组的燃煤成本、启动成本和弃风惩罚；碳排放量则可以用机组出力与对应排放强度的线性关系估算。风电自身的运行成本很低，但把弃风惩罚直接设为固定值，优化器会倾向于让风电满发，即使那会挤占备用容量；更好的做法是把弃风惩罚设置成与旋转备用容量需求挂钩的分段函数。也可以把系统失负荷期望引入目标，但那样会变成三目标问题，求解和分析复杂度会明显上升。

目标函数的数学形态可以写成：

f1(Pg) = sum_t,h (a_h * P_h,t^2 + b_h * P_h,t + c_h * u_h,t + c_su_h * (1-u_h,t-1) * u_h,t) + w_curt * sum_t P_curt,t f2(Pg) = sum_t,h e_h * P_h,t + e_t * sum_t P_curt,t

公式里的 u 是启停状态，e_h 是机组碳排放强度，e_t 是弃风惩罚所对应的“虚拟排放系数”，当系统要把弃风造成的清洁能源浪费折算成碳成本时才会用到这一项。目标函数写完以后，一般会先做单目标求解，得到两个目标的理想点和最高点，用于后续帕累托前沿的归一化。

3.2 旋转备用与功率平衡约束的写法

风电不确定性对调度的最大冲击是功率平衡约束。如果只要求“预测出力下功率平衡”，那就完全没有考虑预测误差。对应地，对每一个风电场景 s，每个时段都需要满足：

sum_h P_h,t,s + P_w,t,s - P_curt,t = L_t

其中 P_w,t,s 是场景 s 下风电出力，L_t 是负荷。这个等式约束必须对每个场景、每个时段成立。常规机组出力和爬坡约束仍然按机组单位描述，备用约束要额外加一组不等式：

P_h,t + R_h,t^up <= P_h^max * u_h,t - eta_w * sigma_w,t

这里的 sigma_w,t 是风电预测误差标准差的函数，eta 是根据置信水平调整的系数，95% 置信区间对应 1.96 倍标准差。由于场景已经带了分布信息，也可以写成累计机会约束：要求场景集合中至少 95% 的场景下备用容量不发生越限，用整数变量记录越限场景，并把越限数量设置为不大于总场景数的 5%。

这种概率备用约束在混合整数规划里可以写成带大 M 形式的线性不等式；在多目标进化算法里则作为一个约束函数计算违反量。注意它是约束函数，不是目标函数。如果一个解在 5% 场景下备用越限，应该通过约束支配机制淘汰，而不是让它在帕累托前沿上形成一个“成本低廉但不可靠”的点。

3.3 用 Python 定义多目标调度问题的骨架

给一个用 pymoo 定义问题的骨架。这个例子缩小到单时段 6 个调度时段，只为了展示决策变量拆分和约束返回的写法；扩展到 24 时段只是维度变化。

import numpy as np from pymoo.core.problem import Problem class MultiPeriodDispatch(Problem): def __init__(self, n_gen=3, n_scenarios=10, n_periods=6): # 决策变量：每台机组每时段出力 + 每个时段弃风量 self.n_gen = n_gen self.n_periods = n_periods n_var = n_gen * n_periods + n_periods self.wind = np.linspace(15, 30, n_periods) # 六个时段的预测风电 self.load = np.linspace(40, 60, n_periods) super().__init__( n_var=n_var, n_obj=2, n_ieq_constr=n_periods, xl=np.concatenate([np.full(n_gen * n_periods, 10), np.zeros(n_periods)]), xu=np.concatenate([np.full(n_gen * n_periods, 80), self.wind]) ) def _evaluate(self, X, out, *args, **kwargs): n_gen, n_periods = self.n_gen, self.n_periods # 把一维染色体拆成 (场景数, 时段数, 机组数) 和弃风量 Pg = X[:, :n_gen * n_periods].reshape(X.shape[0], n_periods, n_gen) Pcurt = X[:, n_gen * n_periods:] # (n_pop, n_periods) # 目标1：发电成本 + 弃风惩罚 cost_coal = np.array([0.04, 0.06, 0.05]) # 每台机组成本系数 f1 = np.sum(Pg * cost_coal, axis=2).sum(axis=1) + 30 * Pcurt.sum(axis=1) # 目标2：碳排放量 emiss = np.array([0.8, 1.0, 0.9]) f2 = np.sum(Pg * emiss, axis=2).sum(axis=1) # 功率平衡：负荷 = 风电 + 机组出力 wind_supply = self.wind - Pcurt balance_resid = self.load - (Pg.sum(axis=2) + wind_supply) out["F"] = np.column_stack([f1, f2]) out["G"] = balance_resid # 等式约束，需转换成立不等式

代码里的 xl 和 xu 分别给出决策变量下界和上界。机组出力下界设成 10MW，上界 80MW；弃风变量的上界是当前时段的预测风电，保证弃风量不可能超过可发电量。cost_coal 和 emiss 数组的顺序必须与机组顺序一致。pymoo 的等式约束需要转换成不等式，返回前要做一次np.abs(balance_resid) - 1e-3 <= 0，否则直接用 balance_resid 会被当成“小于等于 0”的不等式约束，含义完全不同。这个错误很隐蔽，我在移植模型时踩过一次。

4. 用 NSGA-II 求解多目标调度问题与参数调优

模型建好之后，下一步就是让多目标优化算法在决策空间里搜索帕累托前沿。对于这种带非线性、整数变量和场景约束的调度模型，商用 MILP 求解器只能处理一个目标或带权重的多目标求解，不能直接告诉你在成本和碳排之间有哪些不可替代的折中方案。NSGA-II 是这类调度里最常见的启发式算法，它能在一轮计算中给出整条前沿，代价是求解时间不稳定。这里把算法调参说透。

4.1 NSGA-II 的核心流程与关键参数

NSGA-II 的核心是快速非支配排序和拥挤度距离。每一代中，当前解集通过非支配排序分成多个层级，同一层级内部用拥挤度距离来保持解的散布度。选择时用二元锦标赛：只要一个解的排序层级更低，就赢；层级相同则看拥挤度距离更大者。交叉和变异对解集做局部更新，然后父代子代合并共同排序，取最优的部分进入下一代。这个过程保证了帕累托前沿的收敛性和多样性同时被照顾。

调度问题关心的算法参数一般集中在种群大小、迭代次数、交叉变异操作，以及约束处理方法。下面是一组从多次调试里得到的推荐范围：

表中“24时段调度建议”特别要注意决策变量数。同样是 24 时段，3 台机组和 10 台机组的决策空间完全不同，变异概率最好按 1/n_var 设置，否则后期种群会被无规则变异打乱。种群大小和迭代次数可以先跑 10 次取帕累托前沿稳定所需次数，再放宽到该值的 1.5 倍，而不要一上来就设置 500。

4.2 约束支配与可行性修复

调度问题的难点不在目标计算，而在约束多。功率平衡、爬坡、备用、网络安全，这些约束在 NSGA-II 里都以不等式形式进入约束违反量。如果直接对所有约束用加权罚函数，算法前期会把大部分计算力浪费在寻找可行解上。更推荐的方法是约束支配：任意非支配解优先于被支配解；如果两个解互为非支配关系，那么可行解优先于不可行解；如果都不可行，违反量更小的解获胜。

除了排序规则，还可以在编码阶段做可行性修复。例如功率平衡约束，在求出所有机组的初始出力后，把偏差量加到当前时段最便宜或最贵的机组上，使等式满足。这样做会损失一点点目标精度，但能让遗传算子找到的新解通常是可行的。注意修复逻辑不要破坏爬坡约束，否则会在相邻时段之间制造更大的矛盾。

4.3 从帕累托前沿选折中解：模糊隶属度法

NSGA-II 跑完会得到一组解，而不是一个唯一方案。工程调度需要一个结果，所以要从帕累托前沿中选一个“折中解”。不同目标量纲不同，不能简单把目标值直接相加，也不能照搬欧氏距离，因为成本和碳排放的数字尺度差一个数量级。常见做法是先用模糊隶属度把目标归一化，再找综合满意度最好的点。

def pick_compromise(front): # front 是 NSGA-II 输出的目标矩阵，shape (n_solutions, 2) f_min = front.min(axis=0) f_max = front.max(axis=0) f_norm = (front - f_min) / (f_max - f_min + 1e-9) # 平均满意度越小，说明距离最优点越近 satisfaction = f_norm.mean(axis=1) return int(np.argmin(satisfaction)), front[np.argmin(satisfaction)]

f_max - f_min是帕累托前沿在单个目标上的跨度，加 1e-9 是为了防止消除零除异常。平均满意度的含义是把两个目标看成等权重；如果调度员更看重经济性，可以把成本对应的权重调大，例如改成0.7 * f_norm[:,0] + 0.3 * f_norm[:,1]。这种选点方式只对帕累托前沿有效；如果输入点中包含被支配解，归一化区间会被拖偏，所以必须确认这是 NSGA-II 最终输出的非支配解集。

5. 案例验证与进阶排错：把多目标调度推向实时

模型和算法都定了，接下来用一个小算例验证，然后解决工程落地时更常见的问题：为什么实时数据和调度计划对不上。

5.1 小型系统算例与约束违反量检查

假设一个 3 火电机组加单风电场的系统，机组容量都在 80MW 以下，风电场额定容量 50MW。用 30 个风电场景做日前调度，NSGA-II 得到 100 个解，再用模糊隶属度选折中解。对比确定性方案、成本最优方案和折中方案如下：

表格数字是示意，不是标准算例结果。检查约束违反量比盯着目标值更重要，如果备用越限小时数不是 0，需要先查看是哪些场景、哪些时段造成越限，再决定调整备用系数还是削减场景。

5.2 预测误差分布不匹配时的排查

实时运行中如果频繁出现备用调用，先检查事前假定的 Weibull 分布与实测误差是否一致。最简单的方式是用历史预测和实际出力做残差分析，看 90% 分位误差是否与模型设定一致。如果不匹配，不是算法问题，而是不确定性建模参数需要重新估计。更新风速分布参数后，重新生成场景，再跑一次多目标优化，通常就能让备用调用频率回到正常范围。

5.3 与模型预测控制结合的滚动调度

多目标优化调度不应该是“一天只算一次”的黑箱。更可靠的落地方式是把场景法封装成一个滚动求解单元，每个控制周期拿到最新风电预测后，仅执行下一个时段的决策。这一步可以用一个很短的循环表达：

for t in range(0, T - 24, 1): new_scenarios = generate_scenarios(forecast_at_t) front = run_nsga2(new_scenarios) plan = choose_compromise(front) apply_control(plan.iloc[:, t % 24])

这段循环把多目标求解嵌入到模型预测控制主循环里。滚动计算对求解速度有要求，如果 NSGA-II 单次需要一分钟以上，就只能把窗口缩小到 6 小时，或者用 GPU 批量评估种群。调度员刚开始做模型预测控制时，总是想把窗口设成全天，结果刷新频率太低；实际上把控制时域缩短到未来 4~8 小时，再配合每 15 分钟滚动更新，效果通常比一次性求全天好得多。

本文还有配套的精品资源，点击获取